diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C index 192236ab29..e281fc989f 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C @@ -345,6 +345,8 @@ void Foam::AMIInterpolation::agglomerate // targetRestrictAddressing or allRestrict (since the same // for local data). const labelList& elems = map.subMap()[procI]; + const labelList& elemsMap = + map.constructMap()[Pstream::myProcNo()]; labelList& newSubMap = tgtSubMap[procI]; newSubMap.setSize(elems.size()); @@ -353,7 +355,7 @@ void Foam::AMIInterpolation::agglomerate forAll(elems, i) { - label fineElem = elems[i]; + label fineElem = elemsMap[elems[i]]; label coarseElem = allRestrict[fineElem]; if (oldToNew[coarseElem] == -1) { @@ -371,16 +373,26 @@ void Foam::AMIInterpolation::agglomerate // the sending map labelListList tgtConstructMap(Pstream::nProcs()); - labelList tgtCompactMap; // Local constructMap is just identity { tgtConstructMap[Pstream::myProcNo()] = identity(targetCoarseSize); - - tgtCompactMap = targetRestrictAddressing; } - tgtCompactMap.setSize(map.constructSize()); + + labelList tgtCompactMap(map.constructSize()); + + { + const labelList& elemsMap = + map.constructMap()[Pstream::myProcNo()]; + forAll(elemsMap, i) + { + label fineElem = elemsMap[i]; + label coarseElem = allRestrict[fineElem]; + tgtCompactMap[fineElem] = coarseElem; + } + } + label compactI = targetCoarseSize; // Compact data from other processors @@ -439,7 +451,6 @@ void Foam::AMIInterpolation::agglomerate } } - srcAddress.setSize(sourceCoarseSize); srcWeights.setSize(sourceCoarseSize); @@ -1095,7 +1106,7 @@ void Foam::AMIInterpolation::append const TargetPatch& tgtPatch ) { - // create a new interpolation + // Create a new interpolation autoPtr > newPtr ( new AMIInterpolation @@ -1110,7 +1121,145 @@ void Foam::AMIInterpolation::append ) ); - // combine new and current data + // If parallel then combine the mapDistribution and re-index + if (singlePatchProc_ == -1) + { + labelListList& srcSubMap = srcMapPtr_->subMap(); + labelListList& srcConstructMap = srcMapPtr_->constructMap(); + + labelListList& tgtSubMap = tgtMapPtr_->subMap(); + labelListList& tgtConstructMap = tgtMapPtr_->constructMap(); + + labelListList& newSrcSubMap = newPtr->srcMapPtr_->subMap(); + labelListList& newSrcConstructMap = newPtr->srcMapPtr_->constructMap(); + + labelListList& newTgtSubMap = newPtr->tgtMapPtr_->subMap(); + labelListList& newTgtConstructMap = newPtr->tgtMapPtr_->constructMap(); + + // Re-calculate the source indices + { + labelList mapMap(0), newMapMap(0); + forAll(srcSubMap, procI) + { + mapMap.append + ( + identity(srcConstructMap[procI].size()) + + mapMap.size() + newMapMap.size() + ); + newMapMap.append + ( + identity(newSrcConstructMap[procI].size()) + + mapMap.size() + newMapMap.size() + ); + } + + forAll(srcSubMap, procI) + { + forAll(srcConstructMap[procI], srcI) + { + srcConstructMap[procI][srcI] = + mapMap[srcConstructMap[procI][srcI]]; + } + } + + forAll(srcSubMap, procI) + { + forAll(newSrcConstructMap[procI], srcI) + { + newSrcConstructMap[procI][srcI] = + newMapMap[newSrcConstructMap[procI][srcI]]; + } + } + + forAll(tgtAddress_, tgtI) + { + forAll(tgtAddress_[tgtI], tgtJ) + { + tgtAddress_[tgtI][tgtJ] = + mapMap[tgtAddress_[tgtI][tgtJ]]; + } + } + + forAll(newPtr->tgtAddress_, tgtI) + { + forAll(newPtr->tgtAddress_[tgtI], tgtJ) + { + newPtr->tgtAddress_[tgtI][tgtJ] = + newMapMap[newPtr->tgtAddress_[tgtI][tgtJ]]; + } + } + } + + // Re-calculate the target indices + { + labelList mapMap(0), newMapMap(0); + forAll(srcSubMap, procI) + { + mapMap.append + ( + identity(tgtConstructMap[procI].size()) + + mapMap.size() + newMapMap.size() + ); + newMapMap.append + ( + identity(newTgtConstructMap[procI].size()) + + mapMap.size() + newMapMap.size() + ); + } + + forAll(srcSubMap, procI) + { + forAll(tgtConstructMap[procI], tgtI) + { + tgtConstructMap[procI][tgtI] = + mapMap[tgtConstructMap[procI][tgtI]]; + } + } + + forAll(srcSubMap, procI) + { + forAll(newTgtConstructMap[procI], tgtI) + { + newTgtConstructMap[procI][tgtI] = + newMapMap[newTgtConstructMap[procI][tgtI]]; + } + } + + forAll(srcAddress_, srcI) + { + forAll(srcAddress_[srcI], srcJ) + { + srcAddress_[srcI][srcJ] = + mapMap[srcAddress_[srcI][srcJ]]; + } + } + + forAll(newPtr->srcAddress_, srcI) + { + forAll(newPtr->srcAddress_[srcI], srcJ) + { + newPtr->srcAddress_[srcI][srcJ] = + newMapMap[newPtr->srcAddress_[srcI][srcJ]]; + } + } + } + + // Sum the construction sizes + srcMapPtr_->constructSize() += newPtr->srcMapPtr_->constructSize(); + tgtMapPtr_->constructSize() += newPtr->tgtMapPtr_->constructSize(); + + // Combine the maps + forAll(srcSubMap, procI) + { + srcSubMap[procI].append(newSrcSubMap[procI]); + srcConstructMap[procI].append(newSrcConstructMap[procI]); + + tgtSubMap[procI].append(newTgtSubMap[procI]); + tgtConstructMap[procI].append(newTgtConstructMap[procI]); + } + } + + // Combine new and current source data forAll(srcMagSf_, srcFaceI) { srcAddress_[srcFaceI].append(newPtr->srcAddress()[srcFaceI]); @@ -1118,6 +1267,7 @@ void Foam::AMIInterpolation::append srcWeightsSum_[srcFaceI] += newPtr->srcWeightsSum()[srcFaceI]; } + // Combine new and current target data forAll(tgtMagSf_, tgtFaceI) { tgtAddress_[tgtFaceI].append(newPtr->tgtAddress()[tgtFaceI]); diff --git a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C index 59a5359d8d..4a9c9d4f89 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C @@ -44,6 +44,151 @@ namespace Foam // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // +void Foam::cyclicPeriodicAMIPolyPatch::syncTransforms() const +{ + if (owner()) + { + // At this point we guarantee that the transformations have been + // updated. There is one particular case, where if the periodic patch + // locally has zero faces then its transformation will not be set. We + // need to synchronise the transforms over the zero-sized patches as + // well. + // + // We can't put the logic into cyclicPolyPatch as + // processorCyclicPolyPatch uses cyclicPolyPatch functionality. + // Synchronisation will fail because processor-type patches do not exist + // on all processors. + // + // The use in cyclicPeriodicAMI is special; we use the patch even + // though we have no faces in it. Ideally, in future, the transformation + // logic will be abstracted out, and we won't need a periodic patch + // here. Until then, we can just fix the behaviour of the zero-sized + // coupled patches here + + // Get the periodic patch + const coupledPolyPatch& periodicPatch + ( + refCast + ( + boundaryMesh()[periodicPatchID()] + ) + ); + + // If there are any zero-sized periodic patches + if (returnReduce((size() && !periodicPatch.size()), orOp())) + { + if (periodicPatch.separation().size() > 1) + { + FatalErrorIn + ( + "cyclicPeriodicAMIPolyPatch::resetAMI" + "(const AMIPatchToPatchInterpolation::interpolationMethod&" + ") const" + ) << "Periodic patch " << periodicPatchName_ + << " has non-uniform separation vector " + << periodicPatch.separation() + << "This is not allowed inside " << type() + << " patch " << name() + << exit(FatalError); + } + + if (periodicPatch.forwardT().size() > 1) + { + FatalErrorIn + ( + "cyclicPeriodicAMIPolyPatch::resetAMI" + "(const AMIPatchToPatchInterpolation::interpolationMethod&" + ") const" + ) << "Periodic patch " << periodicPatchName_ + << " has non-uniform transformation tensor " + << periodicPatch.forwardT() + << "This is not allowed inside " << type() + << " patch " << name() + << exit(FatalError); + } + + // Note that zero-sized patches will have zero-sized fields for the + // separation vector, forward and reverse transforms. These need + // replacing with the transformations from other processors. + + // Parallel in this context refers to a parallel transformation, + // rather than a rotational transformation + bool isParallel = periodicPatch.parallel(); + reduce(isParallel, orOp()); + + if (isParallel) + { + // Sync a list of separation vectors + List sep(Pstream::nProcs()); + sep[Pstream::myProcNo()] = periodicPatch.separation(); + Pstream::gatherList(sep); + Pstream::scatterList(sep); + + List coll(Pstream::nProcs()); + coll[Pstream::myProcNo()] = periodicPatch.collocated(); + Pstream::gatherList(coll); + Pstream::scatterList(coll); + + // If locally we have zero faces pick the first one that has a + // separation vector + if (!periodicPatch.size()) + { + forAll(sep, procI) + { + if (sep[procI].size()) + { + const_cast + ( + periodicPatch.separation() + ) = sep[procI]; + const_cast + ( + periodicPatch.collocated() + ) = coll[procI]; + + break; + } + } + } + } + else + { + // Sync a list of forward and reverse transforms + List forwardT(Pstream::nProcs()); + forwardT[Pstream::myProcNo()] = periodicPatch.forwardT(); + Pstream::gatherList(forwardT); + Pstream::scatterList(forwardT); + + List reverseT(Pstream::nProcs()); + reverseT[Pstream::myProcNo()] = periodicPatch.reverseT(); + Pstream::gatherList(reverseT); + Pstream::scatterList(reverseT); + + // If locally we have zero faces pick the first one that has a + // transformation vector + if (!periodicPatch.size()) + { + forAll(forwardT, procI) + { + if (forwardT[procI].size()) + { + const_cast + ( + periodicPatch.forwardT() + ) = forwardT[procI]; + const_cast + ( + periodicPatch.reverseT() + ) = reverseT[procI]; + } + } + } + } + } + } +} + + void Foam::cyclicPeriodicAMIPolyPatch::resetAMI ( const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod @@ -60,6 +205,9 @@ void Foam::cyclicPeriodicAMIPolyPatch::resetAMI ) ); + // Synchronise the transforms + syncTransforms(); + // Create copies of both patches' points, transformed to the owner pointField thisPoints0(localPoints()); pointField nbrPoints0(neighbPatch().localPoints()); diff --git a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H index 715cc5429d..ae59b1d5d7 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H @@ -72,6 +72,9 @@ private: // Private Member Functions + //- Synchronise the periodic transformations + void syncTransforms() const; + //- Reset the AMI interpolator virtual void resetAMI (