Go to the documentation of this file.
61 if (mesh_.time().writeTime())
68 mesh_.time().timeName(),
83 Foam::fv::meanVelocityForce::meanVelocityForce
85 const word& sourceName,
86 const word& modelType,
92 Ubar_(coeffs_.get<
vector>(
"Ubar")),
95 flowDir_(Ubar_/
mag(Ubar_)),
96 relaxation_(coeffs_.lookupOrDefault<scalar>(
"relaxation", 1.0)),
99 coeffs_.readEntry(
"fields", fieldNames_);
101 if (fieldNames_.size() != 1)
107 applied_.setSize(fieldNames_.size(),
false);
112 mesh_.time().timePath()/
"uniform"/(name_ +
"Properties")
115 if (propsFile.
good())
117 Info<<
" Reading pressure gradient from file" <<
endl;
119 propsDict.readEntry(
"gradient", gradP0_);
122 Info<<
" Initial pressure gradient = " << gradP0_ <<
nl <<
endl;
133 scalar magUbarAve = 0.0;
138 label celli = cells_[i];
139 scalar volCell = cv[celli];
140 magUbarAve += (flowDir_ &
U[celli])*volCell;
161 scalar volCell = cv[celli];
162 rAUave +=
rAU[celli]*volCell;
184 U.correctBoundaryConditions();
188 Info<<
"Pressure gradient source: uncorrected Ubar = " <<
magUbarAve
189 <<
", pressure gradient = " <<
gradP <<
endl;
205 name_ + fieldNames_[fieldi] +
"Sup",
206 mesh_.time().timeName(),
215 scalar
gradP = gradP0_ + dGradP_;
230 this->addSup(eqn, fieldi);
249 mesh_.time().timeName(),
260 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))
A class for handling words, derived from Foam::string.
Cell-set options abtract base class. Provides a base set of controls, e.g.:
static constexpr const zero Zero
Global zero.
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.
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_
Average velocity.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
messageStream Info
Information stream (uses stdout - output is on the master only)
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,...
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)
defineTypeNameAndDebug(option, 0)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
virtual void constrain(fvMatrix< vector > &eqn, const label fieldi)
Set 1/A coefficient.
virtual scalar magUbarAve(const volVectorField &U) const
Calculate and return the magnitude of the mean velocity.
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....
tmp< volScalarField > rAU
A List with indirect addressing.
bool good() const
Return true if next operation might succeed.
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.