Go to the documentation of this file.
40 namespace regionModels
48 void Foam::regionModels::regionModel::constructMeshObjects()
52 fvMesh* regionMeshPtr =
64 regionMeshPtr->objectRegistry::store();
69 void Foam::regionModels::regionModel::initialise()
73 Pout<<
"regionModel::initialise()" <<
endl;
76 label nBoundaryFaces = 0;
77 DynamicList<label> primaryPatchIDs;
78 DynamicList<label> intCoupledPatchIDs;
79 const polyBoundaryMesh& rbm = regionMesh().boundaryMesh();
83 const polyPatch& regionPatch = rbm[patchi];
84 if (isA<mappedPatchBase>(regionPatch))
88 Pout<<
"found " << mappedWallPolyPatch::typeName
89 <<
" " << regionPatch.
name() <<
endl;
92 intCoupledPatchIDs.append(patchi);
94 nBoundaryFaces += regionPatch.faceCells().size();
96 const mappedPatchBase& mapPatch =
97 refCast<const mappedPatchBase>(regionPatch);
101 primaryMesh_.time().foundObject<polyMesh>
103 mapPatch.sampleRegion()
108 const label primaryPatchi = mapPatch.samplePolyPatch().index();
109 primaryPatchIDs.append(primaryPatchi);
114 primaryPatchIDs_.transfer(primaryPatchIDs);
115 intCoupledPatchIDs_.transfer(intCoupledPatchIDs);
120 <<
"Region model has no mapped boundary conditions - transfer "
121 <<
"between regions will not be possible" <<
endl;
124 if (!outputPropertiesPtr_)
126 const fileName uniformPath(word(
"uniform")/
"regionModels");
128 outputPropertiesPtr_.reset
134 regionName_ +
"OutputProperties",
136 uniformPath/regionName_,
155 if (
const dictionary* dictptr = findDict(modelName_ +
"Coeffs"))
157 coeffs_ <<= *dictptr;
160 infoOutput_.readIfPresent(
"infoOutput", *
this);
174 if (
const dictionary* dictptr =
dict.findDict(modelName_ +
"Coeffs"))
176 coeffs_ <<= *dictptr;
179 infoOutput_.readIfPresent(
"infoOutput",
dict);
191 const label regionPatchi,
192 const label nbrPatchi,
196 label nbrRegionID = interRegionAMINames_.find(nbrRegion.
name());
200 if (nbrRegionID != -1)
202 if (!interRegionAMI_[nbrRegionID].
set(regionPatchi))
204 const polyPatch&
p = regionMesh().boundaryMesh()[regionPatchi];
210 interRegionAMI_[nbrRegionID].set
215 faceAreaWeightAMI::typeName,
221 interRegionAMI_[nbrRegionID][regionPatchi].calculate(
p, nbrP);
226 return interRegionAMI_[nbrRegionID][regionPatchi];
230 label nbrRegionID = interRegionAMINames_.size();
232 interRegionAMINames_.append(nbrRegion.
name());
234 const polyPatch&
p = regionMesh().boundaryMesh()[regionPatchi];
237 label nPatch = regionMesh().boundaryMesh().size();
240 interRegionAMI_.
resize(nbrRegionID + 1);
251 interRegionAMI_[nbrRegionID].set
256 faceAreaWeightAMI::typeName,
262 interRegionAMI_[nbrRegionID][regionPatchi].calculate(
p, nbrP);
266 return interRegionAMI_[nbrRegionID][regionPatchi];
274 const label regionPatchi
277 label nbrPatchi = -1;
287 if (regionPatchi > pbm.size() - 1)
290 <<
"region patch index out of bounds: "
291 <<
"region patch index = " << regionPatchi
292 <<
", maximum index = " << pbm.size() - 1
298 if (!isA<mappedPatchBase>(pp))
301 <<
"Expected a " << mappedPatchBase::typeName
316 refCast<const mappedPatchBase>(nbrPbm[nbrRegionPatchi]);
320 nbrPatchi = nbrRegionPatchi;
327 const polyPatch&
p = regionMesh().boundaryMesh()[regionPatchi];
330 <<
"Unable to find patch pair for local patch "
331 <<
p.name() <<
" and region " << nbrRegion.
name()
341 Foam::regionModels::regionModel::regionModel
344 const word& regionType
351 regionType +
"Properties",
352 mesh.time().constant(),
364 outputPropertiesPtr_(
nullptr),
366 intCoupledPatchIDs_(),
369 interRegionAMINames_(),
374 Foam::regionModels::regionModel::regionModel
377 const word& regionType,
378 const word& modelName,
386 regionType +
"Properties",
387 mesh.time().constant(),
395 active_(get<Switch>(
"active")),
397 modelName_(modelName),
398 coeffs_(subOrEmptyDict(modelName +
"Coeffs")),
399 outputPropertiesPtr_(
nullptr),
401 intCoupledPatchIDs_(),
402 regionName_(
lookup(
"region")),
403 functions_(*
this, subOrEmptyDict(
"functions"))
407 constructMeshObjects();
418 Foam::regionModels::regionModel::regionModel
421 const word& regionType,
422 const word& modelName,
431 regionType +
"Properties",
432 mesh.time().constant(),
444 modelName_(modelName),
445 coeffs_(
dict.subOrEmptyDict(modelName +
"Coeffs")),
446 outputPropertiesPtr_(
nullptr),
448 intCoupledPatchIDs_(),
449 regionName_(
dict.lookup(
"region")),
450 functions_(*
this, subOrEmptyDict(
"functions"))
454 constructMeshObjects();
471 Info<<
"\nEvolving " << modelName_ <<
" for region "
472 << regionMesh().name() <<
endl;
488 if (time_.writeTime())
490 outputProperties().writeObject
502 functions_.preEvolveRegion();
512 functions_.postEvolveRegion();
int debug
Static debugging option.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
void set(List< bool > &bools, const labelRange &range)
Set the specified range 'on' in a boolList.
A class for handling words, derived from Foam::string.
defineTypeNameAndDebug(KirchhoffShell, 0)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
virtual void evolveRegion()
Evolve the region.
virtual void info()
Provide some feedback.
bool read(const char *buf, int32_t &val)
Same as readInt32.
static word timeName(const scalar t, const int precision=precision_)
virtual bool read()
Read object.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
IOobject(const IOobject &)=default
Copy construct.
Ostream & endl(Ostream &os)
Add newline and flush stream.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
static autoPtr< AMIInterpolation > New(const word &modelName, const dictionary &dict, const bool reverseTarget=false)
Selector for dictionary.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
const word & name() const
Ostream & incrIndent(Ostream &os)
Increment the indent level.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
#define forAll(list, i)
Loop across all elements in list.
Base class for region models.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
virtual void postEvolveRegion()
Post-evolve region.
messageStream Info
Information stream (stdout output on master, null elsewhere)
A patch is a list of labels that address the faces in the global face list.
virtual const fileName & name() const
Get the name of the stream.
const fvMesh & regionMesh() const
Return the region mesh database.
The IOstreamOption is a simple container for options an IOstream can normally have.
virtual void preEvolveRegion()
Pre-evolve region.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
virtual bool read()
Read control parameters from dictionary.
label nbrCoupledPatchID(const regionModel &nbrRegion, const label regionPatchi) const
Return the coupled patch ID paired with coupled patch.
Lookup type of boundary radiation properties.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
errorManip< error > abort(error &err)
virtual const AMIPatchToPatchInterpolation & interRegionAMI(const regionModel &nbrRegion, const label regionPatchi, const label nbrPatchi, const bool flip) const
Create or return a new inter-region AMI object.
void resize(const label newLen)
Adjust size of PtrList.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const Time & time_
Reference to the time database.
static int & msgType() noexcept
Message tag of standard messages.
virtual void evolve()
Main driver routing to evolve the region - calls other evolves.
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
#define WarningInFunction
Report a warning using Foam::Warning.
word regionName_
Region name.
const word & samplePatch() const
Patch (only if NEARESTPATCHFACE)