33template<
class ZoneType,
class ZoneMesh>
34void Foam::fvMeshDistribute::reorderZones
36 const wordList& zoneNames,
40 zones.clearAddressing();
43 UPtrList<ZoneType> newZonePtrs(zoneNames.size());
46 auto* zonePtr = zones.get(zonei);
52 const label newIndex = zoneNames.find(zonePtr->name());
53 zonePtr->index() = newIndex;
54 newZonePtrs.set(newIndex, zonePtr);
60 if (!newZonePtrs.get(i))
76 zones.swap(newZonePtrs);
80template<
class GeoField>
85 typename GeoField::value_type,
92 mesh.objectRegistry::lookupClass<GeoField>()
97 const GeoField&
fld = *iter();
98 if (!isA<excludeType>(
fld))
100 Pout<<
"Field:" << iter.key() <<
" internalsize:" <<
fld.size()
108template<
class GeoField>
113 mesh.objectRegistry::lookupClass<GeoField>()
118 const GeoField&
fld = *iter();
120 Pout<<
"Field:" << iter.key() <<
" internalsize:" <<
fld.size()
124 for (
const auto& patchFld :
fld.boundaryField())
126 Pout<<
" " << patchFld.patch().index()
127 <<
' ' << patchFld.patch().
name()
128 <<
' ' << patchFld.
type()
129 <<
' ' << patchFld.size()
136template<
class T,
class Mesh>
137void Foam::fvMeshDistribute::saveBoundaryFields
148 mesh_.objectRegistry::lookupClass<
const fldType>()
151 bflds.resize(flds.size());
156 const fldType&
fld = *iter();
158 bflds.set(i,
fld.boundaryField().clone());
165template<
class T,
class Mesh>
166void Foam::fvMeshDistribute::mapBoundaryFields
168 const mapPolyMesh& map,
169 const PtrList<FieldField<fvsPatchField, T>>& oldBflds
174 const labelList& oldPatchStarts = map.oldPatchStarts();
177 typedef GeometricField<T, fvsPatchField, Mesh> fldType;
179 HashTable<fldType*> flds
181 mesh_.objectRegistry::lookupClass<fldType>()
184 if (flds.size() != oldBflds.size())
194 fldType&
fld = *iter();
195 auto& bfld =
fld.boundaryFieldRef();
197 const FieldField<fvsPatchField, T>& oldBfld = oldBflds[fieldi++];
203 fvsPatchField<T>& patchFld = bfld[patchi];
204 label facei = patchFld.patch().start();
208 label oldFacei =
faceMap[facei++];
211 forAll(oldPatchStarts, oldPatchi)
213 label oldLocalI = oldFacei - oldPatchStarts[oldPatchi];
215 if (oldLocalI >= 0 && oldLocalI < oldBfld[oldPatchi].size())
217 patchFld[i] = oldBfld[oldPatchi][oldLocalI];
227void Foam::fvMeshDistribute::saveInternalFields
229 PtrList<Field<T>>& iflds
232 typedef GeometricField<T, fvsPatchField, surfaceMesh> fldType;
234 HashTable<const fldType*> flds
236 mesh_.objectRegistry::lookupClass<
const fldType>()
239 iflds.resize(flds.size());
245 const fldType&
fld = *iter();
247 iflds.set(i,
fld.primitiveField().clone());
255void Foam::fvMeshDistribute::mapExposedFaces
257 const mapPolyMesh& map,
258 const PtrList<Field<T>>& oldFlds
265 typedef GeometricField<T, fvsPatchField, surfaceMesh> fldType;
267 HashTable<fldType*> flds
269 mesh_.objectRegistry::lookupClass<fldType>()
272 if (flds.size() != oldFlds.size())
284 fldType&
fld = *iter();
285 const bool oriented =
fld.oriented()();
287 typename fldType::Boundary& bfld =
fld.boundaryFieldRef();
289 const Field<T>& oldInternal = oldFlds[fieldI++];
295 fvsPatchField<T>& patchFld = bfld[patchi];
299 const label faceI = patchFld.patch().start()+i;
301 label oldFaceI =
faceMap[faceI];
303 if (oldFaceI < oldInternal.size())
305 patchFld[i] = oldInternal[oldFaceI];
307 if (oriented && map.flipFaceFlux().found(faceI))
309 patchFld[i] = flipOp()(patchFld[i]);
318template<
class GeoField,
class PatchFieldType>
319void Foam::fvMeshDistribute::initPatchFields
321 const typename GeoField::value_type& initVal
326 HashTable<GeoField*> flds
328 mesh_.objectRegistry::lookupClass<GeoField>()
333 GeoField&
fld = *iter();
335 auto& bfld =
fld.boundaryFieldRef();
339 if (isA<PatchFieldType>(bfld[patchi]))
341 bfld[patchi] == initVal;
366template<
class GeoField>
367void Foam::fvMeshDistribute::getFieldNames
370 HashTable<wordList>& allFieldNames,
371 const word& excludeType,
375 wordList& list = allFieldNames(GeoField::typeName);
378 if (!excludeType.empty())
380 const wordList& excludeList = allFieldNames(excludeType);
382 DynamicList<word> newList(list.size());
383 for(
const auto&
name : list)
385 if (!excludeList.found(
name))
387 newList.append(
name);
390 if (newList.size() < list.size())
392 list = std::move(newList);
406 globalIndex::gatherOnly{}
416 if (procNames != list)
419 <<
"When checking for equal " << GeoField::typeName
421 <<
"processor0 has:" << list <<
nl
422 <<
"processor" << proci <<
" has:" << procNames <<
nl
423 << GeoField::typeName
424 <<
" need to be synchronised on all processors."
433template<
class GeoField>
434void Foam::fvMeshDistribute::sendFields
437 const HashTable<wordList>& allFieldNames,
438 const fvMeshSubset& subsetter,
466 for (
const word& fieldName : fieldNames)
470 Pout<<
"Subsetting " << GeoField::typeName
471 <<
" field " << fieldName
472 <<
" for domain:" << domain <<
endl;
477 const GeoField&
fld =
478 subsetter.baseMesh().lookupObject<GeoField>(fieldName);
482 tmp<GeoField> tsubfld = subsetter.interpolate(
fld,
true);
493template<
class GeoField>
494void Foam::fvMeshDistribute::receiveFields
497 const HashTable<wordList>& allFieldNames,
499 PtrList<GeoField>&
fields,
500 const dictionary& allFieldsDict
508 const dictionary& fieldDicts =
509 allFieldsDict.subDict(GeoField::typeName);
514 Pout<<
"Receiving:" << GeoField::typeName
515 <<
" fields:" << fieldNames
516 <<
" from domain:" << domain <<
endl;
519 fields.resize(fieldNames.size());
522 for (
const word& fieldName : fieldNames)
526 Pout<<
"Constructing type:" << GeoField::typeName
527 <<
" field:" << fieldName
528 <<
" from domain:" << domain <<
endl;
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
globalIndex procAddr(aMesh.nFaces())
A field of fields is a PtrList of fields with reference counting.
Generic GeometricField class.
A HashTable similar to std::unordered_map.
static const List< word > & null()
Return a null List.
virtual const fileName & name() const
Get the name of the stream.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
static word timeName(const scalar t, const int precision=precision_)
static bool & parRun() noexcept
Test if this a parallel run.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Type type(bool followLink=true, bool checkGzip=false) const
static void printIntFieldInfo(const fvMesh &)
Print some field info.
static void printFieldInfo(const fvMesh &)
Print some field info.
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
labelRange range() const
Return start/size range of local processor data.
static void gather(const labelUList &offsets, const label comm, const ProcIDsContainer &procIDs, const UList< Type > &fld, List< Type > &allFld, const int tag=UPstream::msgType(), const UPstream::commsTypes=UPstream::commsTypes::nonBlocking)
Collect data in processor order on master (== procIDs[0]).
labelRange subProcs() const noexcept
Range of process indices for addressed sub-offsets (processes)
wordList sortedNames() const
The sorted names of all objects.
splitCell * master() const
@ BEGIN_BLOCK
Begin block [isseparator].
@ END_BLOCK
End block [isseparator].
Mesh data needed to do the Finite Volume discretisation.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
List< word > wordList
A List of words.
List< label > labelList
A List of labels.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
#define forAll(list, i)
Loop across all elements in list.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.