multiphaseVoFMixture.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) 2023 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 "multiphaseVoFMixture.H"
27 #include "correctContactAngle.H"
28 #include "surfaceInterpolate.H"
29 #include "fvcGrad.H"
30 #include "fvcSnGrad.H"
31 #include "fvcDiv.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
38 }
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
44 {
45  scalar level = 0.0;
46  alphas_ == 0.0;
47 
48  forAll(phases_, phasei)
49  {
50  alphas_ += level*phases_[phasei];
51  level += 1.0;
52  }
53 }
54 
55 
56 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
57 
59 (
60  const fvMesh& mesh,
61  const VoFphase::iNew& inew
62 )
63 :
64  VoFMixture(mesh),
66  (
67  IOobject
68  (
69  "phaseProperties",
70  mesh.time().constant(),
71  mesh,
72  IOobject::MUST_READ_IF_MODIFIED,
73  IOobject::NO_WRITE
74  )
75  ),
76 
77  phases_(lookup("phases"), inew),
78 
79  mesh_(mesh),
80 
81  alphas_
82  (
83  IOobject
84  (
85  "alphas",
86  mesh_.time().name(),
87  mesh_,
88  IOobject::NO_READ,
89  IOobject::AUTO_WRITE
90  ),
91  mesh_,
93  ),
94 
95  sigmas_(lookup("sigmas")),
96  dimSigma_(1, 0, -2, 0, 0),
97 
98  deltaN_
99  (
100  "deltaN",
101  1e-8/pow(average(mesh_.V()), 1.0/3.0)
102  )
103 {
104  calcAlphas();
105  alphas_.write();
106 }
107 
108 
109 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
110 
113 {
115  (
117  (
118  "surfaceTensionForce",
119  mesh_,
120  dimensionedScalar(dimensionSet(1, -2, -2, 0, 0), 0)
121  )
122  );
123 
124  surfaceScalarField& stf = tstf.ref();
125 
126  forAll(phases_, phasei)
127  {
128  const VoFphase& alpha1 = phases_[phasei];
129 
130  for (label phasej = phasei+1; phasej<phases_.size(); phasej++)
131  {
132  const VoFphase& alpha2 = phases_[phasej];
133 
134  sigmaTable::const_iterator sigma =
135  sigmas_.find(interfacePair(alpha1, alpha2));
136 
137  if (sigma == sigmas_.end())
138  {
140  << "Cannot find interface " << interfacePair(alpha1, alpha2)
141  << " in list of sigma values"
142  << exit(FatalError);
143  }
144 
145  stf += dimensionedScalar(dimSigma_, sigma())
146  *fvc::interpolate(K(alpha1, alpha2, U))*
147  (
148  fvc::interpolate(alpha2)*fvc::snGrad(alpha1)
149  - fvc::interpolate(alpha1)*fvc::snGrad(alpha2)
150  );
151  }
152  }
153 
154  return tstf;
155 }
156 
157 
159 (
160  const volScalarField& alpha1,
161  const volScalarField& alpha2
162 ) const
163 {
164  /*
165  // Cell gradient of alpha
166  volVectorField gradAlpha =
167  alpha2*fvc::grad(alpha1) - alpha1*fvc::grad(alpha2);
168 
169  // Interpolated face-gradient of alpha
170  surfaceVectorField gradAlphaf = fvc::interpolate(gradAlpha);
171  */
172 
173  surfaceVectorField gradAlphaf
174  (
176  - fvc::interpolate(alpha1)*fvc::interpolate(fvc::grad(alpha2))
177  );
178 
179  // Face unit interface normal
180  return gradAlphaf/(mag(gradAlphaf) + deltaN_);
181 }
182 
183 
185 (
186  const volScalarField& alpha1,
187  const volScalarField& alpha2
188 ) const
189 {
190  // Face unit interface normal flux
191  return nHatfv(alpha1, alpha2) & mesh_.Sf();
192 }
193 
194 
196 (
197  const VoFphase& alpha1,
198  const VoFphase& alpha2,
199  const volVectorField& U
200 ) const
201 {
202  tmp<surfaceVectorField> tnHatfv = nHatfv(alpha1, alpha2);
203 
205  (
206  alpha1,
207  alpha2,
208  U.boundaryField(),
209  deltaN_,
210  tnHatfv.ref().boundaryFieldRef()
211  );
212 
213  // Simple expression for curvature
214  return -fvc::div(tnHatfv & mesh_.Sf());
215 }
216 
217 
220 {
221  tmp<volScalarField> tnearInt
222  (
224  (
225  "nearInterface",
226  mesh_,
228  )
229  );
230 
231  forAll(phases_, phasei)
232  {
233  tnearInt.ref() = max
234  (
235  tnearInt(),
236  pos0(phases_[phasei] - 0.01)*pos0(0.99 - phases_[phasei])
237  );
238  }
239 
240  return tnearInt;
241 }
242 
243 
245 {
246  if (regIOobject::read())
247  {
248  lookup("sigmas") >> sigmas_;
249  return true;
250  }
251  else
252  {
253  return false;
254  }
255 }
256 
257 
258 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
Two-phase VoF mixture.
Definition: VoFMixture.H:52
Return a pointer to a new VoFphase.
Definition: VoFphase.H:77
Single incompressible VoFphase derived from the phase-fraction. Used as part of the multiPhaseMixture...
Definition: VoFphase.H:53
Dimension set for the base types.
Definition: dimensionSet.H:125
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Multiphase VoF mixture with support for interface properties.
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
tmp< volScalarField > K(const VoFphase &alpha1, const VoFphase &alpha2, const volVectorField &U) const
tmp< surfaceScalarField > surfaceTensionForce(const volVectorField &U) const
PtrListDictionary< VoFphase > phases_
Dictionary of phases.
tmp< surfaceScalarField > nHatf(const volScalarField &alpha1, const volScalarField &alpha2) const
tmp< surfaceVectorField > nHatfv(const volScalarField &alpha1, const volScalarField &alpha2) const
multiphaseVoFMixture(const fvMesh &mesh, const VoFphase::iNew &)
Construct from fvMesh.
bool read()
Read base phaseProperties dictionary.
virtual bool write(const bool write=true) const
Write using setting from DB.
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
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Calculate the divergence of the given field.
Calculate the gradient of the given field.
Calculate the snGrad of the given volField.
K
Definition: pEqn.H:75
U
Definition: pEqn.H:72
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
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
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const doubleScalar e
Definition: doubleScalar.H:106
dimensionedScalar pos0(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void correctContactAngle(const volScalarField &alpha1, const volScalarField &alpha2, const volVectorField::Boundary &Ubf, const dimensionedScalar &deltaN, surfaceVectorField::Boundary &nHatbf)
Correct the contact angle for the two volume fraction fields.
dimensioned< scalar > mag(const dimensioned< Type > &)
defineTypeNameAndDebug(combustionModel, 0)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
error FatalError
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.