30#include "alphaContactAngleFvPatchScalarField.H"
60void Foam::multiphaseSystem::calcAlphas()
67 alphas_ += level*
phases()[i];
73void 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_);
383void Foam::multiphaseSystem::correctContactAngle
393 const fvBoundaryMesh&
boundary = mesh_.boundary();
399 isA<reactingMultiphaseEuler::alphaContactAngleFvPatchScalarField>
408 const reactingMultiphaseEuler::
409 alphaContactAngleFvPatchScalarField
416 mesh_.Sf().boundaryField()[patchi]
417 /mesh_.magSf().boundaryField()[patchi]
420 reactingMultiphaseEuler::alphaContactAngleFvPatchScalarField::
421 thetaPropsTable::const_iterator tp =
425 if (tp == acap.thetaProps().end())
428 <<
"Cannot find interface "
430 <<
"\n in table of theta properties for patch "
431 << acap.patch().name()
435 bool matched = (tp.key().first() ==
phase1.
name());
437 scalar theta0 =
degToRad(tp().theta0(matched));
440 scalar uTheta = tp().uTheta();
445 const scalar thetaA =
degToRad(tp().thetaA(matched));
446 const scalar thetaR =
degToRad(tp().thetaR(matched));
454 Uwall -= (AfHatPatch & Uwall)*AfHatPatch;
459 nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch
463 nWall /= (
mag(nWall) + SMALL);
469 theta += (thetaA - thetaR)*
tanh(uwall/uTheta);
483 b2[facei] =
cos(
acos(a12[facei]) - theta[facei]);
491 nHatPatch = a*AfHatPatch +
b*nHatPatch;
493 nHatPatch /= (
mag(nHatPatch) + deltaN_.value());
505 tmp<surfaceVectorField> tnHatfv = nHatfv(
phase1,
phase2);
510 return -
fvc::div(tnHatfv & mesh_.Sf());
537 cAlphas_(
lookup(
"interfaceCompression")),
570 tSurfaceTension.ref().setOriented();
580 cAlphaTable::const_iterator cAlpha(cAlphas_.find(key12));
582 if (cAlpha != cAlphas_.end())
584 tSurfaceTension.ref() +=
594 tSurfaceTension->setOriented();
596 return tSurfaceTension;
603 tmp<volScalarField> tnearInt
637 tmp<volScalarField> trSubDeltaT;
645 List<volScalarField*> alphaPtrs(
phases().size());
646 PtrList<surfaceScalarField> alphaPhiSums(
phases().size());
674 subCycleTime alphaSubCycle
679 !(++alphaSubCycle).
end();
693 if (phase.stationary())
continue;
706 if (phase.stationary())
continue;
708 phase.alphaRhoPhiRef() =
711 phase.clip(SMALL, 1 - SMALL);
CGAL::Exact_predicates_exact_constructions_kernel K
surfaceScalarField phic(mixture.cAlpha() *mag(alphaPhic/mesh.magSf()))
MULES: Multidimensional universal limiter for explicit solution.
const volScalarField & alpha1
const volScalarField & alpha2
DimensionedField< scalar, volMesh > Internal
The internal field type from which this GeometricField is derived.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
GeometricBoundaryField< vector, fvsPatchField, surfaceMesh > Boundary
Type of boundary fields.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
virtual bool enabled() const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Mesh data needed to do the Finite Volume discretisation.
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
const surfaceScalarField & nHatf() const
areaScalarField & surfaceTension()
Return surface tension field.
Incompressible multi-phase mixture with built in solution for the phase fractions with interface comp...
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
const PtrDictionary< phaseModel > & phases() const
Return the phases.
void solve()
Solve for the mixture phase-fractions.
const Time & time() const noexcept
Return time registry.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
virtual tmp< volScalarField > divU() const =0
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
const volVectorField & U() const
const surfaceScalarField & phi() const
const word & name() const
An ordered or unorder pair of phase names. Typically specified as follows.
Class to represent a system of phases and model interfacial transfers between them.
Lookup type of boundary radiation properties.
void setFluxRequired(const word &name) const
Get flux-required for given name, or default.
A class for managing temporary objects.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
surfaceScalarField phir(fvc::flux(UdmModel.Udm()))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the face-flux of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Calculate the snGrad of the given volField.
Calculate the field for explicit evaluation of implicit and explicit sources.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the laplacian of the field.
Calculate the finiteVolume matrix for implicit and explicit sources.
PtrList< surfaceScalarField > alphafs(phases.size())
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
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)
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
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.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
dimensionedScalar det(const dimensionedSphericalTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
To & refCast(From &r)
Reference type cast template function.
const dimensionSet dimless
Dimensionless.
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar tanh(const dimensionedScalar &ds)
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
#define defineRunTimeSelectionTable(baseType, argNames)
Define run-time selection table.
multiphaseSystem::phaseModelList & phases
#define forAll(list, i)
Loop across all elements in list.
surfaceScalarField alphaPhi(phi.name()+alpha1.name(), fvc::flux(phi, alpha1, alphaScheme))
Unit conversion functions.