All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
blendingMethod.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2014-2022 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "blendingMethod.H"
27 #include "phaseSystem.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(blendingMethod, 0);
34  defineRunTimeSelectionTable(blendingMethod, dictionary);
35 }
36 
37 
38 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
39 
41 (
42  const word& name,
43  const dictionary& dict,
44  const Pair<scalar>& bounds,
45  const bool allowNone
46 )
47 {
48  if (dict.found(name) || allowNone)
49  {
50  const token t(dict.lookup(name));
51 
52  if (allowNone && t.isWord() && t.wordToken() == "none")
53  {
54  return NaN;
55  }
56 
57  if (t.isNumber())
58  {
59  forAll(bounds, i)
60  {
61  if (!std::isnan(bounds[i]))
62  {
63  const label s = i == 0 ? -1 : +1;
64 
65  if (s*t.number() > s*bounds[i])
66  {
68  << "Blending parameter " << name << " is "
69  << (i == 0 ? "less" : "greater") << " than "
70  << bounds[i] << exit(FatalError);
71  }
72  }
73  }
74 
75  return t.number();
76  }
77 
79  << "wrong token type - expected Scalar or the word 'none', found "
80  << t.info() << exit(FatalIOError);
81  }
82 
83  return NaN;
84 }
85 
86 
88 (
89  const word& name,
90  const dictionary& dict,
91  const phaseInterface& interface,
92  const Pair<scalar>& bounds,
93  const bool allowNone
94 )
95 {
96  const word name1 = IOobject::groupName(name, interface.phase1().name());
97  const word name2 = IOobject::groupName(name, interface.phase2().name());
98  return
99  Pair<scalar>
100  (
101  readParameter(name1, dict, bounds, allowNone),
102  readParameter(name2, dict, bounds, allowNone)
103  );
104 }
105 
106 
107 bool Foam::blendingMethod::isParameter(const scalar parameter)
108 {
109  return !std::isnan(parameter);
110 }
111 
112 
113 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
114 
116 (
117  const UPtrList<const volScalarField>& alphas,
118  const scalar k
119 ) const
120 {
121  return
123  (
124  name(k),
125  alphas.first().mesh(),
127  );
128 }
129 
130 
132 (
133  const UPtrList<const volScalarField>& alphas,
134  const label set,
135  const bool protect
136 ) const
137 {
138  tmp<volScalarField> talpha = constant(alphas, 0);
139 
140  forAllConstIter(phaseInterface, interface_, iter)
141  {
142  if (0b01 << iter.index() & set)
143  {
144  talpha.ref() +=
145  protect
146  ? max(iter().residualAlpha(), alphas[iter().index()])()
147  : alphas[iter().index()];
148  }
149  }
150 
151  return talpha;
152 }
153 
154 
156 (
157  const UPtrList<const volScalarField>& alphas,
158  const label set,
159  const Pair<scalar>& parameters
160 ) const
161 {
162  tmp<volScalarField> talphaParameter = constant(alphas, 0);
163 
164  forAllConstIter(phaseInterface, interface_, iter)
165  {
166  if (0b01 << iter.index() & set)
167  {
168  talphaParameter.ref() +=
169  max(iter().residualAlpha(), alphas[iter().index()])
170  *parameters[iter.index()];
171  }
172  }
173 
174  return talphaParameter/alpha(alphas, set, true);
175 }
176 
177 
179 (
180  const UPtrList<const volScalarField>& alphas,
181  const label phaseSet,
182  const label systemSet
183 ) const
184 {
185  return
186  systemSet == 0b00
187  ? alpha(alphas, phaseSet, false)
188  : alpha(alphas, phaseSet, false)/alpha(alphas, systemSet, true);
189 }
190 
191 
193 (
194  const UPtrList<const volScalarField>& alphas,
195  const label phaseSet,
196  const label systemSet
197 ) const
198 {
199  label canBeContinuousPhaseSet = 0b00;
200  label canBeContinuousSystemSet = 0b00;
201  forAllConstIter(phaseInterface, interface_, iter)
202  {
203  if (canBeContinuous(iter.index()))
204  {
205  canBeContinuousPhaseSet += 0b01 << iter.index() & phaseSet;
206  canBeContinuousSystemSet += 0b01 << iter.index() & systemSet;
207  }
208  }
209 
210  if (canBeContinuousPhaseSet == 0)
211  {
212  return constant(alphas, 0);
213  }
214 
215  if (canBeContinuousPhaseSet == canBeContinuousSystemSet)
216  {
217  return constant(alphas, 1);
218  }
219 
220  return
221  fContinuous
222  (
223  alphas,
224  canBeContinuousPhaseSet,
225  canBeContinuousSystemSet
226  );
227 }
228 
229 
230 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
231 
233 (
234  const dictionary& dict,
235  const phaseInterface& interface
236 )
237 :
238  interface_(interface)
239 {}
240 
241 
242 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
243 
245 {}
246 
247 
248 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
249 
251 (
252  const UPtrList<const volScalarField>& alphas
253 ) const
254 {
255  return f(alphas, 0b10, 0b11);
256 }
257 
258 
260 (
261  const UPtrList<const volScalarField>& alphas
262 ) const
263 {
264  return f(alphas, 0b01, 0b11);
265 }
266 
267 
269 (
270  const UPtrList<const volScalarField>& alphas
271 ) const
272 {
273  return 1 - f(alphas, 0b11, 0b00);
274 }
275 
276 
277 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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)
Definition: errorManip.H:124
tmp< volScalarField > alpha(const UPtrList< const volScalarField > &alphas, const label set, const bool protect) const
Get the volume fraction of the given set.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
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 > &parameters) 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.
Definition: contiguous.H:49
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.
Definition: complex.C:47
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.
Definition: error.H:318
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.
Definition: PtrList.H:53
Namespace for OpenFOAM.
tmp< volScalarField > fDisplaced(const UPtrList< const volScalarField > &alphas) const
Return the coefficient for when the interface is displaced by a.
IOerror FatalIOError