solidification.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) 2013-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 
26 #include "solidification.H"
28 #include "thermoSingleLayer.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace regionModels
35 {
36 namespace surfaceFilmModels
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 defineTypeNameAndDebug(solidification, 0);
42 
44 (
45  phaseChangeModel,
46  solidification,
47  dictionary
48 );
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
53 (
55  const dictionary& dict
56 )
57 :
58  phaseChangeModel(typeName, film, dict),
59  T0_(coeffDict_.lookup<scalar>("T0")),
60  maxSolidificationFrac_
61  (
62  coeffDict_.lookupOrDefault("maxSolidificationFrac", 0.2)
63  ),
64  maxSolidificationRate_
65  (
67  (
68  "maxSolidificationRate",
69  coeffDict_,
71  great
72  )
73  ),
74  mass_
75  (
76  IOobject
77  (
78  IOobject::modelName("mass", typeName),
79  film.regionMesh().time().timeName(),
80  film.regionMesh(),
83  ),
84  film.regionMesh(),
86  ),
87  thickness_
88  (
89  IOobject
90  (
91  IOobject::modelName("thickness", typeName),
92  film.regionMesh().time().timeName(),
93  film.regionMesh(),
96  ),
97  film.regionMesh(),
99  )
100 {}
101 
102 
103 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
104 
106 {}
107 
108 
109 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
110 
112 (
113  const scalar dt,
114  scalarField& availableMass,
115  scalarField& dMass,
116  scalarField& dEnergy
117 )
118 {
119  const thermoSingleLayer& film = filmType<thermoSingleLayer>();
120 
121  const scalarField& T = film.thermo().T();
122  const scalarField& he = film.thermo().he();
123  const scalarField& coverage = film.coverage();
124 
125  const scalar rateLimiter = min
126  (
128  (
131  ).value()
132  );
133 
134  forAll(coverage, celli)
135  {
136  if (coverage[celli] > 0.5)
137  {
138  if (T[celli] < T0_)
139  {
140  const scalar dm = rateLimiter*availableMass[celli];
141 
142  mass_[celli] += dm;
143  dMass[celli] += dm;
144 
145  // Heat is assumed to be removed by heat-transfer to the wall
146  // so the energy remains unchanged by the phase-change.
147  dEnergy[celli] += dm*he[celli];
148  }
149  }
150  }
151 
152  thickness_ = mass_/film.magSf()/film.rho();
153 }
154 
155 
156 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
157 
158 } // End namespace surfaceFilmModels
159 } // End namespace regionModels
160 } // End namespace Foam
161 
162 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
const dimensionSet dimless
Generic dimensioned Type class.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
Macros for easy insertion into run-time selection tables.
const dimensionSet dimLength
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
Base class for surface film phase change models.
const dimensionSet dimTime
virtual void correctModel(const scalar dt, scalarField &availableMass, scalarField &dMass, scalarField &dEnergy)
Correct.
scalar maxSolidificationFrac_
Solidification limiter.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
const volScalarField & coverage() const
Return the film coverage, 1 = covered, 0 = uncovered [].
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:55
dimensionedScalar maxSolidificationRate_
Solidification limiter.
thermo he()
solidification(surfaceFilmRegionModel &film, const dictionary &dict)
Construct from surface film model.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const dimensionSet dimMass
volScalarField::Internal thickness_
Accumulated solid thickness [m].
volScalarField::Internal mass_
Accumulated solid mass [kg].
scalar T0_
Temperature at which solidification starts.
const volScalarField::Internal & magSf() const
Return the face area magnitudes [m^2].
virtual const volScalarField & T() const =0
Temperature [K].
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
surfaceFilmRegionModel & filmModel_
Reference to the film surface film model.
const volScalarField & rho() const
Return the film density [kg/m^3].
static word modelName(Name name, const word &model)
Return the name of the object within the given model.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
defineTypeNameAndDebug(kinematicSingleLayer, 0)
addToRunTimeSelectionTable(surfaceFilmModel, noFilm, mesh)
Thermodynamic form of single-cell layer surface film model.
Namespace for OpenFOAM.