26 #include "BlendedInterfacialModel.H" 28 #include "surfaceInterpolate.H" 32 template<
class modelType>
33 template<
class GeometricField>
39 typename GeometricField::Boundary& fieldBf =
40 field.boundaryFieldRef();
46 isA<fixedValueFvsPatchScalarField>
48 pair_.phase1().phi().boundaryField()[
patchi]
60 template<
class modelType>
63 const phasePair::dictTable& modelTable,
64 const blendingMethod& blending,
65 const phasePair& pair,
66 const orderedPhasePair& pair1In2,
67 const orderedPhasePair& pair2In1,
68 const bool correctFixedFluxBCs
75 correctFixedFluxBCs_(correctFixedFluxBCs)
77 if (modelTable.found(pair_))
89 if (modelTable.found(pair1In2_))
95 modelTable[pair1In2_],
101 if (modelTable.found(pair2In1_))
107 modelTable[pair2In1_],
117 template<
class modelType>
124 template<
class modelType>
128 tmp<volScalarField>
f1, f2;
130 if (model_.valid() || model1In2_.valid())
132 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
135 if (model_.valid() || model2In1_.valid())
137 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
140 tmp<volScalarField>
x 146 modelType::typeName +
":K",
147 pair_.phase1().mesh().time().timeName(),
148 pair_.phase1().mesh(),
153 pair_.phase1().mesh(),
160 x.ref() += model_->K()*(
f1() - f2());
163 if (model1In2_.valid())
165 x.ref() += model1In2_->K()*(1 -
f1);
168 if (model2In1_.valid())
170 x.ref() += model2In1_->K()*f2;
176 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
179 correctFixedFluxBCs(
x.ref());
186 template<
class modelType>
190 tmp<surfaceScalarField>
f1, f2;
192 if (model_.valid() || model1In2_.valid())
196 blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
200 if (model_.valid() || model2In1_.valid())
204 blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
208 tmp<surfaceScalarField>
x 214 modelType::typeName +
":Kf",
215 pair_.phase1().mesh().time().timeName(),
216 pair_.phase1().mesh(),
221 pair_.phase1().mesh(),
228 x.ref() += model_->Kf()*(
f1() - f2());
231 if (model1In2_.valid())
233 x.ref() += model1In2_->Kf()*(1 -
f1);
236 if (model2In1_.valid())
238 x.ref() += model2In1_->Kf()*f2;
244 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
247 correctFixedFluxBCs(
x.ref());
254 template<
class modelType>
259 tmp<volScalarField>
f1, f2;
261 if (model_.valid() || model1In2_.valid())
263 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
266 if (model_.valid() || model2In1_.valid())
268 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
271 tmp<GeometricField<Type, fvPatchField, volMesh>>
x 273 new GeometricField<Type, fvPatchField, volMesh>
277 modelType::typeName +
":F",
278 pair_.phase1().mesh().time().timeName(),
279 pair_.phase1().mesh(),
284 pair_.phase1().mesh(),
285 dimensioned<Type>(
"zero", modelType::dimF,
Zero)
291 x.ref() += model_->F()*(
f1() - f2());
294 if (model1In2_.valid())
296 x.ref() += model1In2_->F()*(1 -
f1);
299 if (model2In1_.valid())
301 x.ref() -= model2In1_->F()*f2;
307 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
310 correctFixedFluxBCs(
x.ref());
317 template<
class modelType>
321 tmp<surfaceScalarField>
f1, f2;
323 if (model_.valid() || model1In2_.valid())
327 blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
331 if (model_.valid() || model2In1_.valid())
335 blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
339 tmp<surfaceScalarField>
x 345 modelType::typeName +
":Ff",
346 pair_.phase1().mesh().time().timeName(),
347 pair_.phase1().mesh(),
352 pair_.phase1().mesh(),
359 x.ref() += model_->Ff()*(
f1() - f2());
362 if (model1In2_.valid())
364 x.ref() += model1In2_->Ff()*(1 -
f1);
367 if (model2In1_.valid())
369 x.ref() -= model2In1_->Ff()*f2;
375 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
378 correctFixedFluxBCs(
x.ref());
385 template<
class modelType>
389 tmp<volScalarField>
f1, f2;
391 if (model_.valid() || model1In2_.valid())
393 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
396 if (model_.valid() || model2In1_.valid())
398 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
401 tmp<volScalarField>
x 407 modelType::typeName +
":D",
408 pair_.phase1().mesh().time().timeName(),
409 pair_.phase1().mesh(),
414 pair_.phase1().mesh(),
421 x.ref() += model_->D()*(
f1() - f2());
424 if (model1In2_.valid())
426 x.ref() += model1In2_->D()*(1 -
f1);
429 if (model2In1_.valid())
431 x.ref() += model2In1_->D()*f2;
437 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
440 correctFixedFluxBCs(
x.ref());
447 template<
class modelType>
450 const class phaseModel& phase
454 &phase == &(pair_.phase1())
456 : model2In1_.valid();
460 template<
class modelType>
463 const class phaseModel& phase
466 return &phase == &(pair_.phase1()) ? model1In2_ : model2In1_;
tmp< volScalarField > D() const
Return the blended diffusivity.
#define forAll(list, i)
Loop across all elements in list.
~BlendedInterfacialModel()
Destructor.
const modelType & phaseModel(const phaseModel &phase) const
Return the model for the supplied phase.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
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< Type, fvPatchField, volMesh > > F() const
Return the blended force.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
bool hasModel(const phaseModel &phase) const
Return true if a model is specified for the supplied phase.
tmp< surfaceScalarField > Ff() const
Return the face blended force.
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
tmp< volScalarField > K() const
Return the blended force coefficient.
const dimensionSet dimArea(sqr(dimLength))
tmp< surfaceScalarField > Kf() const
Return the face blended force coefficient.