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 144 modelType::typeName +
":K",
145 pair_.phase1().mesh(),
152 x.ref() += model_->K()*(
f1() - f2());
155 if (model1In2_.valid())
157 x.ref() += model1In2_->K()*(1 -
f1);
160 if (model2In1_.valid())
162 x.ref() += model2In1_->K()*f2;
168 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
171 correctFixedFluxBCs(
x.ref());
178 template<
class modelType>
182 tmp<surfaceScalarField>
f1, f2;
184 if (model_.valid() || model1In2_.valid())
188 blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
192 if (model_.valid() || model2In1_.valid())
196 blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
200 tmp<surfaceScalarField>
x 204 modelType::typeName +
":Kf",
205 pair_.phase1().mesh(),
212 x.ref() += model_->Kf()*(
f1() - f2());
215 if (model1In2_.valid())
217 x.ref() += model1In2_->Kf()*(1 -
f1);
220 if (model2In1_.valid())
222 x.ref() += model2In1_->Kf()*f2;
228 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
231 correctFixedFluxBCs(
x.ref());
238 template<
class modelType>
243 tmp<volScalarField>
f1, f2;
245 if (model_.valid() || model1In2_.valid())
247 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
250 if (model_.valid() || model2In1_.valid())
252 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
255 tmp<GeometricField<Type, fvPatchField, volMesh>>
x 257 new GeometricField<Type, fvPatchField, volMesh>
261 modelType::typeName +
":F",
262 pair_.phase1().mesh().time().timeName(),
263 pair_.phase1().mesh(),
268 pair_.phase1().mesh(),
269 dimensioned<Type>(
"zero", modelType::dimF,
Zero)
275 x.ref() += model_->F()*(
f1() - f2());
278 if (model1In2_.valid())
280 x.ref() += model1In2_->F()*(1 -
f1);
283 if (model2In1_.valid())
285 x.ref() -= model2In1_->F()*f2;
291 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
294 correctFixedFluxBCs(
x.ref());
301 template<
class modelType>
305 tmp<surfaceScalarField>
f1, f2;
307 if (model_.valid() || model1In2_.valid())
311 blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
315 if (model_.valid() || model2In1_.valid())
319 blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
323 tmp<surfaceScalarField>
x 327 modelType::typeName +
":Ff",
328 pair_.phase1().mesh(),
335 x.ref() += model_->Ff()*(
f1() - f2());
338 if (model1In2_.valid())
340 x.ref() += model1In2_->Ff()*(1 -
f1);
343 if (model2In1_.valid())
345 x.ref() -= model2In1_->Ff()*f2;
351 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
354 correctFixedFluxBCs(
x.ref());
361 template<
class modelType>
365 tmp<volScalarField>
f1, f2;
367 if (model_.valid() || model1In2_.valid())
369 f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
372 if (model_.valid() || model2In1_.valid())
374 f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
377 tmp<volScalarField>
x 381 modelType::typeName +
":D",
382 pair_.phase1().mesh(),
389 x.ref() += model_->D()*(
f1() - f2());
392 if (model1In2_.valid())
394 x.ref() += model1In2_->D()*(1 -
f1);
397 if (model2In1_.valid())
399 x.ref() += model2In1_->D()*f2;
405 && (model_.valid() || model1In2_.valid() || model2In1_.valid())
408 correctFixedFluxBCs(
x.ref());
415 template<
class modelType>
418 const class phaseModel& phase
422 &phase == &(pair_.phase1())
424 : model2In1_.valid();
428 template<
class modelType>
431 const class phaseModel& phase
434 return &phase == &(pair_.phase1()) ? model1In2_ : model2In1_;
#define forAll(list, i)
Loop across all elements in list.
bool hasModel(const phaseModel &phase) const
Return true if a model is specified for the supplied phase.
const modelType & phaseModel(const phaseModel &phase) const
Return the model for the supplied phase.
~BlendedInterfacialModel()
Destructor.
tmp< volScalarField > K() const
Return the blended force coefficient.
tmp< surfaceScalarField > Ff() const
Return the face blended force.
tmp< surfaceScalarField > Kf() const
Return the face blended force coefficient.
tmp< GeometricField< Type, fvPatchField, volMesh > > F() const
Return the blended force.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
BlendedInterfacialModel(const phaseModel &phase1, const phaseModel &phase2, const blendingMethod &blending, autoPtr< ModelType > model, autoPtr< ModelType > model1In2, autoPtr< ModelType > model2In1, const bool correctFixedFluxBCs=true)
Construct from two phases, blending method and three models.
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.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< volScalarField > D() const
Return the blended diffusivity.
A class for managing temporary objects.
const dimensionSet dimArea(sqr(dimLength))