Go to the documentation of this file.
29 #include "alphaContactAngleFvPatchScalarField.H"
54 const Foam::scalar Foam::multiphaseSystem::convertToRad =
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 = convertToRad*tp().theta0(matched);
430 scalar uTheta = tp().uTheta();
435 scalar thetaA = convertToRad*tp().thetaA(matched);
436 scalar thetaR = convertToRad*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());
518 mesh_.time().timeName(),
527 cAlphas_(
lookup(
"interfaceCompression")),
538 mesh_.setFluxRequired(alphai.name());
566 tSurfaceTension.
ref().setOriented();
576 cAlphaTable::const_iterator cAlpha(cAlphas_.find(key12));
578 if (cAlpha != cAlphas_.end())
580 tSurfaceTension.
ref() +=
590 tSurfaceTension->setOriented();
592 return tSurfaceTension;
675 !(++alphaSubCycle).
end();
689 if (
phase.stationary())
continue;
702 if (
phase.stationary())
continue;
704 phase.alphaRhoPhiRef() =
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
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.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual ~multiphaseSystem()
Destructor.
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
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)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
virtual void solve()
Solve for the 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.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
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: [].
const Time & time() const
Return time.
Calculate the snGrad of the given volField.
static word timeName(const scalar t, const int precision=precision_)
multiphaseSystem(const fvMesh &mesh)
Construct from fvMesh.
Calculate the divergence of the given field.
const word & name() const
Return const reference to name.
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
PtrList< surfaceScalarField > alphafs(phases.size())
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Dimension set for the base types.
tmp< surfaceScalarField > surfaceTension(const phaseModel &phase) const
Return the surface tension force.
void clip(const dimensioned< MinMax< Type >> &range)
Clip the field to be bounded within the specified range.
A class for managing sub-cycling times.
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 (uses stdout - output is on the master only)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
const cellModel * lookup(const word &modelName)
Deprecated(2017-11) equivalent to cellModel::ptr static method.
Calculate the matrix for the laplacian of the field.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
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.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Mesh data needed to do the Finite Volume discretisation.
const T * set(const label i) const
Return const pointer to element (if set) or nullptr.
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.
Calculate the matrix for implicit and explicit sources.
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual tmp< surfaceScalarField > phi() const =0
Return the volumetric flux.
const word & name() const
Return the name of this phase.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr scalar pi(M_PI)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Calculate the face-flux of the given field.
multiphaseSystem::phaseModelList & phases
const dimensionedScalar e
Elementary charge.
virtual tmp< volVectorField > U() const =0
Return the velocity.
Calculate the first temporal derivative.
Class to represent a system of phases and model interfacial transfers between them.
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)
bool valid() const noexcept
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys and label hasher.
MULES: Multidimensional universal limiter for explicit solution.
defineTypeNameAndDebug(combustionModel, 0)
const UPtrList< phaseModel > & phases() const
Return phases.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
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 dimensionedScalar & rho() const
Return const-access to phase1 density.
dimensionedScalar cos(const dimensionedScalar &ds)