meanVelocityForce.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) 2011-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 
26 #include "meanVelocityForce.H"
27 #include "volFields.H"
28 #include "fvMatrices.H"
29 #include "timeIOdictionary.H"
30 #include "IFstream.H"
32 
33 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace fv
38 {
40 
42  (
46  );
47 }
48 }
49 
50 
51 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
52 
53 void Foam::fv::meanVelocityForce::readCoeffs()
54 {
55  UName_ = coeffs().lookupOrDefault<word>("U", "U");
56  Ubar_ = coeffs().lookup<vector>("Ubar");
57  relaxation_ = coeffs().lookupOrDefault<scalar>("relaxation", 1);
58 }
59 
60 
61 void Foam::fv::meanVelocityForce::writeProps
62 (
63  const scalar gradP
64 ) const
65 {
66  // Only write on output time
67  if (mesh().time().writeTime())
68  {
69  timeIOdictionary propsDict
70  (
71  IOobject
72  (
73  name() + "Properties",
74  mesh().time().name(),
75  "uniform",
76  mesh(),
79  )
80  );
81  propsDict.add("gradient", gradP);
82  propsDict.regIOobject::write();
83  }
84 }
85 
86 
87 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
88 
90 (
91  const word& name,
92  const word& modelType,
93  const fvMesh& mesh,
94  const dictionary& dict
95 )
96 :
97  fvConstraint(name, modelType, mesh, dict),
98  set_(mesh, coeffs()),
99  UName_(word::null),
100  Ubar_(vector::uniform(NaN)),
101  relaxation_(NaN),
102  gradP0_(0),
103  dGradP_(0),
104  rAPtr_(nullptr)
105 {
106  readCoeffs();
107 
108  // Read the initial pressure gradient from file if it exists
109  IFstream propsFile
110  (
111  mesh.time().timePath()/"uniform"/(this->name() + "Properties")
112  );
113  if (propsFile.good())
114  {
115  Info<< " Reading pressure gradient from file" << endl;
117  propsDict.lookup("gradient") >> gradP0_;
118  }
119 
120  Info<< " Initial pressure gradient = " << gradP0_ << nl << endl;
121 }
122 
123 
124 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
125 
127 {
128  return wordList(1, UName_);
129 }
130 
131 
132 Foam::scalar Foam::fv::meanVelocityForce::magUbarAve
133 (
134  const volVectorField& U
135 ) const
136 {
137  const labelUList cells = set_.cells();
138  const scalarField& cv = mesh().V();
139 
140  scalar magUbarAve = 0;
141  forAll(cells, i)
142  {
143  const label celli = cells[i];
144  magUbarAve += (normalised(Ubar_) & U[celli])*cv[celli];
145  }
146  reduce(magUbarAve, sumOp<scalar>());
147  magUbarAve /= set_.V();
148 
149  return magUbarAve;
150 }
151 
152 
154 (
155  fvMatrix<vector>& eqn,
156  const word& fieldName
157 ) const
158 {
160  (
161  IOobject
162  (
163  name() + fieldName + "Sup",
164  mesh().time().name(),
165  mesh(),
168  ),
169  mesh(),
171  );
172 
173  const scalar gradP = gradP0_ + dGradP_;
174 
175  UIndirectList<vector>(Su, set_.cells()) = normalised(Ubar_)*gradP;
176 
177  eqn -= Su;
178 
179  if (rAPtr_.empty())
180  {
181  rAPtr_.reset
182  (
183  new volScalarField
184  (
185  IOobject
186  (
187  name() + ":rA",
188  mesh().time().name(),
189  mesh(),
192  false
193  ),
194  1/eqn.A()
195  )
196  );
197  }
198  else
199  {
200  rAPtr_() = 1/eqn.A();
201  }
202 
203  gradP0_ += dGradP_;
204  dGradP_ = 0;
205 
206  return true;
207 }
208 
209 
211 {
212  const scalarField& rAU = rAPtr_();
213 
214  const labelUList cells = set_.cells();
215  const scalarField& cv = mesh().V();
216 
217  // Average rAU over the cell set
218  scalar rAUave = 0;
219  forAll(cells, i)
220  {
221  const label celli = cells[i];
222  rAUave += rAU[celli]*cv[celli];
223  }
224  reduce(rAUave, sumOp<scalar>());
225  rAUave /= set_.V();
226 
227  const scalar magUbarAve = this->magUbarAve(U);
228 
229  // Calculate the pressure gradient increment needed to adjust the average
230  // flow-rate to the desired value
231  dGradP_ = relaxation_*(mag(Ubar_) - magUbarAve)/rAUave;
232 
233  // Apply correction to velocity field
234  forAll(cells, i)
235  {
236  label celli = cells[i];
237  U[celli] += normalised(Ubar_)*rAU[celli]*dGradP_;
238  }
239 
240  const scalar gradP = gradP0_ + dGradP_;
241 
242  Info<< "Pressure gradient source: uncorrected Ubar = " << magUbarAve
243  << ", pressure gradient = " << gradP << endl;
244 
245  writeProps(gradP);
246 
247  return true;
248 }
249 
250 
252 {
253  set_.movePoints();
254  return true;
255 }
256 
257 
259 {
260  set_.topoChange(map);
261 }
262 
263 
265 {
266  set_.mapMesh(map);
267 }
268 
269 
271 {
272  set_.distribute(map);
273 }
274 
275 
277 {
279  {
280  set_.read(coeffs());
281  return true;
282  }
283  else
284  {
285  return false;
286  }
287 }
288 
289 
290 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
Input from file stream.
Definition: IFstream.H:85
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:330
fileName timePath() const
Return current time path.
Definition: Time.H:281
A List with indirect addressing.
Definition: UIndirectList.H:60
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
static const dictionary null
Null dictionary.
Definition: dictionary.H:258
Finite volume options abstract base class.
Definition: fvConstraint.H:57
const dictionary & coeffs() const
Return dictionary.
Definition: fvConstraintI.H:40
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvConstraint.C:166
const fvMesh & mesh() const
Return const access to the mesh database.
Definition: fvConstraintI.H:34
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:808
const dimensionSet & dimensions() const
Definition: fvMatrix.H:302
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:418
Calculates and applies the force necessary to maintain the specified mean velocity.
virtual bool movePoints()
Update for mesh motion.
virtual bool constrain(fvMatrix< vector > &eqn, const word &fieldName) const
Add the momentum source and set the 1/A coefficient.
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
virtual bool read(const dictionary &dict)
Read source dictionary.
meanVelocityForce(const word &sourceName, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from explicit source name and mesh.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
virtual wordList constrainedFields() const
Return the list of fields constrained by the fvConstraint.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
A class for handling words, derived from string.
Definition: word.H:62
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const cellShapeList & cells
volScalarField rAU(1.0/UEqn.A())
U
Definition: pEqn.H:72
addToRunTimeSelectionTable(fvConstraint, bound, dictionary)
defineTypeNameAndDebug(bound, 0)
tmp< VolField< Type > > Su(const VolField< Type > &su, const VolField< Type > &vf)
Definition: fvcSup.C:44
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
List< word > wordList
A List of words.
Definition: fileName.H:54
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
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
messageStream Info
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
const dimensionSet dimVolume
dimensionSet normalised(const dimensionSet &)
Definition: dimensionSet.C:510
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
static const char nl
Definition: Ostream.H:266
labelList fv(nPoints)
dictionary dict
dimensionedVector gradP("gradP", dimensionSet(0, 1, -2, 0, 0), Zero)
IOdictionary propsDict(systemDict("particleTracksDict", args, runTime))