Go to the documentation of this file.
35 #include "surfaceInterpolate.H"
66 limitedPhiAlphas_(phaseModels_.size()),
67 Su_(phaseModels_.size()),
68 Sp_(phaseModels_.size())
71 phases_.setSize(phaseModels_.size());
75 phases_.set(
phasei++, &pm);
78 mesh.solverDict(
"alpha").readEntry(
"cAlphas", cAlphas_);
93 mesh_.time().timeName(),
109 mesh_.time().timeName(),
124 this->alphaTransfer(
Su_,
Sp_);
142 mesh_.time().timeName(),
154 !(++alphaSubCycle).
end();
158 rhoPhiSum += (mesh_.time().deltaT()/totalDeltaT)*rhoPhi_;
170 void Foam::multiphaseSystem::solveAlphas()
177 PtrList<surfaceScalarField> phiAlphaCorrs(phases_.size());
190 if (!phicp.coupled())
199 for (phaseModel&
phase1 : phases_)
208 "phi" +
alpha1.name() +
"Corr",
213 "div(phi," +
alpha1.name() +
')'
220 for (phaseModel&
phase2 : phases_)
228 if (!cAlphas_.found(key12))
231 <<
"Phase compression factor (cAlpha) not found for : "
235 scalar cAlpha = cAlphas_.find(key12)();
255 forAll(phiAlphaCorr.boundaryField(), patchi)
258 phiAlphaCorr.boundaryFieldRef()[patchi];
260 if (!phiAlphaCorrp.coupled())
266 forAll(phiAlphaCorrp, facei)
268 if (phi1p[facei] < 0)
270 phiAlphaCorrp[facei] = alpha1p[facei]*phi1p[facei];
280 for (
const phaseModel& phase : phases_)
296 for (phaseModel& phase : phases_)
307 1.0/mesh_.time().deltaT().value(),
328 mesh_.time().timeName(),
336 for (phaseModel& phase : phases_)
350 fv::EulerDdtScheme<scalar>(mesh_).fvmDdt(
alpha1)
351 + fv::gaussConvectionScheme<scalar>
355 upwind<scalar>(mesh_,
phi)
365 phiAlpha += alpha1Eqn.flux();
379 phase.alphaPhi() = phiAlpha;
391 mesh_.time().timeName(),
401 for (phaseModel& phase : phases_)
410 Info<<
"Phase-sum volume fraction, min, max = "
411 << sumAlpha.weightedAverage(mesh_.V()).value()
412 <<
' ' <<
min(sumAlpha).value()
413 <<
' ' <<
max(sumAlpha).value()
418 for (phaseModel& phase : phases_)
424 <<
alpha.weightedAverage(mesh_.V()).value()
425 <<
" Min(alpha) = " <<
min(
alpha).value()
426 <<
" Max(alpha) = " <<
max(
alpha).value()
466 auto iter = phaseModels_.cbegin();
468 scalar maxVal =
max(iter()->diffNo()).value();
470 for (++iter; iter != phaseModels_.cend(); ++iter)
472 maxVal =
max(maxVal,
max(iter()->diffNo()).value());
475 return maxVal * mesh_.time().deltaT().value();
482 return limitedPhiAlphas_;
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const compressionFluxTable & limitedPhiAlphas() const
Access to compression fluxes for phaes.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
surfaceScalarField::Boundary & phicBf
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)
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
void calculateSuSp()
Calculate Sp and Su.
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 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)
fvsPatchField< scalar > fvsPatchScalarField
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Calculate the snGrad of the given volField.
Calculate the divergence of the given field.
const word & name() const
Return const reference to name.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const phaseModel & phase(const label i) const
Constant access phase model i.
const volScalarField & alpha1
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
surfaceScalarField phic(mixture.cAlpha() *mag(alphaPhic/mesh.magSf()))
SuSpTable Su_
Su phase source terms.
scalar maxDiffNo() const
Maximum diffusion number.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
surfaceScalarField phir(fvc::flux(UdmModel.Udm()))
fvMatrix< scalar > fvScalarMatrix
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
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.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar ddtAlphaMax() const
Access to ddtAlphaMax.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
SuSpTable & Sp()
Access Sp.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
SuSpTable Sp_
Sp phase source terms.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
bool read()
Read base transportProperties dictionary.
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 dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
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.
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.
forAllConstIters(mixture.phases(), phase)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Calculate the gradient of the given field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Calculate the face-flux of the given field.
Calculate the first temporal derivative.
Class to represent a system of phases and model interfacial transfers between them.
Area-weighted average a surfaceField creating a volField.
SuSpTable & Su()
Access Su.
defineTypeNameAndDebug(combustionModel, 0)
multiphaseSystem(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Calculate the matrix for the first temporal derivative.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimless
Dimensionless.
Perform a subCycleTime on a field.