Go to the documentation of this file.
44 List<T> newElems(newToOld.size(), nullValue);
48 const label oldI = newToOld[i];
52 newElems[i] = elems[oldI];
63 const bitSet& isMasterElem,
70 <<
"Number of elements in list " <<
values.size()
71 <<
" does not correspond to number of elements in isMasterElem "
72 << isMasterElem.
size()
81 if (isMasterElem.
test(i))
112 const label nBFaces = mesh_.nBoundaryFaces();
114 if (faceData.
size() != nBFaces || syncedFaceData.
size() != nBFaces)
117 <<
"Boundary faces:" << nBFaces
118 <<
" faceData:" << faceData.
size()
119 <<
" syncedFaceData:" << syncedFaceData.
size()
129 label bFacei = pp.
start() - mesh_.nInternalFaces();
133 const T&
data = faceData[bFacei];
134 const T& syncData = syncedFaceData[bFacei];
136 if (
mag(
data - syncData) > tol)
138 label facei = pp.
start()+i;
144 <<
" fc:" << mesh_.faceCentres()[facei]
145 <<
" patch:" << pp.
name()
146 <<
" faceData:" <<
data
147 <<
" syncedFaceData:" << syncData
148 <<
" diff:" <<
mag(
data - syncData)
175 Pstream::commsTypes::blocking
184 Pstream::commsTypes::blocking
191 for (
const label allPointi : visitOrder)
199 template<
class GeoField>
200 void Foam::meshRefinement::addPatchFields
203 const word& patchFieldType
208 mesh.objectRegistry::lookupClass<GeoField>()
213 GeoField&
fld = *iter();
214 auto& fldBf =
fld.boundaryFieldRef();
216 label sz = fldBf.size();
232 template<
class GeoField>
233 void Foam::meshRefinement::reorderPatchFields
239 HashTable<GeoField*> flds
241 mesh.objectRegistry::lookupClass<GeoField>()
246 iter()->boundaryFieldRef().reorder(oldToNew);
251 template<
class EnumContainer>
254 const EnumContainer& namedEnum,
260 for (
const word& w : words)
262 flags |= namedEnum[w];
273 const bitSet& isMasterEdge,
283 edges.size() != isMasterEdge.
size()
284 || edges.size() != edgeWeights.size()
285 || meshPoints.size() != pointData.size()
289 <<
"Inconsistent sizes for edge or point data:"
290 <<
" isMasterEdge:" << isMasterEdge.
size()
291 <<
" edgeWeights:" << edgeWeights.size()
292 <<
" edges:" << edges.size()
293 <<
" pointData:" << pointData.size()
294 <<
" meshPoints:" << meshPoints.size()
298 sum.setSize(meshPoints.size());
303 if (isMasterEdge.
test(edgeI))
305 const edge&
e = edges[edgeI];
307 scalar eWeight = edgeWeights[edgeI];
312 sum[v0] += eWeight*pointData[v1];
313 sum[v1] += eWeight*pointData[v0];
343 <<
"Entry '" << keyword <<
"' not found in dictionary "
List< label > labelList
A List of labels.
static Type get(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt=keyType::REGEX, const Type &deflt=Zero)
Wrapper around dictionary::get which does not exit.
static void weightedSum(const polyMesh &mesh, const bitSet &isMasterEdge, const labelList &meshPoints, const edgeList &edges, const scalarField &edgeWeights, const Field< Type > &data, Field< Type > &sum)
Helper: weighted sum (over all subset of mesh points) by.
static T gAverage(const bitSet &isMasterElem, const UList< T > &values)
Helper: calculate average.
Calculates points shared by more than two processor patches or cyclic patches.
A class for handling words, derived from Foam::string.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
static constexpr const zero Zero
Global zero (0)
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool test(const label pos) const
Test value at specified position, never auto-vivify entries.
const fvMesh & mesh() const
Reference to mesh.
Mesh consisting of general polyhedral cells.
static void updateList(const labelList &newToOld, const T &nullValue, List< T > &elems)
Helper: reorder list according to map.
#define forAll(list, i)
Loop across all elements in list.
const fileName & name() const
The dictionary name.
static int readFlags(const EnumContainer &namedEnum, const wordList &words)
Helper: convert wordList into bit pattern using provided Enum.
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
void transfer(List< T > &list)
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
static void collectAndPrint(const UList< point > &points, const UList< T > &data)
Print list according to (collected and) sorted coordinate.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
reduce(hasMovingMesh, orOp< bool >())
void testSyncBoundaryFaceList(const scalar mergeDistance, const string &, const UList< T > &, const UList< T > &) const
Compare two lists over all boundary faces.
Mesh data needed to do the Finite Volume discretisation.
errorManip< error > abort(error &err)
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
label start() const
Return start label of this patch in the polyMesh face list.
errorManipArg< error, int > exit(error &err, const int errNo=1)
A HashTable similar to std::unordered_map.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Traits class for primitives.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label size() const noexcept
Number of entries.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
const dimensionedScalar e
Elementary charge.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const polyBoundaryMesh & patches
void size(const label n) noexcept
Override size to be inconsistent with allocated storage.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
vector point
Point is a vector.
const word & name() const
The patch name.
Database for solution data, solver performance and other reduced data.
option
Enumeration for the data type and search/match modes (bitmask)