fixedCoeff.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) 2012-2023 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 "fixedCoeff.H"
28 #include "fvMatrices.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  namespace porosityModels
35  {
38  }
39 }
40 
41 
42 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 void Foam::porosityModels::fixedCoeff::apply
45 (
46  scalarField& Udiag,
47  vectorField& Usource,
48  const scalarField& V,
49  const vectorField& U,
50  const scalar rho
51 ) const
52 {
54 
55  forAll(cells, i)
56  {
57  const label celli = cells[i];
58  const label j = fieldIndex(i);
59  const tensor Cd = rho*(alpha_[j] + beta_[j]*mag(U[celli]));
60  const scalar isoCd = tr(Cd);
61 
62  Udiag[celli] += V[celli]*isoCd;
63  Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
64  }
65 }
66 
67 
68 void Foam::porosityModels::fixedCoeff::apply
69 (
70  tensorField& AU,
71  const vectorField& U,
72  const scalar rho
73 ) const
74 {
75  const labelList& cells = mesh_.cellZones()[zoneName_];
76 
77  forAll(cells, i)
78  {
79  const label celli = cells[i];
80  const label j = fieldIndex(i);
81  const tensor alpha = alpha_[j];
82  const tensor beta = beta_[j];
83 
84  AU[celli] += rho*(alpha + beta*mag(U[celli]));
85  }
86 }
87 
88 
89 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
90 
92 (
93  const word& name,
94  const word& modelType,
95  const fvMesh& mesh,
96  const dictionary& dict,
97  const word& cellZoneName
98 )
99 :
100  porosityModel(name, modelType, mesh, dict, cellZoneName),
101  alphaXYZ_("alpha", dimless/dimTime, coeffs_),
102  betaXYZ_("beta", dimless/dimLength, coeffs_)
103 {
104  adjustNegativeResistance(alphaXYZ_);
105  adjustNegativeResistance(betaXYZ_);
106 
108 }
109 
110 
111 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
112 
114 {}
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
120 {
121  if (coordSys_.R().uniform())
122  {
123  alpha_.setSize(1);
124  beta_.setSize(1);
125 
126  alpha_[0] = Zero;
127  alpha_[0].xx() = alphaXYZ_.value().x();
128  alpha_[0].yy() = alphaXYZ_.value().y();
129  alpha_[0].zz() = alphaXYZ_.value().z();
130  alpha_[0] = coordSys_.R().transform(Zero, alpha_[0]);
131 
132  beta_[0] = Zero;
133  beta_[0].xx() = betaXYZ_.value().x();
134  beta_[0].yy() = betaXYZ_.value().y();
135  beta_[0].zz() = betaXYZ_.value().z();
136  beta_[0] = coordSys_.R().transform(Zero, beta_[0]);
137  }
138  else
139  {
140  const labelList& cells = mesh_.cellZones()[zoneName_];
141 
142  alpha_.setSize(cells.size());
143  beta_.setSize(cells.size());
144 
145  forAll(cells, i)
146  {
147  alpha_[i] = Zero;
148  alpha_[i].xx() = alphaXYZ_.value().x();
149  alpha_[i].yy() = alphaXYZ_.value().y();
150  alpha_[i].zz() = alphaXYZ_.value().z();
151 
152  beta_[i] = Zero;
153  beta_[i].xx() = betaXYZ_.value().x();
154  beta_[i].yy() = betaXYZ_.value().y();
155  beta_[i].zz() = betaXYZ_.value().z();
156  }
157 
158  const coordinateRotation& R = coordSys_.R
159  (
160  UIndirectList<vector>(mesh_.C(), cells)()
161  );
162 
163  alpha_ = R.transform(alpha_);
164  beta_ = R.transform(beta_);
165  }
166 }
167 
168 
170 (
171  const volVectorField& U,
172  const volScalarField& rho,
173  const volScalarField& mu,
174  vectorField& force
175 ) const
176 {
177  scalarField Udiag(U.size(), 0.0);
178  vectorField Usource(U.size(), Zero);
179  const scalarField& V = mesh_.V();
180  scalar rhoRef = coeffs_.lookup<scalar>("rhoRef");
181 
182  apply(Udiag, Usource, V, U, rhoRef);
183 
184  force = Udiag*U - Usource;
185 }
186 
187 
189 (
191 ) const
192 {
193  const vectorField& U = UEqn.psi();
194  const scalarField& V = mesh_.V();
195  scalarField& Udiag = UEqn.diag();
196  vectorField& Usource = UEqn.source();
197 
198  scalar rho = 1.0;
199  if (UEqn.dimensions() == dimForce)
200  {
201  coeffs_.lookup("rhoRef") >> rho;
202  }
203 
204  apply(Udiag, Usource, V, U, rho);
205 }
206 
207 
209 (
210  const fvVectorMatrix& UEqn,
211  volTensorField& AU
212 ) const
213 {
214  const vectorField& U = UEqn.psi();
215 
216  scalar rho = 1.0;
217  if (UEqn.dimensions() == dimForce)
218  {
219  coeffs_.lookup("rhoRef") >> rho;
220  }
221 
222  apply(AU, U, rho);
223 }
224 
225 
227 {
228  os << indent << name_ << endl;
229  dict_.write(os);
230 
231  return true;
232 }
233 
234 
235 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
virtual Ostream & write(const char)=0
Write character.
A List with indirect addressing.
Definition: UIndirectList.H:60
Abstract base class for coordinate rotation.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const cellZoneList & cellZones() const
Return cell zones.
Definition: polyMesh.H:452
Top level model for porosity models.
Definition: porosityModel.H:57
const fvMesh & mesh_
Reference to the mesh database.
Definition: porosityModel.H:75
word zoneName_
Name of cellZone.
Definition: porosityModel.H:84
void adjustNegativeResistance(dimensionedVector &resist)
Adjust negative resistance values to be multiplier of max value.
Definition: porosityModel.C:40
label fieldIndex(const label index) const
Return label index.
Definition: porosityModel.C:66
Fixed coefficient form of porosity model.
Definition: fixedCoeff.H:61
virtual void calcForce(const volVectorField &U, const volScalarField &rho, const volScalarField &mu, vectorField &force) const
Calculate the porosity force.
Definition: fixedCoeff.C:170
bool writeData(Ostream &os) const
Write.
Definition: fixedCoeff.C:226
virtual ~fixedCoeff()
Destructor.
Definition: fixedCoeff.C:113
virtual void correct(fvVectorMatrix &UEqn) const
Add resistance.
Definition: fixedCoeff.C:189
fixedCoeff(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict, const word &cellZoneName)
Definition: fixedCoeff.C:92
virtual void calcTransformModelData()
Transform the model data wrt mesh changes.
Definition: fixedCoeff.C:119
A class for handling words, derived from string.
Definition: word.H:62
fvVectorMatrix & UEqn
Definition: UEqn.H:11
A special matrix type and solver, designed for finite volume solutions of scalar equations.
volScalarField scalarField(fieldObject, mesh)
volVectorField vectorField(fieldObject, mesh)
const cellShapeList & cells
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const dimensionedScalar mu
Atomic mass unit.
addToRunTimeSelectionTable(porosityModel, powerLawLopesdaCosta, mesh)
defineTypeNameAndDebug(powerLawLopesdaCosta, 0)
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
Tensor< scalar > tensor
Tensor of scalars.
Definition: tensor.H:51
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
const dimensionSet dimLength
static const Identity< scalar > I
Definition: Identity.H:93
const dimensionSet dimForce
const dimensionSet dimTime
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
static scalar R(const scalar a, const scalar x)
Definition: invIncGamma.C:102
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:227
dictionary dict