twoPhaseMixture.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) 2011-2024 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 "twoPhaseMixture.H"
27 #include "viscosityModel.H"
28 #include "surfaceInterpolate.H"
29 #include "fvcGrad.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
36 }
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * //
39 
41 Foam::twoPhaseMixture::readPhasePropertiesDict
42 (
43  const objectRegistry& obr
44 )
45 {
46  typeIOobject<IOdictionary> phasePropertiesIO
47  (
48  "phaseProperties",
49  obr.time().constant(),
50  obr,
53  true
54  );
55 
56  if (phasePropertiesIO.headerOk())
57  {
58  return phasePropertiesIO;
59  }
60  else
61  {
62  typeIOobject<IOdictionary> thermophysicalPropertiesIO
63  (
64  "thermophysicalProperties",
65  obr.time().constant(),
66  obr,
69  true
70  );
71 
72  if (thermophysicalPropertiesIO.headerOk())
73  {
74  IOdictionary phasePropertiesDict(thermophysicalPropertiesIO);
75  phasePropertiesDict.rename("phaseProperties");
76  return phasePropertiesDict;
77  }
78  else
79  {
80  typeIOobject<IOdictionary> transportPropertiesIO
81  (
82  "transportProperties",
83  obr.time().constant(),
84  obr,
87  true
88  );
89 
90  if (transportPropertiesIO.headerOk())
91  {
92  IOdictionary phasePropertiesDict(transportPropertiesIO);
93  phasePropertiesDict.rename("phaseProperties");
94 
95  const wordList phases(phasePropertiesDict.lookup("phases"));
96 
97  forAll(phases, i)
98  {
99  IOdictionary phaseDict
100  (
101  IOobject
102  (
104  (
105  physicalProperties::typeName,
106  phases[i]
107  ),
108  obr.time().constant(),
109  obr,
112  true
113  )
114  );
115 
116  phaseDict.merge(phasePropertiesDict.subDict(phases[i]));
117 
118  phaseDict.changeKeyword
119  (
120  "transportModel",
121  viscosityModel::typeName
122  );
123 
124  phaseDict.writeObject
125  (
129  true
130  );
131 
132  phasePropertiesDict.remove(phases[i]);
133  }
134 
135  phasePropertiesDict.writeObject
136  (
140  true
141  );
142 
144  << "Upgrading case by "
145  "converting transportProperties into phaseProperties, "
147  (
148  physicalProperties::typeName,
149  phases[0]
150  )
151  << " and "
153  (
154  physicalProperties::typeName,
155  phases[1]
156  )
157  << nl << endl;
158 
159  return phasePropertiesDict;
160  }
161  else
162  {
163  return phasePropertiesIO;
164  }
165  }
166  }
167 }
168 
169 
170 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
171 
173 :
174  IOdictionary(readPhasePropertiesDict(mesh)),
175 
176  phase1Name_(wordList(lookup("phases"))[0]),
177  phase2Name_(wordList(lookup("phases"))[1]),
178 
179  alpha1_
180  (
181  IOobject
182  (
183  IOobject::groupName("alpha", phase1Name_),
184  mesh.time().name(),
185  mesh,
186  IOobject::MUST_READ,
187  IOobject::AUTO_WRITE
188  ),
189  mesh
190  ),
191 
192  alpha2_
193  (
194  IOobject
195  (
196  IOobject::groupName("alpha", phase2Name_),
197  mesh.time().name(),
198  mesh
199  ),
200  1.0 - alpha1_
201  )
202 {}
203 
204 
205 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
206 
208 {
209  return regIOobject::read();
210 }
211 
212 
215 {
216  const fvMesh& mesh = alpha1_.mesh();
217 
219  (
221  (
222  "A",
223  mesh,
225  )
226  );
228 
229  const surfaceVectorField& Sf = mesh.Sf();
230  const labelUList& own = mesh.owner();
231  const labelUList& nei = mesh.neighbour();
232 
233  const surfaceScalarField alphaf(fvc::interpolate(alpha1_));
234 
235  const volVectorField::Internal gradAlpha(fvc::grad(alpha1_));
237  (
238  gradAlpha/(mag(gradAlpha) + dimensionedScalar(dimless/dimLength, small))
239  );
240 
241  const scalarField& ialpha = alpha1_;
242  const scalarField& ialphaf = alphaf;
243  scalarField sumnSf(mesh.nCells(), 0);
244 
245  forAll(own, facei)
246  {
247  {
248  const scalar nSf(mag(n[own[facei]] & Sf[facei]));
249  A[own[facei]] += nSf*(ialphaf[facei] - ialpha[own[facei]]);
250  sumnSf[own[facei]] += nSf;
251  }
252  {
253  const scalar nSf(mag(n[nei[facei]] & Sf[facei]));
254  A[nei[facei]] += nSf*(ialphaf[facei] - ialpha[nei[facei]]);
255  sumnSf[nei[facei]] += nSf;
256  }
257  }
258 
259  forAll(mesh.boundary(), patchi)
260  {
261  const labelUList& own = mesh.boundary()[patchi].faceCells();
262  const fvsPatchScalarField& palphaf = alphaf.boundaryField()[patchi];
263 
264  forAll(mesh.boundary()[patchi], facei)
265  {
266  const scalar nSf(mag(n[own[facei]] & Sf[facei]));
267  A[own[facei]] += nSf*(palphaf[facei] - ialpha[own[facei]]);
268  sumnSf[own[facei]] += nSf;
269  }
270  }
271 
272  scalarField& a = A.primitiveFieldRef();
273  forAll(a, i)
274  {
275  if (sumnSf[i] > small)
276  {
277  a[i] = 2*mag(a[i])/sumnSf[i];
278  }
279  else
280  {
281  a[i] = 0;
282  }
283  }
284 
285  return tA;
286 }
287 
288 
289 // ************************************************************************* //
static const Foam::dimensionedScalar A("A", Foam::dimPressure, 611.21)
label n
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &, const Field< Type > &)
Return a temporary field constructed from name, mesh,.
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
IOdictionary(const IOobject &io, const word &wantedType)
Construct given an IOobject, supply wanted typeName.
Definition: IOdictionary.C:47
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
@ MUST_READ_IF_MODIFIED
Definition: IOobject.H:119
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:167
static word groupName(Name name, const word &group)
static const versionNumber currentVersion
Current version number.
Definition: IOstream.H:203
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:469
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:857
const surfaceVectorField & Sf() const
Return cell face area vectors.
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:475
const polyMesh & mesh() const
Return reference to polyMesh.
Definition: fvMesh.H:441
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:83
label nCells() const
virtual bool read()
Read object.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Class to represent a mixture of two phases.
twoPhaseMixture(const fvMesh &mesh)
Construct from components.
virtual bool read()=0
Read base phaseProperties dictionary.
tmp< volScalarField::Internal > interfaceFraction() const
Interface fraction in a cell.
Templated form of IOobject providing type information for file reading and header type checking.
Definition: IOobject.H:531
Calculate the gradient of the given field.
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
Namespace for OpenFOAM.
List< word > wordList
A List of words.
Definition: fileName.H:54
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
const dimensionSet dimLength
dimensioned< scalar > mag(const dimensioned< Type > &)
defineTypeNameAndDebug(combustionModel, 0)
static const char nl
Definition: Ostream.H:266
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.