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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2018-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 
29 #include "meanVelocityForce.H"
30 #include "fvMatrices.H"
31 #include "DimensionedField.H"
32 #include "IFstream.H"
34 
35 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace fv
40 {
41  defineTypeNameAndDebug(meanVelocityForce, 0);
42  addToRunTimeSelectionTable(option, meanVelocityForce, dictionary);
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
48 
50 (
51  const scalar gradP
52 ) const
53 {
54  // Only write on output time
55  if (mesh_.time().writeTime())
56  {
58  (
59  IOobject
60  (
61  name_ + "Properties",
62  mesh_.time().timeName(),
63  "uniform",
64  mesh_,
65  IOobject::NO_READ,
66  IOobject::NO_WRITE
67  )
68  );
69  propsDict.add("gradient", gradP);
70  propsDict.regIOobject::write();
71  }
72 }
73 
74 
75 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
76 
78 (
79  const word& sourceName,
80  const word& modelType,
81  const dictionary& dict,
82  const fvMesh& mesh
83 )
84 :
85  fv::cellSetOption(sourceName, modelType, dict, mesh),
86  Ubar_(coeffs_.get<vector>("Ubar")),
87  gradP0_(0.0),
88  dGradP_(0.0),
89  flowDir_(Ubar_/mag(Ubar_)),
90  relaxation_(coeffs_.getOrDefault<scalar>("relaxation", 1)),
91  rAPtr_(nullptr)
92 {
93  coeffs_.readEntry("fields", fieldNames_);
94 
95  if (fieldNames_.size() != 1)
96  {
98  << "settings are:" << fieldNames_ << exit(FatalError);
99  }
100 
101  fv::option::resetApplied();
102 
103  // Read the initial pressure gradient from file if it exists
104  IFstream propsFile
105  (
106  mesh_.time().timePath()/"uniform"/(name_ + "Properties")
107  );
108 
109  if (propsFile.good())
110  {
111  Info<< " Reading pressure gradient from file" << endl;
112  dictionary propsDict(propsFile);
113  propsDict.readEntry("gradient", gradP0_);
114  }
115 
116  Info<< " Initial pressure gradient = " << gradP0_ << nl << endl;
117 }
118 
119 
120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
121 
123 (
124  const volVectorField& U
125 ) const
126 {
127  scalar magUbarAve = 0.0;
128 
129  const scalarField& cv = mesh_.V();
130  forAll(cells_, i)
131  {
132  label celli = cells_[i];
133  scalar volCell = cv[celli];
134  magUbarAve += (flowDir_ & U[celli])*volCell;
135  }
136 
137  reduce(magUbarAve, sumOp<scalar>());
138 
139  magUbarAve /= V_;
140 
141  return magUbarAve;
142 }
143 
144 
146 {
147  const scalarField& rAU = rAPtr_();
148 
149  // Integrate flow variables over cell set
150  scalar rAUave = 0.0;
151  const scalarField& cv = mesh_.V();
152  forAll(cells_, i)
153  {
154  label celli = cells_[i];
155  scalar volCell = cv[celli];
156  rAUave += rAU[celli]*volCell;
157  }
158 
159  // Collect across all processors
160  reduce(rAUave, sumOp<scalar>());
161 
162  // Volume averages
163  rAUave /= V_;
164 
165  scalar magUbarAve = this->magUbarAve(U);
166 
167  // Calculate the pressure gradient increment needed to adjust the average
168  // flow-rate to the desired value
169  dGradP_ = relaxation_*(mag(Ubar_) - magUbarAve)/rAUave;
170 
171  // Apply correction to velocity field
172  forAll(cells_, i)
173  {
174  label celli = cells_[i];
175  U[celli] += flowDir_*rAU[celli]*dGradP_;
176  }
177 
178  U.correctBoundaryConditions();
179 
180  scalar gradP = gradP0_ + dGradP_;
181 
182  Info<< "Pressure gradient source: uncorrected Ubar = " << magUbarAve
183  << ", pressure gradient = " << gradP << endl;
184 
185  writeProps(gradP);
186 }
187 
188 
190 (
191  fvMatrix<vector>& eqn,
192  const label fieldi
193 )
194 {
196  (
197  IOobject
198  (
199  name_ + fieldNames_[fieldi] + "Sup",
200  mesh_.time().timeName(),
201  mesh_,
204  ),
205  mesh_,
207  );
208 
209  scalar gradP = gradP0_ + dGradP_;
210 
211  UIndirectList<vector>(Su, cells_) = flowDir_*gradP;
212 
213  eqn += Su;
214 }
215 
216 
218 (
219  const volScalarField& rho,
220  fvMatrix<vector>& eqn,
221  const label fieldi
222 )
223 {
224  this->addSup(eqn, fieldi);
225 }
226 
227 
229 (
230  fvMatrix<vector>& eqn,
231  const label
232 )
233 {
234  if (!rAPtr_)
235  {
236  rAPtr_.reset
237  (
238  new volScalarField
239  (
240  IOobject
241  (
242  name_ + ":rA",
243  mesh_.time().timeName(),
244  mesh_,
247  ),
248  1.0/eqn.A()
249  )
250  );
251  }
252  else
253  {
254  rAPtr_() = 1.0/eqn.A();
255  }
256 
257  gradP0_ += dGradP_;
258  dGradP_ = 0.0;
259 }
260 
261 
263 {
265 
266  return false;
267 }
268 
269 
270 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::IOdictionary
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:54
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
gradP
volVectorField gradP(fvc::grad(p))
rAU
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::fv::cellSetOption
Intermediate abstract class for handling cell-set options for the derived fvOptions.
Definition: cellSetOption.H:163
Foam::fv::meanVelocityForce::meanVelocityForce
meanVelocityForce(const word &sourceName, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from explicit source name and mesh.
Definition: meanVelocityForce.C:78
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::IFstream
Input from file stream, using an ISstream.
Definition: IFstream.H:53
DimensionedField.H
Foam::fv::cellSetOption::V_
scalar V_
Sum of cell volumes.
Definition: cellSetOption.H:207
Foam::fv::meanVelocityForce::rAPtr_
autoPtr< volScalarField > rAPtr_
Matrix 1/A coefficients field pointer.
Definition: meanVelocityForce.H:172
Foam::fvMatrix::dimensions
const dimensionSet & dimensions() const
Definition: fvMatrix.H:440
Foam::fv::meanVelocityForce::relaxation_
scalar relaxation_
Relaxation factor.
Definition: meanVelocityForce.H:169
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::fv::meanVelocityForce::writeProps
void writeProps(const scalar gradP) const
Write the pressure gradient to file (for restarts etc)
Definition: meanVelocityForce.C:50
rho
rho
Definition: readInitialConditions.H:88
Foam::fv::option::mesh_
const fvMesh & mesh_
Reference to the mesh database.
Definition: fvOption.H:139
Foam::sumOp
Definition: ops.H:213
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::IOstream::good
bool good() const noexcept
True if next operation might succeed.
Definition: IOstream.H:233
Foam::dimensionedVector
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Definition: dimensionedVector.H:50
fvMatrices.H
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Su
zeroField Su
Definition: alphaSuSp.H:1
NotImplemented
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:517
Foam::fv::meanVelocityForce::flowDir_
vector flowDir_
Flow direction.
Definition: meanVelocityForce.H:166
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Definition: PstreamReduceOps.H:51
Foam::fv::meanVelocityForce::Ubar_
vector Ubar_
Desired mean velocity.
Definition: meanVelocityForce.H:157
Foam::Field< scalar >
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::fv::meanVelocityForce::read
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: meanVelocityForce.C:262
Foam::fv::meanVelocityForce::dGradP_
scalar dGradP_
Change in pressure gradient.
Definition: meanVelocityForce.H:163
Foam::fv::cellSetOption::cells_
labelList cells_
Set of cells to apply source to.
Definition: cellSetOption.H:204
propsDict
IOdictionary propsDict(IOobject("particleTrackProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED))
Foam::fv::meanVelocityForce::addSup
virtual void addSup(fvMatrix< vector > &eqn, const label fieldi)
Add explicit contribution to momentum equation.
Definition: meanVelocityForce.C:190
IFstream.H
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
reduce
reduce(hasMovingMesh, orOp< bool >())
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::fvMatrix::A
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:1394
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
fv
labelList fv(nPoints)
U
U
Definition: pEqn.H:72
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::fv::meanVelocityForce::constrain
virtual void constrain(fvMatrix< vector > &eqn, const label fieldi)
Set 1/A coefficient.
Definition: meanVelocityForce.C:229
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::Vector< scalar >
Foam::fv::meanVelocityForce::magUbarAve
virtual scalar magUbarAve(const volVectorField &U) const
Definition: meanVelocityForce.C:123
Foam::fv::meanVelocityForce::gradP0_
scalar gradP0_
Pressure gradient before correction.
Definition: meanVelocityForce.H:160
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::fvMatrix
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvPatchField.H:68
meanVelocityForce.H
Foam::UIndirectList
A List with indirect addressing.
Definition: faMatrix.H:60
Foam::fv::defineTypeNameAndDebug
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
Foam::fv::addToRunTimeSelectionTable
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
Foam::dimVolume
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Foam::GeometricField< vector, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::fv::meanVelocityForce::correct
virtual void correct(volVectorField &U)
Correct the pressure gradient.
Definition: meanVelocityForce.C:145
Foam::fvMesh::V
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
Definition: fvMeshGeometry.C:179