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) 2017-2018 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 "solidification.H"
28 #include "geometricOneField.H"
29 #include "fvMatrices.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  namespace porosityModels
36  {
37  defineTypeNameAndDebug(solidification, 0);
38  addToRunTimeSelectionTable(porosityModel, solidification, mesh);
39  }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
45 Foam::porosityModels::solidification::solidification
46 (
47  const word& name,
48  const word& modelType,
49  const fvMesh& mesh,
50  const dictionary& dict,
51  const word& cellZoneName
52 )
53 :
54  porosityModel(name, modelType, mesh, dict, cellZoneName),
55  TName_(coeffs_.lookupOrDefault<word>("T", "T")),
56  alphaName_(coeffs_.lookupOrDefault<word>("alpha", "none")),
57  rhoName_(coeffs_.lookupOrDefault<word>("rho", "rho")),
58  D_(Function1<scalar>::New("D", coeffs_))
59 {}
60 
61 
62 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
63 
65 {}
66 
67 
68 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
69 
71 {}
72 
73 
75 (
76  const volVectorField& U,
77  const volScalarField& rho,
78  const volScalarField& mu,
79  vectorField& force
80 ) const
81 {
82  scalarField Udiag(U.size(), 0.0);
83  const scalarField& V = mesh_.V();
84 
85  apply(Udiag, V, rho, U);
86 
87  force = Udiag*U;
88 }
89 
90 
92 (
93  fvVectorMatrix& UEqn
94 ) const
95 {
96  const volVectorField& U = UEqn.psi();
97  const scalarField& V = mesh_.V();
98  scalarField& Udiag = UEqn.diag();
99 
100  if (UEqn.dimensions() == dimForce)
101  {
102  const volScalarField& rho = mesh_.lookupObject<volScalarField>
103  (
104  IOobject::groupName(rhoName_, U.group())
105  );
106 
107  apply(Udiag, V, rho, U);
108  }
109  else
110  {
111  apply(Udiag, V, geometricOneField(), U);
112  }
113 }
114 
115 
117 (
118  fvVectorMatrix& UEqn,
119  const volScalarField& rho,
120  const volScalarField& mu
121 ) const
122 {
123  const volVectorField& U = UEqn.psi();
124  const scalarField& V = mesh_.V();
125  scalarField& Udiag = UEqn.diag();
126 
127  apply(Udiag, V, rho, U);
128 }
129 
130 
132 (
133  const fvVectorMatrix& UEqn,
134  volTensorField& AU
135 ) const
136 {
137  const volVectorField& U = UEqn.psi();
138 
139  if (UEqn.dimensions() == dimForce)
140  {
141  const volScalarField& rho = mesh_.lookupObject<volScalarField>
142  (
143  IOobject::groupName(rhoName_, U.group())
144  );
145 
146  apply(AU, rho, U);
147  }
148  else
149  {
150  apply(AU, geometricOneField(), U);
151  }
152 }
153 
154 
156 {
157  os << indent << name_ << endl;
158  dict_.write(os);
159 
160  return true;
161 }
162 
163 
164 // ************************************************************************* //
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:176
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
virtual void calcTransformModelData()
Transform the model data wrt mesh changes.
defineTypeNameAndDebug(powerLawLopesdaCosta, 0)
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
addToRunTimeSelectionTable(porosityModel, powerLawLopesdaCosta, mesh)
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:284
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Macros for easy insertion into run-time selection tables.
virtual void correct(fvVectorMatrix &UEqn) const
Add resistance.
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
static word groupName(Name name, const word &group)
virtual void calcForce(const volVectorField &U, const volScalarField &rho, const volScalarField &mu, vectorField &force) const
Calculate the porosity force.
virtual ~solidification()
Destructor.
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
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const dimensionSet dimForce
U
Definition: pEqn.H:72
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.
virtual Ostream & write(const token &)=0
Write next token to stream.
scalarField & diag()
Definition: lduMatrix.C:186
bool writeData(Ostream &os) const
Write.
Top level model for porosity models.
Definition: porosityModel.H:55
Namespace for OpenFOAM.
const dimensionSet & dimensions() const
Definition: fvMatrix.H:289