Go to the documentation of this file.
30 #include "alphaContactAngleFvPatchScalarField.H"
60 void Foam::multiphaseSystem::calcAlphas()
67 alphas_ += level*
phases()[i];
73 void Foam::multiphaseSystem::solveAlphas()
86 mesh_.time().timeName(),
92 forAll(stationaryPhases(), stationaryPhasei)
94 alphaVoid -= stationaryPhases()[stationaryPhasei];
108 upwind<scalar>(mesh_, phi_).interpolate(phase)
114 PtrList<surfaceScalarField> alphaPhiCorrs(
phases().size());
115 forAll(stationaryPhases(), stationaryPhasei)
117 phaseModel& phase = stationaryPhases()[stationaryPhasei];
125 - upwind<scalar>(mesh_, phi_).flux(phase)
129 forAll(movingPhases(), movingPhasei)
131 phaseModel& phase = movingPhases()[movingPhasei];
151 if (&
phase2 == &phase)
continue;
155 cAlphaTable::const_iterator cAlpha
157 cAlphas_.find(phasePairKey(phase.name(),
phase2.
name()))
160 if (cAlpha != cAlphas_.end())
183 phase.correctInflowOutflow(alphaPhiCorr);
193 min(alphaVoid.primitiveField(), phase.alphaMax())(),
201 forAll(stationaryPhases(), stationaryPhasei)
203 fixedAlphaPhiCorrs.insert(stationaryPhases()[stationaryPhasei].index());
208 forAll(movingPhases(), movingPhasei)
210 phaseModel& phase = movingPhases()[movingPhasei];
214 alphaPhi += upwind<scalar>(mesh_, phi_).flux(phase);
215 phase.correctInflowOutflow(
alphaPhi);
222 mesh_.time().timeName(),
235 if (phase.divU().valid())
241 if (dgdt[celli] > 0.0)
243 Sp[celli] -= dgdt[celli];
244 Su[celli] += dgdt[celli];
246 else if (dgdt[celli] < 0.0)
260 if (&
phase2 == &phase)
continue;
268 if (dgdt2[celli] < 0.0)
278 else if (dgdt2[celli] > 0.0)
280 Sp[celli] -= dgdt2[celli];
303 Info<< phase.name() <<
" fraction, min, max = "
304 << phase.weightedAverage(mesh_.V()).value()
305 <<
' ' <<
min(phase).value()
306 <<
' ' <<
max(phase).value()
315 mesh_.time().timeName(),
321 forAll(movingPhases(), movingPhasei)
323 sumAlphaMoving += movingPhases()[movingPhasei];
326 Info<<
"Phase-sum volume fraction, min, max = "
327 << (sumAlphaMoving + 1 - alphaVoid)().weightedAverage(mesh_.V()).value()
328 <<
' ' <<
min(sumAlphaMoving + 1 - alphaVoid).value()
329 <<
' ' <<
max(sumAlphaMoving + 1 - alphaVoid).value()
333 forAll(movingPhases(), movingPhasei)
335 movingPhases()[movingPhasei] *= alphaVoid/sumAlphaMoving;
362 return gradAlphaf/(
mag(gradAlphaf) + deltaN_);
383 void Foam::multiphaseSystem::correctContactAngle
387 surfaceVectorField::Boundary& nHatb
390 const volScalarField::Boundary& gbf
393 const fvBoundaryMesh&
boundary = mesh_.boundary();
397 if (isA<alphaContactAngleFvPatchScalarField>(gbf[patchi]))
399 const alphaContactAngleFvPatchScalarField& acap =
400 refCast<const alphaContactAngleFvPatchScalarField>(gbf[patchi]);
406 mesh_.Sf().boundaryField()[patchi]
407 /mesh_.magSf().boundaryField()[patchi]
410 alphaContactAngleFvPatchScalarField::thetaPropsTable::
415 if (tp == acap.thetaProps().end())
418 <<
"Cannot find interface "
420 <<
"\n in table of theta properties for patch "
421 << acap.patch().name()
425 bool matched = (tp.key().first() ==
phase1.
name());
427 scalar theta0 =
degToRad(tp().theta0(matched));
430 scalar uTheta = tp().uTheta();
435 const scalar thetaA =
degToRad(tp().thetaA(matched));
436 const scalar thetaR =
degToRad(tp().thetaR(matched));
441 phase1.
U()().boundaryField()[patchi].patchInternalField()
442 -
phase1.
U()().boundaryField()[patchi]
444 Uwall -= (AfHatPatch & Uwall)*AfHatPatch;
449 nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch
453 nWall /= (
mag(nWall) + SMALL);
459 theta += (thetaA - thetaR)*
tanh(uwall/uTheta);
473 b2[facei] =
cos(
acos(a12[facei]) - theta[facei]);
481 nHatPatch = a*AfHatPatch +
b*nHatPatch;
483 nHatPatch /= (
mag(nHatPatch) + deltaN_.value());
495 tmp<surfaceVectorField> tnHatfv = nHatfv(
phase1,
phase2);
500 return -
fvc::div(tnHatfv & mesh_.Sf());
527 cAlphas_(lookup(
"interfaceCompression")),
550 tmp<surfaceScalarField> tSurfaceTension
560 tSurfaceTension.ref().setOriented();
570 cAlphaTable::const_iterator cAlpha(cAlphas_.find(key12));
572 if (cAlpha != cAlphas_.end())
574 tSurfaceTension.ref() +=
584 tSurfaceTension->setOriented();
586 return tSurfaceTension;
593 tmp<volScalarField> tnearInt
627 tmp<volScalarField> trSubDeltaT;
635 List<volScalarField*> alphaPtrs(
phases().size());
636 PtrList<surfaceScalarField> alphaPhiSums(
phases().size());
664 subCycleTime alphaSubCycle
669 !(++alphaSubCycle).
end();
683 if (phase.stationary())
continue;
696 if (phase.stationary())
continue;
698 phase.alphaRhoPhiRef() =
701 phase.clip(SMALL, 1 - SMALL);
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
void solve()
Solve for the mixture phase-fractions.
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)
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
const volScalarField & alpha2
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Calculate the snGrad of the given volField.
static word timeName(const scalar t, const int precision=precision_)
void setFluxRequired(const word &name) const
Get flux-required for given name, or default.
Calculate the divergence of the given field.
const word & name() const
Return const reference to name.
Unit conversion functions.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pos0(const dimensionedScalar &ds)
Calculate the field for explicit evaluation of implicit and explicit sources.
const volScalarField & alpha1
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
tmp< surfaceScalarField > surfaceTension(const phaseModel &phase) const
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
virtual tmp< volScalarField > divU() const =0
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
surfaceScalarField phic(mixture.cAlpha() *mag(alphaPhic/mesh.magSf()))
static tmp< GeometricField< scalar, fvsPatchField, surfaceMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &ds, const word &patchFieldType=fvsPatchField< scalar >::calculatedType())
Return tmp field from name, mesh, dimensions and patch type.
CGAL::Exact_predicates_exact_constructions_kernel K
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
surfaceScalarField phir(fvc::flux(UdmModel.Udm()))
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
dimensionedScalar tanh(const dimensionedScalar &ds)
surfaceScalarField alphaPhi(phi.name()+alpha1.name(), fvc::flux(phi, alpha1, alphaScheme))
messageStream Info
Information stream (stdout output on master, null elsewhere)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
Calculate the matrix for the laplacian of the field.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
PtrList< surfaceScalarField > alphafs(phases.size())
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.
const PtrDictionary< phaseModel > & phases() const
Return the phases.
const surfaceScalarField & phi() const
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
Calculate the matrix for implicit and explicit sources.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const word & name() const
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const volVectorField & U() const
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Calculate the face-flux of the given field.
const dimensionedScalar e
Elementary charge.
Calculate the first temporal derivative.
dimensionedScalar det(const dimensionedSphericalTensor &dt)
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
MULES: Multidimensional universal limiter for explicit solution.
defineTypeNameAndDebug(combustionModel, 0)
const Time & time() const noexcept
Return time registry.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
multiphaseSystem(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
Calculate the matrix for the first temporal derivative.
const dimensionSet dimless
Dimensionless.
multiphaseSystem::phaseModelList & phases
dimensionedScalar cos(const dimensionedScalar &ds)