43 const dictionary& dict,
44 const Pair<scalar>& bounds,
48 if (dict.found(name) || allowNone)
50 const token t(dict.lookup(name));
52 if (allowNone && t.isWord() && t.wordToken() ==
"none")
61 if (!std::isnan(bounds[i]))
63 const label s = i == 0 ? -1 : +1;
65 if (s*t.number() > s*bounds[i])
68 <<
"Blending parameter " << name <<
" is " 69 << (i == 0 ?
"less" :
"greater") <<
" than " 79 <<
"wrong token type - expected Scalar or the word 'none', found " 90 const dictionary& dict,
91 const phaseInterface& interface,
92 const Pair<scalar>& bounds,
101 readParameter(name1, dict, bounds, allowNone),
102 readParameter(name2, dict, bounds, allowNone)
109 return !std::isnan(parameter);
117 const UPtrList<const volScalarField>& alphas,
125 alphas.first().mesh(),
133 const UPtrList<const volScalarField>& alphas,
138 tmp<volScalarField> talpha = constant(alphas, 0);
142 if (0b01 << iter.index() &
set)
146 ?
max(iter().residualAlpha(), alphas[iter().index()])()
147 : alphas[iter().index()];
157 const UPtrList<const volScalarField>& alphas,
159 const Pair<scalar>& parameters
162 tmp<volScalarField> talphaParameter = constant(alphas, 0);
166 if (0b01 << iter.index() &
set)
168 talphaParameter.ref() +=
169 max(iter().residualAlpha(), alphas[iter().index()])
170 *parameters[iter.index()];
174 return talphaParameter/
alpha(alphas,
set,
true);
180 const UPtrList<const volScalarField>& alphas,
181 const label phaseSet,
182 const label systemSet
187 ?
alpha(alphas, phaseSet,
false)
188 : alpha(alphas, phaseSet, false)/alpha(alphas, systemSet, true);
194 const UPtrList<const volScalarField>& alphas,
195 const label phaseSet,
196 const label systemSet
199 label canBeContinuousPhaseSet = 0b00;
200 label canBeContinuousSystemSet = 0b00;
203 if (canBeContinuous(iter.index()))
205 canBeContinuousPhaseSet += 0b01 << iter.index() & phaseSet;
206 canBeContinuousSystemSet += 0b01 << iter.index() & systemSet;
210 if (canBeContinuousPhaseSet == 0)
212 return constant(alphas, 0);
215 if (canBeContinuousPhaseSet == canBeContinuousSystemSet)
217 return constant(alphas, 1);
224 canBeContinuousPhaseSet,
225 canBeContinuousSystemSet
234 const dictionary& dict,
235 const phaseInterface& interface
238 interface_(interface)
252 const UPtrList<const volScalarField>& alphas
255 return f(alphas, 0b10, 0b11);
261 const UPtrList<const volScalarField>& alphas
264 return f(alphas, 0b01, 0b11);
270 const UPtrList<const volScalarField>& alphas
273 return 1 -
f(alphas, 0b11, 0b00);
#define forAll(list, i)
Loop across all elements in list.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static Pair< scalar > readParameters(const word &name, const dictionary &dict, const phaseInterface &interface, const Pair< scalar > &bounds, const bool allowNone)
Read a parameter for each phase in the interface.
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< volScalarField > alpha(const UPtrList< const volScalarField > &alphas, const label set, const bool protect) const
Get the volume fraction of the given set.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
tmp< volScalarField > parameter(const UPtrList< const volScalarField > &alphas, const label set, const Pair< scalar > ¶meters) const
Get a blending parameter averaged for the given set.
tmp< volScalarField > f2DispersedIn1(const UPtrList< const volScalarField > &alphas) const
Return the coefficient for models in which phase 2 is dispersed in.
const dimensionSet dimless
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
tmp< volScalarField > constant(const UPtrList< const volScalarField > &alphas, const scalar k) const
Return a constant field with the given value.
An ordered pair of two objects of type <T> with first() and second() elements.
static word groupName(Name name, const word &group)
static scalar readParameter(const word &name, const dictionary &dict, const Pair< scalar > &bounds, const bool allowNone)
Read a parameter and check it lies within specified bounds.
static bool isParameter(const scalar parameter)
Test if the given scalar is a valid parameter.
virtual ~blendingMethod()
Destructor.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
defineTypeNameAndDebug(combustionModel, 0)
word name(const complex &)
Return a string representation of a complex.
virtual tmp< volScalarField > f(const UPtrList< const volScalarField > &alphas, const label phaseSet, const label systemSet) const
Evaluate the blending function. Filters out phases that cannot.
tmp< volScalarField > f1DispersedIn2(const UPtrList< const volScalarField > &alphas) const
Return the coefficient for models in which phase 1 is dispersed in.
blendingMethod(const dictionary &dict, const phaseInterface &interface)
Construct from a dictionary and an interface.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< volScalarField > x(const UPtrList< const volScalarField > &alphas, const label phaseSet, const label systemSet) const
Return the coordinate of the blending function.
A class for managing temporary objects.
tmp< volScalarField > fDisplaced(const UPtrList< const volScalarField > &alphas) const
Return the coefficient for when the interface is displaced by a.