36Foam::supersonicFreestreamFvPatchVectorField::
37supersonicFreestreamFvPatchVectorField
43 mixedFvPatchVectorField(
p, iF),
46 psiName_(
"thermo:psi"),
52 refValue() = patchInternalField();
58Foam::supersonicFreestreamFvPatchVectorField::
59supersonicFreestreamFvPatchVectorField
66 mixedFvPatchVectorField(
p, iF),
67 TName_(
dict.getOrDefault<
word>(
"T",
"T")),
68 pName_(
dict.getOrDefault<
word>(
"p",
"p")),
69 psiName_(
dict.getOrDefault<
word>(
"psi",
"thermo:psi")),
71 pInf_(
dict.get<scalar>(
"pInf")),
72 TInf_(
dict.get<scalar>(
"TInf")),
73 gamma_(
dict.get<scalar>(
"gamma"))
96 <<
" unphysical pInf specified (pInf <= 0.0)"
97 <<
"\n on patch " << this->patch().name()
98 <<
" of field " << this->internalField().name()
99 <<
" in file " << this->internalField().objectPath()
105Foam::supersonicFreestreamFvPatchVectorField::
106supersonicFreestreamFvPatchVectorField
114 mixedFvPatchVectorField(ptf,
p, iF, mapper),
117 psiName_(ptf.psiName_),
125Foam::supersonicFreestreamFvPatchVectorField::
126supersonicFreestreamFvPatchVectorField
131 mixedFvPatchVectorField(sfspvf),
132 TName_(sfspvf.TName_),
133 pName_(sfspvf.pName_),
134 psiName_(sfspvf.psiName_),
138 gamma_(sfspvf.gamma_)
142Foam::supersonicFreestreamFvPatchVectorField::
143supersonicFreestreamFvPatchVectorField
149 mixedFvPatchVectorField(sfspvf, iF),
150 TName_(sfspvf.TName_),
151 pName_(sfspvf.pName_),
152 psiName_(sfspvf.psiName_),
156 gamma_(sfspvf.gamma_)
164 if (!size() || updated())
180 scalar
R = 1.0/(ppsi[0]*pT[0]);
182 scalar MachInf =
mag(UInf_)/
sqrt(gamma_*
R*TInf_);
187 <<
"\n on patch " << this->patch().name()
188 <<
" of field " << this->internalField().name()
189 <<
" in file " << this->internalField().objectPath()
222 sqrt((gamma_ + 1)/(gamma_ - 1))
223 *
atan(
sqrt((gamma_ - 1)/(gamma_ + 1)*(
sqr(MachInf) - 1)))
232 if (pp[facei] >= pInf_)
239 /(gamma_*
sqr(MachInf))*
mag(Ut[facei])*
log(pp[facei]/pInf_);
241 Up[facei] = Ut[facei] + fpp*nHatInf[facei];
244 scalar Mach =
mag(Up[facei])/
sqrt(gamma_/ppsi[facei]);
250 Up[facei] =
U[facei];
251 valueFraction()[facei] = 0;
262 (2/(gamma_ - 1))*(1 + ((gamma_ - 1)/2)*
sqr(MachInf))
263 *
pow(pp[facei]/pInf_, (1 - gamma_)/gamma_)
273 sqrt((gamma_ + 1)/(gamma_ - 1))
274 *
atan(
sqrt((gamma_ - 1)/(gamma_ + 1)*(
sqr(Mach) - 1)))
277 scalar fpp = (nuMachInf - nuMachf)*
mag(Ut[facei]);
279 Up[facei] = Ut[facei] + fpp*nHatInf[facei];
284 <<
"unphysical subsonic inflow has been generated"
285 <<
"\n on patch " << this->patch().name()
286 <<
" of field " << this->internalField().name()
288 << this->internalField().objectPath()
294 mixedFvPatchVectorField::updateCoeffs();
308 writeEntry(
"value",
os);
#define R(A, B, C, D, E, F, K, M)
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
void size(const label n)
Older name for setAddressableSize.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void operator=(const UList< Type > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Lookup type of boundary radiation properties.
This boundary condition provides a supersonic free-stream condition.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
static const word null
An empty word.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
dimensionedScalar atan(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define forAll(list, i)
Loop across all elements in list.