Go to the documentation of this file.
39 namespace porosityModels
52 const word& modelType,
57 zoneName_(
name +
":porousZone")
62 vector zDir(coeffs.lookup(
"zDir"));
65 scalar searchSpan(coeffs.get<scalar>(
"searchSpan"));
68 const word topSurfaceFileName(coeffs.get<
word>(
"topSurface"));
72 List<wordRe> groundPatches(coeffs.lookup(
"groundPatches"));
87 mesh.time().constant(),
99 surfBounds.min() - searchSpan*zDir, surfBounds.max()
106 label porousCelli = 0;
110 if (searchBounds.contains(
C[celli]))
112 porousCells[porousCelli++] = celli;
116 porousCells.setSize(porousCelli);
121 forAll(porousCells, porousCelli)
123 start[porousCelli] =
C[porousCells[porousCelli]];
124 end[porousCelli] =
start[porousCelli] + searchSpan*zDir;
138 forAll(porousCells, celli)
144 porousCells[porousCelli] = porousCells[celli];
147 (hit.
hitPoint() -
C[porousCells[porousCelli]]) & zDir;
154 porousCells.setSize(porousCelli);
155 zTop.setSize(porousCelli);
160 triSurfaceTools::triangulate
163 mesh.boundaryMesh().patchSet(groundPatches),
169 if (Pstream::parRun())
174 groundSurfaceProcTris[Pstream::myProcNo()] = groundSurface;
175 groundSurfaceProcPoints[Pstream::myProcNo()] = groundSurface.
points();
177 Pstream::gatherList(groundSurfaceProcTris);
178 Pstream::scatterList(groundSurfaceProcTris);
179 Pstream::gatherList(groundSurfaceProcPoints);
180 Pstream::scatterList(groundSurfaceProcPoints);
183 forAll(groundSurfaceProcTris, i)
185 nTris += groundSurfaceProcTris[i].size();
191 forAll(groundSurfaceProcTris, i)
193 forAll(groundSurfaceProcTris[i], j)
195 groundSurfaceTris[trii] = groundSurfaceProcTris[i][j];
196 groundSurfaceTris[trii][0] += offset;
197 groundSurfaceTris[trii][1] += offset;
198 groundSurfaceTris[trii][2] += offset;
201 offset += groundSurfaceProcPoints[i].size();
205 forAll(groundSurfaceProcPoints, i)
207 nPoints += groundSurfaceProcPoints[i].size();
213 forAll(groundSurfaceProcPoints, i)
215 forAll(groundSurfaceProcPoints[i], j)
217 groundSurfacePoints[pointi++] = groundSurfaceProcPoints[i][j];
221 groundSurface =
triSurface(groundSurfaceTris, groundSurfacePoints);
229 start.setSize(porousCelli);
230 end.setSize(porousCelli);
232 forAll(porousCells, porousCelli)
234 start[porousCelli] =
C[porousCells[porousCelli]];
235 end[porousCelli] =
start[porousCelli] - searchSpan*zDir;
242 forAll(porousCells, porousCelli)
248 zBottom[porousCelli] =
249 (
C[porousCells[porousCelli]] - hit.
hitPoint()) & zDir;
257 Sigma_ = SigmaFunc->value(zNorm);
267 zoneID = cellZones.size();
268 cellZones.setSize(
zoneID + 1);
285 <<
"Unable to create porous cellZone " << zoneName_
286 <<
": zone already exists"
292 Foam::porosityModels::powerLawLopesdaCosta::powerLawLopesdaCosta
295 const word& modelType,
298 const word& dummyCellZoneName
308 powerLawLopesdaCostaZone::zoneName_
310 Cd_(coeffs_.get<scalar>(
"Cd")),
311 C1_(coeffs_.get<scalar>(
"C1")),
312 rhoName_(coeffs_.lookupOrDefault<
word>(
"rho",
"rho"))
414 dict_.writeEntry(name_, os);
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const Field< PointType > & points() const
Return reference to global points.
A class for handling words, derived from Foam::string.
const dimensionedScalar mu
Atomic mass unit.
static constexpr const zero Zero
Global zero.
virtual void calcForce(const volVectorField &U, const volScalarField &rho, const volScalarField &mu, vectorField &force) const
Calculate the porosity force.
bool hit() const
Is there a hit.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
IOoject and searching on triSurface.
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...
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Triangulated surface description with patch information.
const Point & hitPoint() const
Return hit point.
word name(const complex &c)
Return string representation of complex.
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.
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 given a 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)
label ListType::const_reference const label start
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.
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.