34 #include "surfaceInterpolate.H" 52 directionalPressureGradientExplicitSource,
65 Foam::fv::directionalPressureGradientExplicitSource::pressureDropModelNames_
67 { pressureDropModel::pVolumetricFlowRateTable,
"volumetricFlowRateTable" },
68 { pressureDropModel::pConstant,
"constant" },
69 { pressureDropModel::pDarcyForchheimer,
"DarcyForchheimer" },
75 void Foam::fv::directionalPressureGradientExplicitSource::initialise()
80 label numFaces = fZone.
size();
91 const label meshFacei = fZone[i];
95 label facePatchId = -1;
103 if (isA<emptyPolyPatch>(
pp))
108 const auto* cpp = isA<coupledPolyPatch>(
pp);
110 if (cpp && !cpp->owner())
120 faceId_[numFaces] =
faceId;
121 facePatchId_[numFaces] = facePatchId;
130 facePatchId_.
resize(numFaces);
134 void Foam::fv::directionalPressureGradientExplicitSource::writeProps
140 if (mesh_.time().writeTime())
146 name_ +
"Properties",
147 mesh_.time().timeName(),
166 const word& sourceName,
167 const word& modelType,
173 model_(pressureDropModelNames_.
get(
"model", coeffs_)),
174 gradP0_(cells_.size(),
Zero),
175 dGradP_(cells_.size(),
Zero),
176 gradPporous_(cells_.size(),
Zero),
177 flowDir_(coeffs_.
get<
vector>(
"flowDir")),
184 faceZoneName_(coeffs_.
get<
word>(
"faceZone")),
185 zoneID_(mesh_.faceZones().findZoneID(faceZoneName_)),
188 relaxationFactor_(coeffs_.getOrDefault<scalar>(
"relaxationFactor", 0.3)),
189 cellFaceMap_(cells_.size(), -1)
193 flowDir_.normalise();
198 <<
"Source can only be applied to a single field. Current " 205 <<
type() <<
" " << this->
name() <<
": " 206 <<
" Unknown face zone name: " << faceZoneName_
213 flowRate_ = interpolationTable<scalar>(
coeffs_);
228 <<
"Did not find mode " << model_ <<
nl 229 <<
"Please set 'model' to one of " 230 << pressureDropModelNames_
242 if (propsFile.good())
244 Info<<
" Reading pressure gradient from file" <<
endl;
246 propsDict.readEntry(
"gradient", gradP0_);
249 Info<<
" Initial pressure gradient = " << gradP0_ <<
nl <<
endl;
270 case pDarcyForchheimer:
274 const auto& turbModel =
282 gradPporous_ = -flowDir_*(D_*
nu + I_*0.5*magUn)*magUn*length_;
286 const auto& turbModel =
297 - flowDir_*(D_*
mu + I_*0.5*
rho*magUn)*magUn*length_;
303 gradPporous_ = -flowDir_*pressureDrop_;
307 case pVolumetricFlowRateTable:
309 scalar volFlowRate = 0;
314 label faceI = faceId_[i];
315 if (facePatchId_[i] != -1)
317 label patchI = facePatchId_[i];
318 totalphi +=
phi.boundaryField()[patchI][faceI];
322 totalphi +=
phi[faceI];
325 reduce(totalphi, sumOp<scalar>());
329 volFlowRate =
mag(totalphi);
333 const auto& turbModel =
341 volFlowRate =
mag(totalphi)/rhoAve;
344 gradPporous_ = -flowDir_*flowRate_(volFlowRate);
349 const faceZone& fZone = mesh_.faceZones()[zoneID_];
351 labelList meshToLocal(mesh_.nCells(), -1);
354 meshToLocal[cells_[i]] = i;
357 labelList faceToCellIndex(fZone.size(), -1);
358 const labelList& mc = fZone.masterCells();
359 const labelList& sc = fZone.slaveCells();
363 label masterCellI = mc[i];
365 if (meshToLocal[masterCellI] != -1 && masterCellI != -1)
367 faceToCellIndex[i] = meshToLocal[masterCellI];
369 else if (meshToLocal[masterCellI] == -1)
372 <<
"Did not find cell " << masterCellI
373 <<
"in cellZone :" << zoneName()
382 const polyBoundaryMesh&
pbm = mesh_.boundaryMesh();
384 FieldField<Field, vector> upwindValues(
pbm.size());
386 forAll(
U.boundaryField(), patchI)
392 upwindValues.set(patchI, pf.patchNeighbourField());
394 else if (!isA<emptyFvPatchScalarField>(pf))
402 label faceI = fZone[i];
403 label
cellId = faceToCellIndex[i];
407 label sourceCellId = sc[i];
408 if (mesh_.isInternalFace(faceI))
410 scalar w = mesh_.magSf()[faceI];
411 UfCells[
cellId] +=
U[sourceCellId]*w;
412 UfCellWeights[
cellId] += w;
414 else if (fZone.flipMap()[i])
416 const label patchI =
pbm.patchID(faceI);
417 label localFaceI =
pbm[patchI].whichFace(faceI);
419 scalar w = mesh_.magSf().boundaryField()[patchI][localFaceI];
421 if (upwindValues.set(patchI))
423 UfCells[
cellId] += upwindValues[patchI][localFaceI]*w;
424 UfCellWeights[
cellId] += w;
430 UfCells /= UfCellWeights;
434 label cellI = cells_[i];
436 const vector Ufnorm(UfCells[i]/(
mag(UfCells[i]) + SMALL));
443 (
D & UfCells[i]) -
U[cellI]
449 Info<<
"Difference mag(U) = " 450 <<
mag(UfCells[i]) -
mag(
U[cellI])
452 Info<<
"Pressure drop in flowDir direction : " 453 << gradPporous_[i] <<
endl;
454 Info<<
"UfCell:= " << UfCells[i] <<
"U : " <<
U[cellI] <<
endl;
457 writeProps(gradP0_ + dGradP_);
463 fvMatrix<vector>& eqn,
469 name_ + fieldNames_[fieldI] +
"Sup",
474 auto&
Su = tSu.ref();
489 this->addSup(eqn, fieldI);
512 invAPtr_() = 1.0/eqn.
A();
539 coeffs.readEntry(
"flowDir", flowDir_);
540 flowDir_.normalise();
542 if (model_ == pConstant)
544 coeffs.readEntry(
"pressureDrop", pressureDrop_);
546 else if (model_ == pDarcyForchheimer)
548 coeffs.readEntry(
"D", D_);
549 coeffs.readEntry(
"I", I_);
550 coeffs.readEntry(
"length", length_);
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
const polyBoundaryMesh & pbm
void size(const label n)
Older name for setAddressableSize.
fvPatchField< vector > fvPatchVectorField
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
wordList fieldNames_
Field names to apply source to - populated by derived models.
void resize(const label len)
Adjust allocated size of list.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
fileName timePath() const
Return current time path = path/timeName.
const fvMesh & mesh_
Reference to the mesh database.
pressureDropModel
Modes of pressure drop.
constexpr char nl
The newline '\n' character (0x0a)
scalarProduct< Type, Type >::type gSumProd(const UList< Type > &f1, const UList< Type > &f2, const label comm)
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
Ostream & endl(Ostream &os)
Add newline and flush stream.
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
Ignore writing from objectRegistry::writeObject()
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
const Time & time() const
Return the top-level database.
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from vector n1 to n2.
tmp< volScalarField > rAU
Macros for easy insertion into run-time selection tables.
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
IOdictionary propsDict(dictIO)
virtual void correct(volVectorField &U)
Correct the pressure gradient.
#define forAll(list, i)
Loop across all elements in list.
volVectorField gradP(fvc::grad(p))
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, IOobjectOption::readOption readOpt=IOobjectOption::MUST_READ) const
Find entry and assign to T val. FatalIOError if it is found and the number of tokens is incorrect...
Templated abstract base class for single-phase incompressible turbulence models.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionSet dimVolume(pow3(dimLength))
static word scopedName(const std::string &scope, const word &name)
Create scope:name or scope_name string.
Type gSum(const FieldField< Field, Type > &f)
const word name_
Source name.
static const word propertiesName
Default name of the turbulence properties dictionary.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
label size() const noexcept
The number of entries in the list.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
const word & name() const noexcept
Return const access to the source name.
label whichPatch(const label meshFacei) const
Return patch index for a given mesh face index. Uses binary search.
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
virtual void writeData(Ostream &os) const
Write the source properties.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
OBJstream os(runTime.globalPath()/outputName)
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
virtual bool read(const dictionary &dict)
Read source dictionary.
const dimensionedScalar mu
Atomic mass unit.
virtual void constrain(fvMatrix< vector > &eqn, const label fieldI)
Set 1/A coefficient.
tmp< volScalarField > A() const
Return the central coefficient.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
directionalPressureGradientExplicitSource(const word &sourceName, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from explicit source name and mesh.
Mesh data needed to do the Finite Volume discretisation.
A List with indirect addressing. Like IndirectList but does not store addressing. ...
void resetApplied()
Resize/reset applied flag list for all fieldNames_ entries.
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, IOobjectOption::registerOption regOpt, const Mesh &mesh, const dimensionSet &dims, const Field< Type > &iField)
Return tmp field (NO_READ, NO_WRITE) from name, mesh, dimensions, copy of internal field...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
A special matrix type and solver, designed for finite volume solutions of scalar equations.
wordList names() const
A list of the zone names.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< vector > vectorField
Specialisation of Field<T> for vector.
Intermediate abstract class for handling cell-set options for the derived fvOptions.
const dimensionedScalar & D
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T, or return the given default value. FatalIOError if it is found and the number of...
List< label > labelList
A List of labels.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
dictionary coeffs_
Dictionary containing source coefficients.
Do not request registration (bool: false)
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
static constexpr const zero Zero
Global zero (0)
virtual void addSup(fvMatrix< vector > &eqn, const label fieldI)
Add explicit contribution to momentum equation.