Go to the documentation of this file.
55 if (mesh_.time().writeTime())
62 mesh_.time().timeName(),
79 const word& sourceName,
80 const word& modelType,
86 Ubar_(coeffs_.get<
vector>(
"Ubar")),
89 flowDir_(Ubar_/
mag(Ubar_)),
90 relaxation_(coeffs_.getOrDefault<scalar>(
"relaxation", 1)),
93 coeffs_.readEntry(
"fields", fieldNames_);
95 if (fieldNames_.size() != 1)
101 fv::option::resetApplied();
106 mesh_.time().timePath()/
"uniform"/(name_ +
"Properties")
109 if (propsFile.
good())
111 Info<<
" Reading pressure gradient from file" <<
endl;
113 propsDict.readEntry(
"gradient", gradP0_);
116 Info<<
" Initial pressure gradient = " << gradP0_ <<
nl <<
endl;
127 scalar magUbarAve = 0.0;
132 label celli = cells_[i];
133 scalar volCell = cv[celli];
134 magUbarAve += (flowDir_ &
U[celli])*volCell;
155 scalar volCell = cv[celli];
156 rAUave +=
rAU[celli]*volCell;
178 U.correctBoundaryConditions();
182 Info<<
"Pressure gradient source: uncorrected Ubar = " <<
magUbarAve
183 <<
", pressure gradient = " <<
gradP <<
endl;
199 name_ + fieldNames_[fieldi] +
"Sup",
200 mesh_.time().timeName(),
209 scalar
gradP = gradP0_ + dGradP_;
224 this->addSup(eqn, fieldi);
243 mesh_.time().timeName(),
254 rAPtr_() = 1.0/eqn.
A();
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
volVectorField gradP(fvc::grad(p))
tmp< volScalarField > rAU
A class for handling words, derived from Foam::string.
Intermediate abstract class for handling cell-set options for the derived fvOptions.
meanVelocityForce(const word &sourceName, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from explicit source name and mesh.
static constexpr const zero Zero
Global zero (0)
Input from file stream, using an ISstream.
scalar V_
Sum of cell volumes.
autoPtr< volScalarField > rAPtr_
Matrix 1/A coefficients field pointer.
const dimensionSet & dimensions() const
scalar relaxation_
Relaxation factor.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void writeProps(const scalar gradP) const
Write the pressure gradient to file (for restarts etc)
const fvMesh & mesh_
Reference to the mesh database.
#define forAll(list, i)
Loop across all elements in list.
bool good() const noexcept
True if next operation might succeed.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
vector flowDir_
Flow direction.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
vector Ubar_
Desired mean velocity.
messageStream Info
Information stream (stdout output on master, null elsewhere)
virtual bool read(const dictionary &dict)
Read source dictionary.
scalar dGradP_
Change in pressure gradient.
labelList cells_
Set of cells to apply source to.
IOdictionary propsDict(IOobject("particleTrackProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED))
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 a number of values (eg,...
reduce(hasMovingMesh, orOp< bool >())
Macros for easy insertion into run-time selection tables.
tmp< volScalarField > A() const
Return the central coefficient.
Mesh data needed to do the Finite Volume discretisation.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual void constrain(fvMatrix< vector > &eqn, const label fieldi)
Set 1/A coefficient.
virtual scalar magUbarAve(const volVectorField &U) const
scalar gradP0_
Pressure gradient before correction.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A special matrix type and solver, designed for finite volume solutions of scalar equations....
A List with indirect addressing.
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
const dimensionSet dimVolume(pow3(dimLength))
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void correct(volVectorField &U)
Correct the pressure gradient.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.