VoFSolidificationMeltingSource.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017 OpenFOAM Foundation
9 Copyright (C) 2020-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
33
34// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
35
36namespace Foam
37{
38 namespace fv
39 {
40 defineTypeNameAndDebug(VoFSolidificationMeltingSource, 0);
41
43 (
44 option,
45 VoFSolidificationMeltingSource,
46 dictionary
47 );
48 }
49}
50
51
52// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53
54void Foam::fv::VoFSolidificationMeltingSource::update()
55{
56 if (curTimeIndex_ == mesh_.time().timeIndex())
57 {
58 return;
59 }
60
61 if (debug)
62 {
63 Info<< type() << ": " << name_
64 << " - updating solid phase fraction" << endl;
65 }
66
67 alphaSolid_.oldTime();
68
69 const twoPhaseMixtureThermo& thermo
70 (
71 mesh_.lookupObject<twoPhaseMixtureThermo>
72 (
74 )
75 );
76
77 const volScalarField& TVoF = thermo.thermo1().T();
78 const volScalarField CpVoF(thermo.thermo1().Cp());
79 const volScalarField& alphaVoF = thermo.alpha1();
80
81 forAll(cells_, i)
82 {
83 const label celli = cells_[i];
84
85 alphaSolid_[celli] = min
86 (
87 relax_*alphaVoF[celli]*alphaSolidT_->value(TVoF[celli])
88 + (1 - relax_)*alphaSolid_[celli],
89 alphaVoF[celli]
90 );
91 }
92
93 alphaSolid_.correctBoundaryConditions();
94
95 curTimeIndex_ = mesh_.time().timeIndex();
96}
97
98
99Foam::word Foam::fv::VoFSolidificationMeltingSource::alphaSolidName() const
100{
101 const twoPhaseMixtureThermo& thermo
102 (
103 mesh_.lookupObject<twoPhaseMixtureThermo>
104 (
106 )
107 );
108
109 const volScalarField& alphaVoF = thermo.alpha1();
110
111 return IOobject::groupName(alphaVoF.name(), "solid");
112}
113
114
115// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
116
118(
119 const word& sourceName,
120 const word& modelType,
121 const dictionary& dict,
122 const fvMesh& mesh
123)
124:
125 fv::cellSetOption(sourceName, modelType, dict, mesh),
126 alphaSolidT_(Function1<scalar>::New("alphaSolidT", coeffs_, &mesh)),
127 L_("L", dimEnergy/dimMass, coeffs_),
128 relax_(coeffs_.getOrDefault("relax", 0.9)),
129 Cu_(coeffs_.getOrDefault<scalar>("Cu", 100000)),
130 q_(coeffs_.getOrDefault<scalar>("q", 0.001)),
131 alphaSolid_
132 (
133 IOobject
134 (
135 alphaSolidName(),
136 mesh.time().timeName(),
137 mesh,
138 IOobject::READ_IF_PRESENT,
139 IOobject::AUTO_WRITE
140 ),
141 mesh,
143 zeroGradientFvPatchScalarField::typeName
144 ),
145 curTimeIndex_(-1)
146{
147 fieldNames_.resize(2);
148 fieldNames_[0] = "U";
149 fieldNames_[1] = "T";
150
151 fv::option::resetApplied();
152}
153
154
155// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156
158(
159 fvMatrix<scalar>& eqn,
160 const label fieldi
161)
162{
163 apply(geometricOneField(), eqn);
164}
165
166
168(
169 const volScalarField& rho,
170 fvMatrix<scalar>& eqn,
171 const label fieldi
172)
173{
174 apply(rho, eqn);
175}
176
177
179(
180 fvMatrix<vector>& eqn,
181 const label fieldi
182)
183{
184 if (debug)
185 {
186 Info<< type() << ": applying source to " << eqn.psi().name() << endl;
187 }
188
189 update();
190
191 scalarField& Sp = eqn.diag();
192 const scalarField& V = mesh_.V();
193
194 forAll(cells_, i)
195 {
196 const label celli = cells_[i];
197 const scalar Vc = V[celli];
198 const scalar alphaFluid = 1 - alphaSolid_[celli];
199
200 const scalar S = Cu_*sqr(1 - alphaFluid)/(pow3(alphaFluid) + q_);
201
202 Sp[celli] -= Vc*S;
203 }
204}
205
206
208(
209 const volScalarField& rho,
210 fvMatrix<vector>& eqn,
211 const label fieldi
212)
213{
214 // Momentum source uses a Boussinesq approximation - redirect
215 addSup(eqn, fieldi);
216}
217
218
219// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
void correctBoundaryConditions()
Correct boundary field.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
label timeIndex() const noexcept
Return current time index.
Definition: TimeStateI.H:37
static const word dictName
Definition: basicThermo.H:256
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Solidification and melting model for VoF simulations.
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Add explicit contribution to enthalpy equation.
labelList cells_
Set of cells to apply source to.
const fvMesh & mesh_
Reference to the mesh database.
Definition: fvOption.H:139
const word name_
Source name.
Definition: fvOption.H:133
const Type & lookupObject(const word &name, const bool recursive=false) const
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
mesh update()
dynamicFvMesh & mesh
zeroField Sp
Definition: alphaSuSp.H:2
word timeName
Definition: getTimeIndex.H:3
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimEnergy
dimensionedScalar pow3(const dimensionedScalar &ds)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
static void apply(bitSet &selection, const Detail::parcelSelection::actionType action, const Predicate &accept, const UList< Type > &list, const AccessOp &aop)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
labelList fv(nPoints)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333