59 template<
class ListType>
63 typename ListType::const_reference t,
67 for (
label i = start; i < l.size(); i++)
111 box =
combine(box,
srcBox(srcFace, srcPoints0, srcPointNormals0));
127 forAll(srcPatch, srcFacei)
133 srcBox(srcPatch, srcPointNormals, srcPointNormals0, srcFacei)
170 forAll(tgtPatch, tgtFacei)
172 box =
combine(box, tgtBox(tgtPatch, tgtFacei));
187 const label srcFacei,
192 auto intersected = []
196 const label otherFacei
199 forAll(faceLocalFaces[facei], i)
201 if (faceLocalFaces[facei][i] == otherFacei)
210 intersected(srcLocalTgtFaces_, srcFacei, tgtFacei)
211 || intersected(tgtLocalSrcFaces_, tgtFacei, srcFacei)
253 label nFaceComplete = 0;
257 const label otherFacei = queue[queuei].second();
262 otherFaceVisited[otherFacei] =
true;
263 label otherPerimeterReached =
false;
264 while (otherCurrentFaces.
size())
267 forAll(otherCurrentFaces, otheri)
269 const label otherFacej = otherCurrentFaces[otheri];
279 isSrc ? facei : otherFacej,
280 isSrc ? otherFacej : facei
300 forAll(otherPatch[otherFacej], otherFacePointj)
302 const label otherPointj =
303 otherLocalFaces[otherFacej][otherFacePointj];
304 forAll(otherPointFaces[otherPointj], otherPointFacej)
306 const label otherFacek =
307 otherPointFaces[otherPointj][otherPointFacej];
308 if (!otherFaceVisited[otherFacek])
310 otherFaceVisited[otherFacek] =
true;
311 otherVisitedFaces.
append(otherFacek);
312 otherNextFaces.
append(otherFacek);
318 forAll(otherPatch[otherFacej], otherFaceEdgej)
320 const label otherEdgej =
321 otherPatch.
faceEdges()[otherFacej][otherFaceEdgej];
328 otherPerimeterReached =
true;
336 otherFaceComplete[otherFacej] == 0
337 && !otherFaceQueued[otherFacej]
340 otherFaceQueued[otherFacej] =
true;
341 otherQueuedFaces.
append(otherFacej);
347 otherCurrentFaces.
transfer(otherNextFaces);
351 forAll(otherVisitedFaces, otherVisitedFacei)
353 otherFaceVisited[otherVisitedFaces[otherVisitedFacei]] =
false;
357 nFaceComplete -= faceComplete[facei] == 2;
358 faceComplete[facei] =
359 max(faceComplete[facei], otherPerimeterReached ? 1 : 2);
360 nFaceComplete += faceComplete[facei] == 2;
364 forAll(otherQueuedFaces, otherQueuedFacei)
366 otherFaceQueued[otherQueuedFaces[otherQueuedFacei]] =
false;
369 return nFaceComplete;
381 if (srcPatch.empty() || tgtPatch.empty())
return;
387 const fileName srcFileName = typeName +
"_srcPatch.vtk";
400 const fileName tgtFileName = typeName +
"_tgtPatch.vtk";
416 auto writeQueue = [&]
418 const label restarti,
419 const label iterationi,
427 typeName +
"_restarti=" +
name(restarti) +
"_iterationi=";
429 word(isSrc ?
"src" :
"tgt") +
"Queue";
433 if (iterationi == 0 || iterationi % 2 == !isSrc)
436 queueFileNamePart0 +
name(iterationi) +
'_'
437 + queueFileNamePart1 +
".vtk";
454 queueFileNamePart0 +
name(iterationi - 1) +
'_'
455 + queueFileNamePart1 +
".vtk",
456 queueFileNamePart0 +
name(iterationi) +
'_'
457 + queueFileNamePart1 +
".vtk"
464 auto writeQueues = [&]
466 const label restarti,
467 const label iterationi,
475 <<
": processing " << srcQueue.size() <<
'/'
476 << srcPatch.size() <<
" source and "
477 << tgtQueue.size() <<
'/' << tgtPatch.size()
478 <<
" target faces " <<
endl;
483 writeQueue(restarti, iterationi,
true, srcQueue);
484 writeQueue(restarti, iterationi,
false, tgtQueue);
514 label nSrcFaceComplete = 0, nTgtFaceComplete = 0;
515 labelList srcFaceComplete(srcPatch.size(), 0);
516 labelList tgtFaceComplete(tgtPatch.size(), 0);
517 boolList srcFaceQueued(srcPatch.size(),
false);
518 boolList tgtFaceQueued(tgtPatch.size(),
false);
519 boolList srcFaceVisited(srcPatch.size(),
false);
520 boolList tgtFaceVisited(tgtPatch.size(),
false);
523 while (srcFacei < srcPatch.size() && srcFacei != -1)
526 srcFaceComplete[srcFacei] = 2;
533 srcBox(srcPatch, srcPointNormals, srcPointNormals0, srcFacei)
536 if (!seedTgtFaces.
empty())
541 <<
" from at source face at "
548 forAll(seedTgtFaces, seedTgtFacei)
550 const label tgtFacei = seedTgtFaces[seedTgtFacei];
557 writeQueues(restarti, 0, srcQueue, tgtQueue);
561 label iterationi = 0;
584 writeQueues(restarti, 2*iterationi + 1, srcQueue, tgtQueue);
587 if (!tgtQueue.
size())
break;
609 writeQueues(restarti, 2*iterationi + 2, srcQueue, tgtQueue);
612 if (!srcQueue.
size())
break;
619 Info<<
indent <<
"Completed " << nSrcFaceComplete <<
'/'
620 << srcPatch.size() <<
" source faces and "
621 << nTgtFaceComplete <<
'/' << tgtPatch.size()
647 srcLocalTgtFaces_.resize(srcPatch.size());
648 forAll(srcLocalTgtFaces_, i)
650 srcLocalTgtFaces_[i].clear();
653 tgtLocalSrcFaces_.resize(tgtPatch.size());
654 forAll(tgtLocalSrcFaces_, i)
656 tgtLocalSrcFaces_[i].clear();
670 boolList oldLocalTgtFaceIsUsed(tgtPatch.size(),
false);
671 forAll(srcLocalTgtFaces_, srcFacei)
673 forAll(srcLocalTgtFaces_[srcFacei], i)
675 oldLocalTgtFaceIsUsed[srcLocalTgtFaces_[srcFacei][i]] =
true;
680 labelList oldToNewLocalTgtFace, newToOldLocalTgtFace;
683 oldLocalTgtFaceIsUsed,
685 oldToNewLocalTgtFace,
690 forAll(srcLocalTgtFaces_, srcFacei)
692 forAll(srcLocalTgtFaces_[srcFacei], i)
694 srcLocalTgtFaces_[srcFacei][i] =
695 oldToNewLocalTgtFace[srcLocalTgtFaces_[srcFacei][i]];
702 localTgtProcFacesPtr_() =
703 List<remote>(localTgtProcFacesPtr_(), newToOldLocalTgtFace);
705 return newToOldLocalTgtFace;
714 localSrcProcFacesPtr_.reset
730 localSrcProcFacesPtr_(),
747 forAll(srcLocalTgtFaces_, srcFacei)
749 nCouples += srcLocalTgtFaces_[srcFacei].size();
751 forAll(tgtLocalSrcFaces_, tgtFacei)
753 nCouples += tgtLocalSrcFaces_[tgtFacei].size();
772 localSrcProcFacesPtr_(nullptr),
773 localTgtProcFacesPtr_(nullptr)
787 const word& patchToPatchType,
791 boolConstructorTable::iterator cstrIter =
792 boolConstructorTablePtr_->find(patchToPatchType);
794 if (cstrIter == boolConstructorTablePtr_->end())
797 <<
"Unknown " << typeName <<
" type "
799 <<
"Valid " << typeName <<
" types are : " <<
endl
800 << boolConstructorTablePtr_->sortedToc()
813 forAll(srcLocalTgtFaces_, srcFacei)
815 result[srcFacei] = !srcLocalTgtFaces_[srcFacei].
empty();
824 forAll(tgtLocalSrcFaces_, tgtFacei)
826 result[tgtFacei] = !tgtLocalSrcFaces_[tgtFacei].
empty();
841 localTgtProcFacesPtr_()
855 localSrcProcFacesPtr_()
875 if (srcTotalSize == 0 || tgtTotalSize == 0)
901 tTgtPoints0Ptr.
ref(),
922 Info<<
indent << typeName <<
": Calculating couplings between "
923 << srcTotalSize <<
" source faces and " << tgtTotalSize
935 if (isSingleProcess())
938 initialise(srcPatch, srcPointNormals, srcPointNormals0, tTgtPatch);
942 srcBox(srcPatch, srcPointNormals, srcPointNormals0);
944 if (srcPatchBox.
overlaps(tgtPatchBox))
973 localTgtProcFacesPtr_.
reset
992 SubList<face>(localTTgtPatchPtr(), localTTgtPatchPtr().size()),
993 localTTgtPatchPtr().
points(),
994 localTTgtPatchPtr().points0()
998 SubList<face>(localTTgtPatchPtr(), localTTgtPatchPtr().size()),
999 localTTgtPatchPtr().
points()
1003 initialise(srcPatch, srcPointNormals, srcPointNormals, localTTgtPatch);
1006 if (localTTgtPatch.size())
1033 localTgtProcFacesPtr_()
1036 distributeSrc(srcPatch);
1039 rDistributeTgt(tgtPatch);
1043 const label nCouples =
1055 Info<<
indent << nCouples <<
" couplings calculated in "
1079 NullObjectRef<vectorField>(),
#define forAll(list, i)
Loop across all elements in list.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
void transfer(List< T > &)
Transfer contents of the argument List into this.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
void clear()
Clear the addressed list, i.e. set the size to zero.
void max(const dimensioned< Type > &)
void min(const dimensioned< Type > &)
A List with indirect addressing.
void append(const T &)
Append an element at the end of the list.
void size(const label)
Override size to be inconsistent with allocated storage.
bool empty() const
Return true if the list is empty (ie, size() is zero).
bool has0() const
Return whether or not old-time geometry is available.
const Field< PointType > & localPoints0() const
Return pointField of old-time points in patch.
A list of faces which address into the list of points.
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
const labelListList & pointFaces() const
Return point-face addressing.
const labelListList & edgeFaces() const
Return edge-face addressing.
const labelListList & faceEdges() const
Return face-edge addressing.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
const Field< PointType > & faceCentres() const
Return face centres for patch.
A List obtained as a section of another List.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
T & first()
Return the first element of the list.
bool empty() const
Return true if the UList is empty (ie, size() is zero)
void reset(T *=nullptr)
If object pointer already set, delete object and set to given.
const point & min() const
Minimum point defining the bounding box.
const point & max() const
Maximum point defining the bounding box.
void inflate(const scalar s)
Inflate box by factor*mag(span) in all dimensions.
Starts timing CPU usage and return elapsed time from start.
double cpuTimeIncrement() const
Return CPU time (in seconds) since last call to cpuTimeIncrement()
A face is a list of labels corresponding to mesh vertices.
A class for handling file names.
Non-pointer based hierarchical recursive searching.
Class to generate coupling geometry between two primitive patches.
virtual treeBoundBox srcBox(const face &srcFace, const pointField &srcPoints, const vectorField &srcPointNormals) const =0
Get the bound box for a source face.
virtual ~patchToPatch()
Destructor.
virtual void rDistributeTgt(const primitiveOldTimePatch &tgtPatch)
Send data that resulted from an intersection between the source.
virtual labelList finaliseLocal(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch)
Finalise the intersection locally. Trims the local target patch.
virtual void distributeSrc(const primitiveOldTimePatch &srcPatch)
Distribute the source patch so that any information needed by.
treeBoundBox tgtBox(const primitiveOldTimePatch &tgtPatch, const label tgtFacei) const
Get the bound box for a target face.
static autoPtr< patchToPatch > New(const word &patchToPatchType, const bool reverse)
Select from name.
List< List< remote > > tgtSrcProcFaces() const
For each target face, the coupled source procs and faces.
PackedBoolList srcCoupled() const
Return a list indicating which source faces are coupled.
bool findOrIntersectFaces(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const label srcFacei, const label tgtFacei)
Intersect two faces.
List< List< remote > > srcTgtProcFaces() const
For each source face, the coupled target procs and faces.
virtual void initialise(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch)
Initialisation.
void intersectPatches(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch)
Intersect the patches.
label intersectPatchQueue(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const bool isSrc, const DynamicList< labelPair > &queue, labelList &faceComplete, DynamicList< labelPair > &otherQueue, const labelList &otherFaceComplete, boolList &otherFaceQueued, boolList &otherFaceVisited)
Intersect a queue of source-target face pairs. Update completion.
void update(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const transformer &tgtToSrc=NullObjectRef< transformer >())
Update addressing and weights for the given patches.
PackedBoolList tgtCoupled() const
Return a list indicating which target faces are coupled.
patchToPatch(const bool reverse)
Construct from components.
virtual label finalise(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const transformer &tgtToSrc)
Finalising.
A class for managing temporary objects without reference counting.
T & ref()
Return non-const reference or generate a fatal error.
Standard boundBox + extra functionality for use in octree.
treeBoundBox extend(const scalar s) const
Return asymetrically extended bounding box, with guaranteed.
bool overlaps(const boundBox &) const
Overlaps other bounding box?
static const treeBoundBox invertedBox
As per boundBox::invertedBox, but with great instead of vGreat.
Encapsulation of data needed to search on PrimitivePatches.
A class for handling words, derived from string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
void write(const fileName &file, const word &title, const bool binary, const PointField &points, const VertexList &vertices, const LineList &lines, const FaceList &faces, const wordList &fieldNames, const boolList &fieldIsPointValues, const UPtrList< const Field< label >> &fieldLabelValues #define FieldTypeValuesConstArg(Type, nullArg))
Write VTK polygonal data to a file. Takes a PtrList of fields of labels and.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Pair< label > labelPair
Label pair.
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const bool)
Return a word representation of a bool.
label findNotIndex(const ListType &l, typename ListType::const_reference t, const label start=0)
labelList second(const UList< labelPair > &p)
vectorField pointField
pointField is a vectorField.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
void reverse(UList< T > &, const label n)
labelList first(const UList< labelPair > &p)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
PrimitiveOldTimePatch< SubList< face >, const pointField & > primitiveOldTimePatch
Addressing for a faceList slice.
List< labelList > labelListList
A List of labelList.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
defineTypeNameAndDebug(combustionModel, 0)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
bool isNull(const T &t)
Return true if t is a reference to the nullObject of type T.
static const label labelMax
treeBoundBox combine(const treeBoundBox &a, const treeBoundBox &b)
Ostream & indent(Ostream &os)
Indent stream.
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.