31 #include "surfaceInterpolate.H" 35 template<
class modelType>
36 template<
class GeometricField>
42 typename GeometricField::Boundary& fieldBf =
43 field.boundaryFieldRef();
45 forAll(pair_.phase1().phi().boundaryField(), patchi)
49 isA<fixedValueFvsPatchScalarField>
51 pair_.phase1().phi().boundaryField()[patchi]
55 fieldBf[patchi] =
Zero;
63 template<
class modelType>
66 const phasePair::dictTable& modelTable,
67 const blendingMethod& blending,
68 const phasePair& pair,
69 const orderedPhasePair& pair1In2,
70 const orderedPhasePair& pair2In1,
71 const bool correctFixedFluxBCs
78 correctFixedFluxBCs_(correctFixedFluxBCs)
80 if (modelTable.found(pair_))
92 if (modelTable.found(pair1In2_))
98 modelTable[pair1In2_],
104 if (modelTable.found(pair2In1_))
110 modelTable[pair2In1_],
120 template<
class modelType>
127 template<
class modelType>
131 tmp<volScalarField> f1, f2;
133 if (model_ || model1In2_)
135 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
138 if (model_ || model2In1_)
140 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
143 tmp<volScalarField>
x 149 modelType::typeName +
":K",
150 pair_.phase1().mesh().time().timeName(),
151 pair_.phase1().mesh(),
154 IOobject::NO_REGISTER
156 pair_.phase1().mesh(),
163 x.ref() += model_->K()*(f1() - f2());
168 x.ref() += model1In2_->K()*(1 - f1);
173 x.ref() += model2In1_->K()*f2;
179 && (model_ || model1In2_ || model2In1_)
182 correctFixedFluxBCs(
x.ref());
189 template<
class modelType>
193 tmp<surfaceScalarField> f1, f2;
195 if (model_ || model1In2_)
199 blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
203 if (model_ || model2In1_)
207 blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
211 tmp<surfaceScalarField>
x 217 modelType::typeName +
":Kf",
218 pair_.phase1().mesh().time().timeName(),
219 pair_.phase1().mesh(),
222 IOobject::NO_REGISTER
224 pair_.phase1().mesh(),
231 x.ref() += model_->Kf()*(f1() - f2());
236 x.ref() += model1In2_->Kf()*(1 - f1);
241 x.ref() += model2In1_->Kf()*f2;
247 && (model_ || model1In2_ || model2In1_)
250 correctFixedFluxBCs(
x.ref());
257 template<
class modelType>
262 tmp<volScalarField> f1, f2;
264 if (model_ || model1In2_)
266 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
269 if (model_ || model2In1_)
271 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
274 auto x = tmp<GeometricField<Type, fvPatchField, volMesh>>
::New 278 modelType::typeName +
":F",
279 pair_.phase1().mesh().time().timeName(),
280 pair_.phase1().mesh(),
283 IOobject::NO_REGISTER
285 pair_.phase1().mesh(),
286 dimensioned<Type>(modelType::dimF,
Zero)
291 x.ref() += model_->F()*(f1() - f2());
296 x.ref() += model1In2_->F()*(1 - f1);
301 x.ref() -= model2In1_->F()*f2;
307 && (model_ || model1In2_ || model2In1_)
310 correctFixedFluxBCs(
x.ref());
317 template<
class modelType>
321 tmp<surfaceScalarField> f1, f2;
323 if (model_ || model1In2_)
327 blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
331 if (model_ || model2In1_)
335 blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
343 modelType::typeName +
":Ff",
344 pair_.phase1().mesh().time().timeName(),
345 pair_.phase1().mesh(),
348 IOobject::NO_REGISTER
350 pair_.phase1().mesh(),
354 x.ref().setOriented();
358 x.ref() += model_->Ff()*(f1() - f2());
363 x.ref() += model1In2_->Ff()*(1 - f1);
368 x.ref() -= model2In1_->Ff()*f2;
374 && (model_ || model1In2_ || model2In1_)
377 correctFixedFluxBCs(
x.ref());
384 template<
class modelType>
388 tmp<volScalarField> f1, f2;
390 if (model_ || model1In2_)
392 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
395 if (model_ || model2In1_)
397 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
400 tmp<volScalarField>
x 406 modelType::typeName +
":D",
407 pair_.phase1().mesh().time().timeName(),
408 pair_.phase1().mesh(),
411 IOobject::NO_REGISTER
413 pair_.phase1().mesh(),
420 x.ref() += model_->D()*(f1() - f2());
425 x.ref() += model1In2_->D()*(1 - f1);
430 x.ref() += model2In1_->D()*f2;
436 && (model_ || model1In2_ || model2In1_)
439 correctFixedFluxBCs(
x.ref());
446 template<
class modelType>
449 const class phaseModel& phase
454 &phase == &(pair_.phase1())
461 template<
class modelType>
467 return &
phase == &(pair_.phase1()) ? *model1In2_ : *model2In1_;
bool hasModel(const phaseModel &phase) const
Return true if a model is specified for the supplied phase.
static tmp< GeometricField< Type, faePatchField, edgeMesh > > interpolate(const GeometricField< Type, faPatchField, areaMesh > &tvf, const edgeScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
~BlendedInterfacialModel()
Destructor.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
tmp< volScalarField > K() const
Return the blended force coefficient.
tmp< surfaceScalarField > Ff() const
Return the face blended force.
tmp< GeometricField< Type, fvPatchField, volMesh > > F() const
Return the blended force.
#define forAll(list, i)
Loop across all elements in list.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const modelType & phaseModel(const phaseModel &phase) const
Return the model for the supplied phase.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< surfaceScalarField > Kf() const
Return the face blended force coefficient.
tmp< volScalarField > D() const
Return the blended diffusivity.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const dimensionSet dimArea(sqr(dimLength))
static constexpr const zero Zero
Global zero (0)