fixedCoeff.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) 2012-2016 OpenFOAM Foundation
9 Copyright (C) 2018-2022 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
30#include "fixedCoeff.H"
31#include "fvMatrices.H"
32#include "pointIndList.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
38 namespace porosityModels
39 {
42 }
43}
44
45
46// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47
49(
50 scalarField& Udiag,
51 vectorField& Usource,
52 const scalarField& V,
53 const vectorField& U,
54 const scalar rho
55) const
56{
57 forAll(cellZoneIDs_, zoneI)
58 {
59 const tensorField& alphaZones = alpha_[zoneI];
60 const tensorField& betaZones = beta_[zoneI];
61
62 const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
63
64 forAll(cells, i)
65 {
66 const label celli = cells[i];
67 const label j = fieldIndex(i);
68 const tensor Cd = rho*(alphaZones[j] + betaZones[j]*mag(U[celli]));
69 const scalar isoCd = tr(Cd);
70
71 Udiag[celli] += V[celli]*isoCd;
72 Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
73 }
74 }
75}
76
77
79(
80 tensorField& AU,
81 const vectorField& U,
82 const scalar rho
83) const
84{
85 forAll(cellZoneIDs_, zoneI)
86 {
87 const tensorField& alphaZones = alpha_[zoneI];
88 const tensorField& betaZones = beta_[zoneI];
89
90 const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
91
92 forAll(cells, i)
93 {
94 const label celli = cells[i];
95 const label j = fieldIndex(i);
96 const tensor alpha = alphaZones[j];
97 const tensor beta = betaZones[j];
98
99 AU[celli] += rho*(alpha + beta*mag(U[celli]));
100 }
101 }
102}
103
104
105// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
106
108(
109 const word& name,
110 const word& modelType,
111 const fvMesh& mesh,
112 const dictionary& dict,
113 const wordRe& cellZoneName
114)
115:
116 porosityModel(name, modelType, mesh, dict, cellZoneName),
117 alphaXYZ_("alpha", dimless/dimTime, coeffs_),
118 betaXYZ_("beta", dimless/dimLength, coeffs_),
119 alpha_(cellZoneIDs_.size()),
120 beta_(cellZoneIDs_.size())
121{
122 adjustNegativeResistance(alphaXYZ_);
123 adjustNegativeResistance(betaXYZ_);
124
126}
127
128
129// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
130
132{
133 // The alpha coefficient as a tensor
134 tensor alphaCoeff(Zero);
135 alphaCoeff.xx() = alphaXYZ_.value().x();
136 alphaCoeff.yy() = alphaXYZ_.value().y();
137 alphaCoeff.zz() = alphaXYZ_.value().z();
138
139 // The beta coefficient as a tensor
140 tensor betaCoeff(Zero);
141 betaCoeff.xx() = betaXYZ_.value().x();
142 betaCoeff.yy() = betaXYZ_.value().y();
143 betaCoeff.zz() = betaXYZ_.value().z();
144
145 if (csys().uniform())
146 {
147 forAll(cellZoneIDs_, zonei)
148 {
149 alpha_[zonei].resize(1);
150 beta_[zonei].resize(1);
151
152 alpha_[zonei] = csys().transform(alphaCoeff);
153 beta_[zonei] = csys().transform(betaCoeff);
154 }
155 }
156 else
157 {
158 forAll(cellZoneIDs_, zonei)
159 {
160 const pointUIndList cc
161 (
162 mesh_.cellCentres(),
163 mesh_.cellZones()[cellZoneIDs_[zonei]]
164 );
165
166 alpha_[zonei] = csys().transform(cc, alphaCoeff);
167 beta_[zonei] = csys().transform(cc, betaCoeff);
168 }
169 }
170}
171
172
174(
175 const volVectorField& U,
176 const volScalarField& rho,
177 const volScalarField& mu,
178 vectorField& force
179) const
180{
181 scalarField Udiag(U.size(), Zero);
182 vectorField Usource(U.size(), Zero);
183 const scalarField& V = mesh_.V();
184 const scalar rhoRef = coeffs_.get<scalar>("rhoRef");
185
186 apply(Udiag, Usource, V, U, rhoRef);
187
188 force = Udiag*U - Usource;
189}
190
191
193(
195) const
196{
197 const vectorField& U = UEqn.psi();
198 const scalarField& V = mesh_.V();
199 scalarField& Udiag = UEqn.diag();
200 vectorField& Usource = UEqn.source();
201
202 scalar rho = 1.0;
203 if (UEqn.dimensions() == dimForce)
204 {
205 coeffs_.readEntry("rhoRef", rho);
206 }
207
208 apply(Udiag, Usource, V, U, rho);
209}
210
211
213(
215 const volScalarField&,
216 const volScalarField&
217) const
218{
219 const vectorField& U = UEqn.psi();
220 const scalarField& V = mesh_.V();
221 scalarField& Udiag = UEqn.diag();
222 vectorField& Usource = UEqn.source();
223
224 scalar rho = 1.0;
225 if (UEqn.dimensions() == dimForce)
226 {
227 coeffs_.readEntry("rhoRef", rho);
228 }
229
230 apply(Udiag, Usource, V, U, rho);
231}
232
233
235(
236 const fvVectorMatrix& UEqn,
238) const
239{
240 const vectorField& U = UEqn.psi();
241
242 scalar rho = 1.0;
243 if (UEqn.dimensions() == dimForce)
244 {
245 coeffs_.readEntry("rhoRef", rho);
246 }
247
248 apply(AU, U, rho);
249}
250
251
253{
254 dict_.writeEntry(name_, os);
255
256 return true;
257}
258
259
260// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
const Cmpt & xx() const
Definition: TensorI.H:153
const Cmpt & zz() const
Definition: TensorI.H:209
const Cmpt & yy() const
Definition: TensorI.H:181
A List with indirect addressing. Like IndirectList but does not store addressing.
Definition: IndirectList.H:79
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type get(const label i) const
Definition: UList.H:528
virtual void apply()=0
Apply bins.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
Definition: polyMesh.H:504
Top level model for porosity models.
Definition: porosityModel.H:60
const fvMesh & mesh_
Reference to the mesh database.
Definition: porosityModel.H:69
void adjustNegativeResistance(dimensionedVector &resist)
Adjust negative resistance values to be multiplier of max value.
Definition: porosityModel.C:43
labelList cellZoneIDs_
Cell zone IDs.
Definition: porosityModel.H:84
label fieldIndex(const label index) const
Return label index.
Fixed coefficient form of porosity model.
Definition: fixedCoeff.H:64
virtual void calcForce(const volVectorField &U, const volScalarField &rho, const volScalarField &mu, vectorField &force) const
Calculate the porosity force.
Definition: fixedCoeff.C:174
bool writeData(Ostream &os) const
Write.
Definition: fixedCoeff.C:252
virtual void calcTransformModelData()
Transform the model data wrt mesh changes.
Definition: fixedCoeff.C:131
Tensor of scalars, i.e. Tensor<scalar>.
A wordRe is a Foam::word, but can contain a regular expression for matching words or strings.
Definition: wordRe.H:83
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
U
Definition: pEqn.H:72
fvVectorMatrix & UEqn
Definition: UEqn.H:13
const volScalarField & mu
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const cellShapeList & cells
Namespace for OpenFOAM.
const dimensionSet dimless
Dimensionless.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
List< label > labelList
A List of labels.
Definition: List.H:66
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
static const Identity< scalar > I
Definition: Identity.H:94
static void apply(bitSet &selection, const Detail::parcelSelection::actionType action, const Predicate &accept, const UList< Type > &list, const AccessOp &aop)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionSet dimForce
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
volScalarField & alpha
dictionary dict
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333