41 mixedFvPatchVectorField(
p, iF,
dict, false),
42 TName_(
dict.lookupOrDefault<
word>(
"T",
"T")),
43 pName_(
dict.lookupOrDefault<
word>(
"p",
"p")),
44 psiName_(
dict.lookupOrDefault<
word>(
"psi",
"psi")),
45 UInf_(
dict.lookup(
"UInf")),
46 pInf_(
dict.lookup<scalar>(
"pInf")),
47 TInf_(
dict.lookup<scalar>(
"TInf")),
48 gamma_(
dict.lookup<scalar>(
"gamma"))
50 if (
dict.found(
"value"))
71 ) <<
" unphysical pInf specified (pInf <= 0.0)"
72 <<
"\n on patch " << this->patch().name()
73 <<
" of field " << this->internalField().name()
74 <<
" in file " << this->internalField().objectPath()
89 mixedFvPatchVectorField(ptf,
p, iF, mapper),
92 psiName_(ptf.psiName_),
107 mixedFvPatchVectorField(sfspvf, iF),
108 TName_(sfspvf.TName_),
109 pName_(sfspvf.pName_),
110 psiName_(sfspvf.psiName_),
114 gamma_(sfspvf.gamma_)
122 if (!size() || updated())
138 scalar
R = 1.0/(ppsi[0]*pT[0]);
140 scalar MachInf =
mag(UInf_)/
sqrt(gamma_*
R*TInf_);
145 <<
"\n on patch " << this->patch().name()
146 <<
" of field " << this->internalField().name()
147 <<
" in file " << this->internalField().objectPath()
180 sqrt((gamma_ + 1)/(gamma_ - 1))
181 *
atan(
sqrt((gamma_ - 1)/(gamma_ + 1)*(
sqr(MachInf) - 1)))
190 if (pp[facei] >= pInf_)
197 /(gamma_*
sqr(MachInf))*
mag(Ut[facei])*
log(pp[facei]/pInf_);
199 Up[facei] = Ut[facei] + fpp*nHatInf[facei];
202 scalar Mach =
mag(Up[facei])/
sqrt(gamma_/ppsi[facei]);
208 Up[facei] =
U[facei];
209 valueFraction()[facei] = 0;
220 (2/(gamma_ - 1))*(1 + ((gamma_ - 1)/2)*
sqr(MachInf))
221 *
pow(pp[facei]/pInf_, (1 - gamma_)/gamma_)
231 sqrt((gamma_ + 1)/(gamma_ - 1))
232 *
atan(
sqrt((gamma_ - 1)/(gamma_ + 1)*(
sqr(Mach) - 1)))
235 scalar fpp = (nuMachInf - nuMachf)*
mag(Ut[facei]);
237 Up[facei] = Ut[facei] + fpp*nHatInf[facei];
242 <<
"unphysical subsonic inflow has been generated"
243 <<
"\n on patch " << this->patch().name()
244 <<
" of field " << this->internalField().name()
246 << this->internalField().objectPath()
252 mixedFvPatchVectorField::updateCoeffs();
259 writeEntryIfDifferent<word>(os,
"T",
"T", TName_);
260 writeEntryIfDifferent<word>(os,
"p",
"p", pName_);
261 writeEntryIfDifferent<word>(os,
"psi",
"psi", psiName_);
#define forAll(list, i)
Loop across all elements in list.
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,...
A list of keyword definitions, which are a keyword followed by any number of values (e....
Foam::fvPatchFieldMapper.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void write(Ostream &) const
Write.
virtual void operator=(const UList< Type > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
This boundary condition provides a supersonic free-stream condition.
virtual void write(Ostream &) const
Write.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
supersonicFreestreamFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
A class for managing temporary objects.
A class for handling words, derived from string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensionedScalar sqrt(const dimensionedScalar &ds)
VolField< scalar > volScalarField
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
static scalar R(const scalar a, const scalar x)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dimensionedScalar atan(const dimensionedScalar &ds)