Go to the documentation of this file.
40 namespace porosityModels
53 const word& modelType,
58 zoneName_(
name +
":porousZone")
66 scalar searchSpan(coeffs.get<scalar>(
"searchSpan"));
69 const word topSurfaceFileName(coeffs.get<
word>(
"topSurface"));
88 mesh.time().constant(),
100 surfBounds.min() - searchSpan*zDir, surfBounds.max()
107 label porousCelli = 0;
111 if (searchBounds.contains(
C[celli]))
113 porousCells[porousCelli++] = celli;
117 porousCells.setSize(porousCelli);
122 forAll(porousCells, porousCelli)
124 start[porousCelli] =
C[porousCells[porousCelli]];
125 end[porousCelli] = start[porousCelli] + searchSpan*zDir;
139 forAll(porousCells, celli)
145 porousCells[porousCelli] = porousCells[celli];
148 (hit.
hitPoint() -
C[porousCells[porousCelli]]) & zDir;
155 porousCells.setSize(porousCelli);
156 zTop.setSize(porousCelli);
161 triSurfaceTools::triangulate
164 mesh.boundaryMesh().patchSet(groundPatches),
170 if (Pstream::parRun())
175 groundSurfaceProcTris[Pstream::myProcNo()] = groundSurface;
176 groundSurfaceProcPoints[Pstream::myProcNo()] = groundSurface.
points();
178 Pstream::gatherList(groundSurfaceProcTris);
179 Pstream::scatterList(groundSurfaceProcTris);
180 Pstream::gatherList(groundSurfaceProcPoints);
181 Pstream::scatterList(groundSurfaceProcPoints);
184 forAll(groundSurfaceProcTris, i)
186 nTris += groundSurfaceProcTris[i].size();
192 forAll(groundSurfaceProcTris, i)
194 forAll(groundSurfaceProcTris[i], j)
196 groundSurfaceTris[trii] = groundSurfaceProcTris[i][j];
197 groundSurfaceTris[trii][0] += offset;
198 groundSurfaceTris[trii][1] += offset;
199 groundSurfaceTris[trii][2] += offset;
202 offset += groundSurfaceProcPoints[i].size();
206 forAll(groundSurfaceProcPoints, i)
208 nPoints += groundSurfaceProcPoints[i].size();
214 forAll(groundSurfaceProcPoints, i)
216 forAll(groundSurfaceProcPoints[i], j)
218 groundSurfacePoints[pointi++] = groundSurfaceProcPoints[i][j];
222 groundSurface =
triSurface(groundSurfaceTris, groundSurfacePoints);
230 start.setSize(porousCelli);
231 end.setSize(porousCelli);
233 forAll(porousCells, porousCelli)
235 start[porousCelli] =
C[porousCells[porousCelli]];
236 end[porousCelli] = start[porousCelli] - searchSpan*zDir;
243 forAll(porousCells, porousCelli)
249 zBottom[porousCelli] =
250 (
C[porousCells[porousCelli]] - hit.
hitPoint()) & zDir;
258 Sigma_ = SigmaFunc->value(zNorm);
268 zoneID = cellZones.size();
269 cellZones.setSize(
zoneID + 1);
286 <<
"Unable to create porous cellZone " << zoneName_
287 <<
": zone already exists"
293 Foam::porosityModels::powerLawLopesdaCosta::powerLawLopesdaCosta
296 const word& modelType,
299 const word& dummyCellZoneName
309 powerLawLopesdaCostaZone::zoneName_
311 Cd_(coeffs_.get<scalar>(
"Cd")),
312 C1_(coeffs_.get<scalar>(
"C1")),
313 rhoName_(coeffs_.getOrDefault<
word>(
"rho",
"rho"))
415 dict_.writeEntry(name_,
os);
const Field< point_type > & points() const noexcept
Return reference to global points.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A class for handling words, derived from Foam::string.
const dimensionedScalar mu
Atomic mass unit.
static constexpr const zero Zero
Global zero (0)
virtual void calcForce(const volVectorField &U, const volScalarField &rho, const volScalarField &mu, vectorField &force) const
Calculate the porosity force.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
IOoject and searching on triSurface.
const point_type & hitPoint() const
Return hit point. Fatal if not hit.
const dimensionSet dimForce
bool writeData(Ostream &os) const
Write.
Helper class to search on triSurface.
powerLawLopesdaCostaZone(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Constructor.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
#define forAll(list, i)
Loop across all elements in list.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
virtual void correct(fvVectorMatrix &UEqn) const
Add resistance.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
bool hit() const noexcept
Is there a hit?
Triangulated surface description with patch information.
defineTypeNameAndDebug(powerLawLopesdaCosta, 0)
virtual void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Find first intersection on segment from start to end.
virtual ~powerLawLopesdaCosta()
Destructor.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const scalarField & Sigma() const
Return the porosity surface area per unit volume zone field.
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
const labelIOList & zoneID
Mesh data needed to do the Finite Volume discretisation.
virtual tmp< pointField > points() const
Get the points that define the surface.
errorManip< error > abort(error &err)
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Top level model for porosity models.
void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &info) const
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static void apply(bitSet &selection, const Detail::parcelSelection::actionType action, const Predicate &accept, const UList< Type > &list, const AccessOp &aop)
A List of wordRe with additional matching capabilities.
virtual void calcTransformModelData()
Transform the model data wrt mesh changes.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
A bounding box defined in terms of min/max extrema points.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
addToRunTimeSelectionTable(porosityModel, powerLawLopesdaCosta, mesh)
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Graphite solid properties.
Variant of the power law porosity model with spatially varying drag coefficient.