42 #define PhiDimensionErrorInFunction(phi) \
43 FatalErrorInFunction \
44 << "Incompatible dimensions for " << phi.name() << ": " \
45 << phi.dimensions() << nl \
46 << "Dimensions should be " << dimMass/dimTime << " or " \
47 << dimVolume/dimTime << exit(FatalError)
53 namespace functionObjects
81 PhiPatchFieldTypes[
patchi] =
82 p.boundaryField()[
patchi].fixesValue()
113 Foam::functionObjects::phaseScalarTransport::alphaPhi()
123 << alphaPhiName_ <<
" was not found, so generating it" <<
endl;
131 tmp<surfaceScalarField> tAlphaPhi
145 const word laplacianScheme =
"laplacian(" + pName_ +
")";
151 auto writeDDt = [&](
const label i)
164 Info<<
type() <<
": Writing " << DDtAlpha.name() <<
endl;
167 if (debug && time_.writeTime())
173 nonOrthogonalSolutionControl& control =
174 mesh_.lookupObjectRef<nonOrthogonalSolutionControl>
176 solutionControl::typeName
182 while (control.correctNonOrthogonal())
191 PhiEqn.solve(pName_);
193 if (control.finalNonOrthogonalIter())
195 alphaPhi += PhiEqn.flux();
204 while (control.correctNonOrthogonal())
213 PhiEqn.solve(pName_);
215 if (control.finalNonOrthogonalIter())
217 alphaPhi += PhiEqn.flux();
227 if (debug && time_.writeTime())
237 Foam::functionObjects::phaseScalarTransport::D
252 const word& nameNoPhase = momentumTransportModel::typeName;
289 fieldName_(
dict.lookup(
"field")),
292 residualAlpha_(rootSmall),
311 <<
"Field \"" << fieldName_ <<
"\" does not have a phase extension "
312 <<
"in its name. If it is associated with \"phaseA\" then it "
313 <<
"should be named \"" << fieldName_ <<
".phaseA\"."
345 phiName_ =
dict.lookupOrDefault<
word>(
"phi",
"phi");
352 pName_ =
dict.lookupOrDefault<
word>(
"p",
"p");
353 schemesField_ =
dict.lookupOrDefault<
word>(
"schemesField", fieldName_);
355 constantD_ =
dict.readIfPresent(
"D", D_);
356 alphaD_ =
dict.lookupOrDefault<scalar>(
"alphaD", 1);
357 alphaDt_ =
dict.lookupOrDefault<scalar>(
"alphaDt", 1);
359 dict.readIfPresent(
"nCorr", nCorr_);
360 dict.readIfPresent(
"residualAlpha", residualAlpha_);
361 writeAlphaField_ =
dict.lookupOrDefault<
bool>(
"writeAlphaField",
true);
369 return wordList{alphaName_, alphaPhiName_, phiName_, pName_};
388 const word divScheme =
389 "div(" + alphaPhi.
name() +
"," + schemesField_ +
")";
390 const word laplacianScheme =
391 "laplacian(" +
D.
name() +
"," + schemesField_ +
")";
394 const scalar relaxCoeff =
395 mesh_.solution().relaxEquation(schemesField_)
396 ? mesh_.solution().equationRelaxationFactor(schemesField_)
408 for (
int i=0; i<=nCorr_; i++)
426 fieldEqn.
relax(relaxCoeff);
428 fieldEqn.
solve(schemesField_);
437 for (
int i=0; i<=nCorr_; i++)
455 fieldEqn.
relax(relaxCoeff);
457 fieldEqn.
solve(schemesField_);
467 if (writeAlphaField_)
469 if (!alphaSPtr_.valid())
478 +
word(toupper(fieldName_[0]))
479 + fieldName_(1, fieldName_.size() - 1),
491 alphaSPtr_() =
alpha*s_;
495 if (alphaSPtr_.valid())
497 alphaSPtr_().clear();
static const Foam::dimensionedScalar D("D", Foam::dimTemperature, 257.14)
#define forAll(list, i)
Loop across all elements in list.
Macros for easy insertion into run-time selection tables.
const dimensionSet & dimensions() const
Return dimensions.
static const char *const typeName
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &)
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
word group() const
Return group (extension part of name)
word member() const
Return member (name without the extension)
const word & name() const
Return name.
static word groupName(Name name, const word &group)
virtual Ostream & write(const char)=0
Write character.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
label size() const
Return the number of elements in the UPtrList.
A list of keyword definitions, which are a keyword followed by any number of values (e....
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
const word & name() const
Return const reference to name.
Abstract base-class for Time/database functionObjects.
const Time & time_
Reference to time.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const fvMesh & mesh_
Reference to the fvMesh.
Evolves a passive scalar transport equation within one phase of a multiphase simulation....
virtual wordList fields() const
Return the list of fields required.
virtual ~phaseScalarTransport()
Destructor.
virtual bool execute()
Solve for the evolution of the field.
virtual bool write()
Do nothing. The field is registered and written automatically.
virtual bool read(const dictionary &)
Read the settings from the given dictionary.
phaseScalarTransport(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
virtual bool read(const dictionary &)
Read optional controls.
Finite volume constraints.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
const fvSchemes & schemes() const
Return the fvSchemes.
void setFluxRequired(const word &name) const
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
A class for managing temporary objects.
A class for handling words, derived from string.
static const word null
An empty word.
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
compressibleMomentumTransportModel momentumTransportModel
const char *const group
Group name for atomic constants.
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
const dimensionSet dimViscosity
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
List< word > wordList
A List of words.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
fvMatrix< scalar > fvScalarMatrix
SurfaceField< scalar > surfaceScalarField
const dimensionSet dimLength
const dimensionSet dimTime
const dimensionSet dimVolume
VolField< scalar > volScalarField
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
#define PhiDimensionErrorInFunction(phi)
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))