solidificationMeltingSource.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2014-2017 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 "fvMatrices.H"
28 #include "basicThermo.H"
33 #include "geometricOneField.H"
34 
35 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  template<>
40  const char* NamedEnum
41  <
43  2
44  >::names[] =
45  {
46  "thermo",
47  "lookup"
48  };
49 
50  namespace fv
51  {
52  defineTypeNameAndDebug(solidificationMeltingSource, 0);
53 
55  (
56  option,
57  solidificationMeltingSource,
59  );
60  }
61 }
62 
65 
66 
67 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
68 
70 Foam::fv::solidificationMeltingSource::Cp() const
71 {
72  switch (mode_)
73  {
74  case mdThermo:
75  {
76  const basicThermo& thermo =
77  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
78 
79  return thermo.Cp();
80  break;
81  }
82  case mdLookup:
83  {
84  if (CpName_ == "CpRef")
85  {
86  scalar CpRef = readScalar(coeffs_.lookup("CpRef"));
87 
88  return tmp<volScalarField>
89  (
90  new volScalarField
91  (
92  IOobject
93  (
94  name_ + ":Cp",
95  mesh_.time().timeName(),
96  mesh_,
99  ),
100  mesh_,
102  (
103  "Cp",
105  CpRef
106  ),
107  extrapolatedCalculatedFvPatchScalarField::typeName
108  )
109  );
110  }
111  else
112  {
113  return mesh_.lookupObject<volScalarField>(CpName_);
114  }
115 
116  break;
117  }
118  default:
119  {
121  << "Unhandled thermo mode: " << thermoModeTypeNames_[mode_]
122  << abort(FatalError);
123  }
124  }
125 
126  return tmp<volScalarField>(nullptr);
127 }
128 
129 
130 Foam::vector Foam::fv::solidificationMeltingSource::g() const
131 {
132  if (mesh_.foundObject<uniformDimensionedVectorField>("g"))
133  {
134  const uniformDimensionedVectorField& value =
135  mesh_.lookupObject<uniformDimensionedVectorField>("g");
136  return value.value();
137  }
138  else
139  {
140  return coeffs_.lookup("g");
141  }
142 }
143 
144 
145 void Foam::fv::solidificationMeltingSource::update(const volScalarField& Cp)
146 {
147  if (curTimeIndex_ == mesh_.time().timeIndex())
148  {
149  return;
150  }
151 
152  if (debug)
153  {
154  Info<< type() << ": " << name_ << " - updating phase indicator" << endl;
155  }
156 
157  // update old time alpha1 field
158  alpha1_.oldTime();
159 
160  const volScalarField& T = mesh_.lookupObject<volScalarField>(TName_);
161 
162  forAll(cells_, i)
163  {
164  label celli = cells_[i];
165 
166  scalar Tc = T[celli];
167  scalar Cpc = Cp[celli];
168  scalar alpha1New = alpha1_[celli] + relax_*Cpc*(Tc - Tmelt_)/L_;
169 
170  alpha1_[celli] = max(0, min(alpha1New, 1));
171  deltaT_[i] = Tc - Tmelt_;
172  }
173 
174  alpha1_.correctBoundaryConditions();
175 
176  curTimeIndex_ = mesh_.time().timeIndex();
177 }
178 
179 
180 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
181 
182 Foam::fv::solidificationMeltingSource::solidificationMeltingSource
183 (
184  const word& sourceName,
185  const word& modelType,
186  const dictionary& dict,
187  const fvMesh& mesh
188 )
189 :
190  cellSetOption(sourceName, modelType, dict, mesh),
191  Tmelt_(readScalar(coeffs_.lookup("Tmelt"))),
192  L_(readScalar(coeffs_.lookup("L"))),
193  relax_(coeffs_.lookupOrDefault("relax", 0.9)),
194  mode_(thermoModeTypeNames_.read(coeffs_.lookup("thermoMode"))),
195  rhoRef_(readScalar(coeffs_.lookup("rhoRef"))),
196  TName_(coeffs_.lookupOrDefault<word>("T", "T")),
197  CpName_(coeffs_.lookupOrDefault<word>("Cp", "Cp")),
198  UName_(coeffs_.lookupOrDefault<word>("U", "U")),
199  phiName_(coeffs_.lookupOrDefault<word>("phi", "phi")),
200  Cu_(coeffs_.lookupOrDefault<scalar>("Cu", 100000)),
201  q_(coeffs_.lookupOrDefault("q", 0.001)),
202  beta_(readScalar(coeffs_.lookup("beta"))),
203  alpha1_
204  (
205  IOobject
206  (
207  name_ + ":alpha1",
208  mesh.time().timeName(),
209  mesh,
212  ),
213  mesh,
214  dimensionedScalar("alpha1", dimless, 0.0),
215  zeroGradientFvPatchScalarField::typeName
216  ),
217  curTimeIndex_(-1),
218  deltaT_(cells_.size(), 0)
219 {
220  fieldNames_.setSize(2);
221  fieldNames_[0] = UName_;
222 
223  switch (mode_)
224  {
225  case mdThermo:
226  {
227  const basicThermo& thermo =
228  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
229 
230  fieldNames_[1] = thermo.he().name();
231  break;
232  }
233  case mdLookup:
234  {
235  fieldNames_[1] = TName_;
236  break;
237  }
238  default:
239  {
241  << "Unhandled thermo mode: " << thermoModeTypeNames_[mode_]
242  << abort(FatalError);
243  }
244  }
245 
246  applied_.setSize(fieldNames_.size(), false);
247 }
248 
249 
250 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
251 
253 (
254  fvMatrix<scalar>& eqn,
255  const label fieldi
256 )
257 {
258  apply(geometricOneField(), eqn);
259 }
260 
261 
263 (
264  const volScalarField& rho,
265  fvMatrix<scalar>& eqn,
266  const label fieldi
267 )
268 {
269  apply(rho, eqn);
270 }
271 
272 
274 (
275  fvMatrix<vector>& eqn,
276  const label fieldi
277 )
278 {
279  if (debug)
280  {
281  Info<< type() << ": applying source to " << eqn.psi().name() << endl;
282  }
283 
284  const volScalarField Cp(this->Cp());
285 
286  update(Cp);
287 
288  vector g = this->g();
289 
290  scalarField& Sp = eqn.diag();
291  vectorField& Su = eqn.source();
292  const scalarField& V = mesh_.V();
293 
294  forAll(cells_, i)
295  {
296  label celli = cells_[i];
297 
298  scalar Vc = V[celli];
299  scalar alpha1c = alpha1_[celli];
300 
301  scalar S = -Cu_*sqr(1.0 - alpha1c)/(pow3(alpha1c) + q_);
302  vector Sb = rhoRef_*g*beta_*deltaT_[i];
303 
304  Sp[celli] += Vc*S;
305  Su[celli] += Vc*Sb;
306  }
307 }
308 
309 
311 (
312  const volScalarField& rho,
313  fvMatrix<vector>& eqn,
314  const label fieldi
315 )
316 {
317  // Momentum source uses a Boussinesq approximation - redirect
318  addSup(eqn, fieldi);
319 }
320 
321 
322 // ************************************************************************* //
defineTypeNameAndDebug(option, 0)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
const word & name() const
Return name.
Definition: IOobject.H:291
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
zeroField Su
Definition: alphaSuSp.H:1
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
static const NamedEnum< thermoMode, 2 > thermoModeTypeNames_
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:284
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:644
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:52
Macros for easy insertion into run-time selection tables.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
const Type & value() const
Return const reference to value.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
const dimensionedVector & g
Field< Type > & source()
Definition: fvMatrix.H:294
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure [J/kg/K].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
const dimensionSet dimEnergy
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Add explicit contribution to enthalpy equation.
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A special matrix type and solver, designed for finite volume solutions of scalar equations.
messageStream Info
Cell-set options abtract base class. Provides a base set of controls, e.g.:
Definition: cellSetOption.H:69
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
scalarField & diag()
Definition: lduMatrix.C:186
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
zeroField Sp
Definition: alphaSuSp.H:2