All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
VoFSolidificationMeltingSource.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) 2017-2021 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 
27 #include "twoPhaseMixtureThermo.H"
28 #include "fvcDdt.H"
31 
32 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  namespace fv
37  {
38  defineTypeNameAndDebug(VoFSolidificationMeltingSource, 0);
39 
41  (
42  fvModel,
43  VoFSolidificationMeltingSource,
44  dictionary
45  );
46  }
47 }
48 
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
52 void Foam::fv::VoFSolidificationMeltingSource::readCoeffs()
53 {
54  alphaSolidT_.reset(Function1<scalar>::New("alphaSolidT", coeffs()).ptr());
55  L_ = dimensionedScalar("L", dimEnergy/dimMass, coeffs());
56  relax_ = coeffs().lookupOrDefault<scalar>("relax", 0.9);
57  Cu_ = coeffs().lookupOrDefault<scalar>("Cu", 100000);
58  q_ = coeffs().lookupOrDefault<scalar>("q", 0.001);
59 }
60 
61 
62 void Foam::fv::VoFSolidificationMeltingSource::update() const
63 {
64  if (curTimeIndex_ == mesh().time().timeIndex())
65  {
66  return;
67  }
68 
69  if (debug)
70  {
71  Info<< type() << ": " << name()
72  << " - updating solid phase fraction" << endl;
73  }
74 
75  alphaSolid_.oldTime();
76 
77  const twoPhaseMixtureThermo& thermo
78  (
79  mesh().lookupObject<twoPhaseMixtureThermo>
80  (
82  )
83  );
84 
85  const volScalarField& TVoF = thermo.thermo1().T();
86  const volScalarField CpVoF(thermo.thermo1().Cp());
87  const volScalarField& alphaVoF = thermo.alpha1();
88 
89  const labelList& cells = set_.cells();
90 
91  forAll(cells, i)
92  {
93  const label celli = cells[i];
94 
95  alphaSolid_[celli] = min
96  (
97  relax_*alphaVoF[celli]*alphaSolidT_->value(TVoF[celli])
98  + (1 - relax_)*alphaSolid_[celli],
99  alphaVoF[celli]
100  );
101  }
102 
103  alphaSolid_.correctBoundaryConditions();
104 
105  curTimeIndex_ = mesh().time().timeIndex();
106 }
107 
108 
109 Foam::word Foam::fv::VoFSolidificationMeltingSource::alphaSolidName() const
110 {
111  const twoPhaseMixtureThermo& thermo
112  (
113  mesh().lookupObject<twoPhaseMixtureThermo>
114  (
116  )
117  );
118 
119  const volScalarField& alphaVoF = thermo.alpha1();
120 
121  return IOobject::groupName(alphaVoF.name(), "solid");
122 }
123 
124 
125 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
126 
128 (
129  const word& name,
130  const word& modelType,
131  const dictionary& dict,
132  const fvMesh& mesh
133 )
134 :
135  fvModel(name, modelType, dict, mesh),
136  set_(coeffs(), mesh),
137  alphaSolidT_(),
138  L_("L", dimEnergy/dimMass, NaN),
139  relax_(NaN),
140  Cu_(NaN),
141  q_(NaN),
142  alphaSolid_
143  (
144  IOobject
145  (
146  alphaSolidName(),
147  mesh.time().timeName(),
148  mesh,
149  IOobject::READ_IF_PRESENT,
150  IOobject::AUTO_WRITE
151  ),
152  mesh,
154  zeroGradientFvPatchScalarField::typeName
155  ),
156  curTimeIndex_(-1)
157 {
158  readCoeffs();
159 }
160 
161 
162 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
163 
165 {
166  return wordList({"U", "T"});
167 }
168 
169 
171 (
172  const volScalarField& rho,
173  fvMatrix<scalar>& eqn,
174  const word& fieldName
175 ) const
176 {
177  if (debug)
178  {
179  Info<< type() << ": applying source to " << eqn.psi().name() << endl;
180  }
181 
182  update();
183 
184  const twoPhaseMixtureThermo& thermo
185  (
186  mesh().lookupObject<twoPhaseMixtureThermo>
187  (
189  )
190  );
191 
192  const volScalarField CpVoF(thermo.thermo1().Cp());
193 
194  if (eqn.psi().dimensions() == dimTemperature)
195  {
196  eqn += L_/CpVoF*(fvc::ddt(rho, alphaSolid_));
197  }
198  else
199  {
200  eqn += L_*(fvc::ddt(rho, alphaSolid_));
201  }
202 }
203 
204 
206 (
207  const volScalarField& rho,
208  fvMatrix<vector>& eqn,
209  const word& fieldName
210 ) const
211 {
212  if (debug)
213  {
214  Info<< type() << ": applying source to " << eqn.psi().name() << endl;
215  }
216 
217  update();
218 
219  scalarField& Sp = eqn.diag();
220  const scalarField& V = mesh().V();
221 
222  const labelList& cells = set_.cells();
223 
224  forAll(cells, i)
225  {
226  const label celli = cells[i];
227  const scalar Vc = V[celli];
228  const scalar alphaFluid = 1 - alphaSolid_[celli];
229 
230  const scalar S = Cu_*sqr(1 - alphaFluid)/(pow3(alphaFluid) + q_);
231 
232  Sp[celli] -= Vc*rho[celli]*S;
233  }
234 }
235 
236 
238 (
239  const mapPolyMesh& mpm
240 )
241 {
242  set_.updateMesh(mpm);
243 }
244 
245 
246 bool Foam::fv::VoFSolidificationMeltingSource::read(const dictionary& dict)
247 {
248  if (fvModel::read(dict))
249  {
250  set_.read(coeffs());
251  readCoeffs();
252  return true;
253  }
254  else
255  {
256  return false;
257  }
258 
259  return false;
260 }
261 
262 
263 // ************************************************************************* //
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
fluidReactionThermo & thermo
Definition: createFields.H:28
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvModel.C:158
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const dimensionSet dimless
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
const fvMesh & mesh() const
Return const access to the mesh database.
Definition: fvModelI.H:34
Macros for easy insertion into run-time selection tables.
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:121
VoFSolidificationMeltingSource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from explicit source name and mesh.
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Calculate the first temporal derivative.
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
tmp< fvMatrix< Type > > S(const Pair< tmp< volScalarField::Internal >> &, const GeometricField< Type, fvPatchField, volMesh > &)
static word groupName(Name name, const word &group)
const dictionary & coeffs() const
Return dictionary.
Definition: fvModelI.H:40
word timeName
Definition: getTimeIndex.H:3
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual bool read(const dictionary &dict)
Read source dictionary.
const dimensionSet dimEnergy
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const dimensionSet dimMass
virtual void updateMesh(const mapPolyMesh &)
Update for mesh changes.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
addToRunTimeSelectionTable(fvConstraint, fixedTemperatureConstraint, dictionary)
List< word > wordList
A List of words.
Definition: fileName.H:54
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label timeIndex
Definition: getTimeIndex.H:4
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual void addSup(const volScalarField &rho, fvMatrix< scalar > &eqn, const word &fieldName) const
Add explicit contribution to compressible enthalpy equation.
messageStream Info
virtual wordList addSupFields() const
Return the list of fields for which the option adds source term.
const dimensionSet dimTemperature
Namespace for OpenFOAM.
static autoPtr< Function1< scalar > > New(const word &name, const dictionary &dict)
Selector.
Definition: Function1New.C:32