30#include "alphaContactAngleFvPatchScalarField.H"
35#include "surfaceInterpolate.H"
46void Foam::multiphaseSystem::calcAlphas()
51 for (
const phaseModel& phase : phases_)
53 alphas_ += level * phase;
59void Foam::multiphaseSystem::solveAlphas()
61 PtrList<surfaceScalarField> alphaPhiCorrs(phases_.size());
64 for (phaseModel& phase : phases_)
85 for (phaseModel&
phase2 : phases_)
89 if (&
phase2 == &phase)
continue;
93 const auto cAlpha = cAlphas_.cfind(interfacePair(phase,
phase2));
105 const word phirScheme
118 phase.correctInflowOutflow(alphaPhiCorr);
122 1.0/mesh_.time().deltaT().value(),
144 mesh_.time().timeName(),
153 for (phaseModel& phase : phases_)
156 alphaPhi += upwind<scalar>(mesh_, phi_).flux(phase);
157 phase.correctInflowOutflow(
alphaPhi);
168 Info<< phase.name() <<
" volume fraction, min, max = "
169 << phase.weightedAverage(mesh_.V()).value()
170 <<
' ' <<
min(phase).value()
171 <<
' ' <<
max(phase).value()
179 Info<<
"Phase-sum volume fraction, min, max = "
180 << sumAlpha.weightedAverage(mesh_.V()).value()
181 <<
' ' <<
min(sumAlpha).value()
182 <<
' ' <<
max(sumAlpha).value()
187 for (phaseModel& phase : phases_)
219 return gradAlphaf/(
mag(gradAlphaf) + deltaN_);
240void Foam::multiphaseSystem::correctContactAngle
250 const fvBoundaryMesh&
boundary = mesh_.boundary();
256 isA<multiphaseEuler::alphaContactAngleFvPatchScalarField>
265 const multiphaseEuler::alphaContactAngleFvPatchScalarField
275 mesh_.Sf().boundaryField()[patchi]
276 /mesh_.magSf().boundaryField()[patchi]
285 <<
"Cannot find interface " << interfacePair(
phase1,
phase2)
286 <<
"\n in table of theta properties for patch "
287 << acap.patch().name()
291 bool matched = (tp.key().first() ==
phase1.
name());
293 const scalar theta0 =
degToRad(tp().theta0(matched));
296 scalar uTheta = tp().uTheta();
301 const scalar thetaA =
degToRad(tp().thetaA(matched));
302 const scalar thetaR =
degToRad(tp().thetaR(matched));
310 Uwall -= (AfHatPatch & Uwall)*AfHatPatch;
315 nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch
319 nWall /= (
mag(nWall) + SMALL);
325 theta += (thetaA - thetaR)*
tanh(uwall/uTheta);
339 b2[facei] =
cos(
acos(a12[facei]) - theta[facei]);
347 nHatPatch = a*AfHatPatch +
b*nHatPatch;
349 nHatPatch /= (
mag(nHatPatch) + deltaN_.value());
361 tmp<surfaceVectorField> tnHatfv = nHatfv(
phase1,
phase2);
366 return -
fvc::div(tnHatfv & mesh_.Sf());
382 "transportProperties",
409 sigmas_(
lookup(
"sigmas")),
410 dimSigma_(1, 0, -2, 0, 0),
411 cAlphas_(
lookup(
"interfaceCompression")),
412 Cvms_(
lookup(
"virtualMass")),
432 *phases_.lookup(iter.key().first()),
433 *phases_.lookup(iter.key().second())
449 if (sigmas_.
found(key) && !cAlphas_.
found(key))
452 <<
"Compression coefficient not specified for phase pair ("
454 <<
") for which a surface tension coefficient is specified"
466 auto iter = phases_.cbegin();
471 for (++iter; iter != phases_.cend(); ++iter)
473 rho += iter()*iter().rho();
483 auto iter = phases_.cbegin();
488 for (++iter; iter != phases_.cend(); ++iter)
490 rho += iter().boundaryField()[patchi]*iter().rho().value();
499 auto iter = phases_.cbegin();
504 for (++iter; iter != phases_.cend(); ++iter)
506 mu += iter()*(iter().rho()*iter().nu());
516 auto iter = phases_.cbegin();
519 iter().boundaryField()[patchi]
520 *(iter().rho().value()*iter().nu().value());
523 for (++iter; iter != phases_.cend(); ++iter)
527 *(iter().rho().value()*iter().nu().value());
530 return tmu/
rho(patchi);
546 mesh_.time().timeName(),
594 mesh_.time().timeName(),
632 tSvm.
ref().boundaryFieldRef();
639 isA<fixedValueFvsPatchScalarField>
645 SvmBf[patchi] =
Zero;
688 isA<fixedValueFvsPatchScalarField>
698 dragCoeffsPtr().
set(iter.key(), Kptr);
701 return dragCoeffsPtr;
718 mesh_.time().timeName(),
736 dmIter.good() && dcIter.good();
746 tdragCoeff.
ref() += *dcIter();
766 mesh_.time().timeName(),
778 tSurfaceTension.
ref().setOriented();
791 tSurfaceTension.
ref() +=
801 return tSurfaceTension;
815 mesh_.time().timeName(),
889 !(++alphaSubCycle).end();
942 lookup(
"sigmas") >> sigmas_;
943 lookup(
"interfaceCompression") >> cAlphas_;
944 lookup(
"virtualMass") >> Cvms_;
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
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
label timeIndex() const
Return the time index of the field.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
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.
A HashTable of pointers to objects of type <T>, with deallocation management of the pointers.
bool set(const Key &key, T *ptr)
Assign a new entry, overwriting existing entries.
typename parent_type::const_iterator const_iterator
A HashTable similar to std::unordered_map.
bool found(const Key &key) const
Return true if hashed entry is found in table.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
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.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
const T * set(const label i) const
label timeIndex() const noexcept
Return current time index.
dimensionedScalar deltaT() const
Return time step.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
void set(T *p) noexcept
Deprecated(2018-02) Identical to reset().
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...
const surfaceScalarField & nHatf() const
areaScalarField & surfaceTension()
Return surface tension field.
const phaseModel & phase2() const
const dimensionedScalar & residualSlip() const
virtual tmp< volScalarField > K(const volScalarField &Ur) const =0
The drag function K used in the momentum eq.
const phaseModel & phase1() const
const dimensionedScalar & residualPhaseFraction() const
Name pair for the interface.
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.
autoPtr< dragCoeffFields > dragCoeffs() const
Return the drag coefficients for all of the interfaces.
tmp< volScalarField > nu() const
Return the mixture laminar viscosity.
tmp< volScalarField > rho() const
Return the mixture density.
tmp< volScalarField > dragCoeff(const phaseModel &phase, const dragCoeffFields &dragCoeffs) const
Return the sum of the drag coefficients for the given phase.
tmp< volVectorField > Svm(const phaseModel &phase) const
Return the virtual-mass source for the given phase.
void solve()
Solve for the mixture phase-fractions.
bool read()
Read base transportProperties dictionary.
const Time & time() const noexcept
Return time registry.
constant condensation/saturation model.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const volVectorField & U() const
const surfaceScalarField & phi() const
const dimensionedScalar & rho() const
const volVectorField & DDtU() const
const word & name() const
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
void correct()
Correct the phase properties.
const dimensionedScalar & rho() const
Return const-access to phase1 density.
bool read(const dictionary &phaseDict)
Read base transportProperties dictionary.
Lookup type of boundary radiation properties.
virtual bool write(const bool valid=true) const
Write using setting from DB.
virtual bool read()
Read object.
A class for managing sub-cycling times.
A class for managing temporary objects.
virtual tmp< volScalarField > Cvm() const
Virtual mass coefficient.
const volScalarField & mu
surfaceScalarField phir(fvc::flux(UdmModel.Udm()))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Area-weighted average a surfaceField creating a volField.
Calculate the divergence of the given field.
Calculate the face-flux of the given field.
Calculate the gradient of the given field.
Calculate the snGrad of the given volField.
#define WarningInFunction
Report a warning using Foam::Warning.
autoPtr< multiphaseSystem::dragCoeffFields > dragCoeffs(fluid.dragCoeffs())
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.
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)
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.
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)
dimensionedScalar cbrt(const dimensionedScalar &ds)
const dimensionSet dimDensity
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
tmp< volScalarField > trho
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
surfaceScalarField alphaPhi(phi.name()+alpha1.name(), fvc::flux(phi, alpha1, alphaScheme))
Unit conversion functions.