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