27 #include "alphaContactAngleFvPatchScalarField.H" 31 #include "surfaceInterpolate.H" 39 const Foam::scalar Foam::multiphaseMixture::convertToRad =
45 void Foam::multiphaseMixture::calcAlphas()
50 forAllIter(PtrDictionary<phase>, phases_, iter)
52 alphas_ += level*iter();
70 "transportProperties",
78 phases_(
lookup(
"phases"), phase::iNew(U, phi)),
123 sigmas_(
lookup(
"sigmas")),
124 dimSigma_(1, 0, -2, 0, 0),
143 tmp<volScalarField> trho = iter()*iter().rho();
146 for (++iter; iter != phases_.end(); ++iter)
148 rho += iter()*iter().rho();
160 tmp<scalarField> trho = iter().boundaryField()[
patchi]*iter().rho().value();
163 for (++iter; iter != phases_.end(); ++iter)
165 rho += iter().boundaryField()[
patchi]*iter().rho().value();
177 tmp<volScalarField> tmu = iter()*iter().rho()*iter().nu();
180 for (++iter; iter != phases_.end(); ++iter)
182 mu += iter()*iter().rho()*iter().nu();
194 tmp<scalarField> tmu =
195 iter().boundaryField()[
patchi]
196 *iter().rho().value()
200 for (++iter; iter != phases_.end(); ++iter)
203 iter().boundaryField()[
patchi]
204 *iter().rho().value()
217 tmp<surfaceScalarField> tmuf =
221 for (++iter; iter != phases_.end(); ++iter)
255 tmp<surfaceScalarField> tstf
261 "surfaceTensionForce",
268 "surfaceTensionForce",
269 dimensionSet(1, -2, -2, 0, 0),
279 const phase& alpha1 = iter1();
284 for (; iter2 != phases_.end(); ++iter2)
286 const phase& alpha2 = iter2();
289 sigmas_.
find(interfacePair(alpha1, alpha2));
291 if (sigma == sigmas_.
end())
294 <<
"Cannot find interface " << interfacePair(alpha1, alpha2)
295 <<
" in list of sigma values" 316 const Time& runTime = mesh_.
time();
322 scalar cAlpha(
readScalar(alphaControls.lookup(
"cAlpha")));
343 !(++alphaSubCycle).
end();
347 rhoPhiSum += (runTime.deltaT()/totalDeltaT)*rhoPhi_;
364 forAllIter(PtrDictionary<phase>, phases_, iter)
393 return gradAlphaf/(
mag(gradAlphaf) + deltaN_);
404 return nHatfv(alpha1, alpha2) & mesh_.
Sf();
414 void Foam::multiphaseMixture::correctContactAngle
418 surfaceVectorField::Boundary& nHatb
421 const volScalarField::Boundary& gbf
422 = alpha1.boundaryField();
428 if (isA<alphaContactAngleFvPatchScalarField>(gbf[patchi]))
430 const alphaContactAngleFvPatchScalarField& acap =
431 refCast<const alphaContactAngleFvPatchScalarField>(gbf[
patchi]);
443 acap.thetaProps().find(interfacePair(alpha1, alpha2));
445 if (tp == acap.thetaProps().end())
448 <<
"Cannot find interface " << interfacePair(alpha1, alpha2)
449 <<
"\n in table of theta properties for patch " 450 << acap.patch().name()
454 bool matched = (tp.key().first() == alpha1.name());
456 scalar theta0 = convertToRad*tp().theta0(matched);
459 scalar uTheta = tp().uTheta();
464 scalar thetaA = convertToRad*tp().thetaA(matched);
465 scalar thetaR = convertToRad*tp().thetaR(matched);
473 Uwall -= (AfHatPatch & Uwall)*AfHatPatch;
478 nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch
482 nWall /= (
mag(nWall) + SMALL);
488 theta += (thetaA - thetaR)*
tanh(uwall/uTheta);
502 b2[facei] =
cos(
acos(a12[facei]) - theta[facei]);
510 nHatPatch = a*AfHatPatch +
b*nHatPatch;
512 nHatPatch /= (
mag(nHatPatch) + deltaN_.
value());
524 tmp<surfaceVectorField> tnHatfv = nHatfv(alpha1, alpha2);
526 correctContactAngle(alpha1, alpha2, tnHatfv.ref().boundaryFieldRef());
536 tmp<volScalarField> tnearInt
553 tnearInt.ref() =
max(tnearInt(),
pos(iter() - 0.01)*
pos(0.99 - iter()));
560 void Foam::multiphaseMixture::solveAlphas
565 static label nSolves=-1;
568 word alphaScheme(
"div(phi,alpha)");
574 PtrList<surfaceScalarField> alphaPhiCorrs(phases_.size());
577 forAllIter(PtrDictionary<phase>, phases_, iter)
579 phase& alpha = iter();
586 "phi" + alpha.name() +
"Corr",
598 forAllIter(PtrDictionary<phase>, phases_, iter2)
600 phase& alpha2 = iter2();
602 if (&alpha2 == &alpha)
continue;
649 forAllIter(PtrDictionary<phase>, phases_, iter)
651 phase& alpha = iter();
654 alphaPhi += upwind<scalar>(mesh_, phi_).
flux(alpha);
665 rhoPhi_ += alphaPhi*alpha.rho();
667 Info<< alpha.name() <<
" volume fraction, min, max = " 668 << alpha.weightedAverage(mesh_.V()).value()
669 <<
' ' <<
min(alpha).value()
670 <<
' ' <<
max(alpha).value()
678 Info<<
"Phase-sum volume fraction, min, max = " 679 << sumAlpha.weightedAverage(mesh_.V()).value()
680 <<
' ' <<
min(sumAlpha).value()
681 <<
' ' <<
max(sumAlpha).value()
694 PtrList<entry> phaseData(
lookup(
"phases"));
697 forAllIter(PtrDictionary<phase>, phases_, iter)
699 readOK &= iter().read(phaseData[phasei++].
dict());
702 lookup(
"sigmas") >> sigmas_;
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
dimensionedScalar tanh(const dimensionedScalar &ds)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
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.
const surfaceVectorField & Sf() const
Return cell face area vectors.
const objectRegistry & db() const
Return the local objectRegistry.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
const double e
Elementary charge.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
IOobject(const word &name, const fileName &instance, const objectRegistry ®istry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
word alpharScheme("div(phirb,alpha)")
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
Ostream & endl(Ostream &os)
Add newline and flush stream.
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.
dimensionedScalar det(const dimensionedSphericalTensor &dt)
label size() const
Return number of elements in list.
GeometricField< vector, fvPatchField, volMesh > volVectorField
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.
dictionary()
Construct top-level dictionary null.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar pos(const dimensionedScalar &ds)
surfaceScalarField alphaPhi(IOobject("alphaPhi", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), phi *fvc::interpolate(alpha1))
dimensionedScalar cos(const dimensionedScalar &ds)
bool read()
Read base transportProperties dictionary.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
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.
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.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
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.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
label readLabel(Istream &is)
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.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
dimensionedScalar deltaT() const
Return time step.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
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.
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)
A class for managing temporary objects.
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
friend class const_iterator
tmp< volScalarField > rho() const
Return the mixture density.
const Time & time() const
Return the top-level database.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.