multiphaseMixture.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) 2011-2016 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 "multiphaseMixture.H"
27 #include "alphaContactAngleFvPatchScalarField.H"
28 #include "Time.H"
29 #include "subCycle.H"
30 #include "MULES.H"
31 #include "surfaceInterpolate.H"
32 #include "fvcGrad.H"
33 #include "fvcSnGrad.H"
34 #include "fvcDiv.H"
35 #include "fvcFlux.H"
36 
37 // * * * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * //
38 
39 const Foam::scalar Foam::multiphaseMixture::convertToRad =
41 
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void Foam::multiphaseMixture::calcAlphas()
46 {
47  scalar level = 0.0;
48  alphas_ == 0.0;
49 
50  forAllIter(PtrDictionary<phase>, phases_, iter)
51  {
52  alphas_ += level*iter();
53  level += 1.0;
54  }
55 }
56 
57 
58 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
59 
61 (
62  const volVectorField& U,
63  const surfaceScalarField& phi
64 )
65 :
67  (
68  IOobject
69  (
70  "transportProperties",
71  U.time().constant(),
72  U.db(),
75  )
76  ),
77 
78  phases_(lookup("phases"), phase::iNew(U, phi)),
79 
80  mesh_(U.mesh()),
81  U_(U),
82  phi_(phi),
83 
84  rhoPhi_
85  (
86  IOobject
87  (
88  "rhoPhi",
89  mesh_.time().timeName(),
90  mesh_,
93  ),
94  mesh_,
95  dimensionedScalar("rhoPhi", dimMass/dimTime, 0.0)
96  ),
97 
98  alphas_
99  (
100  IOobject
101  (
102  "alphas",
103  mesh_.time().timeName(),
104  mesh_,
105  IOobject::NO_READ,
107  ),
108  mesh_,
109  dimensionedScalar("alphas", dimless, 0.0)
110  ),
111 
112  nu_
113  (
114  IOobject
115  (
116  "nu",
117  mesh_.time().timeName(),
118  mesh_
119  ),
120  mu()/rho()
121  ),
122 
123  sigmas_(lookup("sigmas")),
124  dimSigma_(1, 0, -2, 0, 0),
125  deltaN_
126  (
127  "deltaN",
128  1e-8/pow(average(mesh_.V()), 1.0/3.0)
129  )
130 {
131  calcAlphas();
132  alphas_.write();
133 }
134 
135 
136 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
137 
140 {
141  PtrDictionary<phase>::const_iterator iter = phases_.begin();
142 
143  tmp<volScalarField> trho = iter()*iter().rho();
144  volScalarField& rho = trho.ref();
145 
146  for (++iter; iter != phases_.end(); ++iter)
147  {
148  rho += iter()*iter().rho();
149  }
150 
151  return trho;
152 }
153 
154 
156 Foam::multiphaseMixture::rho(const label patchi) const
157 {
158  PtrDictionary<phase>::const_iterator iter = phases_.begin();
159 
160  tmp<scalarField> trho = iter().boundaryField()[patchi]*iter().rho().value();
161  scalarField& rho = trho.ref();
162 
163  for (++iter; iter != phases_.end(); ++iter)
164  {
165  rho += iter().boundaryField()[patchi]*iter().rho().value();
166  }
167 
168  return trho;
169 }
170 
171 
174 {
175  PtrDictionary<phase>::const_iterator iter = phases_.begin();
176 
177  tmp<volScalarField> tmu = iter()*iter().rho()*iter().nu();
178  volScalarField& mu = tmu.ref();
179 
180  for (++iter; iter != phases_.end(); ++iter)
181  {
182  mu += iter()*iter().rho()*iter().nu();
183  }
184 
185  return tmu;
186 }
187 
188 
190 Foam::multiphaseMixture::mu(const label patchi) const
191 {
192  PtrDictionary<phase>::const_iterator iter = phases_.begin();
193 
194  tmp<scalarField> tmu =
195  iter().boundaryField()[patchi]
196  *iter().rho().value()
197  *iter().nu(patchi);
198  scalarField& mu = tmu.ref();
199 
200  for (++iter; iter != phases_.end(); ++iter)
201  {
202  mu +=
203  iter().boundaryField()[patchi]
204  *iter().rho().value()
205  *iter().nu(patchi);
206  }
207 
208  return tmu;
209 }
210 
211 
214 {
215  PtrDictionary<phase>::const_iterator iter = phases_.begin();
216 
217  tmp<surfaceScalarField> tmuf =
218  fvc::interpolate(iter())*iter().rho()*fvc::interpolate(iter().nu());
219  surfaceScalarField& muf = tmuf.ref();
220 
221  for (++iter; iter != phases_.end(); ++iter)
222  {
223  muf +=
224  fvc::interpolate(iter())*iter().rho()*fvc::interpolate(iter().nu());
225  }
226 
227  return tmuf;
228 }
229 
230 
233 {
234  return nu_;
235 }
236 
237 
239 Foam::multiphaseMixture::nu(const label patchi) const
240 {
241  return nu_.boundaryField()[patchi];
242 }
243 
244 
247 {
248  return muf()/fvc::interpolate(rho());
249 }
250 
251 
254 {
255  tmp<surfaceScalarField> tstf
256  (
258  (
259  IOobject
260  (
261  "surfaceTensionForce",
262  mesh_.time().timeName(),
263  mesh_
264  ),
265  mesh_,
267  (
268  "surfaceTensionForce",
269  dimensionSet(1, -2, -2, 0, 0),
270  0.0
271  )
272  )
273  );
274 
275  surfaceScalarField& stf = tstf.ref();
276 
277  forAllConstIter(PtrDictionary<phase>, phases_, iter1)
278  {
279  const phase& alpha1 = iter1();
280 
282  ++iter2;
283 
284  for (; iter2 != phases_.end(); ++iter2)
285  {
286  const phase& alpha2 = iter2();
287 
289  sigmas_.find(interfacePair(alpha1, alpha2));
290 
291  if (sigma == sigmas_.end())
292  {
294  << "Cannot find interface " << interfacePair(alpha1, alpha2)
295  << " in list of sigma values"
296  << exit(FatalError);
297  }
298 
299  stf += dimensionedScalar("sigma", dimSigma_, sigma())
300  *fvc::interpolate(K(alpha1, alpha2))*
301  (
302  fvc::interpolate(alpha2)*fvc::snGrad(alpha1)
303  - fvc::interpolate(alpha1)*fvc::snGrad(alpha2)
304  );
305  }
306  }
307 
308  return tstf;
309 }
310 
311 
313 {
314  correct();
315 
316  const Time& runTime = mesh_.time();
317 
318  volScalarField& alpha = phases_.first();
319 
320  const dictionary& alphaControls = mesh_.solverDict("alpha");
321  label nAlphaSubCycles(readLabel(alphaControls.lookup("nAlphaSubCycles")));
322  scalar cAlpha(readScalar(alphaControls.lookup("cAlpha")));
323 
324  if (nAlphaSubCycles > 1)
325  {
326  surfaceScalarField rhoPhiSum
327  (
328  IOobject
329  (
330  "rhoPhiSum",
331  runTime.timeName(),
332  mesh_
333  ),
334  mesh_,
335  dimensionedScalar("0", rhoPhi_.dimensions(), 0)
336  );
337 
338  dimensionedScalar totalDeltaT = runTime.deltaT();
339 
340  for
341  (
342  subCycle<volScalarField> alphaSubCycle(alpha, nAlphaSubCycles);
343  !(++alphaSubCycle).end();
344  )
345  {
346  solveAlphas(cAlpha);
347  rhoPhiSum += (runTime.deltaT()/totalDeltaT)*rhoPhi_;
348  }
349 
350  rhoPhi_ = rhoPhiSum;
351  }
352  else
353  {
354  solveAlphas(cAlpha);
355  }
356 
357  // Update the mixture kinematic viscosity
358  nu_ = mu()/rho();
359 }
360 
361 
363 {
364  forAllIter(PtrDictionary<phase>, phases_, iter)
365  {
366  iter().correct();
367  }
368 }
369 
370 
371 Foam::tmp<Foam::surfaceVectorField> Foam::multiphaseMixture::nHatfv
372 (
373  const volScalarField& alpha1,
374  const volScalarField& alpha2
375 ) const
376 {
377  /*
378  // Cell gradient of alpha
379  volVectorField gradAlpha =
380  alpha2*fvc::grad(alpha1) - alpha1*fvc::grad(alpha2);
381 
382  // Interpolated face-gradient of alpha
383  surfaceVectorField gradAlphaf = fvc::interpolate(gradAlpha);
384  */
385 
386  surfaceVectorField gradAlphaf
387  (
389  - fvc::interpolate(alpha1)*fvc::interpolate(fvc::grad(alpha2))
390  );
391 
392  // Face unit interface normal
393  return gradAlphaf/(mag(gradAlphaf) + deltaN_);
394 }
395 
396 
397 Foam::tmp<Foam::surfaceScalarField> Foam::multiphaseMixture::nHatf
398 (
399  const volScalarField& alpha1,
400  const volScalarField& alpha2
401 ) const
402 {
403  // Face unit interface normal flux
404  return nHatfv(alpha1, alpha2) & mesh_.Sf();
405 }
406 
407 
408 // Correction for the boundary condition on the unit normal nHat on
409 // walls to produce the correct contact angle.
410 
411 // The dynamic contact angle is calculated from the component of the
412 // velocity on the direction of the interface, parallel to the wall.
413 
414 void Foam::multiphaseMixture::correctContactAngle
415 (
416  const phase& alpha1,
417  const phase& alpha2,
418  surfaceVectorField::Boundary& nHatb
419 ) const
420 {
421  const volScalarField::Boundary& gbf
422  = alpha1.boundaryField();
423 
424  const fvBoundaryMesh& boundary = mesh_.boundary();
425 
426  forAll(boundary, patchi)
427  {
428  if (isA<alphaContactAngleFvPatchScalarField>(gbf[patchi]))
429  {
430  const alphaContactAngleFvPatchScalarField& acap =
431  refCast<const alphaContactAngleFvPatchScalarField>(gbf[patchi]);
432 
433  vectorField& nHatPatch = nHatb[patchi];
434 
435  vectorField AfHatPatch
436  (
437  mesh_.Sf().boundaryField()[patchi]
438  /mesh_.magSf().boundaryField()[patchi]
439  );
440 
442  const_iterator tp =
443  acap.thetaProps().find(interfacePair(alpha1, alpha2));
444 
445  if (tp == acap.thetaProps().end())
446  {
448  << "Cannot find interface " << interfacePair(alpha1, alpha2)
449  << "\n in table of theta properties for patch "
450  << acap.patch().name()
451  << exit(FatalError);
452  }
453 
454  bool matched = (tp.key().first() == alpha1.name());
455 
456  scalar theta0 = convertToRad*tp().theta0(matched);
457  scalarField theta(boundary[patchi].size(), theta0);
458 
459  scalar uTheta = tp().uTheta();
460 
461  // Calculate the dynamic contact angle if required
462  if (uTheta > SMALL)
463  {
464  scalar thetaA = convertToRad*tp().thetaA(matched);
465  scalar thetaR = convertToRad*tp().thetaR(matched);
466 
467  // Calculated the component of the velocity parallel to the wall
468  vectorField Uwall
469  (
470  U_.boundaryField()[patchi].patchInternalField()
471  - U_.boundaryField()[patchi]
472  );
473  Uwall -= (AfHatPatch & Uwall)*AfHatPatch;
474 
475  // Find the direction of the interface parallel to the wall
476  vectorField nWall
477  (
478  nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch
479  );
480 
481  // Normalise nWall
482  nWall /= (mag(nWall) + SMALL);
483 
484  // Calculate Uwall resolved normal to the interface parallel to
485  // the interface
486  scalarField uwall(nWall & Uwall);
487 
488  theta += (thetaA - thetaR)*tanh(uwall/uTheta);
489  }
490 
491 
492  // Reset nHatPatch to correspond to the contact angle
493 
494  scalarField a12(nHatPatch & AfHatPatch);
495 
496  scalarField b1(cos(theta));
497 
498  scalarField b2(nHatPatch.size());
499 
500  forAll(b2, facei)
501  {
502  b2[facei] = cos(acos(a12[facei]) - theta[facei]);
503  }
504 
505  scalarField det(1.0 - a12*a12);
506 
507  scalarField a((b1 - a12*b2)/det);
508  scalarField b((b2 - a12*b1)/det);
509 
510  nHatPatch = a*AfHatPatch + b*nHatPatch;
511 
512  nHatPatch /= (mag(nHatPatch) + deltaN_.value());
513  }
514  }
515 }
516 
517 
518 Foam::tmp<Foam::volScalarField> Foam::multiphaseMixture::K
519 (
520  const phase& alpha1,
521  const phase& alpha2
522 ) const
523 {
524  tmp<surfaceVectorField> tnHatfv = nHatfv(alpha1, alpha2);
525 
526  correctContactAngle(alpha1, alpha2, tnHatfv.ref().boundaryFieldRef());
527 
528  // Simple expression for curvature
529  return -fvc::div(tnHatfv & mesh_.Sf());
530 }
531 
532 
535 {
536  tmp<volScalarField> tnearInt
537  (
538  new volScalarField
539  (
540  IOobject
541  (
542  "nearInterface",
543  mesh_.time().timeName(),
544  mesh_
545  ),
546  mesh_,
547  dimensionedScalar("nearInterface", dimless, 0.0)
548  )
549  );
550 
551  forAllConstIter(PtrDictionary<phase>, phases_, iter)
552  {
553  tnearInt.ref() = max(tnearInt(), pos(iter() - 0.01)*pos(0.99 - iter()));
554  }
555 
556  return tnearInt;
557 }
558 
559 
560 void Foam::multiphaseMixture::solveAlphas
561 (
562  const scalar cAlpha
563 )
564 {
565  static label nSolves=-1;
566  nSolves++;
567 
568  word alphaScheme("div(phi,alpha)");
569  word alpharScheme("div(phirb,alpha)");
570 
571  surfaceScalarField phic(mag(phi_/mesh_.magSf()));
572  phic = min(cAlpha*phic, max(phic));
573 
574  PtrList<surfaceScalarField> alphaPhiCorrs(phases_.size());
575  int phasei = 0;
576 
577  forAllIter(PtrDictionary<phase>, phases_, iter)
578  {
579  phase& alpha = iter();
580 
581  alphaPhiCorrs.set
582  (
583  phasei,
585  (
586  "phi" + alpha.name() + "Corr",
587  fvc::flux
588  (
589  phi_,
590  alpha,
591  alphaScheme
592  )
593  )
594  );
595 
596  surfaceScalarField& alphaPhiCorr = alphaPhiCorrs[phasei];
597 
598  forAllIter(PtrDictionary<phase>, phases_, iter2)
599  {
600  phase& alpha2 = iter2();
601 
602  if (&alpha2 == &alpha) continue;
603 
604  surfaceScalarField phir(phic*nHatf(alpha, alpha2));
605 
606  alphaPhiCorr += fvc::flux
607  (
608  -fvc::flux(-phir, alpha2, alpharScheme),
609  alpha,
611  );
612  }
613 
615  (
616  1.0/mesh_.time().deltaT().value(),
617  geometricOneField(),
618  alpha,
619  phi_,
620  alphaPhiCorr,
621  zeroField(),
622  zeroField(),
623  1,
624  0,
625  true
626  );
627 
628  phasei++;
629  }
630 
631  MULES::limitSum(alphaPhiCorrs);
632 
633  rhoPhi_ = dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0);
634 
635  volScalarField sumAlpha
636  (
637  IOobject
638  (
639  "sumAlpha",
640  mesh_.time().timeName(),
641  mesh_
642  ),
643  mesh_,
644  dimensionedScalar("sumAlpha", dimless, 0)
645  );
646 
647  phasei = 0;
648 
649  forAllIter(PtrDictionary<phase>, phases_, iter)
650  {
651  phase& alpha = iter();
652 
653  surfaceScalarField& alphaPhi = alphaPhiCorrs[phasei];
654  alphaPhi += upwind<scalar>(mesh_, phi_).flux(alpha);
655 
657  (
658  geometricOneField(),
659  alpha,
660  alphaPhi,
661  zeroField(),
662  zeroField()
663  );
664 
665  rhoPhi_ += alphaPhi*alpha.rho();
666 
667  Info<< alpha.name() << " volume fraction, min, max = "
668  << alpha.weightedAverage(mesh_.V()).value()
669  << ' ' << min(alpha).value()
670  << ' ' << max(alpha).value()
671  << endl;
672 
673  sumAlpha += alpha;
674 
675  phasei++;
676  }
677 
678  Info<< "Phase-sum volume fraction, min, max = "
679  << sumAlpha.weightedAverage(mesh_.V()).value()
680  << ' ' << min(sumAlpha).value()
681  << ' ' << max(sumAlpha).value()
682  << endl;
683 
684  calcAlphas();
685 }
686 
687 
689 {
690  if (transportModel::read())
691  {
692  bool readOK = true;
693 
694  PtrList<entry> phaseData(lookup("phases"));
695  label phasei = 0;
696 
697  forAllIter(PtrDictionary<phase>, phases_, iter)
698  {
699  readOK &= iter().read(phaseData[phasei++].dict());
700  }
701 
702  lookup("sigmas") >> sigmas_;
703 
704  return readOK;
705  }
706  else
707  {
708  return false;
709  }
710 }
711 
712 
713 // ************************************************************************* //
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
dimensionedScalar tanh(const dimensionedScalar &ds)
dictionary dict
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
tmp< surfaceScalarField > nuf() const
Return the face-interpolated dynamic laminar viscosity.
tmp< surfaceScalarField > muf() const
Return the face-interpolated dynamic laminar viscosity.
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 surfaceVectorField & Sf() const
Return cell face area vectors.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:221
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
label phasei
Definition: pEqn.H:27
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
const double e
Elementary charge.
Definition: doubleFloat.H:78
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:114
word alpharScheme("div(phirb,alpha)")
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
tmp< volScalarField > trho
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
const Type & value() const
Return const reference to value.
Calculate the snGrad of the given volField.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
dimensionedScalar det(const dimensionedSphericalTensor &dt)
label size() const
Return number of elements in list.
Definition: DLListBaseI.H:77
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
tmp< volScalarField > nu() const
Return the kinematic laminar viscosity.
CGAL::Exact_predicates_exact_constructions_kernel K
multiphaseMixture(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
surfaceScalarField phir(phic *interface.nHatf())
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:138
dictionary()
Construct top-level dictionary null.
Definition: dictionary.C:113
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensionedScalar pos(const dimensionedScalar &ds)
surfaceScalarField alphaPhi(IOobject("alphaPhi", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), phi *fvc::interpolate(alpha1))
dynamicFvMesh & mesh
dimensionedScalar cos(const dimensionedScalar &ds)
bool read()
Read base transportProperties dictionary.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
Definition: MULES.C:51
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
Calculate the gradient of the given field.
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
const Time & time() const
Return time.
Definition: IOobject.C:227
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void correct()
Correct the mixture properties.
Calculate the face-flux of the given field.
friend class const_iterator
Declare friendship with the const_iterator.
Definition: HashTable.H:192
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
word timeName
Definition: getTimeIndex.H:3
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
faceListList boundary(nPatches)
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
label readLabel(Istream &is)
Definition: label.H:64
const dimensionSet & dimensions() const
Return dimensions.
Calculate the divergence of the given field.
tmp< volScalarField > mu() const
Return the dynamic laminar viscosity.
void solve()
Solve for the mixture phase-fractions.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Internal & ref()
Return a reference to the dimensioned internal field.
tmp< surfaceScalarField > surfaceTensionForce() const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
IOdictionary(const IOobject &)
Construct given an IOobject.
Definition: IOdictionary.C:45
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
virtual bool read()=0
Read transportProperties dictionary.
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
MULES: Multidimensional universal limiter for explicit solution.
const dictionary & solverDict(const word &name) const
Return the solver controls dictionary for the given field.
Definition: solution.C:342
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)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:54
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:545
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
tmp< volScalarField > rho() const
Return the mixture density.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
label nAlphaSubCycles(readLabel(alphaControls.lookup("nAlphaSubCycles")))
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451