46 const PtrList<volScalarField>& rAUs,
47 const PtrList<surfaceScalarField>& rAUfs
50 const dictionary& alphaControls = mesh_.solverDict(
"alpha");
59 phaseModel* referencePhasePtr =
nullptr;
70 label solvePhasesi = 0;
75 solvePhases.set(solvePhasesi++, &
movingPhases()[movingPhasei]);
89 if (solvePhases.size() == 1)
91 fluxPhases = solvePhases;
103 PtrList<surfaceScalarField> alphaPhiDbyA0s(
phases().size());
106 const PtrList<surfaceScalarField>
DByAfs(this->
DByAfs(rAUs, rAUfs));
108 forAll(solvePhases, solvePhasei)
110 phaseModel& phase = solvePhases[solvePhasei];
128 mesh_.time().timeName(),
139 bool dilatation =
false;
140 forAll(fluxPhases, fluxPhasei)
142 if (fluxPhases[fluxPhasei].
divU().valid())
151 PtrList<volScalarField::Internal> Sps(
phases().size());
152 PtrList<volScalarField::Internal> Sus(
phases().size());
154 forAll(fluxPhases, fluxPhasei)
156 phaseModel& phase = fluxPhases[fluxPhasei];
158 const label phasei = phase.index();
168 mesh_.time().timeName(),
182 min(alpha, scalar(1))
202 const phaseModel& phase2 =
phases()[phasej];
205 if (&phase2 != &phase)
207 if (phase.divU().valid())
209 dgdt +=
alpha2()*phase.divU()()();
212 if (phase2.divU().valid())
214 dgdt -=
alpha()*phase2.divU()()();
226 Sp[celli] -= dgdt[celli]/
max(1 - alpha[celli], 1
e-4);
227 Su[celli] += dgdt[celli]/
max(1 - alpha[celli], 1
e-4);
229 else if (dgdt[celli] < 0)
231 Sp[celli] += dgdt[celli]/
max(alpha[celli], 1
e-4);
237 tmp<volScalarField> trSubDeltaT;
245 List<volScalarField*> alphaPtrs(
phases().size());
253 subCycle<volScalarField, subCycleFields> alphaSubCycle
258 !(++alphaSubCycle).
end();
263 if (solvePhases.size() > 1)
274 upwind<scalar>(mesh_,
phi_).interpolate(phase)
281 PtrList<surfaceScalarField> alphaPhiCorrs(
phases().size());
283 if (solvePhases.size() > 1)
295 - upwind<scalar>(mesh_,
phi_).flux(phase)
301 forAll(fluxPhases, fluxPhasei)
303 phaseModel& phase = fluxPhases[fluxPhasei];
323 if (&phase2 == &phase)
continue;
329 cAlphas_.
find(phasePairKey(phase.name(), phase2.name()))
340 min(cAlpha()*phic,
max(phic))
341 *
nHatf(alpha, alpha2);
346 "div(phir," + alpha2.name() +
',' + alpha.name() +
')' 357 if (alphaPhiDbyA0s.set(phase.index()))
362 *alphaPhiDbyA0s[phase.index()];
365 phase.correctInflowOutflow(alphaPhiCorr);
375 min(alphaVoid.primitiveField(), phase.alphaMax())(),
381 if (solvePhases.size() > 1)
387 fixedAlphaPhiCorrs.insert
395 forAll(solvePhases, solvePhasei)
397 phaseModel& phase = solvePhases[solvePhasei];
402 phase.correctInflowOutflow(alphaPhi);
413 if (alphaSubCycle.index() == 1)
425 const PtrList<surfaceScalarField>
DByAfs 430 forAll(solvePhases, solvePhasei)
432 phaseModel& phase = solvePhases[solvePhasei];
450 phase.alphaPhiRef() += alphaEqn.flux();
455 forAll(solvePhases, solvePhasei)
457 phaseModel& phase = solvePhases[solvePhasei];
459 Info<< phase.name() <<
" fraction, min, max = " 460 << phase.weightedAverage(mesh_.V()).value()
461 <<
' ' <<
min(phase).value()
462 <<
' ' <<
max(phase).value()
466 if (!referencePhasePtr)
473 mesh_.time().timeName(),
484 Info<<
"Phase-sum volume fraction, min, max = " 486 .weightedAverage(mesh_.V()).value()
501 forAll(solvePhases, solvePhasei)
503 phaseModel& phase = solvePhases[solvePhasei];
509 forAll(solvePhases, solvePhasei)
511 phaseModel& phase = solvePhases[solvePhasei];
513 phase.alphaRhoPhiRef() =
519 if (referencePhasePtr)
521 phaseModel& referencePhase = *referencePhasePtr;
523 referencePhase.alphaPhiRef() =
phi_;
525 forAll(solvePhases, solvePhasei)
527 phaseModel& phase = solvePhases[solvePhasei];
528 referencePhase.alphaPhiRef() -= phase.alphaPhi();
531 referencePhase.correctInflowOutflow(referencePhase.alphaPhiRef());
532 referencePhase.alphaRhoPhiRef() =
534 *referencePhase.alphaPhi();
537 referenceAlpha = alphaVoid;
539 forAll(solvePhases, solvePhasei)
541 referenceAlpha -= solvePhases[solvePhasei];
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)
fvMatrix< scalar > fvScalarMatrix
#define forAll(list, i)
Loop across all elements in list.
const phaseModelPartialList & movingPhases() const
Return the models for phases that are moving.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const word & name() const
Return name.
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &)
Return a temporary field constructed from name, mesh.
tmp< surfaceScalarField > nHatf(const volScalarField &alpha1, const volScalarField &alpha2) const
Normal to interface between two phases dotted with face areas.
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
PtrList< surfaceScalarField > alphafs(phases.size())
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
IOobject(const word &name, const fileName &instance, const objectRegistry ®istry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Calculate the matrix for the laplacian of the field.
Calculate the snGrad of the given volField.
surfaceScalarField phir(fvc::flux(UdmModel.Udm()))
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
const phaseModelPartialList & stationaryPhases() const
Return the models for phases that are stationary.
virtual void solve(const PtrList< volScalarField > &rAUs, const PtrList< surfaceScalarField > &rAUfs)
Solve for the phase fractions.
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
dictionary()
Construct top-level dictionary null.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Calculate the first temporal derivative.
const phaseModelList & phases() const
Return the phase models.
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
void limitSum(UPtrList< scalarField > &phiPsiCorrs)
Calculate the face-flux of the given field.
static word groupName(Name name, const word &group)
friend class const_iterator
Declare friendship with the const_iterator.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Calculate the matrix for the first temporal derivative.
word referencePhaseName_
Name of optional reference phase which is not solved for.
static const word null
An empty word.
UPtrList< phaseModel > phaseModelPartialList
surfaceScalarField alphaPhi(phi.name()+alpha1.name(), fvc::flux(phi, alpha1, alphaScheme))
const fvMesh & mesh_
Reference to the mesh.
Calculate the field for explicit evaluation of implicit and explicit sources.
cAlphaTable cAlphas_
Interface compression coefficients.
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Calculate the divergence of the given field.
bool set(const label) const
Is element set.
tmp< volScalarField > sumAlphaMoving() const
Return the sum of the phase fractions of the moving phases.
const label nAlphaSubCycles(alphaControls.lookup< label >("nAlphaSubCycles"))
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
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 mesh motion flux and convert fluxes from absolute to relative and back.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
virtual PtrList< surfaceScalarField > DByAfs(const PtrList< volScalarField > &rAUs, const PtrList< surfaceScalarField > &rAUfs) const =0
Return the phase diffusivity.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
dimensioned< scalar > mag(const dimensioned< Type > &)
const label nAlphaCorr(alphaControls.lookup< label >("nAlphaCorr"))
MULES: Multidimensional universal limiter for explicit solution.
const doubleScalar e
Elementary charge.
surfaceScalarField phi_
Total volumetric flux.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
virtual bool implicitPhasePressure() const
Returns true if the phase pressure is treated implicitly.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Calculate the matrix for implicit and explicit sources.