Go to the documentation of this file.
42 namespace functionObjects
58 Foam::functionObjects::dsmcFields::dsmcFields
94 word rhoNMeanName =
"rhoNMean";
95 word rhoMMeanName =
"rhoMMean";
96 word momentumMeanName =
"momentumMean";
97 word linearKEMeanName =
"linearKEMean";
98 word internalEMeanName =
"internalEMean";
99 word iDofMeanName =
"iDofMean";
100 word fDMeanName =
"fDMean";
123 if (
min(
mag(rhoNMean)).value() > VSMALL)
125 Log <<
"Calculating dsmcFields." <<
endl;
127 Log <<
" Calculating UMean field." <<
endl;
133 obr_.time().timeName(),
137 momentumMean/rhoMMean
140 Log <<
" Calculating translationalT field." <<
endl;
146 obr_.time().timeName(),
152 *(linearKEMean - 0.5*rhoMMean*(
UMean &
UMean))
155 Log <<
" Calculating internalT field." <<
endl;
161 obr_.time().timeName(),
168 Log <<
" Calculating overallT field." <<
endl;
174 obr_.time().timeName(),
179 *(linearKEMean - 0.5*rhoMMean*(
UMean &
UMean) + internalEMean)
182 Log <<
" Calculating pressure field." <<
endl;
188 obr_.time().timeName(),
195 volScalarField::Boundary& pBf =
p.boundaryFieldRef();
197 forAll(mesh_.boundaryMesh(), i)
201 if (isA<wallPolyPatch>(
patch))
204 fDMean.boundaryField()[i]
210 Log <<
" mag(UMean) max/min : "
214 <<
" translationalT max/min : "
215 <<
max(translationalT).value() <<
" "
216 <<
min(translationalT).value() <<
nl
218 <<
" internalT max/min : "
219 <<
max(internalT).value() <<
" "
220 <<
min(internalT).value() <<
nl
222 <<
" overallT max/min : "
223 <<
max(overallT).value() <<
" "
224 <<
min(overallT).value() <<
nl
227 <<
max(
p).value() <<
" "
232 translationalT.write();
240 Log <<
"dsmcFields written." <<
nl <<
endl;
246 Log <<
"Small value (" <<
min(
mag(rhoNMean))
247 <<
") found in rhoNMean field. "
248 <<
"Not calculating dsmcFields to avoid division by zero."
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
virtual bool execute()
Do nothing.
const dimensionedScalar k
Boltzmann constant.
Different types of constants.
bool read(const char *buf, int32_t &val)
Same as readInt32.
virtual bool write()
Calculate and write the DSMC fields.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
virtual bool read(const dictionary &)
Read the dsmcFields data.
A patch is a list of labels that address the faces in the global face list.
word name(const complex &c)
Return string representation of complex.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
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.
const std::string patch
OpenFOAM patch number as a std::string.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
volVectorField UMean(UMeanHeader, mesh)
defineTypeNameAndDebug(combustionModel, 0)
virtual ~dsmcFields()
Destructor.