BlendedInterfacialModel.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-2018 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 "BlendedInterfacialModel.H"
28 #include "surfaceInterpolate.H"
29 
30 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31 
32 template<class modelType>
33 template<class GeometricField>
35 (
36  GeometricField& field
37 ) const
38 {
39  typename GeometricField::Boundary& fieldBf =
40  field.boundaryFieldRef();
41 
42  forAll(pair_.phase1().phi().boundaryField(), patchi)
43  {
44  if
45  (
46  isA<fixedValueFvsPatchScalarField>
47  (
48  pair_.phase1().phi().boundaryField()[patchi]
49  )
50  )
51  {
52  fieldBf[patchi] = Zero;
53  }
54  }
55 }
56 
57 
58 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
59 
60 template<class modelType>
62 (
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
69 )
70 :
71  pair_(pair),
72  pair1In2_(pair1In2),
73  pair2In1_(pair2In1),
74  blending_(blending),
75  correctFixedFluxBCs_(correctFixedFluxBCs)
76 {
77  if (modelTable.found(pair_))
78  {
79  model_.set
80  (
82  (
83  modelTable[pair_],
84  pair_
85  ).ptr()
86  );
87  }
88 
89  if (modelTable.found(pair1In2_))
90  {
91  model1In2_.set
92  (
94  (
95  modelTable[pair1In2_],
96  pair1In2_
97  ).ptr()
98  );
99  }
100 
101  if (modelTable.found(pair2In1_))
102  {
103  model2In1_.set
104  (
106  (
107  modelTable[pair2In1_],
108  pair2In1_
109  ).ptr()
110  );
111  }
112 }
113 
114 
115 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
116 
117 template<class modelType>
119 {}
120 
121 
122 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
123 
124 template<class modelType>
127 {
128  tmp<volScalarField> f1, f2;
129 
130  if (model_.valid() || model1In2_.valid())
131  {
132  f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
133  }
134 
135  if (model_.valid() || model2In1_.valid())
136  {
137  f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
138  }
139 
140  tmp<volScalarField> x
141  (
142  new volScalarField
143  (
144  IOobject
145  (
146  modelType::typeName + ":K",
147  pair_.phase1().mesh().time().timeName(),
148  pair_.phase1().mesh(),
151  false
152  ),
153  pair_.phase1().mesh(),
154  dimensionedScalar("zero", modelType::dimK, 0)
155  )
156  );
157 
158  if (model_.valid())
159  {
160  x.ref() += model_->K()*(f1() - f2());
161  }
162 
163  if (model1In2_.valid())
164  {
165  x.ref() += model1In2_->K()*(1 - f1);
166  }
167 
168  if (model2In1_.valid())
169  {
170  x.ref() += model2In1_->K()*f2;
171  }
172 
173  if
174  (
175  correctFixedFluxBCs_
176  && (model_.valid() || model1In2_.valid() || model2In1_.valid())
177  )
178  {
179  correctFixedFluxBCs(x.ref());
180  }
181 
182  return x;
183 }
184 
185 
186 template<class modelType>
189 {
190  tmp<surfaceScalarField> f1, f2;
191 
192  if (model_.valid() || model1In2_.valid())
193  {
194  f1 = fvc::interpolate
195  (
196  blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
197  );
198  }
199 
200  if (model_.valid() || model2In1_.valid())
201  {
202  f2 = fvc::interpolate
203  (
204  blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
205  );
206  }
207 
208  tmp<surfaceScalarField> x
209  (
211  (
212  IOobject
213  (
214  modelType::typeName + ":Kf",
215  pair_.phase1().mesh().time().timeName(),
216  pair_.phase1().mesh(),
219  false
220  ),
221  pair_.phase1().mesh(),
222  dimensionedScalar("zero", modelType::dimK, 0)
223  )
224  );
225 
226  if (model_.valid())
227  {
228  x.ref() += model_->Kf()*(f1() - f2());
229  }
230 
231  if (model1In2_.valid())
232  {
233  x.ref() += model1In2_->Kf()*(1 - f1);
234  }
235 
236  if (model2In1_.valid())
237  {
238  x.ref() += model2In1_->Kf()*f2;
239  }
240 
241  if
242  (
243  correctFixedFluxBCs_
244  && (model_.valid() || model1In2_.valid() || model2In1_.valid())
245  )
246  {
247  correctFixedFluxBCs(x.ref());
248  }
249 
250  return x;
251 }
252 
253 
254 template<class modelType>
255 template<class Type>
258 {
259  tmp<volScalarField> f1, f2;
260 
261  if (model_.valid() || model1In2_.valid())
262  {
263  f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
264  }
265 
266  if (model_.valid() || model2In1_.valid())
267  {
268  f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
269  }
270 
271  tmp<GeometricField<Type, fvPatchField, volMesh>> x
272  (
273  new GeometricField<Type, fvPatchField, volMesh>
274  (
275  IOobject
276  (
277  modelType::typeName + ":F",
278  pair_.phase1().mesh().time().timeName(),
279  pair_.phase1().mesh(),
282  false
283  ),
284  pair_.phase1().mesh(),
285  dimensioned<Type>("zero", modelType::dimF, Zero)
286  )
287  );
288 
289  if (model_.valid())
290  {
291  x.ref() += model_->F()*(f1() - f2());
292  }
293 
294  if (model1In2_.valid())
295  {
296  x.ref() += model1In2_->F()*(1 - f1);
297  }
298 
299  if (model2In1_.valid())
300  {
301  x.ref() -= model2In1_->F()*f2; // note : subtraction
302  }
303 
304  if
305  (
306  correctFixedFluxBCs_
307  && (model_.valid() || model1In2_.valid() || model2In1_.valid())
308  )
309  {
310  correctFixedFluxBCs(x.ref());
311  }
312 
313  return x;
314 }
315 
316 
317 template<class modelType>
320 {
321  tmp<surfaceScalarField> f1, f2;
322 
323  if (model_.valid() || model1In2_.valid())
324  {
325  f1 = fvc::interpolate
326  (
327  blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed())
328  );
329  }
330 
331  if (model_.valid() || model2In1_.valid())
332  {
333  f2 = fvc::interpolate
334  (
335  blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed())
336  );
337  }
338 
339  tmp<surfaceScalarField> x
340  (
342  (
343  IOobject
344  (
345  modelType::typeName + ":Ff",
346  pair_.phase1().mesh().time().timeName(),
347  pair_.phase1().mesh(),
350  false
351  ),
352  pair_.phase1().mesh(),
353  dimensionedScalar("zero", modelType::dimF*dimArea, 0)
354  )
355  );
356 
357  if (model_.valid())
358  {
359  x.ref() += model_->Ff()*(f1() - f2());
360  }
361 
362  if (model1In2_.valid())
363  {
364  x.ref() += model1In2_->Ff()*(1 - f1);
365  }
366 
367  if (model2In1_.valid())
368  {
369  x.ref() -= model2In1_->Ff()*f2; // note : subtraction
370  }
371 
372  if
373  (
374  correctFixedFluxBCs_
375  && (model_.valid() || model1In2_.valid() || model2In1_.valid())
376  )
377  {
378  correctFixedFluxBCs(x.ref());
379  }
380 
381  return x;
382 }
383 
384 
385 template<class modelType>
388 {
389  tmp<volScalarField> f1, f2;
390 
391  if (model_.valid() || model1In2_.valid())
392  {
393  f1 = blending_.f1(pair1In2_.dispersed(), pair2In1_.dispersed());
394  }
395 
396  if (model_.valid() || model2In1_.valid())
397  {
398  f2 = blending_.f2(pair1In2_.dispersed(), pair2In1_.dispersed());
399  }
400 
401  tmp<volScalarField> x
402  (
403  new volScalarField
404  (
405  IOobject
406  (
407  modelType::typeName + ":D",
408  pair_.phase1().mesh().time().timeName(),
409  pair_.phase1().mesh(),
412  false
413  ),
414  pair_.phase1().mesh(),
415  dimensionedScalar("zero", modelType::dimD, 0)
416  )
417  );
418 
419  if (model_.valid())
420  {
421  x.ref() += model_->D()*(f1() - f2());
422  }
423 
424  if (model1In2_.valid())
425  {
426  x.ref() += model1In2_->D()*(1 - f1);
427  }
428 
429  if (model2In1_.valid())
430  {
431  x.ref() += model2In1_->D()*f2;
432  }
433 
434  if
435  (
436  correctFixedFluxBCs_
437  && (model_.valid() || model1In2_.valid() || model2In1_.valid())
438  )
439  {
440  correctFixedFluxBCs(x.ref());
441  }
442 
443  return x;
444 }
445 
446 
447 template<class modelType>
449 (
450  const class phaseModel& phase
451 ) const
452 {
453  return
454  &phase == &(pair_.phase1())
455  ? model1In2_.valid()
456  : model2In1_.valid();
457 }
458 
459 
460 template<class modelType>
462 (
463  const class phaseModel& phase
464 ) const
465 {
466  return &phase == &(pair_.phase1()) ? model1In2_ : model2In1_;
467 }
468 
469 
470 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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.
scalar f1
Definition: createFields.H:28
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
static const zero Zero
Definition: zero.H:97
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.
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< volScalarField > D() const
Return the blended diffusivity.
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57