MULESTemplates.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "MULES.H"
27 #include "upwind.H"
28 #include "fvcSurfaceIntegrate.H"
29 #include "localEulerDdtScheme.H"
30 #include "slicedSurfaceFields.H"
31 #include "wedgeFvPatch.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 template<class RdeltaTType, class RhoType, class SpType, class SuType>
37 (
38  const RdeltaTType& rDeltaT,
39  const RhoType& rho,
40  volScalarField& psi,
41  const surfaceScalarField& phiPsi,
42  const SpType& Sp,
43  const SuType& Su
44 )
45 {
46  Info<< "MULES: Solving for " << psi.name() << endl;
47 
48  const fvMesh& mesh = psi.mesh();
49 
50  scalarField& psiIf = psi;
51  const scalarField& psi0 = psi.oldTime();
52 
53  psiIf = 0.0;
54  fvc::surfaceIntegrate(psiIf, phiPsi);
55 
56  if (mesh.moving())
57  {
58  psiIf =
59  (
60  mesh.Vsc0()().field()*rho.oldTime().field()
61  *psi0*rDeltaT/mesh.Vsc()().field()
62  + Su.field()
63  - psiIf
64  )/(rho.field()*rDeltaT - Sp.field());
65  }
66  else
67  {
68  psiIf =
69  (
70  rho.oldTime().field()*psi0*rDeltaT
71  + Su.field()
72  - psiIf
73  )/(rho.field()*rDeltaT - Sp.field());
74  }
75 
77 }
78 
79 
80 template<class RhoType>
82 (
83  const RhoType& rho,
84  volScalarField& psi,
85  const surfaceScalarField& phiPsi
86 )
87 {
88  explicitSolve(rho, psi, phiPsi, zeroField(), zeroField());
89 }
90 
91 
92 template<class RhoType, class SpType, class SuType>
94 (
95  const RhoType& rho,
96  volScalarField& psi,
97  const surfaceScalarField& phiPsi,
98  const SpType& Sp,
99  const SuType& Su
100 )
101 {
102  const fvMesh& mesh = psi.mesh();
103 
104  if (fv::localEulerDdt::enabled(mesh))
105  {
106  const volScalarField& rDeltaT = fv::localEulerDdt::localRDeltaT(mesh);
107  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
108  }
109  else
110  {
111  const scalar rDeltaT = 1.0/mesh.time().deltaTValue();
112  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
113  }
114 }
115 
116 
117 template<class RhoType, class PsiMaxType, class PsiMinType>
119 (
120  const RhoType& rho,
121  volScalarField& psi,
122  const surfaceScalarField& phiBD,
123  surfaceScalarField& phiPsi,
124  const PsiMaxType& psiMax,
125  const PsiMinType& psiMin
126 )
127 {
129  (
130  rho,
131  psi,
132  phiBD,
133  phiPsi,
134  zeroField(),
135  zeroField(),
136  psiMax,
137  psiMin
138  );
139 }
140 
141 
142 template
143 <
144  class RhoType,
145  class SpType,
146  class SuType,
147  class PsiMaxType,
148  class PsiMinType
149 >
151 (
152  const RhoType& rho,
153  volScalarField& psi,
154  const surfaceScalarField& phi,
155  surfaceScalarField& phiPsi,
156  const SpType& Sp,
157  const SuType& Su,
158  const PsiMaxType& psiMax,
159  const PsiMinType& psiMin
160 )
161 {
162  const fvMesh& mesh = psi.mesh();
163 
165 
166  if (fv::localEulerDdt::enabled(mesh))
167  {
168  const volScalarField& rDeltaT = fv::localEulerDdt::localRDeltaT(mesh);
169  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, false);
170  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
171  }
172  else
173  {
174  const scalar rDeltaT = 1.0/mesh.time().deltaTValue();
175  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, false);
176  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
177  }
178 }
179 
180 
181 template
182 <
183  class RdeltaTType,
184  class RhoType,
185  class SpType,
186  class SuType,
187  class PsiMaxType,
188  class PsiMinType
189 >
191 (
192  surfaceScalarField& lambda,
193  const RdeltaTType& rDeltaT,
194  const RhoType& rho,
195  const volScalarField& psi,
196  const surfaceScalarField& phiBD,
197  const surfaceScalarField& phiCorr,
198  const SpType& Sp,
199  const SuType& Su,
200  const PsiMaxType& psiMax,
201  const PsiMinType& psiMin
202 )
203 {
204  const scalarField& psiIf = psi;
205  const volScalarField::Boundary& psiBf = psi.boundaryField();
206 
207  const fvMesh& mesh = psi.mesh();
208 
209  const dictionary& MULEScontrols = mesh.solution().solverDict(psi.name());
210 
211  const label nLimiterIter
212  (
213  MULEScontrols.lookupOrDefault<label>("nLimiterIter", 3)
214  );
215 
216  const scalar smoothLimiter
217  (
218  MULEScontrols.lookupOrDefault<scalar>("smoothLimiter", 0)
219  );
220 
221  const scalar extremaCoeff
222  (
223  MULEScontrols.lookupOrDefault<scalar>("extremaCoeff", 0)
224  );
225 
226  const scalar boundaryExtremaCoeff
227  (
228  MULEScontrols.lookupOrDefault<scalar>
229  (
230  "boundaryExtremaCoeff",
231  extremaCoeff
232  )
233  );
234 
235  const scalar boundaryDeltaExtremaCoeff
236  (
237  max(boundaryExtremaCoeff - extremaCoeff, 0)
238  );
239 
240  const scalarField& psi0 = psi.oldTime();
241 
242  const labelUList& owner = mesh.owner();
243  const labelUList& neighb = mesh.neighbour();
244  tmp<volScalarField::Internal> tVsc = mesh.Vsc();
245  const scalarField& V = tVsc();
246 
247  const scalarField& phiBDIf = phiBD;
248  const surfaceScalarField::Boundary& phiBDBf =
249  phiBD.boundaryField();
250 
251  const scalarField& phiCorrIf = phiCorr;
252  const surfaceScalarField::Boundary& phiCorrBf =
253  phiCorr.boundaryField();
254 
255  scalarField& lambdaIf = lambda;
256  surfaceScalarField::Boundary& lambdaBf = lambda.boundaryFieldRef();
257 
258  scalarField psiMaxn(psiIf.size());
259  scalarField psiMinn(psiIf.size());
260 
261  psiMaxn = psiMin;
262  psiMinn = psiMax;
263 
264  scalarField sumPhiBD(psiIf.size(), 0.0);
265 
266  scalarField sumPhip(psiIf.size(), 0.0);
267  scalarField mSumPhim(psiIf.size(), 0.0);
268 
269  forAll(phiCorrIf, facei)
270  {
271  const label own = owner[facei];
272  const label nei = neighb[facei];
273 
274  psiMaxn[own] = max(psiMaxn[own], psiIf[nei]);
275  psiMinn[own] = min(psiMinn[own], psiIf[nei]);
276 
277  psiMaxn[nei] = max(psiMaxn[nei], psiIf[own]);
278  psiMinn[nei] = min(psiMinn[nei], psiIf[own]);
279 
280  sumPhiBD[own] += phiBDIf[facei];
281  sumPhiBD[nei] -= phiBDIf[facei];
282 
283  const scalar phiCorrf = phiCorrIf[facei];
284 
285  if (phiCorrf > 0)
286  {
287  sumPhip[own] += phiCorrf;
288  mSumPhim[nei] += phiCorrf;
289  }
290  else
291  {
292  mSumPhim[own] -= phiCorrf;
293  sumPhip[nei] -= phiCorrf;
294  }
295  }
296 
297  forAll(phiCorrBf, patchi)
298  {
299  const fvPatchScalarField& psiPf = psiBf[patchi];
300  const scalarField& phiBDPf = phiBDBf[patchi];
301  const scalarField& phiCorrPf = phiCorrBf[patchi];
302 
303  const labelList& pFaceCells = mesh.boundary()[patchi].faceCells();
304 
305  if (psiPf.coupled())
306  {
307  const scalarField psiPNf(psiPf.patchNeighbourField());
308 
309  forAll(phiCorrPf, pFacei)
310  {
311  const label pfCelli = pFaceCells[pFacei];
312 
313  psiMaxn[pfCelli] = max(psiMaxn[pfCelli], psiPNf[pFacei]);
314  psiMinn[pfCelli] = min(psiMinn[pfCelli], psiPNf[pFacei]);
315  }
316  }
317  else if (psiPf.fixesValue())
318  {
319  forAll(phiCorrPf, pFacei)
320  {
321  const label pfCelli = pFaceCells[pFacei];
322 
323  psiMaxn[pfCelli] = max(psiMaxn[pfCelli], psiPf[pFacei]);
324  psiMinn[pfCelli] = min(psiMinn[pfCelli], psiPf[pFacei]);
325  }
326  }
327  else
328  {
329  // Add the optional additional allowed boundary extrema
330  if (boundaryDeltaExtremaCoeff > 0)
331  {
332  forAll(phiCorrPf, pFacei)
333  {
334  const label pfCelli = pFaceCells[pFacei];
335 
336  const scalar extrema =
337  boundaryDeltaExtremaCoeff
338  *(psiMax[pfCelli] - psiMin[pfCelli]);
339 
340  psiMaxn[pfCelli] += extrema;
341  psiMinn[pfCelli] -= extrema;
342  }
343  }
344  }
345 
346  forAll(phiCorrPf, pFacei)
347  {
348  const label pfCelli = pFaceCells[pFacei];
349 
350  sumPhiBD[pfCelli] += phiBDPf[pFacei];
351 
352  const scalar phiCorrf = phiCorrPf[pFacei];
353 
354  if (phiCorrf > 0)
355  {
356  sumPhip[pfCelli] += phiCorrf;
357  }
358  else
359  {
360  mSumPhim[pfCelli] -= phiCorrf;
361  }
362  }
363  }
364 
365  psiMaxn = min(psiMaxn + extremaCoeff*(psiMax - psiMin), psiMax);
366  psiMinn = max(psiMinn - extremaCoeff*(psiMax - psiMin), psiMin);
367 
368  if (smoothLimiter > small)
369  {
370  psiMaxn =
371  min(smoothLimiter*psiIf + (1.0 - smoothLimiter)*psiMaxn, psiMax);
372  psiMinn =
373  max(smoothLimiter*psiIf + (1.0 - smoothLimiter)*psiMinn, psiMin);
374  }
375 
376  if (mesh.moving())
377  {
379 
380  psiMaxn =
381  V
382  *(
383  (rho.field()*rDeltaT - Sp.field())*psiMaxn
384  - Su.field()
385  )
386  - (V0().field()*rDeltaT)*rho.oldTime().field()*psi0
387  + sumPhiBD;
388 
389  psiMinn =
390  V
391  *(
392  Su.field()
393  - (rho.field()*rDeltaT - Sp.field())*psiMinn
394  )
395  + (V0().field()*rDeltaT)*rho.oldTime().field()*psi0
396  - sumPhiBD;
397  }
398  else
399  {
400  psiMaxn =
401  V
402  *(
403  (rho.field()*rDeltaT - Sp.field())*psiMaxn
404  - Su.field()
405  - (rho.oldTime().field()*rDeltaT)*psi0
406  )
407  + sumPhiBD;
408 
409  psiMinn =
410  V
411  *(
412  Su.field()
413  - (rho.field()*rDeltaT - Sp.field())*psiMinn
414  + (rho.oldTime().field()*rDeltaT)*psi0
415  )
416  - sumPhiBD;
417  }
418 
419  scalarField sumlPhip(psiIf.size());
420  scalarField mSumlPhim(psiIf.size());
421 
422  for (int j=0; j<nLimiterIter; j++)
423  {
424  sumlPhip = 0;
425  mSumlPhim = 0;
426 
427  forAll(lambdaIf, facei)
428  {
429  const label own = owner[facei];
430  const label nei = neighb[facei];
431 
432  scalar lambdaPhiCorrf = lambdaIf[facei]*phiCorrIf[facei];
433 
434  if (lambdaPhiCorrf > 0)
435  {
436  sumlPhip[own] += lambdaPhiCorrf;
437  mSumlPhim[nei] += lambdaPhiCorrf;
438  }
439  else
440  {
441  mSumlPhim[own] -= lambdaPhiCorrf;
442  sumlPhip[nei] -= lambdaPhiCorrf;
443  }
444  }
445 
446  forAll(lambdaBf, patchi)
447  {
448  scalarField& lambdaPf = lambdaBf[patchi];
449  const scalarField& phiCorrfPf = phiCorrBf[patchi];
450 
451  const labelList& pFaceCells = mesh.boundary()[patchi].faceCells();
452 
453  forAll(lambdaPf, pFacei)
454  {
455  const label pfCelli = pFaceCells[pFacei];
456  const scalar lambdaPhiCorrf =
457  lambdaPf[pFacei]*phiCorrfPf[pFacei];
458 
459  if (lambdaPhiCorrf > 0)
460  {
461  sumlPhip[pfCelli] += lambdaPhiCorrf;
462  }
463  else
464  {
465  mSumlPhim[pfCelli] -= lambdaPhiCorrf;
466  }
467  }
468  }
469 
470  forAll(sumlPhip, celli)
471  {
472  sumlPhip[celli] =
473  max(min
474  (
475  (sumlPhip[celli] + psiMaxn[celli])
476  /(mSumPhim[celli] + rootVSmall),
477  1.0), 0.0
478  );
479 
480  mSumlPhim[celli] =
481  max(min
482  (
483  (mSumlPhim[celli] + psiMinn[celli])
484  /(sumPhip[celli] + rootVSmall),
485  1.0), 0.0
486  );
487  }
488 
489  const scalarField& lambdam = sumlPhip;
490  const scalarField& lambdap = mSumlPhim;
491 
492  forAll(lambdaIf, facei)
493  {
494  if (phiCorrIf[facei] > 0)
495  {
496  lambdaIf[facei] = min
497  (
498  lambdaIf[facei],
499  min(lambdap[owner[facei]], lambdam[neighb[facei]])
500  );
501  }
502  else
503  {
504  lambdaIf[facei] = min
505  (
506  lambdaIf[facei],
507  min(lambdam[owner[facei]], lambdap[neighb[facei]])
508  );
509  }
510  }
511 
512  forAll(lambdaBf, patchi)
513  {
514  fvsPatchScalarField& lambdaPf = lambdaBf[patchi];
515  const scalarField& phiCorrfPf = phiCorrBf[patchi];
516  const fvPatchScalarField& psiPf = psiBf[patchi];
517 
518  if (isA<wedgeFvPatch>(mesh.boundary()[patchi]))
519  {
520  lambdaPf = 0;
521  }
522  else if (psiPf.coupled())
523  {
524  const labelList& pFaceCells =
525  mesh.boundary()[patchi].faceCells();
526 
527  forAll(lambdaPf, pFacei)
528  {
529  const label pfCelli = pFaceCells[pFacei];
530 
531  if (phiCorrfPf[pFacei] > 0)
532  {
533  lambdaPf[pFacei] =
534  min(lambdaPf[pFacei], lambdap[pfCelli]);
535  }
536  else
537  {
538  lambdaPf[pFacei] =
539  min(lambdaPf[pFacei], lambdam[pfCelli]);
540  }
541  }
542  }
543  }
544 
545  // Take minimum value of limiter across coupled patches
546  surfaceScalarField::Boundary lambdaNbrBf
547  (
548  surfaceScalarField::Internal::null(),
549  lambdaBf.boundaryNeighbourField()
550  );
551  forAll(lambdaBf, patchi)
552  {
553  fvsPatchScalarField& lambdaPf = lambdaBf[patchi];
554  const fvsPatchScalarField& lambdaNbrPf = lambdaNbrBf[patchi];
555  if (lambdaPf.coupled())
556  {
557  lambdaPf = min(lambdaPf, lambdaNbrPf);
558  }
559  }
560  }
561 }
562 
563 
564 template
565 <
566  class RdeltaTType,
567  class RhoType,
568  class SpType,
569  class SuType,
570  class PsiMaxType,
571  class PsiMinType
572 >
574 (
575  const RdeltaTType& rDeltaT,
576  const RhoType& rho,
577  const volScalarField& psi,
578  const surfaceScalarField& phi,
579  surfaceScalarField& phiPsi,
580  const SpType& Sp,
581  const SuType& Su,
582  const PsiMaxType& psiMax,
583  const PsiMinType& psiMin,
584  const bool returnCorr
585 )
586 {
587  const fvMesh& mesh = psi.mesh();
588 
589  surfaceScalarField phiBD(upwind<scalar>(psi.mesh(), phi).flux(psi));
590 
592  const surfaceScalarField::Boundary& phiPsiBf = phiPsi.boundaryField();
593 
594  forAll(phiBDBf, patchi)
595  {
596  fvsPatchScalarField& phiBDPf = phiBDBf[patchi];
597 
598  if (!phiBDPf.coupled())
599  {
600  phiBDPf = phiPsiBf[patchi];
601  }
602  }
603 
604  surfaceScalarField& phiCorr = phiPsi;
605  phiCorr -= phiBD;
606 
608  (
609  IOobject
610  (
611  "lambda",
612  mesh.time().timeName(),
613  mesh,
614  IOobject::NO_READ,
615  IOobject::NO_WRITE,
616  false
617  ),
618  mesh,
620  );
621 
622  limiter
623  (
624  lambda,
625  rDeltaT,
626  rho,
627  psi,
628  phiBD,
629  phiCorr,
630  Sp,
631  Su,
632  psiMax,
633  psiMin
634  );
635 
636  if (returnCorr)
637  {
638  phiCorr *= lambda;
639  }
640  else
641  {
642  phiPsi = phiBD + lambda*phiCorr;
643  }
644 }
645 
646 
647 template
648 <
649  class RhoType,
650  class SpType,
651  class SuType,
652  class PsiMaxType,
653  class PsiMinType
654 >
656 (
657  const RhoType& rho,
658  const volScalarField& psi,
659  const surfaceScalarField& phi,
660  surfaceScalarField& phiPsi,
661  const SpType& Sp,
662  const SuType& Su,
663  const PsiMaxType& psiMax,
664  const PsiMinType& psiMin,
665  const bool rtnCorr
666 )
667 {
668  const fvMesh& mesh = psi.mesh();
669 
670  if (fv::localEulerDdt::enabled(mesh))
671  {
672  const volScalarField& rDeltaT = fv::localEulerDdt::localRDeltaT(mesh);
673  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, rtnCorr);
674  }
675  else
676  {
677  const scalar rDeltaT = 1.0/mesh.time().deltaTValue();
678  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, rtnCorr);
679  }
680 }
681 
682 
683 template<template<class> class AlphaList, template<class> class PhiList>
685 (
686  const AlphaList<const volScalarField>& alphas,
687  PhiList<surfaceScalarField>& phiPsis,
688  const surfaceScalarField& phi
689 )
690 {
691  PtrList<surfaceScalarField> alphaPhiUDs(phiPsis.size());
692 
693  forAll(phiPsis, phasei)
694  {
695  alphaPhiUDs.set
696  (
697  phasei,
698  upwind<scalar>(phi.mesh(), phi).flux(alphas[phasei])
699  );
700 
701  phiPsis[phasei] -= alphaPhiUDs[phasei];
702  }
703 
704  {
705  UPtrList<scalarField> phiPsisInternal(phiPsis.size());
706 
707  forAll(phiPsisInternal, phasei)
708  {
709  phiPsisInternal.set(phasei, &phiPsis[phasei]);
710  }
711 
712  limitSum(phiPsisInternal);
713  }
714 
715  const surfaceScalarField::Boundary& phibf = phi.boundaryField();
716 
717  forAll(phibf, patchi)
718  {
719  if (phibf[patchi].coupled())
720  {
721  UPtrList<scalarField> phiPsisPatch(phiPsis.size());
722 
723  forAll(phiPsisPatch, phasei)
724  {
725  phiPsisPatch.set
726  (
727  phasei,
728  &phiPsis[phasei].boundaryFieldRef()[patchi]
729  );
730  }
731 
732  limitSum(phiPsisPatch);
733  }
734  }
735 
736  forAll(phiPsis, phasei)
737  {
738  phiPsis[phasei] += alphaPhiUDs[phasei];
739  }
740 }
741 
742 
743 // ************************************************************************* //
void limit(const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phi, surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin, const bool returnCorr)
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField...
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const word & name() const
Return name.
Definition: IOobject.H:315
void surfaceIntegrate(Field< Type > &ivf, const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
bool moving() const
Is mesh moving.
Definition: polyMesh.H:525
label phasei
Definition: pEqn.H:27
virtual bool fixesValue() const
Return true if this patch field fixes a value.
Definition: fvPatchField.H:310
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
dimensionedScalar lambda(viscosity->lookup("lambda"))
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const fvSolution & solution() const
Return the fvSchemes.
Definition: fvMesh.C:1683
const dimensionSet dimless
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
const dictionary & solverDict(const word &name) const
Return the solver controls dictionary for the given field.
Definition: solution.C:335
fvMesh & mesh
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:423
tmp< DimensionedField< scalar, volMesh > > Vsc() const
Return sub-cycle cell volumes.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
Definition: MULES.C:30
virtual bool coupled() const
Return true if this patch field is coupled.
A class representing the concept of a field of 0 used to avoid unnecessary manipulations for objects ...
Definition: zeroField.H:50
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:34
const volScalarField & psi
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:323
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
phi
Definition: correctPhi.H:3
void limiter(surfaceScalarField &lambda, const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phiBD, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin)
tmp< GeometricBoundaryField > boundaryNeighbourField() const
Return BoundaryField of the values on the other side of couples.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:54
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
tmp< DimensionedField< scalar, volMesh > > Vsc0() const
Return sub-cycle old-time cell volumes.
const Mesh & mesh() const
Return mesh.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:417
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
MULES: Multidimensional universal limiter for explicit solution.
rDeltaTY field()
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
virtual tmp< Field< Type > > patchNeighbourField(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking) const
Return patchField on the opposite patch of a coupled patch.
Definition: fvPatchField.H:433
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:800