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-2020 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 #include "syncTools.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 template<class RdeltaTType, class RhoType, class SpType, class SuType>
38 (
39  const RdeltaTType& rDeltaT,
40  const RhoType& rho,
41  volScalarField& psi,
42  const surfaceScalarField& phiPsi,
43  const SpType& Sp,
44  const SuType& Su
45 )
46 {
47  Info<< "MULES: Solving for " << psi.name() << endl;
48 
49  const fvMesh& mesh = psi.mesh();
50 
51  scalarField& psiIf = psi;
52  const scalarField& psi0 = psi.oldTime();
53 
54  psiIf = 0.0;
55  fvc::surfaceIntegrate(psiIf, phiPsi);
56 
57  if (mesh.moving())
58  {
59  psiIf =
60  (
61  mesh.Vsc0()().field()*rho.oldTime().field()
62  *psi0*rDeltaT/mesh.Vsc()().field()
63  + Su.field()
64  - psiIf
65  )/(rho.field()*rDeltaT - Sp.field());
66  }
67  else
68  {
69  psiIf =
70  (
71  rho.oldTime().field()*psi0*rDeltaT
72  + Su.field()
73  - psiIf
74  )/(rho.field()*rDeltaT - Sp.field());
75  }
76 
78 }
79 
80 
81 template<class RhoType>
83 (
84  const RhoType& rho,
85  volScalarField& psi,
86  const surfaceScalarField& phiPsi
87 )
88 {
89  explicitSolve(rho, psi, phiPsi, zeroField(), zeroField());
90 }
91 
92 
93 template<class RhoType, class SpType, class SuType>
95 (
96  const RhoType& rho,
97  volScalarField& psi,
98  const surfaceScalarField& phiPsi,
99  const SpType& Sp,
100  const SuType& Su
101 )
102 {
103  const fvMesh& mesh = psi.mesh();
104 
105  if (fv::localEulerDdt::enabled(mesh))
106  {
107  const volScalarField& rDeltaT = fv::localEulerDdt::localRDeltaT(mesh);
108  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
109  }
110  else
111  {
112  const scalar rDeltaT = 1.0/mesh.time().deltaTValue();
113  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
114  }
115 }
116 
117 
118 template<class RhoType, class PsiMaxType, class PsiMinType>
120 (
121  const RhoType& rho,
122  volScalarField& psi,
123  const surfaceScalarField& phiBD,
124  surfaceScalarField& phiPsi,
125  const PsiMaxType& psiMax,
126  const PsiMinType& psiMin
127 )
128 {
130  (
131  rho,
132  psi,
133  phiBD,
134  phiPsi,
135  zeroField(),
136  zeroField(),
137  psiMax,
138  psiMin
139  );
140 }
141 
142 
143 template
144 <
145  class RhoType,
146  class SpType,
147  class SuType,
148  class PsiMaxType,
149  class PsiMinType
150 >
152 (
153  const RhoType& rho,
154  volScalarField& psi,
155  const surfaceScalarField& phi,
156  surfaceScalarField& phiPsi,
157  const SpType& Sp,
158  const SuType& Su,
159  const PsiMaxType& psiMax,
160  const PsiMinType& psiMin
161 )
162 {
163  const fvMesh& mesh = psi.mesh();
164 
166 
167  if (fv::localEulerDdt::enabled(mesh))
168  {
169  const volScalarField& rDeltaT = fv::localEulerDdt::localRDeltaT(mesh);
170  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, false);
171  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
172  }
173  else
174  {
175  const scalar rDeltaT = 1.0/mesh.time().deltaTValue();
176  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, false);
177  explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su);
178  }
179 }
180 
181 
182 template
183 <
184  class RdeltaTType,
185  class RhoType,
186  class SpType,
187  class SuType,
188  class PsiMaxType,
189  class PsiMinType
190 >
192 (
193  scalarField& allLambda,
194  const RdeltaTType& rDeltaT,
195  const RhoType& rho,
196  const volScalarField& psi,
197  const surfaceScalarField& phiBD,
198  const surfaceScalarField& phiCorr,
199  const SpType& Sp,
200  const SuType& Su,
201  const PsiMaxType& psiMax,
202  const PsiMinType& psiMin
203 )
204 {
205  const scalarField& psiIf = psi;
206  const volScalarField::Boundary& psiBf = psi.boundaryField();
207 
208  const fvMesh& mesh = psi.mesh();
209 
210  const dictionary& MULEScontrols = mesh.solverDict(psi.name());
211 
212  const label nLimiterIter
213  (
214  MULEScontrols.lookupOrDefault<label>("nLimiterIter", 3)
215  );
216 
217  const scalar smoothLimiter
218  (
219  MULEScontrols.lookupOrDefault<scalar>("smoothLimiter", 0)
220  );
221 
222  const scalar extremaCoeff
223  (
224  MULEScontrols.lookupOrDefault<scalar>("extremaCoeff", 0)
225  );
226 
227  const scalar boundaryExtremaCoeff
228  (
229  MULEScontrols.lookupOrDefault<scalar>
230  (
231  "boundaryExtremaCoeff",
232  extremaCoeff
233  )
234  );
235 
236  const scalar boundaryDeltaExtremaCoeff
237  (
238  max(boundaryExtremaCoeff - extremaCoeff, 0)
239  );
240 
241  const scalarField& psi0 = psi.oldTime();
242 
243  const labelUList& owner = mesh.owner();
244  const labelUList& neighb = mesh.neighbour();
245  tmp<volScalarField::Internal> tVsc = mesh.Vsc();
246  const scalarField& V = tVsc();
247 
248  const scalarField& phiBDIf = phiBD;
249  const surfaceScalarField::Boundary& phiBDBf =
250  phiBD.boundaryField();
251 
252  const scalarField& phiCorrIf = phiCorr;
253  const surfaceScalarField::Boundary& phiCorrBf =
254  phiCorr.boundaryField();
255 
257  (
258  IOobject
259  (
260  "lambda",
261  mesh.time().timeName(),
262  mesh,
263  IOobject::NO_READ,
264  IOobject::NO_WRITE,
265  false
266  ),
267  mesh,
268  dimless,
269  allLambda,
270  false // Use slices for the couples
271  );
272 
273  scalarField& lambdaIf = lambda;
274  surfaceScalarField::Boundary& lambdaBf = lambda.boundaryFieldRef();
275 
276  scalarField psiMaxn(psiIf.size());
277  scalarField psiMinn(psiIf.size());
278 
279  psiMaxn = psiMin;
280  psiMinn = psiMax;
281 
282  scalarField sumPhiBD(psiIf.size(), 0.0);
283 
284  scalarField sumPhip(psiIf.size(), 0.0);
285  scalarField mSumPhim(psiIf.size(), 0.0);
286 
287  forAll(phiCorrIf, facei)
288  {
289  const label own = owner[facei];
290  const label nei = neighb[facei];
291 
292  psiMaxn[own] = max(psiMaxn[own], psiIf[nei]);
293  psiMinn[own] = min(psiMinn[own], psiIf[nei]);
294 
295  psiMaxn[nei] = max(psiMaxn[nei], psiIf[own]);
296  psiMinn[nei] = min(psiMinn[nei], psiIf[own]);
297 
298  sumPhiBD[own] += phiBDIf[facei];
299  sumPhiBD[nei] -= phiBDIf[facei];
300 
301  const scalar phiCorrf = phiCorrIf[facei];
302 
303  if (phiCorrf > 0)
304  {
305  sumPhip[own] += phiCorrf;
306  mSumPhim[nei] += phiCorrf;
307  }
308  else
309  {
310  mSumPhim[own] -= phiCorrf;
311  sumPhip[nei] -= phiCorrf;
312  }
313  }
314 
315  forAll(phiCorrBf, patchi)
316  {
317  const fvPatchScalarField& psiPf = psiBf[patchi];
318  const scalarField& phiBDPf = phiBDBf[patchi];
319  const scalarField& phiCorrPf = phiCorrBf[patchi];
320 
321  const labelList& pFaceCells = mesh.boundary()[patchi].faceCells();
322 
323  if (psiPf.coupled())
324  {
325  const scalarField psiPNf(psiPf.patchNeighbourField());
326 
327  forAll(phiCorrPf, pFacei)
328  {
329  const label pfCelli = pFaceCells[pFacei];
330 
331  psiMaxn[pfCelli] = max(psiMaxn[pfCelli], psiPNf[pFacei]);
332  psiMinn[pfCelli] = min(psiMinn[pfCelli], psiPNf[pFacei]);
333  }
334  }
335  else if (psiPf.fixesValue())
336  {
337  forAll(phiCorrPf, pFacei)
338  {
339  const label pfCelli = pFaceCells[pFacei];
340 
341  psiMaxn[pfCelli] = max(psiMaxn[pfCelli], psiPf[pFacei]);
342  psiMinn[pfCelli] = min(psiMinn[pfCelli], psiPf[pFacei]);
343  }
344  }
345  else
346  {
347  // Add the optional additional allowed boundary extrema
348  if (boundaryDeltaExtremaCoeff > 0)
349  {
350  forAll(phiCorrPf, pFacei)
351  {
352  const label pfCelli = pFaceCells[pFacei];
353 
354  const scalar extrema =
355  boundaryDeltaExtremaCoeff
356  *(psiMax[pfCelli] - psiMin[pfCelli]);
357 
358  psiMaxn[pfCelli] += extrema;
359  psiMinn[pfCelli] -= extrema;
360  }
361  }
362  }
363 
364  forAll(phiCorrPf, pFacei)
365  {
366  const label pfCelli = pFaceCells[pFacei];
367 
368  sumPhiBD[pfCelli] += phiBDPf[pFacei];
369 
370  const scalar phiCorrf = phiCorrPf[pFacei];
371 
372  if (phiCorrf > 0)
373  {
374  sumPhip[pfCelli] += phiCorrf;
375  }
376  else
377  {
378  mSumPhim[pfCelli] -= phiCorrf;
379  }
380  }
381  }
382 
383  psiMaxn = min(psiMaxn + extremaCoeff*(psiMax - psiMin), psiMax);
384  psiMinn = max(psiMinn - extremaCoeff*(psiMax - psiMin), psiMin);
385 
386  if (smoothLimiter > small)
387  {
388  psiMaxn =
389  min(smoothLimiter*psiIf + (1.0 - smoothLimiter)*psiMaxn, psiMax);
390  psiMinn =
391  max(smoothLimiter*psiIf + (1.0 - smoothLimiter)*psiMinn, psiMin);
392  }
393 
394  if (mesh.moving())
395  {
397 
398  psiMaxn =
399  V
400  *(
401  (rho.field()*rDeltaT - Sp.field())*psiMaxn
402  - Su.field()
403  )
404  - (V0().field()*rDeltaT)*rho.oldTime().field()*psi0
405  + sumPhiBD;
406 
407  psiMinn =
408  V
409  *(
410  Su.field()
411  - (rho.field()*rDeltaT - Sp.field())*psiMinn
412  )
413  + (V0().field()*rDeltaT)*rho.oldTime().field()*psi0
414  - sumPhiBD;
415  }
416  else
417  {
418  psiMaxn =
419  V
420  *(
421  (rho.field()*rDeltaT - Sp.field())*psiMaxn
422  - Su.field()
423  - (rho.oldTime().field()*rDeltaT)*psi0
424  )
425  + sumPhiBD;
426 
427  psiMinn =
428  V
429  *(
430  Su.field()
431  - (rho.field()*rDeltaT - Sp.field())*psiMinn
432  + (rho.oldTime().field()*rDeltaT)*psi0
433  )
434  - sumPhiBD;
435  }
436 
437  scalarField sumlPhip(psiIf.size());
438  scalarField mSumlPhim(psiIf.size());
439 
440  for (int j=0; j<nLimiterIter; j++)
441  {
442  sumlPhip = 0;
443  mSumlPhim = 0;
444 
445  forAll(lambdaIf, facei)
446  {
447  const label own = owner[facei];
448  const label nei = neighb[facei];
449 
450  scalar lambdaPhiCorrf = lambdaIf[facei]*phiCorrIf[facei];
451 
452  if (lambdaPhiCorrf > 0)
453  {
454  sumlPhip[own] += lambdaPhiCorrf;
455  mSumlPhim[nei] += lambdaPhiCorrf;
456  }
457  else
458  {
459  mSumlPhim[own] -= lambdaPhiCorrf;
460  sumlPhip[nei] -= lambdaPhiCorrf;
461  }
462  }
463 
464  forAll(lambdaBf, patchi)
465  {
466  scalarField& lambdaPf = lambdaBf[patchi];
467  const scalarField& phiCorrfPf = phiCorrBf[patchi];
468 
469  const labelList& pFaceCells = mesh.boundary()[patchi].faceCells();
470 
471  forAll(lambdaPf, pFacei)
472  {
473  const label pfCelli = pFaceCells[pFacei];
474  const scalar lambdaPhiCorrf =
475  lambdaPf[pFacei]*phiCorrfPf[pFacei];
476 
477  if (lambdaPhiCorrf > 0)
478  {
479  sumlPhip[pfCelli] += lambdaPhiCorrf;
480  }
481  else
482  {
483  mSumlPhim[pfCelli] -= lambdaPhiCorrf;
484  }
485  }
486  }
487 
488  forAll(sumlPhip, celli)
489  {
490  sumlPhip[celli] =
491  max(min
492  (
493  (sumlPhip[celli] + psiMaxn[celli])
494  /(mSumPhim[celli] + rootVSmall),
495  1.0), 0.0
496  );
497 
498  mSumlPhim[celli] =
499  max(min
500  (
501  (mSumlPhim[celli] + psiMinn[celli])
502  /(sumPhip[celli] + rootVSmall),
503  1.0), 0.0
504  );
505  }
506 
507  const scalarField& lambdam = sumlPhip;
508  const scalarField& lambdap = mSumlPhim;
509 
510  forAll(lambdaIf, facei)
511  {
512  if (phiCorrIf[facei] > 0)
513  {
514  lambdaIf[facei] = min
515  (
516  lambdaIf[facei],
517  min(lambdap[owner[facei]], lambdam[neighb[facei]])
518  );
519  }
520  else
521  {
522  lambdaIf[facei] = min
523  (
524  lambdaIf[facei],
525  min(lambdam[owner[facei]], lambdap[neighb[facei]])
526  );
527  }
528  }
529 
530  forAll(lambdaBf, patchi)
531  {
532  fvsPatchScalarField& lambdaPf = lambdaBf[patchi];
533  const scalarField& phiCorrfPf = phiCorrBf[patchi];
534  const fvPatchScalarField& psiPf = psiBf[patchi];
535 
536  if (isA<wedgeFvPatch>(mesh.boundary()[patchi]))
537  {
538  lambdaPf = 0;
539  }
540  else if (psiPf.coupled())
541  {
542  const labelList& pFaceCells =
543  mesh.boundary()[patchi].faceCells();
544 
545  forAll(lambdaPf, pFacei)
546  {
547  const label pfCelli = pFaceCells[pFacei];
548 
549  if (phiCorrfPf[pFacei] > 0)
550  {
551  lambdaPf[pFacei] =
552  min(lambdaPf[pFacei], lambdap[pfCelli]);
553  }
554  else
555  {
556  lambdaPf[pFacei] =
557  min(lambdaPf[pFacei], lambdam[pfCelli]);
558  }
559  }
560  }
561  }
562 
563  syncTools::syncFaceList(mesh, allLambda, minEqOp<scalar>());
564  }
565 }
566 
567 
568 template
569 <
570  class RdeltaTType,
571  class RhoType,
572  class SpType,
573  class SuType,
574  class PsiMaxType,
575  class PsiMinType
576 >
578 (
579  const RdeltaTType& rDeltaT,
580  const RhoType& rho,
581  const volScalarField& psi,
582  const surfaceScalarField& phi,
583  surfaceScalarField& phiPsi,
584  const SpType& Sp,
585  const SuType& Su,
586  const PsiMaxType& psiMax,
587  const PsiMinType& psiMin,
588  const bool returnCorr
589 )
590 {
591  const fvMesh& mesh = psi.mesh();
592 
593  surfaceScalarField phiBD(upwind<scalar>(psi.mesh(), phi).flux(psi));
594 
595  surfaceScalarField::Boundary& phiBDBf = phiBD.boundaryFieldRef();
596  const surfaceScalarField::Boundary& phiPsiBf = phiPsi.boundaryField();
597 
598  forAll(phiBDBf, patchi)
599  {
600  fvsPatchScalarField& phiBDPf = phiBDBf[patchi];
601 
602  if (!phiBDPf.coupled())
603  {
604  phiBDPf = phiPsiBf[patchi];
605  }
606  }
607 
608  surfaceScalarField& phiCorr = phiPsi;
609  phiCorr -= phiBD;
610 
611  scalarField allLambda(mesh.nFaces(), 1.0);
612 
614  (
615  IOobject
616  (
617  "lambda",
618  mesh.time().timeName(),
619  mesh,
620  IOobject::NO_READ,
621  IOobject::NO_WRITE,
622  false
623  ),
624  mesh,
625  dimless,
626  allLambda,
627  false // Use slices for the couples
628  );
629 
630  limiter
631  (
632  allLambda,
633  rDeltaT,
634  rho,
635  psi,
636  phiBD,
637  phiCorr,
638  Sp,
639  Su,
640  psiMax,
641  psiMin
642  );
643 
644  if (returnCorr)
645  {
646  phiCorr *= lambda;
647  }
648  else
649  {
650  phiPsi = phiBD + lambda*phiCorr;
651  }
652 }
653 
654 
655 template
656 <
657  class RhoType,
658  class SpType,
659  class SuType,
660  class PsiMaxType,
661  class PsiMinType
662 >
664 (
665  const RhoType& rho,
666  const volScalarField& psi,
667  const surfaceScalarField& phi,
668  surfaceScalarField& phiPsi,
669  const SpType& Sp,
670  const SuType& Su,
671  const PsiMaxType& psiMax,
672  const PsiMinType& psiMin,
673  const bool rtnCorr
674 )
675 {
676  const fvMesh& mesh = psi.mesh();
677 
678  if (fv::localEulerDdt::enabled(mesh))
679  {
680  const volScalarField& rDeltaT = fv::localEulerDdt::localRDeltaT(mesh);
681  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, rtnCorr);
682  }
683  else
684  {
685  const scalar rDeltaT = 1.0/mesh.time().deltaTValue();
686  limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, rtnCorr);
687  }
688 }
689 
690 
691 template<template<class> class AlphaList, template<class> class PhiList>
693 (
694  const AlphaList<const volScalarField>& alphas,
695  PhiList<surfaceScalarField>& phiPsis,
696  const surfaceScalarField& phi
697 )
698 {
699  PtrList<surfaceScalarField> alphaPhiUDs(phiPsis.size());
700 
701  forAll(phiPsis, phasei)
702  {
703  alphaPhiUDs.set
704  (
705  phasei,
706  upwind<scalar>(phi.mesh(), phi).flux(alphas[phasei])
707  );
708 
709  phiPsis[phasei] -= alphaPhiUDs[phasei];
710  }
711 
712  {
713  UPtrList<scalarField> phiPsisInternal(phiPsis.size());
714 
715  forAll(phiPsisInternal, phasei)
716  {
717  phiPsisInternal.set(phasei, &phiPsis[phasei]);
718  }
719 
720  limitSum(phiPsisInternal);
721  }
722 
723  const surfaceScalarField::Boundary& phibf = phi.boundaryField();
724 
725  forAll(phibf, patchi)
726  {
727  if (phibf[patchi].coupled())
728  {
729  UPtrList<scalarField> phiPsisPatch(phiPsis.size());
730 
731  forAll(phiPsisPatch, phasei)
732  {
733  phiPsisPatch.set
734  (
735  phasei,
736  &phiPsis[phasei].boundaryFieldRef()[patchi]
737  );
738  }
739 
740  limitSum(phiPsisPatch);
741  }
742  }
743 
744  forAll(phiPsis, phasei)
745  {
746  phiPsis[phasei] += alphaPhiUDs[phasei];
747  }
748 }
749 
750 
751 // ************************************************************************* //
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.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
Definition: IOobject.H:303
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
void surfaceIntegrate(Field< Type > &ivf, const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
bool moving() const
Is mesh moving.
Definition: polyMesh.H:512
label phasei
Definition: pEqn.H:27
virtual bool fixesValue() const
Return true if this patch field fixes a value.
Definition: fvPatchField.H:296
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
label nFaces() const
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
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
const dimensionSet dimless
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
const dictionary & solverDict(const word &name) const
Return the solver controls dictionary for the given field.
Definition: solution.C:342
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:284
tmp< DimensionedField< scalar, volMesh > > Vsc() const
Return sub-cycle cell volumes.
dynamicFvMesh & mesh
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
void limiter(scalarField &allLambda, 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)
virtual tmp< Field< Type > > patchNeighbourField() const
Return patchField on the opposite patch of a coupled patch.
Definition: fvPatchField.H:411
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:309
phi
Definition: correctPhi.H:3
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-cycl old-time cell volumes.
const Mesh & mesh() const
Return mesh.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:278
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.
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
Upwind differencing scheme class.
Definition: upwind.H:51
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
MULES: Multidimensional universal limiter for explicit solution.
const volScalarField & psi
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
Specialisation of GeometricField which holds slices of given complete fields in a form that they act ...
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:92
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:540