meanVelocityForce.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 | Copyright (C) 2011-2015 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 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
53 (
54  const scalar gradP
55 ) const
56 {
57  // Only write on output time
58  if (mesh_.time().outputTime())
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_(NULL)
95 {
96  coeffs_.lookup("fieldNames") >> fieldNames_;
97 
98  if (fieldNames_.size() != 1)
99  {
101  (
102  "Foam::fv::meanVelocityForce::"
103  "meanVelocityForce"
104  "("
105  "const word&, "
106  "const word&, "
107  "const dictionary&, "
108  "const fvMesh&"
109  ")"
110  ) << "Source can only be applied to a single field. Current "
111  << "settings are:" << fieldNames_ << exit(FatalError);
112  }
113 
114  applied_.setSize(fieldNames_.size(), false);
115 
116  // Read the initial pressure gradient from file if it exists
117  IFstream propsFile
118  (
119  mesh_.time().timePath()/"uniform"/(name_ + "Properties")
120  );
121 
122  if (propsFile.good())
123  {
124  Info<< " Reading pressure gradient from file" << endl;
126  propsDict.lookup("gradient") >> gradP0_;
127  }
128 
129  Info<< " Initial pressure gradient = " << gradP0_ << nl << endl;
130 }
131 
132 
133 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
134 
136 (
137  const volVectorField& U
138 ) const
139 {
140  scalar magUbarAve = 0.0;
141 
142  const scalarField& cv = mesh_.V();
143  forAll(cells_, i)
144  {
145  label cellI = cells_[i];
146  scalar volCell = cv[cellI];
147  magUbarAve += (flowDir_ & U[cellI])*volCell;
148  }
149 
150  reduce(magUbarAve, sumOp<scalar>());
151 
152  magUbarAve /= V_;
153 
154  return magUbarAve;
155 }
156 
157 
159 {
160  const scalarField& rAU = rAPtr_().internalField();
161 
162  // Integrate flow variables over cell set
163  scalar rAUave = 0.0;
164  const scalarField& cv = mesh_.V();
165  forAll(cells_, i)
166  {
167  label cellI = cells_[i];
168  scalar volCell = cv[cellI];
169  rAUave += rAU[cellI]*volCell;
170  }
171 
172  // Collect across all processors
173  reduce(rAUave, sumOp<scalar>());
174 
175  // Volume averages
176  rAUave /= V_;
177 
178  scalar magUbarAve = this->magUbarAve(U);
179 
180  // Calculate the pressure gradient increment needed to adjust the average
181  // flow-rate to the desired value
182  dGradP_ = relaxation_*(mag(Ubar_) - magUbarAve)/rAUave;
183 
184  // Apply correction to velocity field
185  forAll(cells_, i)
186  {
187  label cellI = cells_[i];
188  U[cellI] += flowDir_*rAU[cellI]*dGradP_;
189  }
190 
191  scalar gradP = gradP0_ + dGradP_;
192 
193  Info<< "Pressure gradient source: uncorrected Ubar = " << magUbarAve
194  << ", pressure gradient = " << gradP << endl;
195 
196  writeProps(gradP);
197 }
198 
199 
201 (
202  fvMatrix<vector>& eqn,
203  const label fieldI
204 )
205 {
207  (
208  IOobject
209  (
210  name_ + fieldNames_[fieldI] + "Sup",
211  mesh_.time().timeName(),
212  mesh_,
215  ),
216  mesh_,
218  );
219 
220  scalar gradP = gradP0_ + dGradP_;
221 
222  UIndirectList<vector>(Su, cells_) = flowDir_*gradP;
223 
224  eqn += Su;
225 }
226 
227 
229 (
230  const volScalarField& rho,
231  fvMatrix<vector>& eqn,
232  const label fieldI
233 )
234 {
235  this->addSup(eqn, fieldI);
236 }
237 
238 
240 (
241  fvMatrix<vector>& eqn,
242  const label
243 )
244 {
245  if (rAPtr_.empty())
246  {
247  rAPtr_.reset
248  (
249  new volScalarField
250  (
251  IOobject
252  (
253  name_ + ":rA",
254  mesh_.time().timeName(),
255  mesh_,
258  ),
259  1.0/eqn.A()
260  )
261  );
262  }
263  else
264  {
265  rAPtr_() = 1.0/eqn.A();
266  }
267 
268  gradP0_ += dGradP_;
269  dGradP_ = 0.0;
270 }
271 
272 
273 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
defineTypeNameAndDebug(cellSetOption, 0)
Cell-set options abtract base class. Provides a base set of controls, e.g.
Definition: cellSetOption.H:71
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionSet & dimensions() const
Definition: fvMatrix.H:286
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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:68
messageStream Info
virtual void addSup(fvMatrix< vector > &eqn, const label fieldI)
Add explicit contribution to momentum equation.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Input from file stream.
Definition: IFstream.H:81
Namespace for OpenFOAM.
virtual scalar magUbarAve(const volVectorField &U) const
Calculate and return the magnitude of the mean velocity.
static const dictionary null
Null dictionary.
Definition: dictionary.H:193
tmp< GeometricField< Type, fvPatchField, volMesh > > Su(const GeometricField< Type, fvPatchField, volMesh > &su, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:44
virtual void constrain(fvMatrix< vector > &eqn, const label fieldI)
Set 1/A coefficient.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
A List with indirect addressing.
Definition: fvMatrix.H:106
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
IOdictionary propsDict(IOobject( "particleTrackProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED ))
#define forAll(list, i)
Definition: UList.H:421
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
Macros for easy insertion into run-time selection tables.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
virtual void correct(volVectorField &U)
Correct the pressure gradient.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:739
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:740
volScalarField rAU(1.0/UEqn.A())
void writeProps(const scalar gradP) const
Write the pressure gradient to file (for restarts etc)
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
labelList fv(nPoints)
static const Vector zero
Definition: Vector.H:80
A special matrix type and solver, designed for finite volume solutions of scalar equations.