41 fixedValueFvPatchScalarField(
p, iF,
dict),
48 supplyTotalTemperature_
52 plenumVolume_(
dict.lookup<scalar>(
"plenumVolume",
dimVolume)),
61 inletDischargeCoefficient_
65 timeScale_(
dict.lookupOrDefault<scalar>(
"timeScale",
dimTime, 0.0)),
66 phiName_(
dict.lookupOrDefault<
word>(
"phi",
"phi")),
67 UName_(
dict.lookupOrDefault<
word>(
"U",
"U"))
69 if (
dict.found(
"rho"))
85 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
88 supplyMassFlowRate_(ptf.supplyMassFlowRate_),
89 supplyTotalTemperature_(ptf.supplyTotalTemperature_),
90 plenumVolume_(ptf.plenumVolume_),
91 plenumDensity_(ptf.plenumDensity_),
92 plenumTemperature_(ptf.plenumTemperature_),
95 inletAreaRatio_(ptf.inletAreaRatio_),
96 inletDischargeCoefficient_(ptf.inletDischargeCoefficient_),
97 timeScale_(ptf.timeScale_),
98 phiName_(ptf.phiName_),
109 fixedValueFvPatchScalarField(tppsf, iF),
110 gamma_(tppsf.gamma_),
112 supplyMassFlowRate_(tppsf.supplyMassFlowRate_),
113 supplyTotalTemperature_(tppsf.supplyTotalTemperature_),
114 plenumVolume_(tppsf.plenumVolume_),
115 plenumDensity_(tppsf.plenumDensity_),
116 plenumTemperature_(tppsf.plenumTemperature_),
118 hasRho_(tppsf.hasRho_),
119 inletAreaRatio_(tppsf.inletAreaRatio_),
120 inletDischargeCoefficient_(tppsf.inletDischargeCoefficient_),
121 timeScale_(tppsf.timeScale_),
122 phiName_(tppsf.phiName_),
141 internalField().name()
149 const scalar dt = db().time().deltaTValue();
153 if (timeIndex_ != db().time().
timeIndex())
155 timeIndex_ = db().time().timeIndex();
156 plenumDensityOld_ = plenumDensity_;
157 plenumTemperatureOld_ = plenumTemperature_;
161 scalar massFlowRate(1.0);
166 massFlowRate = -
gSum(rho_*phi);
171 <<
"The density must be specified when using a volumetric flux."
177 phi.internalField().dimensions()
184 <<
"The density must be not specified when using a mass flux."
189 massFlowRate = -
gSum(phi);
195 <<
"dimensions of phi are not correct"
196 <<
"\n on patch " << patch().name()
197 <<
" of field " << internalField().name()
198 <<
" in file " << internalField().objectPath() <<
nl
203 const scalar cv = R_/(gamma_ - 1),
cp = R_*gamma_/(gamma_ - 1);
208 + (dt/plenumVolume_)*(supplyMassFlowRate_ - massFlowRate);
210 plenumTemperatureOld_
211 + (dt/(plenumDensity_*cv*plenumVolume_))
214 *(
cp*supplyTotalTemperature_ - cv*plenumTemperature_)
215 - massFlowRate*R_*plenumTemperature_
217 const scalar plenumPressure = plenumDensity_*R_*plenumTemperature_;
225 1.0 - (gamma_ - 1.0)*U_e/(2.0*gamma_*R_*plenumTemperature_)
231 (1.0 - r)/(r*r*inletDischargeCoefficient_*inletDischargeCoefficient_)
243 neg(phi)*t*plenumPressure +
pos0(phi)*
max(
p, plenumPressure)
247 const scalar oneByFraction = timeScale_/dt;
248 const scalar fraction = oneByFraction < 1.0 ? 1.0 : 1.0/oneByFraction;
251 operator==((1.0 - fraction)*p_old + fraction*p_new);
252 fixedValueFvPatchScalarField::updateCoeffs();
261 writeEntry(os,
"supplyMassFlowRate", supplyMassFlowRate_);
262 writeEntry(os,
"supplyTotalTemperature", supplyTotalTemperature_);
263 writeEntry(os,
"plenumVolume", plenumVolume_);
264 writeEntry(os,
"plenumDensity", plenumDensity_);
265 writeEntry(os,
"plenumTemperature", plenumTemperature_);
270 writeEntry(os,
"inletAreaRatio", inletAreaRatio_);
274 "inletDischargeCoefficient",
275 inletDischargeCoefficient_
277 writeEntryIfDifferent<scalar>(os,
"timeScale", 0.0, timeScale_);
278 writeEntryIfDifferent<word>(os,
"phi",
"phi", phiName_);
279 writeEntryIfDifferent<word>(os,
"U",
"U", UName_);
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...
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
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....
Abstract base class for field mapping.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void write(Ostream &) const
Write.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
This boundary condition provides a plenum pressure inlet condition. This condition creates a zero-dim...
virtual void write(Ostream &) const
Write.
plenumPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for handling words, derived from string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
errorManipArg< error, int > exit(error &err, const int errNo=1)
Type gSum(const FieldField< Field, Type > &f)
VolField< vector > volVectorField
dimensionedScalar pos0(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
const dimensionSet dimMassFlux
const dimensionSet dimVolumetricFlux
const dimensionSet dimGasConstant
const dimensionSet dimless
SurfaceField< scalar > surfaceScalarField
const dimensionSet dimTemperature
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimTime
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimDensity
const dimensionSet dimVolume
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar neg(const dimensionedScalar &ds)
const dimensionSet dimMass
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const unitConversion unitFraction