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-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 
26 #include "meanVelocityForce.H"
27 #include "fvMatrices.H"
28 #include "DimensionedField.H"
29 #include "IFstream.H"
31 
32 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace fv
37 {
38  defineTypeNameAndDebug(meanVelocityForce, 0);
39 
41  (
42  option,
43  meanVelocityForce,
44  dictionary
45  );
46 }
47 }
48 
49 
50 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
51 
53 (
54  const scalar gradP
55 ) const
56 {
57  // Only write on output time
58  if (mesh_.time().writeTime())
59  {
61  (
62  IOobject
63  (
64  name_ + "Properties",
65  mesh_.time().timeName(),
66  "uniform",
67  mesh_,
70  )
71  );
72  propsDict.add("gradient", gradP);
73  propsDict.regIOobject::write();
74  }
75 }
76 
77 
78 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
79 
80 Foam::fv::meanVelocityForce::meanVelocityForce
81 (
82  const word& sourceName,
83  const word& modelType,
84  const dictionary& dict,
85  const fvMesh& mesh
86 )
87 :
88  cellSetOption(sourceName, modelType, dict, mesh),
89  Ubar_(coeffs_.lookup("Ubar")),
90  gradP0_(0.0),
91  dGradP_(0.0),
92  flowDir_(Ubar_/mag(Ubar_)),
93  relaxation_(coeffs_.lookupOrDefault<scalar>("relaxation", 1.0)),
94  rAPtr_(nullptr)
95 {
96  coeffs_.lookup("fields") >> fieldNames_;
97 
98  if (fieldNames_.size() != 1)
99  {
101  << "settings are:" << fieldNames_ << exit(FatalError);
102  }
103 
104  applied_.setSize(fieldNames_.size(), false);
105 
106  // Read the initial pressure gradient from file if it exists
107  IFstream propsFile
108  (
109  mesh_.time().timePath()/"uniform"/(name_ + "Properties")
110  );
111 
112  if (propsFile.good())
113  {
114  Info<< " Reading pressure gradient from file" << endl;
116  propsDict.lookup("gradient") >> gradP0_;
117  }
118 
119  Info<< " Initial pressure gradient = " << gradP0_ << nl << endl;
120 }
121 
122 
123 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
124 
126 (
127  const volVectorField& U
128 ) const
129 {
130  scalar magUbarAve = 0.0;
131 
132  const scalarField& cv = mesh_.V();
133  forAll(cells_, i)
134  {
135  label celli = cells_[i];
136  scalar volCell = cv[celli];
137  magUbarAve += (flowDir_ & U[celli])*volCell;
138  }
139 
140  reduce(magUbarAve, sumOp<scalar>());
141 
142  magUbarAve /= V_;
143 
144  return magUbarAve;
145 }
146 
147 
149 {
150  const scalarField& rAU = rAPtr_();
151 
152  // Integrate flow variables over cell set
153  scalar rAUave = 0.0;
154  const scalarField& cv = mesh_.V();
155  forAll(cells_, i)
156  {
157  label celli = cells_[i];
158  scalar volCell = cv[celli];
159  rAUave += rAU[celli]*volCell;
160  }
161 
162  // Collect across all processors
163  reduce(rAUave, sumOp<scalar>());
164 
165  // Volume averages
166  rAUave /= V_;
167 
168  scalar magUbarAve = this->magUbarAve(U);
169 
170  // Calculate the pressure gradient increment needed to adjust the average
171  // flow-rate to the desired value
172  dGradP_ = relaxation_*(mag(Ubar_) - magUbarAve)/rAUave;
173 
174  // Apply correction to velocity field
175  forAll(cells_, i)
176  {
177  label celli = cells_[i];
178  U[celli] += flowDir_*rAU[celli]*dGradP_;
179  }
180 
181  scalar gradP = gradP0_ + dGradP_;
182 
183  Info<< "Pressure gradient source: uncorrected Ubar = " << magUbarAve
184  << ", pressure gradient = " << gradP << endl;
185 
186  writeProps(gradP);
187 }
188 
189 
191 (
192  fvMatrix<vector>& eqn,
193  const label fieldi
194 )
195 {
197  (
198  IOobject
199  (
200  name_ + fieldNames_[fieldi] + "Sup",
201  mesh_.time().timeName(),
202  mesh_,
205  ),
206  mesh_,
208  );
209 
210  scalar gradP = gradP0_ + dGradP_;
211 
212  UIndirectList<vector>(Su, cells_) = flowDir_*gradP;
213 
214  eqn += Su;
215 }
216 
217 
219 (
220  const volScalarField& rho,
221  fvMatrix<vector>& eqn,
222  const label fieldi
223 )
224 {
225  this->addSup(eqn, fieldi);
226 }
227 
228 
230 (
231  fvMatrix<vector>& eqn,
232  const label
233 )
234 {
235  if (rAPtr_.empty())
236  {
237  rAPtr_.reset
238  (
239  new volScalarField
240  (
241  IOobject
242  (
243  name_ + ":rA",
244  mesh_.time().timeName(),
245  mesh_,
248  ),
249  1.0/eqn.A()
250  )
251  );
252  }
253  else
254  {
255  rAPtr_() = 1.0/eqn.A();
256  }
257 
258  gradP0_ += dGradP_;
259  dGradP_ = 0.0;
260 }
261 
262 
263 // ************************************************************************* //
defineTypeNameAndDebug(option, 0)
virtual scalar magUbarAve(const volVectorField &U) const
Calculate and return the magnitude of the mean velocity.
virtual void correct(volVectorField &U)
Correct the pressure gradient.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
zeroField Su
Definition: alphaSuSp.H:1
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static const dictionary null
Null dictionary.
Definition: dictionary.H:202
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
IOdictionary propsDict(IOobject("particleTrackProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED))
Macros for easy insertion into run-time selection tables.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:814
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
void writeProps(const scalar gradP) const
Write the pressure gradient to file (for restarts etc)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
virtual void addSup(fvMatrix< vector > &eqn, const label fieldi)
Add explicit contribution to momentum equation.
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
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
static const zero Zero
Definition: zero.H:97
static const char nl
Definition: Ostream.H:265
Input from file stream.
Definition: IFstream.H:81
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:737
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A List with indirect addressing.
Definition: fvMatrix.H:106
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A special matrix type and solver, designed for finite volume solutions of scalar equations.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
Cell-set options abtract base class. Provides a base set of controls, e.g.:
Definition: cellSetOption.H:69
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
virtual void constrain(fvMatrix< vector > &eqn, const label fieldi)
Set 1/A coefficient.
const dimensionSet & dimensions() const
Definition: fvMatrix.H:289