49 const label startSeedI,
56 for (label i = startSeedI; i < srcCellIDs.
size(); i++)
58 label srcI = srcCellIDs[i];
62 const point& srcCc = srcCcs[srcI];
69 tgtSeedI = hit.
index();
76 <<
"Unable to find nearest target cell"
77 <<
" for source cell " << srcI
78 <<
" with centre " << srcCc
88 Pout<<
"could not find starting seed" <<
endl;
101 const label srcSeedI,
102 const label tgtSeedI,
115 label srcCelli = srcSeedI;
116 label tgtCelli = tgtSeedI;
121 findNearestCell(src_, tgt_, srcCelli, tgtCelli);
124 srcToTgt[srcCelli].
append(tgtCelli);
125 tgtToSrc[tgtCelli].
append(srcCelli);
128 mapFlag[srcCelli] =
false;
131 V_ += srcVc[srcCelli];
143 while (srcCelli >= 0);
151 forAll(tgtToSrc, targetCelli)
153 if (tgtToSrc[targetCelli].size() > 1)
155 const vector& tgtC = tgtCc[targetCelli];
159 label srcCelli = srcCells[0];
160 scalar d =
magSqr(tgtC - srcCc[srcCelli]);
162 for (label i = 1; i < srcCells.
size(); i++)
164 label srcI = srcCells[i];
165 scalar dNew =
magSqr(tgtC - srcCc[srcI]);
174 srcCells.
append(srcCelli);
180 forAll(tgtToSrc, tgtCelli)
182 if (tgtToSrc[tgtCelli].empty())
184 label srcCelli = findMappedSrcCell(tgtCelli, tgtToSrc);
186 findNearestCell(tgt_, src_, tgtCelli, srcCelli);
188 tgtToSrc[tgtCelli].
append(srcCelli);
193 forAll(srcToTgtCellAddr, i)
195 srcToTgtCellWght[i] =
scalarList(srcToTgt[i].size(), srcVc[i]);
196 srcToTgtCellAddr[i].
transfer(srcToTgt[i]);
199 forAll(tgtToSrcCellAddr, i)
201 tgtToSrcCellWght[i] =
scalarList(tgtToSrc[i].size(), tgtVc[i]);
202 tgtToSrcCellAddr[i].
transfer(tgtToSrc[i]);
218 const vector& p1 = Cc1[cell1];
229 label c2 = cells2.
remove();
232 scalar dTest =
magSqr(Cc2[c2] - p1);
237 appendNbrCells(cell2, mesh2, visitedCells, cells2);
240 }
while (cells2.
size() > 0);
253 const labelList& srcNbr = src_.cellCells()[srcCelli];
258 label celli = srcNbr[i];
266 for (label i = startSeedI; i < srcCellIDs.
size(); i++)
268 label celli = srcCellIDs[i];
276 (void)findInitialSeeds
289 const label tgtCelli,
296 testCells.
append(tgtCelli);
301 label tgtI = testCells.
remove();
303 if (!visitedCells.
found(tgtI))
305 visitedCells.
append(tgtI);
307 if (tgtToSrc[tgtI].size())
309 return tgtToSrc[tgtI][0];
313 const labelList& nbrCells = tgt_.cellCells()[tgtI];
315 for (
const label nbrCelli : nbrCells)
317 if (!visitedCells.
found(nbrCelli))
319 testCells.
append(nbrCelli);
324 }
while (testCells.
size());
378 boolList mapFlag(src_.nCells(),
false);
384 label startSeedI = 0;
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
void clear() noexcept
Clear the addressed list, i.e. set the size to zero.
T remove()
Remove and return the last element. Fatal on an empty list.
void append(const T &val)
Copy append an element to the end of this list.
void transfer(List< T > &list)
void append(const T &val)
Append an element at the end of the list.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
label index() const noexcept
Return the hit index.
bool hit() const noexcept
Is there a hit?
bool found(const T &val, label pos=0) const
True if the value if found in the list.
void size(const label n)
Older name for setAddressableSize.
Map nearest mesh-to-mesh interpolation class.
virtual bool findInitialSeeds(const labelList &srcCellIDs, const boolList &mapFlag, const label startSeedI, label &srcSeedI, label &tgtSeedI) const
Find indices of overlapping cells in src and tgt meshes - returns.
virtual void findNearestCell(const polyMesh &mesh1, const polyMesh &mesh2, const label cell1, label &cell2) const
Find the nearest cell on mesh2 for cell1 on mesh1.
virtual void setNextNearestCells(label &startSeedI, label &srcCelli, label &tgtCelli, boolList &mapFlag, const labelList &srcCellIDs) const
Set the next cells for the marching front algorithm.
virtual label findMappedSrcCell(const label tgtCelli, const List< DynamicList< label > > &tgtToSrc) const
Find a source cell mapped to target cell tgtCelli.
virtual void calculateAddressing(labelListList &srcToTgtCellAddr, scalarListList &srcToTgtCellWght, labelListList &tgtToSrcCellAddr, scalarListList &tgtToSrcCellWght, const label srcSeedI, const label tgtSeedI, const labelList &srcCellIDs, boolList &mapFlag, label &startSeedI)
Calculate the mesh-to-mesh addressing and weights.
virtual void calculate(labelListList &srcToTgtAddr, scalarListList &srcToTgtWght, pointListList &srcToTgtVec, labelListList &tgtToSrcAddr, scalarListList &tgtToSrcWght, pointListList &tgtToSrcVec)
Calculate addressing and weights and optionally offset vectors.
virtual ~mapNearestMethod()
Destructor.
Base class for mesh-to-mesh calculation methods.
const polyMesh & tgt_
Reference to the target mesh.
const polyMesh & src_
Reference to the source mesh.
Mesh consisting of general polyhedral cells.
const indexedOctree< treeDataCell > & cellTree() const
Return the cell search tree.
const vectorField & cellCentres() const
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
List< scalar > scalarList
A List of scalars.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
UIndirectList< bool > boolUIndList
UIndirectList of bools.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
#define forAll(list, i)
Loop across all elements in list.