multiphaseSystem.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2015 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 "multiphaseSystem.H"
27 #include "alphaContactAngleFvPatchScalarField.H"
28 
29 #include "MULES.H"
30 #include "subCycle.H"
31 
32 #include "fvcDdt.H"
33 #include "fvcDiv.H"
34 #include "fvcSnGrad.H"
35 #include "fvcFlux.H"
36 #include "fvcMeshPhi.H"
37 #include "fvcSup.H"
38 
39 #include "fvmDdt.H"
40 #include "fvmLaplacian.H"
41 #include "fvmSup.H"
42 
43 // * * * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47  defineTypeNameAndDebug(multiphaseSystem, 0);
48  defineRunTimeSelectionTable(multiphaseSystem, dictionary);
49 }
50 
51 const Foam::scalar Foam::multiphaseSystem::convertToRad =
53 
54 
55 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
56 
57 void Foam::multiphaseSystem::calcAlphas()
58 {
59  scalar level = 0.0;
60  alphas_ == 0.0;
61 
62  forAll(phases(), i)
63  {
64  alphas_ += level*phases()[i];
65  level += 1.0;
66  }
67 
68  alphas_.correctBoundaryConditions();
69 }
70 
71 
72 void Foam::multiphaseSystem::solveAlphas()
73 {
74  bool LTS = fv::localEulerDdt::enabled(mesh_);
75 
76  PtrList<surfaceScalarField> alphaPhiCorrs(phases().size());
77  forAll(phases(), phasei)
78  {
79  phaseModel& phase = phases()[phasei];
80  volScalarField& alpha1 = phase;
81 
82  phase.alphaPhi() =
83  dimensionedScalar("0", dimensionSet(0, 3, -1, 0, 0), 0);
84 
85  alphaPhiCorrs.set
86  (
87  phasei,
89  (
90  "phi" + alpha1.name() + "Corr",
91  fvc::flux
92  (
93  phi_,
94  phase,
95  "div(phi," + alpha1.name() + ')'
96  )
97  )
98  );
99 
100  surfaceScalarField& alphaPhiCorr = alphaPhiCorrs[phasei];
101 
102  forAll(phases(), phasej)
103  {
104  phaseModel& phase2 = phases()[phasej];
105  volScalarField& alpha2 = phase2;
106 
107  if (&phase2 == &phase) continue;
108 
109  surfaceScalarField phir(phase.phi() - phase2.phi());
110 
111  cAlphaTable::const_iterator cAlpha
112  (
113  cAlphas_.find(phasePairKey(phase.name(), phase2.name()))
114  );
115 
116  if (cAlpha != cAlphas_.end())
117  {
119  (
120  (mag(phi_) + mag(phir))/mesh_.magSf()
121  );
122 
123  phir += min(cAlpha()*phic, max(phic))*nHatf(phase, phase2);
124  }
125 
126  word phirScheme
127  (
128  "div(phir," + alpha2.name() + ',' + alpha1.name() + ')'
129  );
130 
131  alphaPhiCorr += fvc::flux
132  (
133  -fvc::flux(-phir, phase2, phirScheme),
134  phase,
135  phirScheme
136  );
137  }
138 
139  // Ensure that the flux at inflow BCs is preserved
140  forAll(alphaPhiCorr.boundaryField(), patchi)
141  {
142  fvsPatchScalarField& alphaPhiCorrp =
143  alphaPhiCorr.boundaryField()[patchi];
144 
145  if (!alphaPhiCorrp.coupled())
146  {
147  const scalarField& phi1p = phase.phi().boundaryField()[patchi];
148  const scalarField& alpha1p = alpha1.boundaryField()[patchi];
149 
150  forAll(alphaPhiCorrp, facei)
151  {
152  if (phi1p[facei] < 0)
153  {
154  alphaPhiCorrp[facei] = alpha1p[facei]*phi1p[facei];
155  }
156  }
157  }
158  }
159 
160  if (LTS)
161  {
163  (
165  geometricOneField(),
166  phase,
167  phi_,
168  alphaPhiCorr,
169  zeroField(),
170  zeroField(),
171  phase.alphaMax(),
172  0,
173  true
174  );
175  }
176  else
177  {
178  const scalar rDeltaT = 1.0/mesh_.time().deltaTValue();
179 
181  (
182  rDeltaT,
183  geometricOneField(),
184  phase,
185  phi_,
186  alphaPhiCorr,
187  zeroField(),
188  zeroField(),
189  phase.alphaMax(),
190  0,
191  true
192  );
193  }
194  }
195 
196  MULES::limitSum(alphaPhiCorrs);
197 
198  volScalarField sumAlpha
199  (
200  IOobject
201  (
202  "sumAlpha",
203  mesh_.time().timeName(),
204  mesh_
205  ),
206  mesh_,
207  dimensionedScalar("sumAlpha", dimless, 0)
208  );
209 
210 
211  volScalarField divU(fvc::div(fvc::absolute(phi_, phases().first().U())));
212 
213  forAll(phases(), phasei)
214  {
215  phaseModel& phase = phases()[phasei];
216  volScalarField& alpha = phase;
217 
218  surfaceScalarField& alphaPhic = alphaPhiCorrs[phasei];
219  alphaPhic += upwind<scalar>(mesh_, phi_).flux(phase);
220 
222  (
223  IOobject
224  (
225  "Sp",
226  mesh_.time().timeName(),
227  mesh_
228  ),
229  mesh_,
230  dimensionedScalar("Sp", divU.dimensions(), 0.0)
231  );
232 
234  (
235  IOobject
236  (
237  "Su",
238  mesh_.time().timeName(),
239  mesh_
240  ),
241  // Divergence term is handled explicitly to be
242  // consistent with the explicit transport solution
243  divU*min(alpha, scalar(1))
244  );
245 
246  if (phase.divU().valid())
247  {
248  const scalarField& dgdt = phase.divU()();
249 
250  forAll(dgdt, celli)
251  {
252  if (dgdt[celli] > 0.0)
253  {
254  Sp[celli] -= dgdt[celli];
255  Su[celli] += dgdt[celli];
256  }
257  else if (dgdt[celli] < 0.0)
258  {
259  Sp[celli] +=
260  dgdt[celli]
261  *(1.0 - alpha[celli])/max(alpha[celli], 1e-4);
262  }
263  }
264  }
265 
266  forAll(phases(), phasej)
267  {
268  const phaseModel& phase2 = phases()[phasej];
269  const volScalarField& alpha2 = phase2;
270 
271  if (&phase2 == &phase) continue;
272 
273  if (phase2.divU().valid())
274  {
275  const scalarField& dgdt2 = phase2.divU()();
276 
277  forAll(dgdt2, celli)
278  {
279  if (dgdt2[celli] < 0.0)
280  {
281  Sp[celli] +=
282  dgdt2[celli]
283  *(1.0 - alpha2[celli])/max(alpha2[celli], 1e-4);
284 
285  Su[celli] -=
286  dgdt2[celli]
287  *alpha[celli]/max(alpha2[celli], 1e-4);
288  }
289  else if (dgdt2[celli] > 0.0)
290  {
291  Sp[celli] -= dgdt2[celli];
292  }
293  }
294  }
295  }
296 
298  (
299  geometricOneField(),
300  alpha,
301  alphaPhic,
302  Sp,
303  Su
304  );
305 
306  phase.alphaPhi() += alphaPhic;
307 
308  Info<< phase.name() << " volume fraction, min, max = "
309  << phase.weightedAverage(mesh_.V()).value()
310  << ' ' << min(phase).value()
311  << ' ' << max(phase).value()
312  << endl;
313 
314  sumAlpha += phase;
315  }
316 
317  Info<< "Phase-sum volume fraction, min, max = "
318  << sumAlpha.weightedAverage(mesh_.V()).value()
319  << ' ' << min(sumAlpha).value()
320  << ' ' << max(sumAlpha).value()
321  << endl;
322 }
323 
324 
325 Foam::tmp<Foam::surfaceVectorField> Foam::multiphaseSystem::nHatfv
326 (
327  const volScalarField& alpha1,
328  const volScalarField& alpha2
329 ) const
330 {
331  /*
332  // Cell gradient of alpha
333  volVectorField gradAlpha =
334  alpha2*fvc::grad(alpha1) - alpha1*fvc::grad(alpha2);
335 
336  // Interpolated face-gradient of alpha
337  surfaceVectorField gradAlphaf = fvc::interpolate(gradAlpha);
338  */
339 
340  surfaceVectorField gradAlphaf
341  (
343  - fvc::interpolate(alpha1)*fvc::interpolate(fvc::grad(alpha2))
344  );
345 
346  // Face unit interface normal
347  return gradAlphaf/(mag(gradAlphaf) + deltaN_);
348 }
349 
350 
351 Foam::tmp<Foam::surfaceScalarField> Foam::multiphaseSystem::nHatf
352 (
353  const volScalarField& alpha1,
354  const volScalarField& alpha2
355 ) const
356 {
357  // Face unit interface normal flux
358  return nHatfv(alpha1, alpha2) & mesh_.Sf();
359 }
360 
361 
362 // Correction for the boundary condition on the unit normal nHat on
363 // walls to produce the correct contact angle.
364 
365 // The dynamic contact angle is calculated from the component of the
366 // velocity on the direction of the interface, parallel to the wall.
367 
368 void Foam::multiphaseSystem::correctContactAngle
369 (
370  const phaseModel& phase1,
371  const phaseModel& phase2,
372  surfaceVectorField::GeometricBoundaryField& nHatb
373 ) const
374 {
375  const volScalarField::GeometricBoundaryField& gbf
376  = phase1.boundaryField();
377 
378  const fvBoundaryMesh& boundary = mesh_.boundary();
379 
380  forAll(boundary, patchi)
381  {
382  if (isA<alphaContactAngleFvPatchScalarField>(gbf[patchi]))
383  {
384  const alphaContactAngleFvPatchScalarField& acap =
385  refCast<const alphaContactAngleFvPatchScalarField>(gbf[patchi]);
386 
387  vectorField& nHatPatch = nHatb[patchi];
388 
389  vectorField AfHatPatch
390  (
391  mesh_.Sf().boundaryField()[patchi]
392  /mesh_.magSf().boundaryField()[patchi]
393  );
394 
395  alphaContactAngleFvPatchScalarField::thetaPropsTable::
396  const_iterator tp =
397  acap.thetaProps()
398  .find(phasePairKey(phase1.name(), phase2.name()));
399 
400  if (tp == acap.thetaProps().end())
401  {
403  << "Cannot find interface "
404  << phasePairKey(phase1.name(), phase2.name())
405  << "\n in table of theta properties for patch "
406  << acap.patch().name()
407  << exit(FatalError);
408  }
409 
410  bool matched = (tp.key().first() == phase1.name());
411 
412  scalar theta0 = convertToRad*tp().theta0(matched);
413  scalarField theta(boundary[patchi].size(), theta0);
414 
415  scalar uTheta = tp().uTheta();
416 
417  // Calculate the dynamic contact angle if required
418  if (uTheta > SMALL)
419  {
420  scalar thetaA = convertToRad*tp().thetaA(matched);
421  scalar thetaR = convertToRad*tp().thetaR(matched);
422 
423  // Calculated the component of the velocity parallel to the wall
424  vectorField Uwall
425  (
426  phase1.U()().boundaryField()[patchi].patchInternalField()
427  - phase1.U()().boundaryField()[patchi]
428  );
429  Uwall -= (AfHatPatch & Uwall)*AfHatPatch;
430 
431  // Find the direction of the interface parallel to the wall
432  vectorField nWall
433  (
434  nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch
435  );
436 
437  // Normalise nWall
438  nWall /= (mag(nWall) + SMALL);
439 
440  // Calculate Uwall resolved normal to the interface parallel to
441  // the interface
442  scalarField uwall(nWall & Uwall);
443 
444  theta += (thetaA - thetaR)*tanh(uwall/uTheta);
445  }
446 
447 
448  // Reset nHatPatch to correspond to the contact angle
449 
450  scalarField a12(nHatPatch & AfHatPatch);
451 
452  scalarField b1(cos(theta));
453 
454  scalarField b2(nHatPatch.size());
455 
456  forAll(b2, facei)
457  {
458  b2[facei] = cos(acos(a12[facei]) - theta[facei]);
459  }
460 
461  scalarField det(1.0 - a12*a12);
462 
463  scalarField a((b1 - a12*b2)/det);
464  scalarField b((b2 - a12*b1)/det);
465 
466  nHatPatch = a*AfHatPatch + b*nHatPatch;
467 
468  nHatPatch /= (mag(nHatPatch) + deltaN_.value());
469  }
470  }
471 }
472 
473 
474 Foam::tmp<Foam::volScalarField> Foam::multiphaseSystem::K
475 (
476  const phaseModel& phase1,
477  const phaseModel& phase2
478 ) const
479 {
480  tmp<surfaceVectorField> tnHatfv = nHatfv(phase1, phase2);
481 
482  correctContactAngle(phase1, phase2, tnHatfv().boundaryField());
483 
484  // Simple expression for curvature
485  return -fvc::div(tnHatfv & mesh_.Sf());
486 }
487 
488 
489 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
490 
492 (
493  const fvMesh& mesh
494 )
495 :
496  phaseSystem(mesh),
497 
498  alphas_
499  (
500  IOobject
501  (
502  "alphas",
503  mesh_.time().timeName(),
504  mesh,
505  IOobject::NO_READ,
506  IOobject::AUTO_WRITE
507  ),
508  mesh,
509  dimensionedScalar("alphas", dimless, 0.0),
510  zeroGradientFvPatchScalarField::typeName
511  ),
512 
513  cAlphas_(lookup("interfaceCompression")),
514 
515  deltaN_
516  (
517  "deltaN",
518  1e-8/pow(average(mesh_.V()), 1.0/3.0)
519  )
520 {
521  forAll(phases(), phasei)
522  {
523  volScalarField& alphai = phases()[phasei];
524  mesh_.setFluxRequired(alphai.name());
525  }
526 }
527 
528 
529 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
530 
532 {}
533 
534 
535 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
536 
538 (
539  const phaseModel& phase1
540 ) const
541 {
542  tmp<surfaceScalarField> tSurfaceTension
543  (
545  (
546  IOobject
547  (
548  "surfaceTension",
549  mesh_.time().timeName(),
550  mesh_
551  ),
552  mesh_,
554  (
555  "surfaceTension",
556  dimensionSet(1, -2, -2, 0, 0),
557  0
558  )
559  )
560  );
561 
562  forAll(phases(), phasej)
563  {
564  const phaseModel& phase2 = phases()[phasej];
565 
566  if (&phase2 != &phase1)
567  {
568  phasePairKey key12(phase1.name(), phase2.name());
569 
570  cAlphaTable::const_iterator cAlpha(cAlphas_.find(key12));
571 
572  if (cAlpha != cAlphas_.end())
573  {
574  tSurfaceTension() +=
575  fvc::interpolate(sigma(key12)*K(phase1, phase2))
576  *(
577  fvc::interpolate(phase2)*fvc::snGrad(phase1)
578  - fvc::interpolate(phase1)*fvc::snGrad(phase2)
579  );
580  }
581  }
582  }
583 
584  return tSurfaceTension;
585 }
586 
587 
590 {
591  tmp<volScalarField> tnearInt
592  (
593  new volScalarField
594  (
595  IOobject
596  (
597  "nearInterface",
598  mesh_.time().timeName(),
599  mesh_
600  ),
601  mesh_,
602  dimensionedScalar("nearInterface", dimless, 0.0)
603  )
604  );
605 
606  forAll(phases(), phasei)
607  {
608  tnearInt() = max
609  (
610  tnearInt(),
611  pos(phases()[phasei] - 0.01)*pos(0.99 - phases()[phasei])
612  );
613  }
614 
615  return tnearInt;
616 }
617 
618 
620 {
621  const Time& runTime = mesh_.time();
622 
623  const dictionary& alphaControls = mesh_.solverDict("alpha");
624  label nAlphaSubCycles(readLabel(alphaControls.lookup("nAlphaSubCycles")));
625 
626  bool LTS = fv::localEulerDdt::enabled(mesh_);
627 
628  if (nAlphaSubCycles > 1)
629  {
630  tmp<volScalarField> trSubDeltaT;
631 
632  if (LTS)
633  {
634  trSubDeltaT =
636  }
637 
638  PtrList<volScalarField> alpha0s(phases().size());
639  PtrList<surfaceScalarField> alphaPhiSums(phases().size());
640 
641  forAll(phases(), phasei)
642  {
643  phaseModel& phase = phases()[phasei];
644  volScalarField& alpha = phase;
645 
646  alpha0s.set
647  (
648  phasei,
649  new volScalarField(alpha.oldTime())
650  );
651 
652  alphaPhiSums.set
653  (
654  phasei,
656  (
657  IOobject
658  (
659  "phiSum" + alpha.name(),
660  runTime.timeName(),
661  mesh_
662  ),
663  mesh_,
664  dimensionedScalar("0", dimensionSet(0, 3, -1, 0, 0), 0)
665  )
666  );
667  }
668 
669  for
670  (
671  subCycleTime alphaSubCycle
672  (
673  const_cast<Time&>(runTime),
675  );
676  !(++alphaSubCycle).end();
677  )
678  {
679  solveAlphas();
680 
681  forAll(phases(), phasei)
682  {
683  alphaPhiSums[phasei] += phases()[phasei].alphaPhi();
684  }
685  }
686 
687  forAll(phases(), phasei)
688  {
689  phaseModel& phase = phases()[phasei];
690  volScalarField& alpha = phase;
691 
692  phase.alphaPhi() = alphaPhiSums[phasei]/nAlphaSubCycles;
693 
694  // Correct the time index of the field
695  // to correspond to the global time
696  alpha.timeIndex() = runTime.timeIndex();
697 
698  // Reset the old-time field value
699  alpha.oldTime() = alpha0s[phasei];
700  alpha.oldTime().timeIndex() = runTime.timeIndex();
701  }
702  }
703  else
704  {
705  solveAlphas();
706  }
707 
708  forAll(phases(), phasei)
709  {
710  phaseModel& phase = phases()[phasei];
711  phase.alphaRhoPhi() = fvc::interpolate(phase.rho())*phase.alphaPhi();
712  }
713 
714  calcAlphas();
715 }
716 
717 
718 // ************************************************************************* //
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf, const surfaceScalarField &faceFlux, Istream &schemeData)
phaseModel & phase1
Definition: createFields.H:12
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
void solve()
Solve for the mixture phase-fractions.
const phaseModel & phase() const
Return the phase.
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
phaseModel & phase2
Definition: createFields.H:13
faceListList boundary(nPatches)
Calculate the snGrad of the given volField.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux(const surfaceScalarField &phi, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcFlux.C:45
volScalarField divU(fvc::div(fvc::absolute(phi, U)))
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:28
Calculate the matrix for implicit and explicit sources.
dimensioned< scalar > mag(const dimensioned< Type > &)
label phasei
Definition: pEqn.H:37
CGAL::Exact_predicates_exact_constructions_kernel K
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
virtual ~multiphaseSystem()
Destructor.
surfaceScalarField phir(phic *interface.nHatf())
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
tmp< GeometricField< Type, fvPatchField, volMesh > > Sp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:67
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
const surfaceScalarField & alphaPhi() const
Definition: phaseModel.H:200
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Calulate the matrix for the first temporal derivative.
messageStream Info
dynamicFvMesh & mesh
multiphaseSystem(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
const dictionary & alphaControls
Definition: alphaControls.H:1
MULES: Multidimensional universal limiter for explicit solution.
Namespace for OpenFOAM.
fvsPatchField< scalar > fvsPatchScalarField
label readLabel(Istream &is)
Definition: label.H:64
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
tmp< GeometricField< Type, fvPatchField, volMesh > > Su(const GeometricField< Type, fvPatchField, volMesh > &su, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:44
const double e
Elementary charge.
Definition: doubleFloat.H:78
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
stressControl lookup("compactNormalStress") >> compactNormalStress
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:187
Calculate the first temporal derivative.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
Calculate the field for explicit evaluation of implicit and explicit sources.
static const volScalarField & localRDeltaT(const fvMesh &mesh)
Return the reciprocal of the local time-step.
Definition: localEulerDdt.C:45
dimensionedScalar det(const dimensionedSphericalTensor &dt)
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
Definition: MULES.C:51
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
Definition: localEulerDdt.C:57
Field< vector > vectorField
Specialisation of Field<T> for vector.
#define forAll(list, i)
Definition: UList.H:421
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
dimensionedScalar cos(const dimensionedScalar &ds)
label patchi
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const dimensionSet & dimensions() const
Return dimensions.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
tmp< surfaceScalarField > surfaceTension(const phaseModel &phase) const
void limit(const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phi, surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su, const scalar psiMax, const scalar psiMin, const bool returnCorr)
Calculate the matrix for the laplacian of the field.
label nAlphaSubCycles(readLabel(alphaControls.lookup("nAlphaSubCycles")))
Calculate the divergence of the given field.
error FatalError
DimensionedField< scalar, volMesh > DimensionedInternalField
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
multiphaseSystem::phaseModelList & phases
Definition: createFields.H:11
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
dimensionedScalar pos(const dimensionedScalar &ds)
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
Definition: localEulerDdt.C:36
U
Definition: pEqn.H:82
const Type & value() const
Return const reference to value.
defineTypeNameAndDebug(combustionModel, 0)
word timeName
Definition: getTimeIndex.H:3
Calculate the face-flux of the given field.