54 int main(
int argc,
char *argv[])
60 "Name of the pressure field"
66 "Initialise U boundary conditions"
72 "Write the velocity potential field"
78 "Calculate and write the pressure field"
89 #include "createFields.H"
93 Info<<
nl <<
"Calculating potential flow" <<
endl;
97 runTime.functionObjects().start();
112 PhiEqn.setReference(PhiRefCell, PhiRefValue);
117 phi -= PhiEqn.flux();
121 Info<<
"Continuity error = "
122 <<
mag(
fvc::div(phi))().weightedAverage(mesh.V()).value()
126 U.correctBoundaryConditions();
128 Info<<
"Interpolated velocity error = "
148 Info<<
nl <<
"Calculating approximate pressure field" <<
endl;
151 scalar pRefValue = 0.0;
185 pEqn.setReference(pRefCell, pRefValue);
192 runTime.functionObjects().end();
194 Info<<
"ExecutionTime = " << runTime.elapsedCpuTime() <<
" s"
195 <<
" ClockTime = " << runTime.elapsedClockTime() <<
" s"
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
Generic GeometricField class.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &phi) const
Return the given relative flux absolute within the MRF region.
void makeRelative(volVectorField &U) const
Make the given absolute velocity relative within the MRF region.
static void addOption(const word &opt, const string ¶m="", const string &usage="")
Add to an option to validOptions with usage information.
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
bool optionFound(const word &opt) const
Return true if the named option is found.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Non-orthogonal solution control class. Provides non-orthogonal-loop control methods.
virtual bool write(const bool write=true) const
Write using setting from DB.
int main(int argc, char *argv[])
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic,...
Calculate the face-flux of the given field.
Reconstruct volField from a face flux field.
Calculate the matrix for the laplacian of the field.
const dimensionedScalar F
Faraday constant: default SI units: [C/kmol].
tmp< SurfaceField< typename innerProduct< vector, Type >::type > > flux(const VolField< Type > &vf)
Return the face-flux field obtained from the given volVectorField.
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionSet dimless
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
bool setRefCell(const volScalarField &field, const volScalarField &fieldRef, const dictionary &dict, label &refCelli, scalar &refValue, const bool forceReference=false)
If the field fieldRef needs referencing find the reference cell nearest.
const dictionary & potentialFlow(mesh.solution().dict().subDict("potentialFlow"))
Foam::argList args(argc, argv)