diff --git a/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.C index beb923b8990b14762f9d25957689e994c524bad9..96415c00e56adaaa1645d264230e02b901037a02 100644 --- a/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.C +++ b/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.C @@ -226,10 +226,9 @@ void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix { const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch(); - // note: only applying coupled contribution + // Note: only applying coupled contribution - const labelUList& nbrFaceCellsCoupled = - cpp.neighbPatch().faceCells(); + const labelUList& nbrFaceCellsCoupled = cpp.neighbPatch().faceCells(); solveScalarField pnf(psiInternal, nbrFaceCellsCoupled); @@ -254,7 +253,7 @@ void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix { const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch(); - // note: only applying coupled contribution + // Note: only applying coupled contribution const labelUList& nbrFaceCellsCoupled = cpp.neighbPatch().faceCells(); @@ -277,7 +276,7 @@ void Foam::cyclicACMIFvPatchField<Type>::manipulateMatrix { const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask(); - // nothing to be done by the AMI, but re-direct to non-overlap patch + // Nothing to be done by the AMI, but re-direct to non-overlap patch // with non-overlap patch weights const fvPatchField<Type>& npf = nonOverlapPatchField(); diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C index 23a12796ad1d9a2c682e398edeb62f53fae10fda..a66ae7d06b98271b7e39a4d4c45ba5f61e995801 100644 --- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C +++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C @@ -46,6 +46,7 @@ namespace Foam void Foam::cyclicACMIFvPatch::resetPatchAreas(const fvPatch& fvp) const { const_cast<vectorField&>(fvp.Sf()) = fvp.patch().faceAreas(); + const_cast<vectorField&>(fvp.Cf()) = fvp.patch().faceCentres(); const_cast<scalarField&>(fvp.magSf()) = mag(fvp.patch().faceAreas()); DebugPout @@ -117,7 +118,7 @@ Foam::tmp<Foam::vectorField> Foam::cyclicACMIFvPatch::delta() const vectorField nbrPatchD(interpolate(nbrPatch.coupledFvPatch::delta())); - tmp<vectorField> tpdv(new vectorField(patchD.size())); + auto tpdv = tmp<vectorField>::New(patchD.size()); vectorField& pdv = tpdv.ref(); // do the transformation if necessary @@ -204,6 +205,9 @@ void Foam::cyclicACMIFvPatch::movePoints() scalarField& phiNonOverlapp = meshPhiBf[nonOverlapPatch.patch().index()]; + const auto& localFaces = cyclicACMIPolyPatch_.localFaces(); + const auto& localPoints = cyclicACMIPolyPatch_.localPoints(); + forAll(phip, facei) { if (newSrcAddr[facei].empty()) @@ -214,12 +218,12 @@ void Foam::cyclicACMIFvPatch::movePoints() else { // Scale the mesh flux according to the area fraction - const face& fAMI = cyclicACMIPolyPatch_.localFaces()[facei]; + const face& fAMI = localFaces[facei]; // Note: using raw point locations to calculate the geometric // area - faces areas are currently scaled (decoupled from // mesh points) - const scalar geomArea = fAMI.mag(cyclicACMIPolyPatch_.localPoints()); + const scalar geomArea = fAMI.mag(localPoints); phip[facei] *= magSf()[facei]/geomArea; } } @@ -234,6 +238,9 @@ void Foam::cyclicACMIFvPatch::movePoints() scalarField& nbrPhiNonOverlapp = meshPhiBf[nbrNonOverlapPatch.patch().index()]; + const auto& nbrLocalFaces = nbrACMI.patch().localFaces(); + const auto& nbrLocalPoints = nbrACMI.patch().localPoints(); + forAll(nbrPhip, facei) { if (newTgtAddr[facei].empty()) @@ -242,12 +249,12 @@ void Foam::cyclicACMIFvPatch::movePoints() } else { - const face& fAMI = nbrACMI.patch().localFaces()[facei]; + const face& fAMI = nbrLocalFaces[facei]; // Note: using raw point locations to calculate the geometric // area - faces areas are currently scaled (decoupled from // mesh points) - const scalar geomArea = fAMI.mag(nbrACMI.patch().localPoints()); + const scalar geomArea = fAMI.mag(nbrLocalPoints); nbrPhip[facei] *= nbrACMI.magSf()[facei]/geomArea; } } diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C index a9c03e340c24d43d29e28aea662dc733344e10f2..baba50933c8bfff3931532bcf37acc6a202ee667 100644 --- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C +++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2019-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -123,7 +124,7 @@ Foam::tmp<Foam::vectorField> Foam::cyclicAMIFvPatch::delta() const const vectorField& nbrPatchD = tnbrPatchD(); - tmp<vectorField> tpdv(new vectorField(patchD.size())); + auto tpdv = tmp<vectorField>::New(patchD.size()); vectorField& pdv = tpdv.ref(); // do the transformation if necessary diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H index 4245dd4be83c232b14262ca5840853d8c5b3ac28..4d9daf3852c25c6908101cce73ee5997b6249d33 100644 --- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H +++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H @@ -159,8 +159,8 @@ public: } //- Return true if this patch is coupled. This is equivalent - // to the coupledPolyPatch::coupled() if parallel running or - // both sides present, false otherwise + //- to the coupledPolyPatch::coupled() if parallel running or + //- both sides present, false otherwise virtual bool coupled() const; //- Return delta (P to N) vectors across coupled patch @@ -190,7 +190,7 @@ public: // Interface transfer functions //- Return the values of the given internal data adjacent to - // the interface as a field + //- the interface as a field virtual tmp<labelField> interfaceInternalField ( const labelUList& internalData diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C index eaf1ebeb5a013419e4b53036835b96f3bf20ce4d..f01f8570da98bf03348f4edf895e0affa672234c 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2015-2018 OpenCFD Ltd. + Copyright (C) 2015-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -27,81 +27,104 @@ License \*---------------------------------------------------------------------------*/ #include "AMIInterpolation.H" -#include "AMIMethod.H" #include "meshTools.H" #include "mapDistribute.H" #include "flipOp.H" #include "profiling.H" - -#define DEBUGAMI(msg){Pout<< "[" << __FILE__ << ":" << __LINE__ << "]: " << #msg << "=" << msg << endl;} +#include "triPointRef.H" +#include "OFstream.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -const Foam::Enum -< - typename Foam::AMIInterpolation<SourcePatch, TargetPatch>:: - interpolationMethod -> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolationMethodNames_ -({ - { interpolationMethod::imNearestFace, "nearestFaceAMI" }, - { interpolationMethod::imFaceAreaWeight, "faceAreaWeightAMI" }, - { interpolationMethod::imPartialFaceAreaWeight, "partialFaceAreaWeightAMI" } -}); - - -template<class SourcePatch, class TargetPatch> -bool Foam::AMIInterpolation<SourcePatch, TargetPatch>::cacheIntersections_ = - false; - - -template<class SourcePatch, class TargetPatch> -template<class Patch> -Foam::tmp<Foam::scalarField> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::patchMagSf +namespace Foam +{ + defineTypeNameAndDebug(AMIInterpolation, 0); + defineRunTimeSelectionTable(AMIInterpolation, dict); + defineRunTimeSelectionTable(AMIInterpolation, component); +} + +bool Foam::AMIInterpolation::cacheIntersections_ = false; + + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +Foam::autoPtr<Foam::indexedOctree<Foam::AMIInterpolation::treeType>> +Foam::AMIInterpolation::createTree ( - const Patch& patch, - const faceAreaIntersect::triangulationMode triMode -) + const primitivePatch& patch +) const { - auto tResult = tmp<scalarField>::New(patch.size(), Zero); - scalarField& result = tResult.ref(); + treeBoundBox bb(patch.points(), patch.meshPoints()); + bb.inflate(0.01); + + return autoPtr<indexedOctree<treeType>>::New + ( + treeType + ( + false, + patch, + indexedOctree<treeType>::perturbTol() + ), + bb, // overall search domain + 8, // maxLevel + 10, // leaf size + 3.0 // duplicity + ); +} - const pointField& patchPoints = patch.localPoints(); - faceList patchFaceTris; +Foam::label Foam::AMIInterpolation::calcDistribution +( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch +) const +{ + label proci = 0; - forAll(result, patchFacei) + if (Pstream::parRun()) { - faceAreaIntersect::triangulate - ( - patch.localFaces()[patchFacei], - patchPoints, - triMode, - patchFaceTris - ); + labelList facesPresentOnProc(Pstream::nProcs(), Zero); + if ((srcPatch.size() > 0) || (tgtPatch.size() > 0)) + { + facesPresentOnProc[Pstream::myProcNo()] = 1; + } + else + { + facesPresentOnProc[Pstream::myProcNo()] = 0; + } + + Pstream::gatherList(facesPresentOnProc); + Pstream::scatterList(facesPresentOnProc); - forAll(patchFaceTris, i) + label nHaveFaces = sum(facesPresentOnProc); + + if (nHaveFaces > 1) { - result[patchFacei] += - triPointRef - ( - patchPoints[patchFaceTris[i][0]], - patchPoints[patchFaceTris[i][1]], - patchPoints[patchFaceTris[i][2]] - ).mag(); + proci = -1; + if (debug) + { + InfoInFunction + << "AMI split across multiple processors" << endl; + } + } + else if (nHaveFaces == 1) + { + proci = facesPresentOnProc.find(1); + if (debug) + { + InfoInFunction + << "AMI local to processor" << proci << endl; + } } } - return tResult; -} + // Either not parallel or no faces on any processor + return proci; +} -// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::projectPointsToSurface +void Foam::AMIInterpolation::projectPointsToSurface ( const searchableSurface& surf, pointField& pts @@ -141,8 +164,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::projectPointsToSurface } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights +void Foam::AMIInterpolation::normaliseWeights ( const scalarList& patchAreas, const word& patchName, @@ -170,7 +192,6 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights scalar s = sum(w); scalar t = s/denom; - if (conformal) { denom = s; @@ -182,10 +203,9 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights } wghtSum[facei] = t; - if (t < lowWeightTol) { - nLowWeight++; + ++nLowWeight; } } else @@ -194,7 +214,6 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights } } - if (output) { const label nFace = returnReduce(wght.size(), sumOp<label>()); @@ -223,8 +242,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::agglomerate +void Foam::AMIInterpolation::agglomerate ( const autoPtr<mapDistribute>& targetMapPtr, const scalarList& fineSrcMagSf, @@ -494,10 +512,10 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::agglomerate forAll(elems, i) { - label elemi = elems[i]; - label coarseElemi = targetRestrictAddressing[elemi]; + const label elemi = elems[i]; + const label coarseElemi = targetRestrictAddressing[elemi]; - label index = newElems.find(coarseElemi); + const label index = newElems.find(coarseElemi); if (index == -1) { newElems.append(coarseElemi); @@ -526,216 +544,87 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::agglomerate } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::constructFromSurface -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const autoPtr<searchableSurface>& surfPtr -) -{ - if (surfPtr.valid()) - { - // Create new patches for source and target - pointField srcPoints(srcPatch.points()); - SourcePatch srcPatch0 - ( - SubList<face> - ( - srcPatch, - srcPatch.size(), - 0 - ), - srcPoints - ); - - if (debug) - { - OFstream os("amiSrcPoints.obj"); - for (const point& pt : srcPoints) - { - meshTools::writeOBJ(os, pt); - } - } - - pointField tgtPoints(tgtPatch.points()); - TargetPatch tgtPatch0 - ( - SubList<face> - ( - tgtPatch, - tgtPatch.size(), - 0 - ), - tgtPoints - ); - - if (debug) - { - OFstream os("amiTgtPoints.obj"); - for (const point& pt : tgtPoints) - { - meshTools::writeOBJ(os, pt); - } - } - - - // Map source and target patches onto projection surface - projectPointsToSurface(surfPtr(), srcPoints); - projectPointsToSurface(surfPtr(), tgtPoints); - - - // Calculate AMI interpolation - update(srcPatch0, tgtPatch0); - } - else - { - update(srcPatch, tgtPatch); - } -} - - // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool requireMatch, - const interpolationMethod& method, - const scalar lowWeightCorrection, - const bool reverseTarget -) -: - AMIInterpolation - ( - srcPatch, - tgtPatch, - nullptr, - triMode, - requireMatch, - method, - lowWeightCorrection, - reverseTarget - ) -{} - - -template<class SourcePatch, class TargetPatch> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation +Foam::AMIInterpolation::AMIInterpolation ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool requireMatch, - const word& methodName, - const scalar lowWeightCorrection, + const dictionary& dict, const bool reverseTarget ) : - AMIInterpolation - ( - srcPatch, - tgtPatch, - nullptr, - triMode, - requireMatch, - interpolationMethodNames_.get(methodName), - lowWeightCorrection, - reverseTarget - ) + requireMatch_(dict.getOrDefault("requireMatch", true)), + reverseTarget_(dict.getOrDefault("reverseTarget", reverseTarget)), + lowWeightCorrection_(dict.getOrDefault<scalar>("lowWeightCorrection", -1)), + singlePatchProc_(-999), + srcMagSf_(), + srcAddress_(), + srcWeights_(), + srcWeightsSum_(), + srcCentroids_(), + srcMapPtr_(nullptr), + tgtMagSf_(), + tgtAddress_(), + tgtWeights_(), + tgtWeightsSum_(), + tgtCentroids_(), + tgtMapPtr_(nullptr), + upToDate_(false) {} -template<class SourcePatch, class TargetPatch> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation +Foam::AMIInterpolation::AMIInterpolation ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const autoPtr<searchableSurface>& surfPtr, - const faceAreaIntersect::triangulationMode& triMode, const bool requireMatch, - const interpolationMethod& method, - const scalar lowWeightCorrection, - const bool reverseTarget + const bool reverseTarget, + const scalar lowWeightCorrection ) : - methodName_(interpolationMethodNames_[method]), - reverseTarget_(reverseTarget), requireMatch_(requireMatch), - singlePatchProc_(-999), + reverseTarget_(reverseTarget), lowWeightCorrection_(lowWeightCorrection), + singlePatchProc_(-999), srcMagSf_(), srcAddress_(), srcWeights_(), srcWeightsSum_(), srcCentroids_(), + srcPatchPts_(), + srcMapPtr_(nullptr), tgtMagSf_(), tgtAddress_(), tgtWeights_(), tgtWeightsSum_(), tgtCentroids_(), - triMode_(triMode), - srcMapPtr_(nullptr), - tgtMapPtr_(nullptr) -{ - constructFromSurface(srcPatch, tgtPatch, surfPtr); -} - - -template<class SourcePatch, class TargetPatch> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const autoPtr<searchableSurface>& surfPtr, - const faceAreaIntersect::triangulationMode& triMode, - const bool requireMatch, - const word& methodName, - const scalar lowWeightCorrection, - const bool reverseTarget -) -: - AMIInterpolation - ( - srcPatch, - tgtPatch, - surfPtr, - triMode, - requireMatch, - interpolationMethodNames_.get(methodName), - lowWeightCorrection, - reverseTarget - ) + tgtPatchPts_(), + tgtMapPtr_(nullptr), + upToDate_(false) {} -template<class SourcePatch, class TargetPatch> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation +Foam::AMIInterpolation::AMIInterpolation ( - const AMIInterpolation<SourcePatch, TargetPatch>& fineAMI, + const AMIInterpolation& fineAMI, const labelList& sourceRestrictAddressing, const labelList& targetRestrictAddressing ) : - methodName_(fineAMI.methodName_), - reverseTarget_(fineAMI.reverseTarget_), requireMatch_(fineAMI.requireMatch_), - singlePatchProc_(fineAMI.singlePatchProc_), + reverseTarget_(fineAMI.reverseTarget_), lowWeightCorrection_(-1.0), + singlePatchProc_(fineAMI.singlePatchProc_), srcMagSf_(), srcAddress_(), srcWeights_(), srcWeightsSum_(), + srcPatchPts_(), + srcMapPtr_(nullptr), tgtMagSf_(), tgtAddress_(), tgtWeights_(), tgtWeightsSum_(), - triMode_(fineAMI.triMode_), - srcMapPtr_(nullptr), - tgtMapPtr_(nullptr) + tgtPatchPts_(), + tgtMapPtr_(nullptr), + upToDate_(false) { label sourceCoarseSize = ( @@ -817,16 +706,67 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation } +Foam::AMIInterpolation::AMIInterpolation(const AMIInterpolation& ami) +: + requireMatch_(ami.requireMatch_), + reverseTarget_(ami.reverseTarget_), + lowWeightCorrection_(ami.lowWeightCorrection_), + singlePatchProc_(ami.singlePatchProc_), + srcMagSf_(ami.srcMagSf_), + srcAddress_(ami.srcAddress_), + srcWeights_(ami.srcWeights_), + srcWeightsSum_(ami.srcWeightsSum_), + srcCentroids_(ami.srcCentroids_), + srcMapPtr_(nullptr), + tgtMagSf_(ami.tgtMagSf_), + tgtAddress_(ami.tgtAddress_), + tgtWeights_(ami.tgtWeights_), + tgtWeightsSum_(ami.tgtWeightsSum_), + tgtCentroids_(ami.tgtCentroids_), + tgtMapPtr_(nullptr), + upToDate_(false) +{} + + // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update +bool Foam::AMIInterpolation::calculate ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr ) { - addProfiling(ami, "AMIInterpolation::update"); + if (upToDate_) + { + return false; + } + + addProfiling(ami, "AMIInterpolation::calculate"); + + if (surfPtr) + { + srcPatchPts_ = srcPatch.points(); + projectPointsToSurface(surfPtr(), srcPatchPts_); + tsrcPatch0_ = tmpNrc<primitivePatch>::New + ( + SubList<face>(srcPatch), + srcPatchPts_ + ); + + tgtPatchPts_ = tgtPatch.points(); + projectPointsToSurface(surfPtr(), tgtPatchPts_); + ttgtPatch0_ = tmpNrc<primitivePatch>::New + ( + SubList<face>(tgtPatch), + tgtPatchPts_ + ); + } + else + { + tsrcPatch0_.cref(srcPatch); + ttgtPatch0_.cref(tgtPatch); + } label srcTotalSize = returnReduce(srcPatch.size(), sumOp<label>()); label tgtTotalSize = returnReduce(tgtPatch.size(), sumOp<label>()); @@ -836,7 +776,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update DebugInfo<< "AMI: no source faces present - no addressing constructed" << endl; - return; + return false; } Info<< indent @@ -845,53 +785,20 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update << tgtTotalSize << " target faces" << endl; - - // Calculate AMI interpolation - autoPtr<AMIMethod<SourcePatch, TargetPatch>> AMIPtr - ( - AMIMethod<SourcePatch, TargetPatch>::New - ( - methodName_, - srcPatch, - tgtPatch, - triMode_, - reverseTarget_, - requireMatch_ && (lowWeightCorrection_ < 0) - ) - ); - - AMIPtr->calculate - ( - srcAddress_, - srcWeights_, - srcCentroids_, - tgtAddress_, - tgtWeights_, - srcMagSf_, - tgtMagSf_, - srcMapPtr_, - tgtMapPtr_ - ); - - AMIPtr->normaliseWeights(true, *this); - - singlePatchProc_ = AMIPtr->distributed() ? -1 : 0; + singlePatchProc_ = calcDistribution(srcPatch, tgtPatch); if (debug) { - Info<< "AMIInterpolation : Constructed addressing and weights" << nl - << " triMode :" - << faceAreaIntersect::triangulationModeNames_[triMode_] << nl + Info<< "AMIInterpolation:" << nl << " singlePatchProc:" << singlePatchProc_ << nl - << " srcMagSf :" << gSum(srcMagSf_) << nl - << " tgtMagSf :" << gSum(tgtMagSf_) << nl << endl; } + + return true; } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update +void Foam::AMIInterpolation::reset ( autoPtr<mapDistribute>&& srcToTgtMap, autoPtr<mapDistribute>&& tgtToSrcMap, @@ -923,32 +830,25 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update srcMapPtr_ = srcToTgtMap; tgtMapPtr_ = tgtToSrcMap; + + upToDate_ = true; } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::append +void Foam::AMIInterpolation::append ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch ) { addProfiling(ami, "AMIInterpolation::append"); // Create a new interpolation - auto newPtr = autoPtr<AMIInterpolation<SourcePatch, TargetPatch>>::New - ( - srcPatch, - tgtPatch, - triMode_, - requireMatch_, - methodName_, - lowWeightCorrection_, - reverseTarget_ - ); + auto newPtr = clone(); + newPtr->calculate(srcPatch, tgtPatch); // If parallel then combine the mapDistribution and re-index - if (singlePatchProc_ == -1) + if (distributed()) { labelListList& srcSubMap = srcMapPtr_->subMap(); labelListList& srcConstructMap = srcMapPtr_->constructMap(); @@ -1007,8 +907,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::append { forAll(tgtAddress_[tgti], tgtj) { - tgtAddress_[tgti][tgtj] = - mapMap[tgtAddress_[tgti][tgtj]]; + tgtAddress_[tgti][tgtj] = mapMap[tgtAddress_[tgti][tgtj]]; } } @@ -1115,8 +1014,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::append } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights +void Foam::AMIInterpolation::normaliseWeights ( const bool conformal, const bool output @@ -1148,313 +1046,10 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights } -template<class SourcePatch, class TargetPatch> -template<class Type, class CombineOp> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToTarget -( - const UList<Type>& fld, - const CombineOp& cop, - List<Type>& result, - const UList<Type>& defaultValues -) const -{ - addProfiling(ami, "AMIInterpolation::interpolateToTarget"); - - if (fld.size() != srcAddress_.size()) - { - FatalErrorInFunction - << "Supplied field size is not equal to source patch size" << nl - << " source patch = " << srcAddress_.size() << nl - << " target patch = " << tgtAddress_.size() << nl - << " supplied field = " << fld.size() - << abort(FatalError); - } - - if (lowWeightCorrection_ > 0) - { - if (defaultValues.size() != tgtAddress_.size()) - { - FatalErrorInFunction - << "Employing default values when sum of weights falls below " - << lowWeightCorrection_ - << " but supplied default field size is not equal to target " - << "patch size" << nl - << " default values = " << defaultValues.size() << nl - << " target patch = " << tgtAddress_.size() << nl - << abort(FatalError); - } - } - - result.setSize(tgtAddress_.size()); - - if (singlePatchProc_ == -1) - { - const mapDistribute& map = srcMapPtr_(); - - List<Type> work(fld); - map.distribute(work); - - forAll(result, facei) - { - if (tgtWeightsSum_[facei] < lowWeightCorrection_) - { - result[facei] = defaultValues[facei]; - } - else - { - const labelList& faces = tgtAddress_[facei]; - const scalarList& weights = tgtWeights_[facei]; - - forAll(faces, i) - { - cop(result[facei], facei, work[faces[i]], weights[i]); - } - } - } - } - else - { - forAll(result, facei) - { - if (tgtWeightsSum_[facei] < lowWeightCorrection_) - { - result[facei] = defaultValues[facei]; - } - else - { - const labelList& faces = tgtAddress_[facei]; - const scalarList& weights = tgtWeights_[facei]; - - forAll(faces, i) - { - cop(result[facei], facei, fld[faces[i]], weights[i]); - } - } - } - } -} - - -template<class SourcePatch, class TargetPatch> -template<class Type, class CombineOp> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToSource +Foam::label Foam::AMIInterpolation::srcPointFace ( - const UList<Type>& fld, - const CombineOp& cop, - List<Type>& result, - const UList<Type>& defaultValues -) const -{ - addProfiling(ami, "AMIInterpolation::interpolateToSource"); - - if (fld.size() != tgtAddress_.size()) - { - FatalErrorInFunction - << "Supplied field size is not equal to target patch size" << nl - << " source patch = " << srcAddress_.size() << nl - << " target patch = " << tgtAddress_.size() << nl - << " supplied field = " << fld.size() - << abort(FatalError); - } - - if (lowWeightCorrection_ > 0) - { - if (defaultValues.size() != srcAddress_.size()) - { - FatalErrorInFunction - << "Employing default values when sum of weights falls below " - << lowWeightCorrection_ - << " but supplied default field size is not equal to source " - << "patch size" << nl - << " default values = " << defaultValues.size() << nl - << " source patch = " << srcAddress_.size() << nl - << abort(FatalError); - } - } - - result.setSize(srcAddress_.size()); - - if (singlePatchProc_ == -1) - { - const mapDistribute& map = tgtMapPtr_(); - - List<Type> work(fld); - map.distribute(work); - - forAll(result, facei) - { - if (srcWeightsSum_[facei] < lowWeightCorrection_) - { - result[facei] = defaultValues[facei]; - } - else - { - const labelList& faces = srcAddress_[facei]; - const scalarList& weights = srcWeights_[facei]; - - forAll(faces, i) - { - cop(result[facei], facei, work[faces[i]], weights[i]); - } - } - } - } - else - { - forAll(result, facei) - { - if (srcWeightsSum_[facei] < lowWeightCorrection_) - { - result[facei] = defaultValues[facei]; - } - else - { - const labelList& faces = srcAddress_[facei]; - const scalarList& weights = srcWeights_[facei]; - - forAll(faces, i) - { - cop(result[facei], facei, fld[faces[i]], weights[i]); - } - } - } - } -} - - -template<class SourcePatch, class TargetPatch> -template<class Type, class CombineOp> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToSource -( - const Field<Type>& fld, - const CombineOp& cop, - const UList<Type>& defaultValues -) const -{ - auto tresult = tmp<Field<Type>>::New(srcAddress_.size(), Zero); - - interpolateToSource - ( - fld, - multiplyWeightedOp<Type, CombineOp>(cop), - tresult.ref(), - defaultValues - ); - - return tresult; -} - - -template<class SourcePatch, class TargetPatch> -template<class Type, class CombineOp> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToSource -( - const tmp<Field<Type>>& tFld, - const CombineOp& cop, - const UList<Type>& defaultValues -) const -{ - return interpolateToSource(tFld(), cop, defaultValues); -} - - -template<class SourcePatch, class TargetPatch> -template<class Type, class CombineOp> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToTarget -( - const Field<Type>& fld, - const CombineOp& cop, - const UList<Type>& defaultValues -) const -{ - auto tresult = tmp<Field<Type>>::New(tgtAddress_.size(), Zero); - - interpolateToTarget - ( - fld, - multiplyWeightedOp<Type, CombineOp>(cop), - tresult.ref(), - defaultValues - ); - - return tresult; -} - - -template<class SourcePatch, class TargetPatch> -template<class Type, class CombineOp> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToTarget -( - const tmp<Field<Type>>& tFld, - const CombineOp& cop, - const UList<Type>& defaultValues -) const -{ - return interpolateToTarget(tFld(), cop, defaultValues); -} - - -template<class SourcePatch, class TargetPatch> -template<class Type> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToSource -( - const Field<Type>& fld, - const UList<Type>& defaultValues -) const -{ - return interpolateToSource(fld, plusEqOp<Type>(), defaultValues); -} - - -template<class SourcePatch, class TargetPatch> -template<class Type> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToSource -( - const tmp<Field<Type>>& tFld, - const UList<Type>& defaultValues -) const -{ - return interpolateToSource(tFld(), plusEqOp<Type>(), defaultValues); -} - - -template<class SourcePatch, class TargetPatch> -template<class Type> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToTarget -( - const Field<Type>& fld, - const UList<Type>& defaultValues -) const -{ - return interpolateToTarget(fld, plusEqOp<Type>(), defaultValues); -} - - -template<class SourcePatch, class TargetPatch> -template<class Type> -Foam::tmp<Foam::Field<Type>> -Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToTarget -( - const tmp<Field<Type>>& tFld, - const UList<Type>& defaultValues -) const -{ - return interpolateToTarget(tFld(), plusEqOp<Type>(), defaultValues); -} - - -template<class SourcePatch, class TargetPatch> -Foam::label Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcPointFace -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const vector& n, const label tgtFacei, point& tgtPoint @@ -1498,11 +1093,10 @@ const } -template<class SourcePatch, class TargetPatch> -Foam::label Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtPointFace +Foam::label Foam::AMIInterpolation::tgtPointFace ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const vector& n, const label srcFacei, point& srcPoint @@ -1546,11 +1140,7 @@ const } -template<class SourcePatch, class TargetPatch> -bool Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkSymmetricWeights -( - const bool log -) const +bool Foam::AMIInterpolation::checkSymmetricWeights(const bool log) const { if (Pstream::parRun() && (singlePatchProc_ == -1)) { @@ -1613,11 +1203,10 @@ bool Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkSymmetricWeights } -template<class SourcePatch, class TargetPatch> -void Foam::AMIInterpolation<SourcePatch, TargetPatch>::writeFaceConnectivity +void Foam::AMIInterpolation::writeFaceConnectivity ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const labelListList& srcAddress ) const @@ -1646,4 +1235,20 @@ const } +void Foam::AMIInterpolation::write(Ostream& os) const +{ + os.writeEntry("AMIMethod", type()); + + if (reverseTarget_) + { + os.writeEntry("flipNormals", reverseTarget_); + } + + if (lowWeightCorrection_ > 0) + { + os.writeEntry("lowWeightCorrection", lowWeightCorrection_); + } +} + + // ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H index 82848433c58d145c6ba35893dde50b84ee6db518..8ce226bbb4f392b1fcf3bb66491c069a0548397b 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2016-2018 OpenCFD Ltd. + Copyright (C) 2016-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -62,75 +62,52 @@ SourceFiles #include "ops.H" #include "Enum.H" #include "pointList.H" +#include "indexedOctree.H" +#include "treeDataPrimitivePatch.H" + +#include "runTimeSelectionTables.H" + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // namespace Foam { -/*---------------------------------------------------------------------------*\ - Class AMIInterpolationName Declaration -\*---------------------------------------------------------------------------*/ - -TemplateName(AMIInterpolation); - - /*---------------------------------------------------------------------------*\ Class AMIInterpolation Declaration \*---------------------------------------------------------------------------*/ -template<class SourcePatch, class TargetPatch> class AMIInterpolation -: - public AMIInterpolationName { public: // Public data types - //- Enumeration specifying interpolation method - enum interpolationMethod - { - imNearestFace, - imFaceAreaWeight, - imPartialFaceAreaWeight - }; - - static const Enum<interpolationMethod> interpolationMethodNames_; - static bool cacheIntersections_; - //- Calculate the patch face magnitudes for the given tri-mode - template<class Patch> - static tmp<scalarField> patchMagSf - ( - const Patch& patch, - const faceAreaIntersect::triangulationMode triMode - ); +protected: -private: + //- Local typedef to octree tree-type + typedef treeDataPrimitivePatch<primitivePatch> treeType; - // Private data + // Protected data - //- Interpolation method - const word methodName_; + //- Flag to indicate that the two patches must be matched/an overlap + //- exists between them + bool requireMatch_; //- Flag to indicate that the two patches are co-directional and //- that the orientation of the target patch should be reversed const bool reverseTarget_; - //- Flag to indicate that the two patches must be matched/an overlap - //- exists between them - const bool requireMatch_; + //- Threshold weight below which interpolation is deactivated + scalar lowWeightCorrection_; //- Index of processor that holds all of both sides. -1 in all other //- cases label singlePatchProc_; - //- Threshold weight below which interpolation is deactivated - scalar lowWeightCorrection_; - // Source patch @@ -149,6 +126,16 @@ private: //- Centroid of target faces per source face pointListList srcCentroids_; + //- Source patch points if input points are manipulated, e.g. + //- projected + pointField srcPatchPts_; + + //- Source patch using manipulated input points + tmpNrc<primitivePatch> tsrcPatch0_; + + //- Source map pointer - parallel running only + autoPtr<mapDistribute> srcMapPtr_; + // Target patch @@ -167,21 +154,21 @@ private: //- Centroid of source faces per target face pointListList tgtCentroids_; + //- Target patch points if input points are manipulated, e.g. + //- projected + pointField tgtPatchPts_; - //- Face triangulation mode - const faceAreaIntersect::triangulationMode triMode_; + //- Target patch using manipulated input points + tmpNrc<primitivePatch> ttgtPatch0_; - //- Source map pointer - parallel running only - autoPtr<mapDistribute> srcMapPtr_; + //- Target map pointer - parallel running only + autoPtr<mapDistribute> tgtMapPtr_; - //- Target map pointer - parallel running only - autoPtr<mapDistribute> tgtMapPtr_; + //- Up-to-date flag + bool upToDate_; - // Private Member Functions - - //- No copy construct - AMIInterpolation(const AMIInterpolation&) = delete; + // Protected Member Functions //- No copy assignment void operator=(const AMIInterpolation&) = delete; @@ -189,6 +176,19 @@ private: // Initialisation + //- Reset the octree for the patch face search + autoPtr<indexedOctree<treeType>> createTree + ( + const primitivePatch& patch + ) const; + + //- Calculate if patches are on multiple processors + label calcDistribution + ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch + ) const; + //- Project points to surface void projectPointsToSurface ( @@ -197,6 +197,15 @@ private: ) const; + // Access + + //- Return the orginal src patch with optionally updated points + inline const primitivePatch& srcPatch0() const; + + //- Return the orginal tgt patch with optionally updated points + inline const primitivePatch& tgtPatch0() const; + + // Evaluation //- Normalise the (area) weights - suppresses numerical error in @@ -236,96 +245,127 @@ private: autoPtr<mapDistribute>& tgtMap ); - void constructFromSurface - ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const autoPtr<searchableSurface>& surfPtr - ); public: - // Constructors + //- Runtime type information + TypeName("AMIInterpolation"); - //- Construct from components - AMIInterpolation + // Selection tables + + //- Selection table for dictionary construction + declareRunTimeSelectionTable ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool requireMatch = true, - const interpolationMethod& method = imFaceAreaWeight, - const scalar lowWeightCorrection = -1, + autoPtr, + AMIInterpolation, + dict, + ( + const dictionary& dict, + const bool reverseTarget + ), + ( + dict, + reverseTarget + ) + ); + + //- Selection table for component-wise construction + declareRunTimeSelectionTable + ( + autoPtr, + AMIInterpolation, + component, + ( + const bool requireMatch, + const bool reverseTarget, + const scalar lowWeightCorrection + ), + ( + requireMatch, + reverseTarget, + lowWeightCorrection + ) + ); + + //- Selector for dictionary + static autoPtr<AMIInterpolation> New + ( + const word& modelName, + const dictionary& dict, const bool reverseTarget = false ); - //- Construct from components - AMIInterpolation + //- Selector for components + static autoPtr<AMIInterpolation> New ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const word& modelName, const bool requireMatch = true, - const word& methodName = - interpolationMethodNames_[imFaceAreaWeight], - const scalar lowWeightCorrection = -1, - const bool reverseTarget = false + const bool reverseTarget = false, + const scalar lowWeightCorrection = -1 ); - //- Construct from components, with projection surface + + // Constructors + + //- Construct from dictionary AMIInterpolation ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const autoPtr<searchableSurface>& surf, - const faceAreaIntersect::triangulationMode& triMode, - const bool requireMatch = true, - const interpolationMethod& method = imFaceAreaWeight, - const scalar lowWeightCorrection = -1, + const dictionary& dict, const bool reverseTarget = false ); - //- Construct from components, with projection surface + //- Construct from components AMIInterpolation ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const autoPtr<searchableSurface>& surf, - const faceAreaIntersect::triangulationMode& triMode, const bool requireMatch = true, - const word& methodName = - interpolationMethodNames_[imFaceAreaWeight], - const scalar lowWeightCorrection = -1, - const bool reverseTarget = false + const bool reverseTarget = false, + const scalar lowWeightCorrection = -1 ); //- Construct from agglomeration of AMIInterpolation. Agglomeration - // passed in as new coarse size and addressing from fine from coarse + //- passed in as new coarse size and addressing from fine from coarse AMIInterpolation ( - const AMIInterpolation<SourcePatch, TargetPatch>& fineAMI, + const AMIInterpolation& fineAMI, const labelList& sourceRestrictAddressing, const labelList& neighbourRestrictAddressing ); + //- Construct as copy + AMIInterpolation(const AMIInterpolation& ami); - //- Destructor - ~AMIInterpolation() = default; + //- Construct and return a clone + virtual autoPtr<AMIInterpolation> clone() const + { + return autoPtr<AMIInterpolation>::New(*this); + } - // Typedef to SourcePatch type this AMIInterpolation is instantiated on - typedef SourcePatch sourcePatchType; - // Typedef to TargetPatch type this AMIInterpolation is instantiated on - typedef TargetPatch targetPatchType; + //- Destructor + virtual ~AMIInterpolation() = default; // Member Functions // Access - //- Set to -1, or the processor holding all faces (both sides) of - //- the AMI - inline label singlePatchProc() const; + //- Access to the up-to-date flag + inline bool upToDate() const; + + //- Access to the up-to-date flag + inline bool& upToDate(); + + //- Access to the distributed flag + inline bool distributed() const; + + //- Access to the requireMatch flag + inline bool requireMatch() const; + + //- Access to the requireMatch flag + inline bool setRequireMatch(const bool flag); + + //- Access to the reverseTarget flag + inline bool reverseTarget() const; //- Threshold weight below which interpolation is deactivated inline scalar lowWeightCorrection() const; @@ -333,6 +373,10 @@ public: //- Return true if employing a 'lowWeightCorrection' inline bool applyLowWeightCorrection() const; + //- Set to -1, or the processor holding all faces (both sides) of + //- the AMI + inline label singlePatchProc() const; + // Source patch @@ -410,14 +454,16 @@ public: // Manipulation - //- Update addressing and weights - void update + //- Update addressing, weights and (optional) centroids + virtual bool calculate ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr = nullptr ); - void update + //- Set the maps, addresses and weights from an external source + void reset ( autoPtr<mapDistribute>&& srcToTgtMap, autoPtr<mapDistribute>&& tgtToSrcMap, @@ -427,17 +473,11 @@ public: scalarListList&& tgtWeights ); - void setAreas(const scalarList& srcMagSf, const scalarList& tgtMagSf) - { - srcMagSf_ = srcMagSf; - tgtMagSf_ = tgtMagSf; - } - //- Append additional addressing and weights void append ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch ); //- Normalise the weights @@ -545,8 +585,8 @@ public: //- Return source patch face index of point on target patch face label srcPointFace ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const vector& n, const label tgtFacei, point& tgtPoint @@ -556,8 +596,8 @@ public: //- Return target patch face index of point on source patch face label tgtPointFace ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const vector& n, const label srcFacei, point& srcPoint @@ -574,10 +614,16 @@ public: //- Write face connectivity as OBJ file void writeFaceConnectivity ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const labelListList& srcAddress ) const; + + + // I-O + + //- Write + virtual void write(Ostream& os) const; }; @@ -592,7 +638,7 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #ifdef NoRepository - #include "AMIInterpolation.C" + #include "AMIInterpolationTemplates.C" #endif // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H index 027a504521d7bf285941b0e4fe8a414ac7f9a72b..22823a23ec6d8d57f334b1731d8dba53d92506f9 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2016 OpenCFD Ltd. + Copyright (C) 2016-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -26,186 +26,203 @@ License \*---------------------------------------------------------------------------*/ -template<class SourcePatch, class TargetPatch> -inline Foam::label -Foam::AMIInterpolation<SourcePatch, TargetPatch>::singlePatchProc() const +inline const Foam::primitivePatch& Foam::AMIInterpolation::srcPatch0() const { - return singlePatchProc_; + if (!tsrcPatch0_.valid()) + { + FatalErrorInFunction + << "tsrcPatch0Ptr_ not set" + << abort(FatalError); + } + + return tsrcPatch0_(); +} + + +inline const Foam::primitivePatch& Foam::AMIInterpolation::tgtPatch0() const +{ + + if (!ttgtPatch0_.valid()) + { + FatalErrorInFunction + << "ttgtPatch0Ptr_ not set" + << abort(FatalError); + } + + return ttgtPatch0_(); } -template<class SourcePatch, class TargetPatch> -inline Foam::scalar -Foam::AMIInterpolation<SourcePatch, TargetPatch>::lowWeightCorrection() const +inline bool Foam::AMIInterpolation::upToDate() const +{ + return upToDate_; +} + + +inline bool& Foam::AMIInterpolation::upToDate() +{ + return upToDate_; +} + + +inline bool Foam::AMIInterpolation::distributed() const +{ + return singlePatchProc_ == -1; +} + + +inline bool Foam::AMIInterpolation::requireMatch() const +{ + return requireMatch_ && lowWeightCorrection_ < 0; +} + + +inline bool Foam::AMIInterpolation::setRequireMatch(const bool flag) +{ + requireMatch_ = flag; + return requireMatch_; +} + + +inline bool Foam::AMIInterpolation::reverseTarget() const +{ + return reverseTarget_; +} + + +inline Foam::scalar Foam::AMIInterpolation::lowWeightCorrection() const { return lowWeightCorrection_; } -template<class SourcePatch, class TargetPatch> -inline bool -Foam::AMIInterpolation<SourcePatch, TargetPatch>:: -applyLowWeightCorrection() const +inline bool Foam::AMIInterpolation::applyLowWeightCorrection() const { return lowWeightCorrection_ > 0; } -template<class SourcePatch, class TargetPatch> -inline const Foam::List<Foam::scalar>& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcMagSf() const +inline Foam::label Foam::AMIInterpolation::singlePatchProc() const +{ + return singlePatchProc_; +} + + +inline const Foam::List<Foam::scalar>& Foam::AMIInterpolation::srcMagSf() const { return srcMagSf_; } -template<class SourcePatch, class TargetPatch> -inline Foam::List<Foam::scalar>& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcMagSf() +inline Foam::List<Foam::scalar>& Foam::AMIInterpolation::srcMagSf() { return srcMagSf_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::labelListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcAddress() const +inline const Foam::labelListList& Foam::AMIInterpolation::srcAddress() const { return srcAddress_; } -template<class SourcePatch, class TargetPatch> -inline Foam::labelListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcAddress() +inline Foam::labelListList& Foam::AMIInterpolation::srcAddress() { return srcAddress_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::scalarListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeights() const +inline const Foam::scalarListList& Foam::AMIInterpolation::srcWeights() const { return srcWeights_; } -template<class SourcePatch, class TargetPatch> -inline Foam::scalarListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeights() +inline Foam::scalarListList& Foam::AMIInterpolation::srcWeights() { return srcWeights_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::scalarField& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeightsSum() const +inline const Foam::scalarField& Foam::AMIInterpolation::srcWeightsSum() const { return srcWeightsSum_; } -template<class SourcePatch, class TargetPatch> -inline Foam::scalarField& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeightsSum() +inline Foam::scalarField& Foam::AMIInterpolation::srcWeightsSum() { return srcWeightsSum_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::pointListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcCentroids() const +inline const Foam::pointListList& Foam::AMIInterpolation::srcCentroids() const { return srcCentroids_; } -template<class SourcePatch, class TargetPatch> -inline Foam::pointListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcCentroids() +inline Foam::pointListList& Foam::AMIInterpolation::srcCentroids() { return srcCentroids_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::mapDistribute& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcMap() const +inline const Foam::mapDistribute& Foam::AMIInterpolation::srcMap() const { return *srcMapPtr_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::List<Foam::scalar>& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtMagSf() const +inline const Foam::List<Foam::scalar>& Foam::AMIInterpolation::tgtMagSf() const { return tgtMagSf_; } -template<class SourcePatch, class TargetPatch> -inline Foam::List<Foam::scalar>& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtMagSf() +inline Foam::List<Foam::scalar>& Foam::AMIInterpolation::tgtMagSf() { return tgtMagSf_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::labelListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtAddress() const +inline const Foam::labelListList& Foam::AMIInterpolation::tgtAddress() const { return tgtAddress_; } -template<class SourcePatch, class TargetPatch> -inline Foam::labelListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtAddress() +inline Foam::labelListList& Foam::AMIInterpolation::tgtAddress() { return tgtAddress_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::scalarListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtWeights() const +inline const Foam::scalarListList& Foam::AMIInterpolation::tgtWeights() const { return tgtWeights_; } -template<class SourcePatch, class TargetPatch> -inline Foam::scalarListList& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtWeights() +inline Foam::scalarListList& Foam::AMIInterpolation::tgtWeights() { return tgtWeights_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::scalarField& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtWeightsSum() const +inline const Foam::scalarField& Foam::AMIInterpolation::tgtWeightsSum() const { return tgtWeightsSum_; } -template<class SourcePatch, class TargetPatch> -inline Foam::scalarField& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtWeightsSum() +inline Foam::scalarField& Foam::AMIInterpolation::tgtWeightsSum() { return tgtWeightsSum_; } -template<class SourcePatch, class TargetPatch> -inline const Foam::mapDistribute& -Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtMap() const +inline const Foam::mapDistribute& Foam::AMIInterpolation::tgtMap() const { return *tgtMapPtr_; } diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationName.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationName.C deleted file mode 100644 index b09ca433fbef4a6a4f8da8d4083a7bb8ec255750..0000000000000000000000000000000000000000 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationName.C +++ /dev/null @@ -1,38 +0,0 @@ -/*---------------------------------------------------------------------------*\ - ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox - \\ / O peration | - \\ / A nd | www.openfoam.com - \\/ M anipulation | -------------------------------------------------------------------------------- - Copyright (C) 2011-2012 OpenFOAM Foundation -------------------------------------------------------------------------------- -License - This file is part of OpenFOAM. - - OpenFOAM is free software: you can redistribute it and/or modify it - under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - OpenFOAM is distributed in the hope that it will be useful, but WITHOUT - ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License - for more details. - - You should have received a copy of the GNU General Public License - along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. - -\*---------------------------------------------------------------------------*/ - -#include "AMIInterpolation.H" - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -namespace Foam -{ -defineTypeNameAndDebug(AMIInterpolationName, 0); -} - - -// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationNew.C similarity index 59% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C rename to src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationNew.C index 5ed107dc23e47de7f3d35a44d06c9d8e44e8eb37..6bc4ba2c2e673869b8afae534ea2ecd5a212a9d2 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationNew.C @@ -5,8 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2019 OpenCFD Ltd. + Copyright (C) 2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -26,46 +25,66 @@ License \*---------------------------------------------------------------------------*/ +#include "AMIInterpolation.H" + // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -Foam::autoPtr<Foam::AMIMethod<SourcePatch, TargetPatch>> -Foam::AMIMethod<SourcePatch, TargetPatch>::New +Foam::autoPtr<Foam::AMIInterpolation> Foam::AMIInterpolation::New ( - const word& methodName, - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const word& modelName, + const dictionary& dict, + const bool reverseTarget +) +{ + DebugInfo << "Selecting model " << modelName << endl; + + auto cstrIter = dictConstructorTablePtr_->cfind(modelName); + + if (!cstrIter.found()) + { + FatalErrorInLookup + ( + typeName, + modelName, + *dictConstructorTablePtr_ + ) << exit(FatalError); + } + + return autoPtr<AMIInterpolation>(cstrIter()(dict, reverseTarget)); +} + + +Foam::autoPtr<Foam::AMIInterpolation> Foam::AMIInterpolation::New +( + const word& modelName, + const bool requireMatch, const bool reverseTarget, - const bool requireMatch + const scalar lowWeightCorrection ) { - DebugInfo << "Selecting AMIMethod " << methodName << endl; + DebugInfo << "Selecting model " << modelName << endl; - auto cstrIter = componentsConstructorTablePtr_->cfind(methodName); + auto cstrIter = componentConstructorTablePtr_->cfind(modelName); if (!cstrIter.found()) { FatalErrorInLookup ( - "AMIMethod", - methodName, - *componentsConstructorTablePtr_ + typeName, + modelName, + *componentConstructorTablePtr_ ) << exit(FatalError); } - return autoPtr<AMIMethod<SourcePatch, TargetPatch>> + return autoPtr<AMIInterpolation> ( cstrIter() ( - srcPatch, - tgtPatch, - triMode, + requireMatch, reverseTarget, - requireMatch + lowWeightCorrection ) ); } - // ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationTemplates.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationTemplates.C new file mode 100644 index 0000000000000000000000000000000000000000..418bf48cc863fe6a977cdd81b54a4ab92f7ec37c --- /dev/null +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationTemplates.C @@ -0,0 +1,318 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2011-2017 OpenFOAM Foundation + Copyright (C) 2015-2018 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "profiling.H" +#include "mapDistribute.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +template<class Type, class CombineOp> +void Foam::AMIInterpolation::interpolateToTarget +( + const UList<Type>& fld, + const CombineOp& cop, + List<Type>& result, + const UList<Type>& defaultValues +) const +{ + addProfiling(ami, "AMIInterpolation::interpolateToTarget"); + + if (fld.size() != srcAddress_.size()) + { + FatalErrorInFunction + << "Supplied field size is not equal to source patch size" << nl + << " source patch = " << srcAddress_.size() << nl + << " target patch = " << tgtAddress_.size() << nl + << " supplied field = " << fld.size() + << abort(FatalError); + } + + if (lowWeightCorrection_ > 0) + { + if (defaultValues.size() != tgtAddress_.size()) + { + FatalErrorInFunction + << "Employing default values when sum of weights falls below " + << lowWeightCorrection_ + << " but supplied default field size is not equal to target " + << "patch size" << nl + << " default values = " << defaultValues.size() << nl + << " target patch = " << tgtAddress_.size() << nl + << abort(FatalError); + } + } + + result.setSize(tgtAddress_.size()); + + if (distributed()) + { + const mapDistribute& map = srcMapPtr_(); + + List<Type> work(fld); + map.distribute(work); + + forAll(result, facei) + { + if (tgtWeightsSum_[facei] < lowWeightCorrection_) + { + result[facei] = defaultValues[facei]; + } + else + { + const labelList& faces = tgtAddress_[facei]; + const scalarList& weights = tgtWeights_[facei]; + + forAll(faces, i) + { + cop(result[facei], facei, work[faces[i]], weights[i]); + } + } + } + } + else + { + forAll(result, facei) + { + if (tgtWeightsSum_[facei] < lowWeightCorrection_) + { + result[facei] = defaultValues[facei]; + } + else + { + const labelList& faces = tgtAddress_[facei]; + const scalarList& weights = tgtWeights_[facei]; + + forAll(faces, i) + { + cop(result[facei], facei, fld[faces[i]], weights[i]); + } + } + } + } +} + + +template<class Type, class CombineOp> +void Foam::AMIInterpolation::interpolateToSource +( + const UList<Type>& fld, + const CombineOp& cop, + List<Type>& result, + const UList<Type>& defaultValues +) const +{ + addProfiling(ami, "AMIInterpolation::interpolateToSource"); + + if (fld.size() != tgtAddress_.size()) + { + FatalErrorInFunction + << "Supplied field size is not equal to target patch size" << nl + << " source patch = " << srcAddress_.size() << nl + << " target patch = " << tgtAddress_.size() << nl + << " supplied field = " << fld.size() + << abort(FatalError); + } + + if (lowWeightCorrection_ > 0) + { + if (defaultValues.size() != srcAddress_.size()) + { + FatalErrorInFunction + << "Employing default values when sum of weights falls below " + << lowWeightCorrection_ + << " but supplied default field size is not equal to source " + << "patch size" << nl + << " default values = " << defaultValues.size() << nl + << " source patch = " << srcAddress_.size() << nl + << abort(FatalError); + } + } + + result.setSize(srcAddress_.size()); + + if (distributed()) + { + const mapDistribute& map = tgtMapPtr_(); + + List<Type> work(fld); + map.distribute(work); + + forAll(result, facei) + { + if (srcWeightsSum_[facei] < lowWeightCorrection_) + { + result[facei] = defaultValues[facei]; + } + else + { + const labelList& faces = srcAddress_[facei]; + const scalarList& weights = srcWeights_[facei]; + + forAll(faces, i) + { + cop(result[facei], facei, work[faces[i]], weights[i]); + } + } + } + } + else + { + forAll(result, facei) + { + if (srcWeightsSum_[facei] < lowWeightCorrection_) + { + result[facei] = defaultValues[facei]; + } + else + { + const labelList& faces = srcAddress_[facei]; + const scalarList& weights = srcWeights_[facei]; + + forAll(faces, i) + { + cop(result[facei], facei, fld[faces[i]], weights[i]); + } + } + } + } +} + + +template<class Type, class CombineOp> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToSource +( + const Field<Type>& fld, + const CombineOp& cop, + const UList<Type>& defaultValues +) const +{ + auto tresult = tmp<Field<Type>>::New(srcAddress_.size(), Zero); + + interpolateToSource + ( + fld, + multiplyWeightedOp<Type, CombineOp>(cop), + tresult.ref(), + defaultValues + ); + + return tresult; +} + + +template<class Type, class CombineOp> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToSource +( + const tmp<Field<Type>>& tFld, + const CombineOp& cop, + const UList<Type>& defaultValues +) const +{ + return interpolateToSource(tFld(), cop, defaultValues); +} + + +template<class Type, class CombineOp> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToTarget +( + const Field<Type>& fld, + const CombineOp& cop, + const UList<Type>& defaultValues +) const +{ + auto tresult = tmp<Field<Type>>::New(tgtAddress_.size(), Zero); + + interpolateToTarget + ( + fld, + multiplyWeightedOp<Type, CombineOp>(cop), + tresult.ref(), + defaultValues + ); + + return tresult; +} + + +template<class Type, class CombineOp> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToTarget +( + const tmp<Field<Type>>& tFld, + const CombineOp& cop, + const UList<Type>& defaultValues +) const +{ + return interpolateToTarget(tFld(), cop, defaultValues); +} + + +template<class Type> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToSource +( + const Field<Type>& fld, + const UList<Type>& defaultValues +) const +{ + return interpolateToSource(fld, plusEqOp<Type>(), defaultValues); +} + + +template<class Type> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToSource +( + const tmp<Field<Type>>& tFld, + const UList<Type>& defaultValues +) const +{ + return interpolateToSource(tFld(), plusEqOp<Type>(), defaultValues); +} + + +template<class Type> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToTarget +( + const Field<Type>& fld, + const UList<Type>& defaultValues +) const +{ + return interpolateToTarget(fld, plusEqOp<Type>(), defaultValues); +} + + +template<class Type> +Foam::tmp<Foam::Field<Type>> Foam::AMIInterpolation::interpolateToTarget +( + const tmp<Field<Type>>& tFld, + const UList<Type>& defaultValues +) const +{ + return interpolateToTarget(tFld(), plusEqOp<Type>(), defaultValues); +} + + +// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H deleted file mode 100644 index 62efb5d48684080300a7c8d7ce7315740201dec0..0000000000000000000000000000000000000000 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H +++ /dev/null @@ -1,380 +0,0 @@ -/*---------------------------------------------------------------------------*\ - ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox - \\ / O peration | - \\ / A nd | www.openfoam.com - \\/ M anipulation | -------------------------------------------------------------------------------- - Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2016-2018 OpenCFD Ltd. -------------------------------------------------------------------------------- -License - This file is part of OpenFOAM. - - OpenFOAM is free software: you can redistribute it and/or modify it - under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - OpenFOAM is distributed in the hope that it will be useful, but WITHOUT - ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License - for more details. - - You should have received a copy of the GNU General Public License - along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. - -Class - Foam::AMIMethod - -Description - Base class for Arbitrary Mesh Interface (AMI) methods - -SourceFiles - AMIMethod.C - -\*---------------------------------------------------------------------------*/ - -#ifndef AMIMethod_H -#define AMIMethod_H - -#include "className.H" -#include "DynamicList.H" -#include "faceAreaIntersect.H" -#include "indexedOctree.H" -#include "treeDataPrimitivePatch.H" -#include "treeBoundBoxList.H" -#include "runTimeSelectionTables.H" -#include "pointList.H" - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -namespace Foam -{ - -template<class SourcePatch, class TargetPatch> class AMIInterpolation; - -/*---------------------------------------------------------------------------*\ - Class AMIMethod Declaration -\*---------------------------------------------------------------------------*/ - -template<class SourcePatch, class TargetPatch> -class AMIMethod -{ - -private: - - // Private Member Functions - - //- No copy construct - AMIMethod(const AMIMethod&) = delete; - - //- No copy assignment - void operator=(const AMIMethod&) = delete; - - // Parallel operations - - //- Calculate if patches are on multiple processors - label calcDistribution - ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch - ) const; - - label calcOverlappingProcs - ( - const List<treeBoundBoxList>& procBb, - const treeBoundBox& bb, - boolList& overlaps - ) const; - - void distributePatches - ( - const mapDistribute& map, - const TargetPatch& pp, - const globalIndex& gi, - List<faceList>& faces, - List<pointField>& points, - List<labelList>& tgtFaceIDs - ) const; - - void distributeAndMergePatches - ( - const mapDistribute& map, - const TargetPatch& tgtPatch, - const globalIndex& gi, - faceList& tgtFaces, - pointField& tgtPoints, - labelList& tgtFaceIDs - ) const; - - autoPtr<mapDistribute> calcProcMap - ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch - ) const; - - -protected: - - //- Local typedef to octree tree-type - typedef treeDataPrimitivePatch<TargetPatch> treeType; - - - // Protected data - - //- Reference to source patch - const SourcePatch& srcPatch0_; - - //- Reference to target patch - const TargetPatch& tgtPatch0_; - - //- Demand-driven extended target mesh (distributed parallel usage) - autoPtr<TargetPatch> extendedTgtPatchPtr_; - - //- Flag to indicate that the two patches are co-directional and - //- that the orientation of the target patch should be reversed - const bool reverseTarget_; - - //- Flag to indicate that the two patches must be matched/an overlap - //- exists between them - const bool requireMatch_; - - //- Labels of faces that are not overlapped by any target faces - //- (should be empty for correct functioning for fully covered AMIs) - labelList srcNonOverlap_; - - //- Octree used to find face seeds - autoPtr<indexedOctree<treeType>> treePtr_; - - //- Face triangulation mode - const faceAreaIntersect::triangulationMode triMode_; - - //- Label of processor containing all meshes - //- Note: set to -1 if distributed - label singleMeshProc_; - - mutable faceList newTgtFaces_; - mutable pointField newTgtPoints_; - - // Protected Member Functions - - // Helper functions - - //- Create a map that extends tgtPatch so that it covers srcPatch - autoPtr<TargetPatch> createExtendedTgtPatch - ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const globalIndex& globalTgtFaces, - autoPtr<mapDistribute>& mapPtr, - labelList& extendedTgtFaceIDs - ) const; - - //- Check AMI patch coupling - void checkPatches() const; - - //- Initialise and return true if all ok - bool initialise - ( - labelListList& srcAddress, - scalarListList& srcWeights, - labelListList& tgtAddress, - scalarListList& tgtWeights, - label& srcFacei, - label& tgtFacei - ); - - //- Write triangle intersection to OBJ file - void writeIntersectionOBJ - ( - const scalar area, - const face& f1, - const face& f2, - const pointField& f1Points, - const pointField& f2Points - ) const; - - - // Common AMI method functions - - //- Reset the octree for the target patch face search - template<class PatchType> - autoPtr<indexedOctree<treeDataPrimitivePatch<PatchType>>> - createTree - ( - const PatchType& patch - ) const; - - label findTargetFace - ( - const label srcFacei, - const UList<label>& excludeFaces = UList<label>::null(), - const label srcFacePti = -1 - ) const; - - //- Add faces neighbouring facei to the ID list - void appendNbrFaces - ( - const label facei, - const TargetPatch& patch, - const DynamicList<label>& visitedFaces, - DynamicList<label>& faceIDs - ) const; - - //- Helper function to decompose a patch - template<class PatchType> - void triangulatePatch - ( - const PatchType& patch, - List<DynamicList<face>>& tris, - List<scalar>& magSf - ) const; - - -public: - - //- Runtime type information - TypeName("AMIMethod"); - - //- Declare runtime constructor selection table - declareRunTimeSelectionTable - ( - autoPtr, - AMIMethod, - components, - ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool reverseTarget, - const bool requireMatch - ), - ( - srcPatch, - tgtPatch, - triMode, - reverseTarget, - requireMatch - ) - ); - - //- Construct from components - AMIMethod - ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool reverseTarget, - const bool requireMatch - ); - - //- Selector - static autoPtr<AMIMethod> New - ( - const word& methodName, - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool reverseTarget, - const bool requireMatch - ); - - - //- Destructor - virtual ~AMIMethod() = default; - - - // Member Functions - - // Access - - //- Return const access to the source patch - inline const SourcePatch& srcPatch() const; - - //- Return const access to the target patch - inline const TargetPatch& tgtPatch() const; - - //- Return true if the patches are split across multiple processors - inline bool distributed() const; - - //- Labels of faces that are not overlapped by any target faces - // Note: this should be empty for correct functioning - inline const labelList& srcNonOverlap() const; - - //- Flag to indicate that interpolation patches are conformal - virtual bool conformal() const; - - - // Manipulation - - //- Update addressing and weights - virtual bool calculate - ( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei = -1, - label tgtFacei = -1 - ) = 0; - - //- Normalise the weight. Can optionally subset addressing - //- (e.g. for mapNearest) - virtual void normaliseWeights - ( - const bool verbose, - AMIInterpolation<SourcePatch, TargetPatch>& inter - ) = 0; -}; - - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -} // End namespace Foam - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -#define makeAMIMethod(AMIType) \ - \ - typedef AMIMethod<AMIType::sourcePatchType,AMIType::targetPatchType> \ - AMIMethod##AMIType; \ - \ - defineNamedTemplateTypeNameAndDebug(AMIMethod##AMIType, 0); \ - defineTemplateRunTimeSelectionTable(AMIMethod##AMIType, components); - - -#define makeAMIMethodType(AMIType, Method) \ - \ - typedef Method<AMIType::sourcePatchType,AMIType::targetPatchType> \ - Method##AMIType; \ - \ - defineNamedTemplateTypeNameAndDebug(Method##AMIType, 0); \ - \ - AMIMethod<AMIType::sourcePatchType,AMIType::targetPatchType>:: \ - addcomponentsConstructorToTable<Method##AMIType> \ - add##Method##AMIType##ConstructorToTable_; - - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -#include "AMIMethodI.H" - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -#ifdef NoRepository - #include "AMIMethod.C" - #include "AMIMethodNew.C" - #include "AMIMethodParallelOps.C" -#endif - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -#endif - -// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/nearestFaceAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/nearestFaceAMI.C deleted file mode 100644 index 0aa96a35c4469d14ca93efcc4a1ec474b0e12bc2..0000000000000000000000000000000000000000 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/nearestFaceAMI.C +++ /dev/null @@ -1,361 +0,0 @@ -/*---------------------------------------------------------------------------*\ - ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox - \\ / O peration | - \\ / A nd | www.openfoam.com - \\/ M anipulation | -------------------------------------------------------------------------------- - Copyright (C) 2020 OpenCFD Ltd. -------------------------------------------------------------------------------- -License - This file is part of OpenFOAM. - - OpenFOAM is free software: you can redistribute it and/or modify it - under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - OpenFOAM is distributed in the hope that it will be useful, but WITHOUT - ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License - for more details. - - You should have received a copy of the GNU General Public License - along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. - -\*---------------------------------------------------------------------------*/ - -#include "nearestFaceAMI.H" - -// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // - -template<class SourcePatch, class TargetPatch> -Foam::autoPtr<Foam::mapDistribute> -Foam::nearestFaceAMI<SourcePatch, TargetPatch>::calcFaceMap -( - const List<nearestAndDist>& localInfo, - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch -) const -{ - // Generate the list of processor bounding boxes - List<boundBox> procBbs(Pstream::nProcs()); - procBbs[Pstream::myProcNo()] = - boundBox(srcPatch.points(), srcPatch.meshPoints(), true); - Pstream::gatherList(procBbs); - Pstream::scatterList(procBbs); - - // Identify which of my local tgt faces intersect with each processor bb - // within the current match's search distance - const pointField& tgtCcs = tgtPatch.faceCentres(); - List<DynamicList<label>> dynSendMap(Pstream::nProcs()); - - forAll(localInfo, tgtFacei) - { - const scalar r2 = localInfo[tgtFacei].second(); - - // Construct local bounding box to test against processor bb - forAll(procBbs, proci) - { - if (proci != Pstream::myProcNo()) - { - if (procBbs[proci].overlaps(tgtCcs[tgtFacei], r2)) - { - dynSendMap[proci].append(tgtFacei); - } - } - } - } - - // Convert dynamicList to labelList - labelListList sendMap(Pstream::nProcs()); - forAll(sendMap, proci) - { - dynSendMap[proci].shrink(); - sendMap[proci].transfer(dynSendMap[proci]); - - if (debug) - { - Pout<< "send map - to proc " << proci << " sending " - << sendMap[proci].size() << " elements" << endl; - } - } - - return autoPtr<mapDistribute>::New(std::move(sendMap)); -} - - -template<class SourcePatch, class TargetPatch> -Foam::autoPtr<Foam::mapDistribute> -Foam::nearestFaceAMI<SourcePatch, TargetPatch>::calcDistributed -( - const SourcePatch& src, - const TargetPatch& tgt, - labelListList& srcToTgtAddr, - scalarListList& srcToTgtWght -) const -{ - const auto tgtTreePtr = this->createTree(tgt); - const auto& tgtTree = tgtTreePtr(); - - // Create global indexing for each patch - globalIndex globalTgtCells(src.size()); - - // First pass: determine local matches - - // Identify local nearest matches - pointField srcCcs(src.faceCentres()); - - List<nearestAndDist> localInfo(src.size()); - forAll(srcCcs, srcCelli) - { - const point& srcCc = srcCcs[srcCelli]; - - pointIndexHit& test = localInfo[srcCelli].first(); - test = tgtTree.findNearest(srcCc, GREAT); - - if (test.hit()) - { - // With a search radius2 of GREAT all cells should receive a hit - localInfo[srcCelli].second() = magSqr(srcCc - test.hitPoint()); - test.setIndex(globalTgtCells.toGlobal(test.index())); - } - } - - // Second pass: determine remote matches - - autoPtr<mapDistribute> mapPtr = calcFaceMap(localInfo, src, tgt); - mapDistribute& map = mapPtr(); - - List<nearestAndDist> remoteInfo(localInfo); - map.distribute(remoteInfo); - - // Note: re-using srcCcs - map.distribute(srcCcs); - - // Test remote target cells against local source cells - nearestAndDist testInfo; - pointIndexHit& test = testInfo.first(); - forAll(remoteInfo, i) - { - test = tgtTree.findNearest(srcCcs[i], remoteInfo[i].second()); - if (test.hit()) - { - testInfo.first().setIndex - ( - globalTgtCells.toGlobal(test.index()) - ); - testInfo.second() = magSqr(test.hitPoint() - srcCcs[i]); - nearestEqOp()(remoteInfo[i], testInfo); - } - } - - // Send back to originating processor. Choose best if sent to multiple - // processors. Note that afterwards all unused entries have the unique - // value nearestZero (distance < 0). This is used later on to see if - // the sample was sent to any processor. - const nearestAndDist nearestZero(pointIndexHit(), -GREAT); - mapDistributeBase::distribute - ( - Pstream::commsTypes::nonBlocking, - List<labelPair>(0), - src.size(), - map.constructMap(), - map.constructHasFlip(), - map.subMap(), - map.subHasFlip(), - remoteInfo, - nearestEqOp(), - noOp(), // no flipping - nearestZero - ); - - // Third pass: combine local and remote info and filter out any - // connections that are further away than threshold distance squared - srcToTgtAddr.setSize(src.size()); - srcToTgtWght.setSize(src.size()); - forAll(srcToTgtAddr, srcCelli) - { - nearestEqOp()(localInfo[srcCelli], remoteInfo[srcCelli]); - if (localInfo[srcCelli].second() < maxDistance2_) - { - const label tgtCelli = localInfo[srcCelli].first().index(); - srcToTgtAddr[srcCelli] = labelList(1, tgtCelli); - srcToTgtWght[srcCelli] = scalarList(1, 1.0); - } - } - - List<Map<label>> cMap; - return autoPtr<mapDistribute>::New(globalTgtCells, srcToTgtAddr, cMap); -} - - -// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // - -template<class SourcePatch, class TargetPatch> -Foam::nearestFaceAMI<SourcePatch, TargetPatch>::nearestFaceAMI -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool reverseTarget, - const bool requireMatch -) -: - AMIMethod<SourcePatch, TargetPatch> - ( - srcPatch, - tgtPatch, - triMode, - reverseTarget, - requireMatch - ), - maxDistance2_(GREAT) -{} - - -// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // - -template<class SourcePatch, class TargetPatch> -bool Foam::nearestFaceAMI<SourcePatch, TargetPatch>::calculate -( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei, - label tgtFacei -) -{ -bool symmetric_ = true; - - if (this->distributed()) - { - tgtMapPtr = - calcDistributed - ( - this->srcPatch0_, - this->tgtPatch0_, - srcAddress, - srcWeights - ); - - if (symmetric_) - { - srcMapPtr = - calcDistributed - ( - this->tgtPatch0_, - this->srcPatch0_, - tgtAddress, - tgtWeights - ); - } - } - else - { - srcAddress.setSize(this->srcPatch0_.size()); - srcWeights.setSize(this->srcPatch0_.size()); - - if (symmetric_) - { - tgtAddress.setSize(this->tgtPatch0_.size()); - tgtWeights.setSize(this->tgtPatch0_.size()); - } - - const pointField& srcCcs = this->srcPatch0_.faceCentres(); - const pointField& tgtCcs = this->tgtPatch0_.faceCentres(); - - const auto tgtTreePtr = this->createTree(this->tgtPatch0_); - const auto& tgtTree = tgtTreePtr(); - - forAll(srcCcs, srcFacei) - { - const point& srcCc = srcCcs[srcFacei]; - const pointIndexHit hit = tgtTree.findNearest(srcCc, GREAT); - - if - ( - hit.hit() - && (magSqr(srcCc - tgtCcs[hit.index()]) < maxDistance2_) - ) - { - label tgtFacei = hit.index(); - srcAddress[srcFacei] = labelList(1, tgtFacei); - srcWeights[srcFacei] = scalarList(1, 1.0); - - if (symmetric_) - { - tgtAddress[tgtFacei] = labelList(1, srcFacei); - tgtWeights[tgtFacei] = scalarList(1, 1.0); - } - } - else - { - if (debug) - { - WarningInFunction - << "Unable to find target face for source face " - << srcFacei << endl; - } - } - } - - if (symmetric_) - { - const auto srcTreePtr = this->createTree(this->srcPatch0_); - const auto& srcTree = srcTreePtr(); - - // Check that all source cells have connections and populate any - // missing entries - forAll(tgtWeights, tgtCelli) - { - if (tgtAddress[tgtCelli].empty()) - { - const point& tgtCc = tgtCcs[tgtCelli]; - pointIndexHit hit = srcTree.findNearest(tgtCc, GREAT); - - if - ( - hit.hit() - && (magSqr(tgtCc - srcCcs[hit.index()]) < maxDistance2_) - ) - { - tgtAddress[tgtCelli] = labelList(1, hit.index()); - tgtWeights[tgtCelli] = scalarList(1, 1.0); - } - } - else - { - if (debug) - { - WarningInFunction - << "Unable to find source face for target face " - << tgtCelli << endl; - } - } - } - } - } - - return true; -} - - -template<class SourcePatch, class TargetPatch> -void Foam::nearestFaceAMI<SourcePatch, TargetPatch>::normaliseWeights -( - const bool verbose, - AMIInterpolation<SourcePatch, TargetPatch>& inter -) -{ - // Do nothing - weights already 1-to-1 and normalised -} - - -// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C deleted file mode 100644 index 67c21724613da2da09aa5b5e04d544a06e452cf3..0000000000000000000000000000000000000000 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C +++ /dev/null @@ -1,46 +0,0 @@ -/*---------------------------------------------------------------------------*\ - ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox - \\ / O peration | - \\ / A nd | www.openfoam.com - \\/ M anipulation | -------------------------------------------------------------------------------- - Copyright (C) 2013 OpenFOAM Foundation -------------------------------------------------------------------------------- -License - This file is part of OpenFOAM. - - OpenFOAM is free software: you can redistribute it and/or modify it - under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - OpenFOAM is distributed in the hope that it will be useful, but WITHOUT - ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License - for more details. - - You should have received a copy of the GNU General Public License - along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. - -\*---------------------------------------------------------------------------*/ - -#include "AMIPatchToPatchInterpolation.H" -#include "AMIMethod.H" -#include "nearestFaceAMI.H" -#include "faceAreaWeightAMI.H" -#include "partialFaceAreaWeightAMI.H" - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -namespace Foam -{ - makeAMIMethod(AMIPatchToPatchInterpolation); - - makeAMIMethodType(AMIPatchToPatchInterpolation, nearestFaceAMI); - makeAMIMethodType(AMIPatchToPatchInterpolation, faceAreaWeightAMI); - makeAMIMethodType(AMIPatchToPatchInterpolation, partialFaceAreaWeightAMI); -} - - -// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.H index c11f82dc3a0f5295149eebd85bfa8e5f8bd7c062..43a8a3e5cc5023531d6d1a1078dbf3e7994cadd7 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.H @@ -30,14 +30,12 @@ License #define AMIPatchToPatchInterpolation_H #include "AMIInterpolation.H" -#include "primitivePatch.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // namespace Foam { - typedef AMIInterpolation<primitivePatch, primitivePatch> - AMIPatchToPatchInterpolation; + typedef AMIInterpolation AMIPatchToPatchInterpolation; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.C similarity index 63% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C rename to src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.C index 5b3d1bc483bdc3fd8ce07a02e4fc52377822dee8..dc41c132ea60a990f89c16e796022146944c596d 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.C @@ -26,17 +26,23 @@ License \*---------------------------------------------------------------------------*/ -#include "AMIMethod.H" +#include "advancingFrontAMI.H" #include "meshTools.H" #include "mapDistribute.H" #include "unitConversion.H" #include "findNearestMaskedOp.H" +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(advancingFrontAMI, 0); +} + // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -void Foam::AMIMethod<SourcePatch, TargetPatch>::checkPatches() const +void Foam::advancingFrontAMI::checkPatches() const { const auto& src = srcPatch(); const auto& tgt = tgtPatch(); @@ -75,101 +81,68 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::checkPatches() const } -template<class SourcePatch, class TargetPatch> -Foam::autoPtr<TargetPatch> -Foam::AMIMethod<SourcePatch, TargetPatch>::createExtendedTgtPatch -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const globalIndex& globalTgtFaces, - autoPtr<mapDistribute>& mapPtr, - labelList& extendedTgtFaceIDs -) const +void Foam::advancingFrontAMI::createExtendedTgtPatch() { - // Create a representation of the src mesh that is extended to overlap the - // tgt mesh - // Create processor map of extended cells. This map gets (possibly // remote) cells from the src mesh such that they (together) cover // all of tgt - mapPtr = calcProcMap(srcPatch, tgtPatch); - const mapDistribute& map = mapPtr(); - - // Create new target patch that fully encompasses source patch - - // Faces and points - // faceList newTgtFaces_; - // pointField newTgtPoints_; + extendedTgtMapPtr_.reset(calcProcMap(srcPatch0(), tgtPatch0())); + const mapDistribute& map = extendedTgtMapPtr_(); - // Original faces from tgtPatch (in globalIndexing since might be - // remote) + // Original faces from tgtPatch + // Note: in globalIndexing since might be remote + globalIndex globalTgtFaces(tgtPatch0().size()); distributeAndMergePatches ( map, - tgtPatch, + tgtPatch0(), globalTgtFaces, - newTgtFaces_, - newTgtPoints_, - extendedTgtFaceIDs + extendedTgtFaces_, + extendedTgtPoints_, + extendedTgtFaceIDs_ ); - return autoPtr<TargetPatch>::New + // Create a representation of the tgt patch that is extended to overlap + // the src patch + extendedTgtPatchPtr_.reset ( - SubList<face> + autoPtr<primitivePatch>::New ( - newTgtFaces_, - newTgtFaces_.size() - ), - newTgtPoints_ + SubList<face>(extendedTgtFaces_), + extendedTgtPoints_ + ) ); } -template<class SourcePatch, class TargetPatch> -bool Foam::AMIMethod<SourcePatch, TargetPatch>::initialise -( - labelListList& srcAddress, - scalarListList& srcWeights, - labelListList& tgtAddress, - scalarListList& tgtWeights, - label& srcFacei, - label& tgtFacei -) +bool Foam::advancingFrontAMI::initialiseWalk(label& srcFacei, label& tgtFacei) { - const auto& src = srcPatch(); - const auto& tgt = tgtPatch(); - - checkPatches(); + const auto& src = this->srcPatch(); + const auto& tgt = this->tgtPatch(); - // Set initial sizes for weights and addressing - must be done even if - // returns false below - srcAddress.setSize(src.size()); - srcWeights.setSize(srcAddress.size()); - tgtAddress.setSize(tgt.size()); - tgtWeights.setSize(tgtAddress.size()); + bool foundFace = false; // Check that patch sizes are valid if (!src.size()) { - return false; + return foundFace; } else if (!tgt.size()) { WarningInFunction << src.size() << " source faces but no target faces" << endl; - return false; + return foundFace; } // Reset the octree - treePtr_.reset(createTree<TargetPatch>(tgtPatch())); + treePtr_.reset(createTree(tgt)); // Find initial face match using brute force/octree search if ((srcFacei == -1) || (tgtFacei == -1)) { srcFacei = 0; tgtFacei = 0; - bool foundFace = false; forAll(src, facei) { tgtFacei = findTargetFace(facei); @@ -190,7 +163,7 @@ bool Foam::AMIMethod<SourcePatch, TargetPatch>::initialise << abort(FatalError); } - return false; + return foundFace; } } @@ -203,8 +176,7 @@ bool Foam::AMIMethod<SourcePatch, TargetPatch>::initialise } -template<class SourcePatch, class TargetPatch> -void Foam::AMIMethod<SourcePatch, TargetPatch>::writeIntersectionOBJ +void Foam::advancingFrontAMI::writeIntersectionOBJ ( const scalar area, const face& f1, @@ -255,37 +227,7 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::writeIntersectionOBJ } -template<class SourcePatch, class TargetPatch> -template<class PatchType> -Foam::autoPtr<Foam::indexedOctree<Foam::treeDataPrimitivePatch<PatchType>>> -Foam::AMIMethod<SourcePatch, TargetPatch>::createTree -( - const PatchType& patch -) const -{ - typedef treeDataPrimitivePatch<PatchType> PatchTreeType; - - treeBoundBox bb(patch.points(), patch.meshPoints()); - bb.inflate(0.01); - - return autoPtr<indexedOctree<PatchTreeType>>::New - ( - PatchTreeType - ( - false, - patch, - indexedOctree<PatchTreeType>::perturbTol() - ), - bb, // overall search domain - 8, // maxLevel - 10, // leaf size - 3.0 // duplicity - ); -} - - -template<class SourcePatch, class TargetPatch> -Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::findTargetFace +Foam::label Foam::advancingFrontAMI::findTargetFace ( const label srcFacei, const UList<label>& excludeFaces, @@ -299,9 +241,9 @@ Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::findTargetFace const pointField& srcPts = src.points(); const face& srcFace = src[srcFacei]; - findNearestMaskedOp<TargetPatch> fnOp(*treePtr_, excludeFaces); + findNearestMaskedOp<primitivePatch> fnOp(*treePtr_, excludeFaces); - boundBox bb(srcPts, srcFace, false); + const boundBox bb(srcPts, srcFace, false); const point srcPt = srcFacePti == -1 ? bb.centre() : srcPts[srcFace[srcFacePti]]; @@ -325,11 +267,10 @@ Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::findTargetFace } -template<class SourcePatch, class TargetPatch> -void Foam::AMIMethod<SourcePatch, TargetPatch>::appendNbrFaces +void Foam::advancingFrontAMI::appendNbrFaces ( const label facei, - const TargetPatch& patch, + const primitivePatch& patch, const DynamicList<label>& visitedFaces, DynamicList<label>& faceIDs ) const @@ -341,30 +282,8 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::appendNbrFaces // Filter out faces already visited from face neighbours for (const label nbrFacei : nbrFaces) { - bool valid = true; - for (const label visitedFacei : visitedFaces) - { - if (nbrFacei == visitedFacei) - { - valid = false; - break; - } - } - - if (valid) - { - for (const label testFacei : faceIDs) - { - if (nbrFacei == testFacei) - { - valid = false; - break; - } - } - } - // Prevent addition of face if it is not on the same plane-ish - if (valid) + if (!visitedFaces.found(nbrFacei) && !faceIDs.found(nbrFacei)) { const vector& n1 = patch.faceNormals()[facei]; const vector& n2 = patch.faceNormals()[nbrFacei]; @@ -380,11 +299,9 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::appendNbrFaces } -template<class SourcePatch, class TargetPatch> -template<class PatchType> -void Foam::AMIMethod<SourcePatch, TargetPatch>::triangulatePatch +void Foam::advancingFrontAMI::triangulatePatch ( - const PatchType& patch, + const primitivePatch& patch, List<DynamicList<face>>& tris, List<scalar>& magSf ) const @@ -396,6 +313,8 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::triangulatePatch // Using methods that index into existing points forAll(patch, facei) { + tris[facei].clear(); + switch (triMode_) { case faceAreaIntersect::tmFan: @@ -428,34 +347,114 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::triangulatePatch // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -Foam::AMIMethod<SourcePatch, TargetPatch>::AMIMethod +Foam::advancingFrontAMI::advancingFrontAMI +( + const dictionary& dict, + const bool reverseTarget +) +: + AMIInterpolation(dict, reverseTarget), + srcTris_(), + tgtTris_(), + extendedTgtPatchPtr_(nullptr), + extendedTgtFaces_(), + extendedTgtPoints_(), + extendedTgtFaceIDs_(), + extendedTgtMapPtr_(nullptr), + srcNonOverlap_(), + triMode_ + ( + faceAreaIntersect::triangulationModeNames_.getOrDefault + ( + "triMode", + dict, + faceAreaIntersect::tmMesh + ) + ) +{} + + +Foam::advancingFrontAMI::advancingFrontAMI ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const bool requireMatch, const bool reverseTarget, - const bool requireMatch + const scalar lowWeightCorrection, + const faceAreaIntersect::triangulationMode triMode ) : - srcPatch0_(srcPatch), - tgtPatch0_(tgtPatch), + AMIInterpolation(requireMatch, reverseTarget, lowWeightCorrection), + srcTris_(), + tgtTris_(), extendedTgtPatchPtr_(nullptr), - reverseTarget_(reverseTarget), - requireMatch_(requireMatch), + extendedTgtFaces_(), + extendedTgtPoints_(), + extendedTgtFaceIDs_(), + extendedTgtMapPtr_(nullptr), srcNonOverlap_(), - triMode_(triMode), - singleMeshProc_(calcDistribution(srcPatch, tgtPatch)) -{ - // Note: setting srcMagSf and tgtMagSf to 1 by default for 1-to-1 methods - // - others will need to overwrite as necessary -} + triMode_(triMode) +{} + + +Foam::advancingFrontAMI::advancingFrontAMI(const advancingFrontAMI& ami) +: + AMIInterpolation(ami), + srcTris_(), + tgtTris_(), + extendedTgtPatchPtr_(nullptr), + extendedTgtFaces_(), + extendedTgtPoints_(), + extendedTgtFaceIDs_(), + extendedTgtMapPtr_(nullptr), + srcNonOverlap_(), + triMode_(ami.triMode_) +{} // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -bool Foam::AMIMethod<SourcePatch, TargetPatch>::conformal() const +bool Foam::advancingFrontAMI::calculate +( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr +) +{ + if (AMIInterpolation::calculate(srcPatch, tgtPatch, surfPtr)) + { + // Create a representation of the target patch that covers the source patch + if (distributed()) + { + createExtendedTgtPatch(); + } + + const auto& src = this->srcPatch(); + const auto& tgt = this->tgtPatch(); + + // Initialise area magnitudes + srcMagSf_.setSize(src.size(), 1.0); + tgtMagSf_.setSize(tgt.size(), 1.0); + + // Source and target patch triangulations + triangulatePatch(src, srcTris_, srcMagSf_); + triangulatePatch(tgt, tgtTris_, tgtMagSf_); + + checkPatches(); + + // Set initial sizes for weights and addressing - must be done even if + // returns false below + srcAddress_.setSize(src.size()); + srcWeights_.setSize(src.size()); + tgtAddress_.setSize(tgt.size()); + tgtWeights_.setSize(tgt.size()); + + return true; + } + + return false; +} + + +bool Foam::advancingFrontAMI::conformal() const { return true; } diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.H new file mode 100644 index 0000000000000000000000000000000000000000..84031daa733f3e53df6d86d6bdaca30a458ee3db --- /dev/null +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.H @@ -0,0 +1,271 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2013-2016 OpenFOAM Foundation + Copyright (C) 2016-2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. + +Class + Foam::advancingFrontAMI + +Description + Base class for Arbitrary Mesh Interface (AMI) methods + +SourceFiles + advancingFrontAMI.C + +\*---------------------------------------------------------------------------*/ + +#ifndef advancingFrontAMI_H +#define advancingFrontAMI_H + +#include "className.H" +#include "DynamicList.H" +#include "faceAreaIntersect.H" +#include "pointList.H" +#include "AMIInterpolation.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class advancingFrontAMI Declaration +\*---------------------------------------------------------------------------*/ + +class advancingFrontAMI +: + public AMIInterpolation +{ + +private: + + // Private Member Functions + + //- No copy assignment + void operator=(const advancingFrontAMI&) = delete; + + + // Parallel operations + + label calcOverlappingProcs + ( + const List<treeBoundBoxList>& procBb, + const treeBoundBox& bb, + boolList& overlaps + ) const; + + void distributePatches + ( + const mapDistribute& map, + const primitivePatch& pp, + const globalIndex& gi, + List<faceList>& faces, + List<pointField>& points, + List<labelList>& tgtFaceIDs + ) const; + + void distributeAndMergePatches + ( + const mapDistribute& map, + const primitivePatch& tgtPatch, + const globalIndex& gi, + faceList& tgtFaces, + pointField& tgtPoints, + labelList& tgtFaceIDs + ) const; + + autoPtr<mapDistribute> calcProcMap + ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch + ) const; + + +protected: + + // Protected data + + //- Storage for src-side triangle decomposition + List<DynamicList<face>> srcTris_; + + //- Storage for tgt-side triangle decomposition + List<DynamicList<face>> tgtTris_; + + //- Demand-driven extended target mesh (distributed parallel usage) + autoPtr<primitivePatch> extendedTgtPatchPtr_; + + //- Extended patch faces + faceList extendedTgtFaces_; + + //- Extended patch points + pointField extendedTgtPoints_; + + //- Extended patch face IDs + labelList extendedTgtFaceIDs_; + + //- Extended patch map + autoPtr<mapDistribute> extendedTgtMapPtr_; + + //- Labels of faces that are not overlapped by any target faces + //- (should be empty for correct functioning for fully covered AMIs) + labelList srcNonOverlap_; + + //- Octree used to find face seeds + autoPtr<indexedOctree<treeType>> treePtr_; + + //- Face triangulation mode + const faceAreaIntersect::triangulationMode triMode_; + + + // Protected Member Functions + + // Helper functions + + //- Create a map that extends tgtPatch so that it covers srcPatch + void createExtendedTgtPatch(); + + //- Check AMI patch coupling + void checkPatches() const; + + virtual bool calculate + ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr = nullptr + ); + + //- Initialise walk and return true if all ok + bool initialiseWalk + ( + label& srcFacei, + label& tgtFacei + ); + + //- Write triangle intersection to OBJ file + void writeIntersectionOBJ + ( + const scalar area, + const face& f1, + const face& f2, + const pointField& f1Points, + const pointField& f2Points + ) const; + + + // Common AMI method functions + + label findTargetFace + ( + const label srcFacei, + const UList<label>& excludeFaces = UList<label>::null(), + const label srcFacePti = -1 + ) const; + + //- Add faces neighbouring facei to the ID list + void appendNbrFaces + ( + const label facei, + const primitivePatch& patch, + const DynamicList<label>& visitedFaces, + DynamicList<label>& faceIDs + ) const; + + //- Helper function to decompose a patch + void triangulatePatch + ( + const primitivePatch& patch, + List<DynamicList<face>>& tris, + List<scalar>& magSf + ) const; + + +public: + + //- Runtime type information + TypeName("advancingFrontAMI"); + + // Constructors + + //- Construct from components + advancingFrontAMI + ( + const dictionary& dict, + const bool reverseTarget + ); + + //- Construct from components + advancingFrontAMI + ( + const bool requireMatch = true, + const bool reverseTarget = false, + const scalar lowWeightCorrection = -1, + const faceAreaIntersect::triangulationMode triMode = + faceAreaIntersect::tmMesh + ); + + //- Construct as copy + advancingFrontAMI(const advancingFrontAMI& ami); + + //- Construct and return a clone + virtual autoPtr<AMIInterpolation> clone() const + { + return autoPtr<AMIInterpolation>(new advancingFrontAMI(*this)); + } + + + //- Destructor + virtual ~advancingFrontAMI() = default; + + + // Member Functions + + //- Return const access to the source patch + inline const primitivePatch& srcPatch() const; + + //- Return const access to the target patch + inline const primitivePatch& tgtPatch() const; + + //- Labels of faces that are not overlapped by any target faces + // Note: this should be empty for correct functioning + inline const labelList& srcNonOverlap() const; + + //- Flag to indicate that interpolation patches are conformal, i.e. + //- should fully cover each other + virtual bool conformal() const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "advancingFrontAMII.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMII.H similarity index 68% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H rename to src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMII.H index 8c740cdcb1c5ce837f2c21b2dbf0d44b2ff1058f..e0441345383976a94c5f1b2ddf26ff89b6ef8ab2 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMII.H @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2013 OpenFOAM Foundation + Copyright (C) 2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -25,37 +25,38 @@ License \*---------------------------------------------------------------------------*/ -template<class SourcePatch, class TargetPatch> -inline const SourcePatch& -Foam::AMIMethod<SourcePatch, TargetPatch>::srcPatch() const +inline const Foam::primitivePatch& Foam::advancingFrontAMI::srcPatch() const { - return srcPatch0_; + if (!tsrcPatch0_.valid()) + { + FatalErrorInFunction + << "tsrcPatch0_ not allocated" + << abort(FatalError); + } + + return tsrcPatch0_(); } -template<class SourcePatch, class TargetPatch> -inline const TargetPatch& -Foam::AMIMethod<SourcePatch, TargetPatch>::tgtPatch() const +inline const Foam::primitivePatch& Foam::advancingFrontAMI::tgtPatch() const { if (extendedTgtPatchPtr_) { return extendedTgtPatchPtr_(); } - return tgtPatch0_; -} - + if (!ttgtPatch0_.valid()) + { + FatalErrorInFunction + << "srcPatch0Ptr not allocated" + << abort(FatalError); + } -template<class SourcePatch, class TargetPatch> -bool Foam::AMIMethod<SourcePatch, TargetPatch>::distributed() const -{ - return singleMeshProc_ == -1; + return ttgtPatch0_(); } -template<class SourcePatch, class TargetPatch> -inline const Foam::labelList& -Foam::AMIMethod<SourcePatch, TargetPatch>::srcNonOverlap() const +inline const Foam::labelList& Foam::advancingFrontAMI::srcNonOverlap() const { return srcNonOverlap_; } diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodParallelOps.C b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMIParallelOps.C similarity index 85% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodParallelOps.C rename to src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMIParallelOps.C index 1066943f293fe26b9cacea3288998abbc22162c4..e67d9e6e400cde827065477af3a849195c9e272e 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodParallelOps.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/advancingFrontAMI/advancingFrontAMIParallelOps.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2018-2020 OpenCFD Ltd. + Copyright (C) 2018 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -25,60 +25,15 @@ License along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. \*---------------------------------------------------------------------------*/ + +#include "advancingFrontAMI.H" #include "mergePoints.H" #include "mapDistribute.H" #include "AABBTree.H" // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::calcDistribution -( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch -) const -{ - label proci = 0; - - if (Pstream::parRun()) - { - labelList facesPresentOnProc(Pstream::nProcs(), Zero); - if ((srcPatch.size() > 0) || (tgtPatch.size() > 0)) - { - facesPresentOnProc[Pstream::myProcNo()] = 1; - } - else - { - facesPresentOnProc[Pstream::myProcNo()] = 0; - } - - Pstream::gatherList(facesPresentOnProc); - Pstream::scatterList(facesPresentOnProc); - - label nHaveFaces = sum(facesPresentOnProc); - - if (nHaveFaces > 1) - { - proci = -1; - DebugInFunction - << "AMI split across multiple processors" << endl; - } - else if (nHaveFaces == 1) - { - proci = facesPresentOnProc.find(1); - DebugInFunction - << "AMI local to processor" << proci << endl; - } - } - - - // Either not parallel or no faces on any processor - return proci; -} - - -template<class SourcePatch, class TargetPatch> -Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::calcOverlappingProcs +Foam::label Foam::advancingFrontAMI::calcOverlappingProcs ( const List<treeBoundBoxList>& procBb, const treeBoundBox& bb, @@ -109,11 +64,10 @@ Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::calcOverlappingProcs } -template<class SourcePatch, class TargetPatch> -void Foam::AMIMethod<SourcePatch, TargetPatch>::distributePatches +void Foam::advancingFrontAMI::distributePatches ( const mapDistribute& map, - const TargetPatch& pp, + const primitivePatch& pp, const globalIndex& gi, List<faceList>& faces, List<pointField>& points, @@ -194,12 +148,10 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::distributePatches } -template<class SourcePatch, class TargetPatch> -void Foam::AMIMethod<SourcePatch, TargetPatch>:: -distributeAndMergePatches +void Foam::advancingFrontAMI::distributeAndMergePatches ( const mapDistribute& map, - const TargetPatch& tgtPatch, + const primitivePatch& tgtPatch, const globalIndex& gi, faceList& tgtFaces, pointField& tgtPoints, @@ -308,12 +260,10 @@ distributeAndMergePatches } -template<class SourcePatch, class TargetPatch> -Foam::autoPtr<Foam::mapDistribute> -Foam::AMIMethod<SourcePatch, TargetPatch>::calcProcMap +Foam::autoPtr<Foam::mapDistribute> Foam::advancingFrontAMI::calcProcMap ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch ) const { // Get decomposition of patch diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.C similarity index 58% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C rename to src/meshTools/AMIInterpolation/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.C index 8eba124a3a9de5870034939ad7f72cecb2d583b5..2be76626bdf415b70290fe4d9f7989a207e40794 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.C @@ -29,41 +29,20 @@ License #include "faceAreaWeightAMI.H" #include "profiling.H" #include "OBJstream.H" +#include "addToRunTimeSelectionTable.H" -// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * // +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::initGeom -( - const globalIndex& globalTgtFaces, - labelList& extendedTgtFaceIDs -) +namespace Foam { - // Create a representation of the target patch that covers the source patch - if (this->distributed()) - { - this->extendedTgtPatchPtr_ = - this->createExtendedTgtPatch - ( - this->srcPatch0_, - this->tgtPatch0_, - globalTgtFaces, - extendedTgtMapPtr_, - extendedTgtFaceIDs - ); - } - - const auto& src = this->srcPatch(); - const auto& tgt = this->tgtPatch(); - - // Initialise area magnitudes - srcMagSf_.setSize(src.size(), 1.0); - tgtMagSf_.setSize(tgt.size(), 1.0); + defineTypeNameAndDebug(faceAreaWeightAMI, 0); + addToRunTimeSelectionTable(AMIInterpolation, faceAreaWeightAMI, dict); + addToRunTimeSelectionTable(AMIInterpolation, faceAreaWeightAMI, component); +} - // Source and target patch triangulations - this->triangulatePatch(src, srcTris_, srcMagSf_); - this->triangulatePatch(tgt, tgtTris_, tgtMagSf_); +// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * // +/* if (debug) { static label nAMI = 0; @@ -99,11 +78,10 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::initGeom ++nAMI; } -} +*/ -template<class SourcePatch, class TargetPatch> -void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing +void Foam::faceAreaWeightAMI::calcAddressing ( List<DynamicList<label>>& srcAddr, List<DynamicList<scalar>>& srcWght, @@ -116,27 +94,28 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing { addProfiling(ami, "faceAreaWeightAMI::calcAddressing"); - // construct weights and addressing + // Construct weights and addressing // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ label nFacesRemaining = srcAddr.size(); - // list of tgt face neighbour faces + // List of tgt face neighbour faces DynamicList<label> nbrFaces(10); - // list of faces currently visited for srcFacei to avoid multiple hits + // List of faces currently visited for srcFacei to avoid multiple hits DynamicList<label> visitedFaces(10); - // list to keep track of tgt faces used to seed src faces + // List to keep track of tgt faces used to seed src faces labelList seedFaces(nFacesRemaining, -1); seedFaces[srcFacei] = tgtFacei; - // list to keep track of whether src face can be mapped + // List to keep track of whether src face can be mapped bitSet mapFlag(nFacesRemaining, true); - // reset starting seed + // Reset starting seed label startSeedi = 0; + bool continueWalk = true; DynamicList<label> nonOverlapFaces; do { @@ -161,34 +140,30 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing mapFlag.unset(srcFacei); - nFacesRemaining--; - if (!faceProcessed) { nonOverlapFaces.append(srcFacei); } - // choose new src face from current src face neighbour - if (nFacesRemaining > 0) - { - setNextFaces - ( - startSeedi, - srcFacei, - tgtFacei, - mapFlag, - seedFaces, - visitedFaces - ); - } - } while (nFacesRemaining > 0); + // Choose new src face from current src face neighbour + continueWalk = setNextFaces + ( + startSeedi, + srcFacei, + tgtFacei, + mapFlag, + seedFaces, + visitedFaces, + requireMatch_ && (lowWeightCorrection_ < 0) + // pass in nonOverlapFaces for failed tree search? + ); + } while (continueWalk); - this->srcNonOverlap_.transfer(nonOverlapFaces); + srcNonOverlap_.transfer(nonOverlapFaces); } -template<class SourcePatch, class TargetPatch> -bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace +bool Foam::faceAreaWeightAMI::processSourceFace ( const label srcFacei, const label tgtStartFacei, @@ -213,15 +188,11 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace return false; } + const auto& tgtPatch = this->tgtPatch(); + // append initial target face and neighbours nbrFaces.append(tgtStartFacei); - this->appendNbrFaces - ( - tgtStartFacei, - this->tgtPatch(), - visitedFaces, - nbrFaces - ); + appendNbrFaces(tgtStartFacei, tgtPatch, visitedFaces, nbrFaces); bool faceProcessed = false; @@ -238,11 +209,7 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace calcInterArea(srcFacei, tgtFacei, interArea, interCentroid); // store when intersection fractional area > tolerance - if - ( - interArea/this->srcMagSf_[srcFacei] - > faceAreaIntersect::tolerance() - ) + if (interArea/srcMagSf_[srcFacei] > faceAreaIntersect::tolerance()) { srcAddr[srcFacei].append(tgtFacei); srcWght[srcFacei].append(interArea); @@ -251,13 +218,7 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace tgtAddr[tgtFacei].append(srcFacei); tgtWght[tgtFacei].append(interArea); - this->appendNbrFaces - ( - tgtFacei, - this->tgtPatch(), - visitedFaces, - nbrFaces - ); + appendNbrFaces(tgtFacei, tgtPatch, visitedFaces, nbrFaces); faceProcessed = true; @@ -275,8 +236,7 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace } -template<class SourcePatch, class TargetPatch> -void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces +bool Foam::faceAreaWeightAMI::setNextFaces ( label& startSeedi, label& srcFacei, @@ -289,12 +249,18 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces { addProfiling(ami, "faceAreaWeightAMI::setNextFaces"); + if (mapFlag.count() == 0) + { + // No more faces to map + return false; + } + const labelList& srcNbrFaces = this->srcPatch().faceFaces()[srcFacei]; - // initialise tgtFacei + // Initialise tgtFacei tgtFacei = -1; - // set possible seeds for later use + // Set possible seeds for later use bool valuesSet = false; for (label faceS: srcNbrFaces) { @@ -303,7 +269,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces for (label faceT : visitedFaces) { const scalar threshold = - this->srcMagSf_[faceS]*faceAreaIntersect::tolerance(); + srcMagSf_[faceS]*faceAreaIntersect::tolerance(); // Store when intersection fractional area > threshold if (overlaps(faceS, faceT, threshold)) @@ -321,73 +287,73 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces } } - // set next src and tgt faces if not set above if (valuesSet) { - return; + return true; } - else + + // Set next src and tgt faces if not set above + // - try to use existing seed + label facei = startSeedi; + if (!mapFlag.test(startSeedi)) { - // try to use existing seed - label facei = startSeedi; - if (!mapFlag.test(startSeedi)) - { - facei = mapFlag.find_next(facei); - } - const label startSeedi0 = facei; + facei = mapFlag.find_next(facei); + } + const label startSeedi0 = facei; - bool foundNextSeed = false; - while (facei != -1) + bool foundNextSeed = false; + while (facei != -1) + { + if (!foundNextSeed) { - if (!foundNextSeed) - { - startSeedi = facei; - foundNextSeed = true; - } - - if (seedFaces[facei] != -1) - { - srcFacei = facei; - tgtFacei = seedFaces[facei]; - - return; - } - - facei = mapFlag.find_next(facei); + startSeedi = facei; + foundNextSeed = true; } - // perform new search to find match - if (debug) + if (seedFaces[facei] != -1) { - Pout<< "Advancing front stalled: searching for new " - << "target face" << endl; + srcFacei = facei; + tgtFacei = seedFaces[facei]; + + return true; } - facei = startSeedi0; - while (facei != -1) - { - srcFacei = facei; - tgtFacei = this->findTargetFace(srcFacei, visitedFaces); + facei = mapFlag.find_next(facei); + } - if (tgtFacei >= 0) - { - return; - } + // Perform new search to find match + if (debug) + { + Pout<< "Advancing front stalled: searching for new " + << "target face" << endl; + } - facei = mapFlag.find_next(facei); - } + facei = startSeedi0; + while (facei != -1) + { + srcFacei = facei; + tgtFacei = findTargetFace(srcFacei, visitedFaces); - if (errorOnNotFound) + if (tgtFacei >= 0) { - FatalErrorInFunction - << "Unable to set source and target faces" << abort(FatalError); + return true; } + + facei = mapFlag.find_next(facei); } + + if (errorOnNotFound) + { + FatalErrorInFunction + << "Unable to set target face for source face " << srcFacei + << abort(FatalError); + } + + return false; } -template<class SourcePatch, class TargetPatch> -void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcInterArea +void Foam::faceAreaWeightAMI::calcInterArea ( const label srcFacei, const label tgtFacei, @@ -397,46 +363,48 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcInterArea { addProfiling(ami, "faceAreaWeightAMI::interArea"); - const pointField& srcPoints = this->srcPatch().points(); - const pointField& tgtPoints = this->tgtPatch().points(); - - // references to candidate faces - const face& src = this->srcPatch()[srcFacei]; - const face& tgt = this->tgtPatch()[tgtFacei]; - - // quick reject if either face has zero area + // Quick reject if either face has zero area if ( - (this->srcMagSf_[srcFacei] < ROOTVSMALL) - || (this->tgtMagSf_[tgtFacei] < ROOTVSMALL) + (srcMagSf_[srcFacei] < ROOTVSMALL) + || (tgtMagSf_[tgtFacei] < ROOTVSMALL) ) { return; } - // create intersection object + const auto& srcPatch = this->srcPatch(); + const auto& tgtPatch = this->tgtPatch(); + + const pointField& srcPoints = srcPatch.points(); + const pointField& tgtPoints = tgtPatch.points(); + + // Create intersection object faceAreaIntersect inter ( srcPoints, tgtPoints, - this->srcTris_[srcFacei], - this->tgtTris_[tgtFacei], - this->reverseTarget_, - AMIInterpolation<SourcePatch, TargetPatch>::cacheIntersections_ + srcTris_[srcFacei], + tgtTris_[tgtFacei], + reverseTarget_, + AMIInterpolation::cacheIntersections_ ); - // crude resultant norm - vector n(-this->srcPatch().faceNormals()[srcFacei]); - if (this->reverseTarget_) + // Crude resultant norm + vector n(-srcPatch.faceNormals()[srcFacei]); + if (reverseTarget_) { - n -= this->tgtPatch().faceNormals()[tgtFacei]; + n -= tgtPatch.faceNormals()[tgtFacei]; } else { - n += this->tgtPatch().faceNormals()[tgtFacei]; + n += tgtPatch.faceNormals()[tgtFacei]; } scalar magN = mag(n); + const face& src = srcPatch[srcFacei]; + const face& tgt = tgtPatch[tgtFacei]; + if (magN > ROOTVSMALL) { inter.calc(src, tgt, n/magN, area, centroid); @@ -451,11 +419,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcInterArea << endl; } - if - ( - AMIInterpolation<SourcePatch, TargetPatch>::cacheIntersections_ - && debug - ) + if (AMIInterpolation::cacheIntersections_ && debug) { static OBJstream tris("intersectionTris.obj"); const auto& triPts = inter.triangles(); @@ -467,58 +431,59 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcInterArea if ((debug > 1) && (area > 0)) { - this->writeIntersectionOBJ(area, src, tgt, srcPoints, tgtPoints); + writeIntersectionOBJ(area, src, tgt, srcPoints, tgtPoints); } } -template<class SourcePatch, class TargetPatch> -bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::overlaps +bool Foam::faceAreaWeightAMI::overlaps ( const label srcFacei, const label tgtFacei, const scalar threshold ) const { - const pointField& srcPoints = this->srcPatch().points(); - const pointField& tgtPoints = this->tgtPatch().points(); - - // references to candidate faces - const face& src = this->srcPatch()[srcFacei]; - const face& tgt = this->tgtPatch()[tgtFacei]; - - // quick reject if either face has zero area + // Quick reject if either face has zero area if ( - (this->srcMagSf_[srcFacei] < ROOTVSMALL) - || (this->tgtMagSf_[tgtFacei] < ROOTVSMALL) + (srcMagSf_[srcFacei] < ROOTVSMALL) + || (tgtMagSf_[tgtFacei] < ROOTVSMALL) ) { return false; } + const auto& srcPatch = this->srcPatch(); + const auto& tgtPatch = this->tgtPatch(); + + const pointField& srcPoints = srcPatch.points(); + const pointField& tgtPoints = tgtPatch.points(); + faceAreaIntersect inter ( srcPoints, tgtPoints, - this->srcTris_[srcFacei], - this->tgtTris_[tgtFacei], - this->reverseTarget_, - AMIInterpolation<SourcePatch, TargetPatch>::cacheIntersections_ + srcTris_[srcFacei], + tgtTris_[tgtFacei], + reverseTarget_, + AMIInterpolation::cacheIntersections_ ); - // crude resultant norm - vector n(-this->srcPatch().faceNormals()[srcFacei]); - if (this->reverseTarget_) + // Crude resultant norm + vector n(-srcPatch.faceNormals()[srcFacei]); + if (reverseTarget_) { - n -= this->tgtPatch().faceNormals()[tgtFacei]; + n -= tgtPatch.faceNormals()[tgtFacei]; } else { - n += this->tgtPatch().faceNormals()[tgtFacei]; + n += tgtPatch.faceNormals()[tgtFacei]; } scalar magN = mag(n); + const face& src = srcPatch[srcFacei]; + const face& tgt = tgtPatch[tgtFacei]; + if (magN > ROOTVSMALL) { return inter.overlaps(src, tgt, n/magN, threshold); @@ -537,9 +502,7 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::overlaps } -template<class SourcePatch, class TargetPatch> -void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>:: -restartUncoveredSourceFace +void Foam::faceAreaWeightAMI::restartUncoveredSourceFace ( List<DynamicList<label>>& srcAddr, List<DynamicList<scalar>>& srcWght, @@ -555,27 +518,29 @@ restartUncoveredSourceFace label nBelowMinWeight = 0; const scalar minWeight = 0.95; - // list of tgt face neighbour faces + // List of tgt face neighbour faces DynamicList<label> nbrFaces(10); - // list of faces currently visited for srcFacei to avoid multiple hits + // List of faces currently visited for srcFacei to avoid multiple hits DynamicList<label> visitedFaces(10); + const auto& srcPatch = this->srcPatch(); + forAll(srcWght, srcFacei) { const scalar s = sum(srcWght[srcFacei]); - const scalar t = s/this->srcMagSf_[srcFacei]; + const scalar t = s/srcMagSf_[srcFacei]; if (t < minWeight) { ++nBelowMinWeight; - const face& f = this->srcPatch()[srcFacei]; + const face& f = srcPatch[srcFacei]; forAll(f, fpi) { const label tgtFacei = - this->findTargetFace(srcFacei, srcAddr[srcFacei], fpi); + findTargetFace(srcFacei, srcAddr[srcFacei], fpi); if (tgtFacei != -1) { @@ -613,77 +578,74 @@ restartUncoveredSourceFace // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::faceAreaWeightAMI +Foam::faceAreaWeightAMI::faceAreaWeightAMI +( + const dictionary& dict, + const bool reverseTarget +) +: + advancingFrontAMI(dict, reverseTarget), + restartUncoveredSourceFace_ + ( + dict.getOrDefault("restartUncoveredSourceFace", true) + ) +{} + + +Foam::faceAreaWeightAMI::faceAreaWeightAMI ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, - const bool reverseTarget, const bool requireMatch, + const bool reverseTarget, + const scalar lowWeightCorrection, + const faceAreaIntersect::triangulationMode triMode, const bool restartUncoveredSourceFace ) : - AMIMethod<SourcePatch, TargetPatch> + advancingFrontAMI ( - srcPatch, - tgtPatch, - triMode, + requireMatch, reverseTarget, - requireMatch + lowWeightCorrection, + triMode ), - restartUncoveredSourceFace_(restartUncoveredSourceFace), - srcMagSf_(), - tgtMagSf_(), - srcTris_(), - tgtTris_(), - extendedTgtMapPtr_(nullptr) + restartUncoveredSourceFace_(restartUncoveredSourceFace) +{} + + +Foam::faceAreaWeightAMI::faceAreaWeightAMI(const faceAreaWeightAMI& ami) +: + advancingFrontAMI(ami), + restartUncoveredSourceFace_(ami.restartUncoveredSourceFace_) {} // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate +bool Foam::faceAreaWeightAMI::calculate ( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei, - label tgtFacei + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr ) { + if (upToDate_) + { + return false; + } + addProfiling(ami, "faceAreaWeightAMI::calculate"); - // Create global indexing for each patch - globalIndex globalSrcFaces(this->srcPatch0_.size()); - globalIndex globalTgtFaces(this->tgtPatch0_.size()); + advancingFrontAMI::calculate(srcPatch, tgtPatch, surfPtr); - // Initialise the geometry - labelList extendedTgtFaceIDs; - initGeom(globalTgtFaces, extendedTgtFaceIDs); + label srcFacei = 0; + label tgtFacei = 0; - bool ok = - this->initialise - ( - srcAddress, - srcWeights, - tgtAddress, - tgtWeights, - srcFacei, - tgtFacei - ); + bool ok = initialiseWalk(srcFacei, tgtFacei); - srcCentroids.setSize(srcAddress.size()); + srcCentroids_.setSize(srcAddress_.size()); const auto& src = this->srcPatch(); - const auto& tgt = this->tgtPatch(); + const auto& tgt = this->tgtPatch(); // might be the extended patch! // Temporary storage for addressing and weights List<DynamicList<label>> srcAddr(src.size()); @@ -705,15 +667,15 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate tgtFacei ); - if (debug && !this->srcNonOverlap_.empty()) + if (debug && !srcNonOverlap_.empty()) { - Pout<< " AMI: " << this->srcNonOverlap_.size() + Pout<< " AMI: " << srcNonOverlap_.size() << " non-overlap faces identified" << endl; } // Check for badly covered faces - if (restartUncoveredSourceFace_) + if (restartUncoveredSourceFace_) // && requireMatch_??? { restartUncoveredSourceFace ( @@ -726,36 +688,38 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate } } - - // Set the patch face areas - srcMagSf = std::move(srcMagSf_); - tgtMagSf = std::move(tgtMagSf_); - // Transfer data to persistent storage forAll(srcAddr, i) { - srcAddress[i].transfer(srcAddr[i]); - srcWeights[i].transfer(srcWght[i]); - srcCentroids[i].transfer(srcCtr[i]); + srcAddress_[i].transfer(srcAddr[i]); + srcWeights_[i].transfer(srcWght[i]); + srcCentroids_[i].transfer(srcCtr[i]); } forAll(tgtAddr, i) { - tgtAddress[i].transfer(tgtAddr[i]); - tgtWeights[i].transfer(tgtWght[i]); + tgtAddress_[i].transfer(tgtAddr[i]); + tgtWeights_[i].transfer(tgtWght[i]); } - if (this->distributed()) + if (distributed()) { - for (labelList& addressing : srcAddress) + const primitivePatch& srcPatch0 = this->srcPatch0(); + const primitivePatch& tgtPatch0 = this->tgtPatch0(); + + // Create global indexing for each original patch + globalIndex globalSrcFaces(srcPatch0.size()); + globalIndex globalTgtFaces(tgtPatch0.size()); + + for (labelList& addressing : srcAddress_) { for (label& addr : addressing) { - addr = extendedTgtFaceIDs[addr]; + addr = extendedTgtFaceIDs_[addr]; } } - for (labelList& addressing : tgtAddress) + for (labelList& addressing : tgtAddress_) { globalSrcFaces.inplaceToGlobal(addressing); } @@ -767,12 +731,12 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate ( Pstream::commsTypes::nonBlocking, List<labelPair>(), - this->tgtPatch0_.size(), + tgtPatch0.size(), extendedTgtMapPtr_->constructMap(), false, // has flip extendedTgtMapPtr_->subMap(), false, // has flip - tgtAddress, + tgtAddress_, ListOps::appendEqOp<label>(), flipOp(), // flip operation labelList() @@ -782,43 +746,55 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate ( Pstream::commsTypes::nonBlocking, List<labelPair>(), - this->tgtPatch0_.size(), + tgtPatch0.size(), extendedTgtMapPtr_->constructMap(), false, extendedTgtMapPtr_->subMap(), false, - tgtWeights, + tgtWeights_, ListOps::appendEqOp<scalar>(), flipOp(), scalarList() ); // Note: using patch face areas calculated by the AMI method - extendedTgtMapPtr_->reverseDistribute - ( - this->tgtPatch0_.size(), - tgtMagSf - ); + extendedTgtMapPtr_->reverseDistribute(tgtPatch0.size(), tgtMagSf_); // Cache maps and reset addresses List<Map<label>> cMapSrc; - srcMapPtr.reset(new mapDistribute(globalSrcFaces, tgtAddress, cMapSrc)); + srcMapPtr_.reset + ( + new mapDistribute(globalSrcFaces, tgtAddress_, cMapSrc) + ); + List<Map<label>> cMapTgt; - tgtMapPtr.reset(new mapDistribute(globalTgtFaces, srcAddress, cMapTgt)); + tgtMapPtr_.reset + ( + new mapDistribute(globalTgtFaces, srcAddress_, cMapTgt) + ); } - return true; + // Convert the weights from areas to normalised values + normaliseWeights(conformal(), true); + + upToDate_ = true; + + return upToDate_; } -template<class SourcePatch, class TargetPatch> -void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::normaliseWeights -( - const bool verbose, - AMIInterpolation<SourcePatch, TargetPatch>& inter -) +void Foam::faceAreaWeightAMI::write(Ostream& os) const { - inter.normaliseWeights(this->conformal(), verbose); + advancingFrontAMI::write(os); + + if (restartUncoveredSourceFace_) + { + os.writeEntry + ( + "restartUncoveredSourceFace", + restartUncoveredSourceFace_ + ); + } } diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.H similarity index 74% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H rename to src/meshTools/AMIInterpolation/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.H index 04fe248b78f74d38139b2f37e873c102d0a22349..818d26ebea8aec4310209e3d46261a90384a3f4d 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.H @@ -38,7 +38,7 @@ SourceFiles #ifndef faceAreaWeightAMI_H #define faceAreaWeightAMI_H -#include "AMIMethod.H" +#include "advancingFrontAMI.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -49,10 +49,9 @@ namespace Foam Class faceAreaWeightAMI Declaration \*---------------------------------------------------------------------------*/ -template<class SourcePatch, class TargetPatch> class faceAreaWeightAMI : - public AMIMethod<SourcePatch, TargetPatch> + public advancingFrontAMI { private: @@ -61,35 +60,20 @@ private: //- Flag to restart uncovered source faces const bool restartUncoveredSourceFace_; - //- Source face areas - List<scalar> srcMagSf_; - - //- Target face areas - List<scalar> tgtMagSf_; - - //- Storage for src-side triangle decomposition - List<DynamicList<face>> srcTris_; - - //- Storage for tgt-side triangle decomposition - List<DynamicList<face>> tgtTris_; protected: - autoPtr<mapDistribute> extendedTgtMapPtr_; - - // Protected Member Functions - //- No copy construct - faceAreaWeightAMI(const faceAreaWeightAMI&) = delete; - //- No copy assignment void operator=(const faceAreaWeightAMI&) = delete; //- Initialise the geometry void initGeom ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, const globalIndex& globalTgtFaces, labelList& extendedTgtFaceIDs ); @@ -135,7 +119,7 @@ protected: ); //- Set the source and target seed faces - virtual void setNextFaces + virtual bool setNextFaces ( label& startSeedi, label& srcFacei, @@ -175,17 +159,33 @@ public: // Constructors + //- Construct from dictionary + faceAreaWeightAMI + ( + const dictionary& dict, + const bool reverseTarget = false + ); + //- Construct from components faceAreaWeightAMI ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const bool requireMatch, const bool reverseTarget = false, - const bool requireMatch = true, + const scalar lowWeightCorrection = -1, + const faceAreaIntersect::triangulationMode triMode = + faceAreaIntersect::tmMesh, const bool restartUncoveredSourceFace = true ); + //- Construct as copy + faceAreaWeightAMI(const faceAreaWeightAMI& ami); + + //- Construct and return a clone + virtual autoPtr<AMIInterpolation> clone() const + { + return autoPtr<AMIInterpolation>(new faceAreaWeightAMI(*this)); + } + //- Destructor virtual ~faceAreaWeightAMI() = default; @@ -193,31 +193,16 @@ public: // Member Functions - // Manipulation - - //- Update addressing, weights and centroids - virtual bool calculate - ( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei = -1, - label tgtFacei = -1 - ); + //- Update addressing, weights and (optional) centroids + virtual bool calculate + ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr = nullptr + ); - //- Normalise the weight. Can optionally subset addressing - //- (e.g. for mapNearest) - virtual void normaliseWeights - ( - const bool verbose, - AMIInterpolation<SourcePatch, TargetPatch>& inter - ); + //- Write + virtual void write(Ostream& os) const; }; @@ -227,12 +212,6 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -#ifdef NoRepository - #include "faceAreaWeightAMI.C" -#endif - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - #endif // ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.C new file mode 100644 index 0000000000000000000000000000000000000000..9108ddbca89f3131065442b65e65c527ad2ffc65 --- /dev/null +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.C @@ -0,0 +1,408 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "nearestFaceAMI.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(nearestFaceAMI, 0); + addToRunTimeSelectionTable(AMIInterpolation, nearestFaceAMI, dict); + addToRunTimeSelectionTable(AMIInterpolation, nearestFaceAMI, component); +} + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +Foam::autoPtr<Foam::mapDistribute> Foam::nearestFaceAMI::calcFaceMap +( + const List<nearestAndDist>& localInfo, + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch +) const +{ + // Generate the list of processor bounding boxes for tgtPatch + List<boundBox> procBbs(Pstream::nProcs()); + procBbs[Pstream::myProcNo()] = + boundBox(tgtPatch.points(), tgtPatch.meshPoints(), true); + Pstream::gatherList(procBbs); + Pstream::scatterList(procBbs); + + // Identify which of my local src faces intersect with each processor + // tgtPatch bb within the current match's search distance + const pointField& srcCcs = srcPatch.faceCentres(); + List<DynamicList<label>> dynSendMap(Pstream::nProcs()); + + forAll(localInfo, srcFacei) + { + // Test local srcPatch face centres against remote processor tgtPatch bb + // using distance from initial pass + + const scalar r2 = localInfo[srcFacei].second(); + + forAll(procBbs, proci) + { + if (proci != Pstream::myProcNo()) + { + if (procBbs[proci].overlaps(srcCcs[srcFacei], r2)) + { + dynSendMap[proci].append(srcFacei); + } + } + } + } + + // Convert dynamicList to labelList + labelListList sendMap(Pstream::nProcs()); + forAll(sendMap, proci) + { + dynSendMap[proci].shrink(); + sendMap[proci].transfer(dynSendMap[proci]); + + if (debug) + { + Pout<< "send map - to proc " << proci << " sending " + << sendMap[proci].size() << " elements" << endl; + } + } + + return autoPtr<mapDistribute>::New(std::move(sendMap)); +} + + +Foam::autoPtr<Foam::mapDistribute> Foam::nearestFaceAMI::calcDistributed +( + const primitivePatch& src, + const primitivePatch& tgt, + labelListList& srcToTgtAddr, + scalarListList& srcToTgtWght +) const +{ + autoPtr<indexedOctree<treeType>> tgtTreePtr; + if (tgt.size()) + { + tgtTreePtr = this->createTree(tgt); + } + + // Create global indexing for tgtPatch + globalIndex globalTgtCells(tgt.size()); + + + // First pass + // ========== + // For each srcPatch face, determine local match on tgtPatch + + // Identify local nearest matches + pointField srcCcs(src.faceCentres()); + + List<nearestAndDist> localInfo(src.size()); + if (tgt.size()) + { + const auto& tgtTree = tgtTreePtr(); + + forAll(srcCcs, srcCelli) + { + const point& srcCc = srcCcs[srcCelli]; + + pointIndexHit& test = localInfo[srcCelli].first(); + test = tgtTree.findNearest(srcCc, GREAT); + + if (test.hit()) + { + // With a search radius2 of GREAT all cells should receive a hit + localInfo[srcCelli].second() = magSqr(srcCc - test.hitPoint()); + test.setIndex(globalTgtCells.toGlobal(test.index())); + } + } + } + else + { + // No local tgtPatch faces - initialise nearest distance for all + // srcPatch faces to GREAT so that they [should be] set by remote + // tgtPatch faces + for (auto& info : localInfo) + { + info.second() = GREAT; + } + } + + + // Second pass + // =========== + // Determine remote matches + + // Map returns labels of src patch faces to send to each proc + autoPtr<mapDistribute> mapPtr = calcFaceMap(localInfo, src, tgt); + mapDistribute& map = mapPtr(); + + List<nearestAndDist> remoteInfo(localInfo); + map.distribute(remoteInfo); + + // Note: re-using srcCcs + map.distribute(srcCcs); + + if (tgt.size()) + { + const auto& tgtTree = tgtTreePtr(); + + // Test remote srcPatch faces against local tgtPatch faces + nearestAndDist testInfo; + pointIndexHit& test = testInfo.first(); + forAll(remoteInfo, i) + { + test = tgtTree.findNearest(srcCcs[i], remoteInfo[i].second()); + if (test.hit()) + { + test.setIndex(globalTgtCells.toGlobal(test.index())); + testInfo.second() = magSqr(test.hitPoint() - srcCcs[i]); + nearestEqOp()(remoteInfo[i], testInfo); + } + } + } + + // Send back to originating processor. Choose best if sent to multiple + // processors. Note that afterwards all unused entries have the unique + // value nearestZero (distance < 0). This is used later on to see if + // the sample was sent to any processor. + const nearestAndDist nearestZero(pointIndexHit(), -GREAT); + mapDistributeBase::distribute + ( + Pstream::commsTypes::nonBlocking, + List<labelPair>(0), + src.size(), + map.constructMap(), + map.constructHasFlip(), + map.subMap(), + map.subHasFlip(), + remoteInfo, + nearestEqOp(), + noOp(), // no flipping + nearestZero + ); + + + // Third pass + // ========== + // Combine local and remote info and filter out any connections that are + // further away than threshold distance squared + + srcToTgtAddr.setSize(src.size()); + srcToTgtWght.setSize(src.size()); + forAll(srcToTgtAddr, srcFacei) + { + nearestEqOp()(localInfo[srcFacei], remoteInfo[srcFacei]); + if (localInfo[srcFacei].second() < maxDistance2_) + { + const label tgtFacei = localInfo[srcFacei].first().index(); + srcToTgtAddr[srcFacei] = labelList(1, tgtFacei); + srcToTgtWght[srcFacei] = scalarList(1, 1.0); + } + } + + List<Map<label>> cMap; + return autoPtr<mapDistribute>::New(globalTgtCells, srcToTgtAddr, cMap); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::nearestFaceAMI::nearestFaceAMI +( + const dictionary& dict, + const bool reverseTarget +) +: + AMIInterpolation(dict, reverseTarget), + maxDistance2_(dict.getOrDefault<scalar>("maxDistance2", GREAT)) +{} + + +Foam::nearestFaceAMI::nearestFaceAMI +( + const bool requireMatch, + const bool reverseTarget, + const scalar lowWeightCorrection +) +: + AMIInterpolation(requireMatch, reverseTarget, lowWeightCorrection), + maxDistance2_(GREAT) +{} + + +Foam::nearestFaceAMI::nearestFaceAMI(const nearestFaceAMI& ami) +: + AMIInterpolation(ami), + maxDistance2_(ami.maxDistance2_) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +bool Foam::nearestFaceAMI::calculate +( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr +) +{ + if (upToDate_) + { + return false; + } + + AMIInterpolation::calculate(srcPatch, tgtPatch, surfPtr); + + const auto& src = this->srcPatch0(); + const auto& tgt = this->tgtPatch0(); + + // Set face area magnitudes + srcMagSf_ = mag(src.faceAreas()); + tgtMagSf_ = mag(tgt.faceAreas()); + + // TODO: implement symmetric calculation controls; assume yes for now + bool symmetric_ = true; + + if (this->distributed()) + { + tgtMapPtr_ = + calcDistributed + ( + src, + tgt, + srcAddress_, + srcWeights_ + ); + + if (symmetric_) + { + srcMapPtr_ = + calcDistributed + ( + tgt, + src, + tgtAddress_, + tgtWeights_ + ); + } + } + else + { + srcAddress_.setSize(src.size()); + srcWeights_.setSize(src.size()); + + if (symmetric_) + { + tgtAddress_.setSize(tgt.size()); + tgtWeights_.setSize(tgt.size()); + } + + const pointField& srcCcs = src.faceCentres(); + const pointField& tgtCcs = tgt.faceCentres(); + + const auto tgtTreePtr = this->createTree(tgtPatch); + const auto& tgtTree = tgtTreePtr(); + + forAll(srcCcs, srcFacei) + { + const point& srcCc = srcCcs[srcFacei]; + const pointIndexHit hit = tgtTree.findNearest(srcCc, GREAT); + + if + ( + hit.hit() + && (magSqr(srcCc - tgtCcs[hit.index()]) < maxDistance2_) + ) + { + label tgtFacei = hit.index(); + srcAddress_[srcFacei] = labelList(1, tgtFacei); + srcWeights_[srcFacei] = scalarList(1, 1.0); + + if (symmetric_) + { + tgtAddress_[tgtFacei] = labelList(1, srcFacei); + tgtWeights_[tgtFacei] = scalarList(1, 1.0); + } + } + else + { + if (debug) + { + WarningInFunction + << "Unable to find target face for source face " + << srcFacei << endl; + } + } + } + + if (symmetric_) + { + const auto srcTreePtr = this->createTree(srcPatch); + const auto& srcTree = srcTreePtr(); + + // Check that all source cells have connections and populate any + // missing entries + forAll(tgtWeights_, tgtCelli) + { + if (tgtAddress_[tgtCelli].empty()) + { + const point& tgtCc = tgtCcs[tgtCelli]; + pointIndexHit hit = srcTree.findNearest(tgtCc, GREAT); + + if + ( + hit.hit() + && (magSqr(tgtCc - srcCcs[hit.index()]) < maxDistance2_) + ) + { + tgtAddress_[tgtCelli] = labelList(1, hit.index()); + tgtWeights_[tgtCelli] = scalarList(1, 1.0); + } + } + else + { + if (debug) + { + WarningInFunction + << "Unable to find source face for target face " + << tgtCelli << endl; + } + } + } + } + } + + srcWeightsSum_.setSize(srcWeights_.size(), 1); + tgtWeightsSum_.setSize(tgtWeights_.size(), 1); + + upToDate_ = true; + + return upToDate_; +} + + +// ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/nearestFaceAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.H similarity index 68% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/nearestFaceAMI.H rename to src/meshTools/AMIInterpolation/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.H index 78bc71d63693c9c25444c7e9e6ee934ac4d9630d..85938ecf4d817422cf1aee4646756949ccefd527 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/nearestFaceAMI.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.H @@ -37,7 +37,7 @@ SourceFiles #ifndef nearestFaceAMI_H #define nearestFaceAMI_H -#include "AMIMethod.H" +#include "AMIInterpolation.H" #include "pointIndexHit.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -49,10 +49,9 @@ namespace Foam Class nearestFaceAMI Declaration \*---------------------------------------------------------------------------*/ -template<class SourcePatch, class TargetPatch> class nearestFaceAMI : - public AMIMethod<SourcePatch, TargetPatch> + public AMIInterpolation { public: @@ -91,23 +90,20 @@ private: // Private Member Functions - //- No copy construct - nearestFaceAMI(const nearestFaceAMI&) = delete; - //- No copy assignment void operator=(const nearestFaceAMI&) = delete; autoPtr<mapDistribute> calcFaceMap ( const List<nearestAndDist>& localInfo, - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch ) const; autoPtr<mapDistribute> calcDistributed ( - const SourcePatch& src, - const TargetPatch& tgt, + const primitivePatch& src, + const primitivePatch& tgt, labelListList& srcToTgtAddr, scalarListList& srcToTgtWght ) const; @@ -121,16 +117,30 @@ public: // Constructors + //- Construct from dictionary + nearestFaceAMI + ( + const dictionary& dict, + const bool reverseTarget = false + ); + //- Construct from components nearestFaceAMI ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const bool requireMatch = true, const bool reverseTarget = false, - const bool requireMatch = true + const scalar lowWeightCorrection = -1 ); + //- Construct as copy + nearestFaceAMI(const nearestFaceAMI& ami); + + //- Construct and return a clone + virtual autoPtr<AMIInterpolation> clone() const + { + return autoPtr<AMIInterpolation>(new nearestFaceAMI(*this)); + } + //- Destructor virtual ~nearestFaceAMI() = default; @@ -138,31 +148,13 @@ public: // Member Functions - // Manipulation - - //- Update addressing and weights - virtual bool calculate - ( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei = -1, - label tgtFacei = -1 - ); - - //- Normalise the weight. Can optionally subset addressing - //- (e.g. for mapNearest) - virtual void normaliseWeights - ( - const bool verbose, - AMIInterpolation<SourcePatch, TargetPatch>& inter - ); + //- Update addressing and weights + virtual bool calculate + ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr = nullptr + ); }; @@ -172,12 +164,6 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -#ifdef NoRepository - #include "nearestFaceAMI.C" -#endif - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - #endif // ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C similarity index 55% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C rename to src/meshTools/AMIInterpolation/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C index 4e15ead5600e8a655e80c5a427f9babd0e4452b9..61b4197f2b0b5354f02bfa7c3dada538d4a46555 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation + Copyright (C) 2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -26,11 +27,30 @@ License \*---------------------------------------------------------------------------*/ #include "partialFaceAreaWeightAMI.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(partialFaceAreaWeightAMI, 0); + addToRunTimeSelectionTable + ( + AMIInterpolation, + partialFaceAreaWeightAMI, + dict + ); + addToRunTimeSelectionTable + ( + AMIInterpolation, + partialFaceAreaWeightAMI, + component + ); +} // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces +bool Foam::partialFaceAreaWeightAMI::setNextFaces ( label& startSeedi, label& srcFacei, @@ -41,7 +61,7 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces const bool errorOnNotFound ) const { - faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces + return faceAreaWeightAMI::setNextFaces ( startSeedi, srcFacei, @@ -56,85 +76,75 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>:: -partialFaceAreaWeightAMI +Foam::partialFaceAreaWeightAMI::partialFaceAreaWeightAMI +( + const dictionary& dict, + const bool reverseTarget +) +: + faceAreaWeightAMI(dict, reverseTarget) +{} + + +Foam::partialFaceAreaWeightAMI::partialFaceAreaWeightAMI ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const bool requireMatch, const bool reverseTarget, - const bool requireMatch + const scalar lowWeightCorrection, + const faceAreaIntersect::triangulationMode triMode, + const bool restartUncoveredSourceFace ) : - faceAreaWeightAMI<SourcePatch, TargetPatch> + faceAreaWeightAMI ( - srcPatch, - tgtPatch, - triMode, - reverseTarget, requireMatch, - true // false // Not performing restart on low weights - valid for partial match + reverseTarget, + lowWeightCorrection, + triMode, + restartUncoveredSourceFace ) {} +Foam::partialFaceAreaWeightAMI::partialFaceAreaWeightAMI +( + const partialFaceAreaWeightAMI& ami +) +: + faceAreaWeightAMI(ami) +{} + + // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -template<class SourcePatch, class TargetPatch> -bool Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::conformal() const +bool Foam::partialFaceAreaWeightAMI::conformal() const { return false; } -template<class SourcePatch, class TargetPatch> -bool Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::calculate +bool Foam::partialFaceAreaWeightAMI::calculate ( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei, - label tgtFacei + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr ) { - bool ok = - faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate - ( - srcAddress, - srcWeights, - srcCentroids, - tgtAddress, - tgtWeights, - srcMagSf, - tgtMagSf, - srcMapPtr, - tgtMapPtr, - srcFacei, - tgtFacei - ); - - if (ok) + if (faceAreaWeightAMI::calculate(srcPatch, tgtPatch, surfPtr)) { - if (this->distributed()) + if (distributed()) { - scalarList newTgtMagSf(std::move(tgtMagSf)); + scalarList newTgtMagSf(std::move(tgtMagSf_)); - // Assign default sizes. Override selected values with - // calculated values. This is to support ACMI - // where some of the target faces are never used (so never get sent - // over and hence never assigned to) - tgtMagSf = this->tgtPatch0_.magFaceAreas(); + // Assign default sizes. Override selected values with calculated + // values. This is to support ACMI where some of the target faces + // are never used (so never get sent over and hence never assigned + // to) + tgtMagSf_ = tgtPatch0().magFaceAreas(); for (const labelList& smap : this->extendedTgtMapPtr_->subMap()) { - UIndirectList<scalar>(tgtMagSf, smap) = + UIndirectList<scalar>(tgtMagSf_, smap) = UIndirectList<scalar>(newTgtMagSf, smap); } } diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H similarity index 70% rename from src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H rename to src/meshTools/AMIInterpolation/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H index ef46251bd6f00007c3d526158c20f875e4f895ff..f03f8d155dd4281eb33555170e4b0b1328f48206 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2016 OpenCFD Ltd. + Copyright (C) 2016-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -49,19 +49,15 @@ namespace Foam Class partialFaceAreaWeightAMI Declaration \*---------------------------------------------------------------------------*/ -template<class SourcePatch, class TargetPatch> class partialFaceAreaWeightAMI : - public faceAreaWeightAMI<SourcePatch, TargetPatch> + public faceAreaWeightAMI { private: // Private Member Functions - //- No copy construct - partialFaceAreaWeightAMI(const partialFaceAreaWeightAMI&) = delete; - //- No copy assignment void operator=(const partialFaceAreaWeightAMI&) = delete; @@ -71,7 +67,7 @@ protected: // Protected Member Functions //- Set the source and target seed faces - virtual void setNextFaces + virtual bool setNextFaces ( label& startSeedi, label& srcFacei, @@ -91,16 +87,34 @@ public: // Constructors + //- Construct from dictionary + partialFaceAreaWeightAMI + ( + const dictionary& dict, + const bool reverseTarget = false + ); + //- Construct from components partialFaceAreaWeightAMI ( - const SourcePatch& srcPatch, - const TargetPatch& tgtPatch, - const faceAreaIntersect::triangulationMode& triMode, + const bool requireMatch = false, const bool reverseTarget = false, - const bool requireMatch = true + const scalar lowWeightCorrection = -1, + const faceAreaIntersect::triangulationMode triMode = + faceAreaIntersect::tmMesh, + const bool restartUncoveredSourceFace = true ); + //- Construct as copy + partialFaceAreaWeightAMI(const partialFaceAreaWeightAMI& ami); + + //- Construct and return a clone + virtual autoPtr<AMIInterpolation> clone() const + { + return + autoPtr<AMIInterpolation>(new partialFaceAreaWeightAMI(*this)); + } + //- Destructor virtual ~partialFaceAreaWeightAMI() = default; @@ -108,29 +122,16 @@ public: // Member Functions - // Access - - //- Flag to indicate that interpolation patches are conformal - virtual bool conformal() const; - - - // Manipulation - - //- Update addressing and weights - virtual bool calculate - ( - labelListList& srcAddress, - scalarListList& srcWeights, - pointListList& srcCentroids, - labelListList& tgtAddress, - scalarListList& tgtWeights, - scalarList& srcMagSf, - scalarList& tgtMagSf, - autoPtr<mapDistribute>& srcMapPtr, - autoPtr<mapDistribute>& tgtMapPtr, - label srcFacei = -1, - label tgtFacei = -1 - ); + //- Flag to indicate that interpolation patches are conformal + virtual bool conformal() const; + + //- Update addressing, weights and (optional) centroids + virtual bool calculate + ( + const primitivePatch& srcPatch, + const primitivePatch& tgtPatch, + const autoPtr<searchableSurface>& surfPtr = nullptr + ); }; @@ -140,12 +141,6 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -#ifdef NoRepository - #include "partialFaceAreaWeightAMI.C" -#endif - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - #endif // ************************************************************************* // diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H index 916d40e36f7cf3d02504b8ac97eedc30f5260b12..d79327de31ddc44d3b331a2e0dcbd2212af0f8e8 100644 --- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H +++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H @@ -46,6 +46,7 @@ SourceFiles #include "face.H" #include "triPoints.H" #include "Enum.H" +#include "searchableSurface.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C index cbba50fb341b641aeee43b308a25d338107b6c38..4e994a6f236185043af6457fa7b607da60537c64 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2017-2019 OpenCFD Ltd. + Copyright (C) 2017-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -72,20 +72,13 @@ void Foam::cyclicACMIPolyPatch::reportCoverage } -void Foam::cyclicACMIPolyPatch::resetAMI -( - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod -) const +void Foam::cyclicACMIPolyPatch::resetAMI() const { - resetAMI(boundaryMesh().mesh().points(), AMIMethod); + resetAMI(boundaryMesh().mesh().points()); } -void Foam::cyclicACMIPolyPatch::resetAMI -( - const UList<point>& points, - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod -) const +void Foam::cyclicACMIPolyPatch::resetAMI(const UList<point>& points) const { if (!owner()) { @@ -94,52 +87,40 @@ void Foam::cyclicACMIPolyPatch::resetAMI const polyPatch& nonOverlapPatch = this->nonOverlapPatch(); - if (debug) - { - Pout<< "cyclicACMIPolyPatch::resetAMI : recalculating weights" - << " for " << name() << " and " << nonOverlapPatch.name() - << endl; - } + DebugPout + << "cyclicACMIPolyPatch::resetAMI : recalculating weights" + << " for " << name() << " and " << nonOverlapPatch.name() + << endl; + const polyMesh& mesh = boundaryMesh().mesh(); -WarningInFunction<< "DEACTIVATED clearGeom()" << endl; -// if (boundaryMesh().mesh().hasCellCentres()) - if (0 && boundaryMesh().mesh().hasCellCentres()) + if (!createAMIFaces_ && mesh.hasCellCentres()) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::resetAMI : clearing cellCentres" - << " for " << name() << " and " << nonOverlapPatch.name() - << endl; - } + DebugPout + << "cyclicACMIPolyPatch::resetAMI : clearing cellCentres" + << " for " << name() << " and " << nonOverlapPatch.name() + << endl; - //WarningInFunction - // << "The mesh already has cellCentres calculated when" - // << " resetting ACMI " << name() << "." << endl - // << "This is a problem since ACMI adapts the face areas" - // << " (to close cells) so this has" << endl - // << "to be done before cell centre calculation." << endl - // << "This can happen if e.g. the cyclicACMI is after" - // << " any processor patches in the boundary." << endl; - const_cast<polyMesh&> - ( - boundaryMesh().mesh() - ).primitiveMesh::clearGeom(); + WarningInFunction + << "The mesh already has cellCentres calculated when" + << " resetting ACMI " << name() << "." << nl + << "This is a problem since ACMI adapts the face areas" + << " (to close cells) so this has" << nl + << "to be done before cell centre calculation." << nl + << "This can happen if e.g. the cyclicACMI is after" + << " any processor patches in the boundary." << endl; + const_cast<polyMesh&>(mesh).primitiveMesh::clearGeom(); } // Trigger re-building of faceAreas - (void)boundaryMesh().mesh().faceAreas(); + (void)mesh.faceAreas(); // Calculate the AMI using partial face-area-weighted. This leaves // the weights as fractions of local areas (sum(weights) = 1 means // face is fully covered) - cyclicAMIPolyPatch::resetAMI - ( - points, - AMIPatchToPatchInterpolation::imPartialFaceAreaWeight - ); + cyclicAMIPolyPatch::resetAMI(points); const AMIPatchToPatchInterpolation& AMI = this->AMI(); @@ -196,15 +177,17 @@ void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas // Primitive patch face areas have been cleared/reset based on the raw // points - need to reset to avoid double-accounting of face areas - DebugPout - << "rescaling non-overlap patch areas" << endl; - const scalar maxTol = scalar(1) - tolerance_; const scalarField& mask = acmipp.mask(); const polyPatch& nonOverlapPatch = acmipp.nonOverlapPatch(); vectorField::subField noSf = nonOverlapPatch.faceAreas(); + DebugPout + << "rescaling non-overlap patch areas for: " << nonOverlapPatch.name() + << endl; + + if (mask.size() != noSf.size()) { WarningInFunction @@ -232,7 +215,7 @@ void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas // updates, e.g. by the calls to cyclicAMIPolyPatch's setTopology and // initMovePoints DebugPout - << "scaling coupled patch areas" << endl; + << "scaling coupled patch areas for: " << acmipp.name() << endl; // Scale the coupled patch face areas vectorField::subField Sf = acmipp.faceAreas(); @@ -266,26 +249,26 @@ void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas void Foam::cyclicACMIPolyPatch::initGeometry(PstreamBuffers& pBufs) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::initGeometry : " << name() << endl; - } + DebugPout << "cyclicACMIPolyPatch::initGeometry : " << name() << endl; // Note: calculates transformation and triggers face centre calculation cyclicAMIPolyPatch::initGeometry(pBufs); - // On start-up there are no topological updates so scale the face areas - // - Note: resetAMI called by cyclicAMIPolyPatch::initGeometry + // Initialise the AMI early to make sure we adapt the face areas before the + // cell centre calculation gets triggered. + if (!createAMIFaces_ && canResetAMI()) + { + resetAMI(); + } + scalePatchFaceAreas(); } void Foam::cyclicACMIPolyPatch::calcGeometry(PstreamBuffers& pBufs) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::calcGeometry : " << name() << endl; - } + DebugPout << "cyclicACMIPolyPatch::calcGeometry : " << name() << endl; + cyclicAMIPolyPatch::calcGeometry(pBufs); } @@ -296,10 +279,7 @@ void Foam::cyclicACMIPolyPatch::initMovePoints const pointField& p ) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::initMovePoints : " << name() << endl; - } + DebugPout<< "cyclicACMIPolyPatch::initMovePoints : " << name() << endl; // Note: calculates transformation and triggers face centre calculation // - Note: resetAMI called by cyclicAMIPolyPatch::initMovePoints @@ -315,10 +295,7 @@ void Foam::cyclicACMIPolyPatch::movePoints const pointField& p ) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::movePoints : " << name() << endl; - } + DebugPout << "cyclicACMIPolyPatch::movePoints : " << name() << endl; // When topology is changing, this will scale the duplicate AMI faces cyclicAMIPolyPatch::movePoints(pBufs, p); @@ -327,30 +304,24 @@ void Foam::cyclicACMIPolyPatch::movePoints void Foam::cyclicACMIPolyPatch::initUpdateMesh(PstreamBuffers& pBufs) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::initUpdateMesh : " << name() << endl; - } + DebugPout << "cyclicACMIPolyPatch::initUpdateMesh : " << name() << endl; + cyclicAMIPolyPatch::initUpdateMesh(pBufs); } void Foam::cyclicACMIPolyPatch::updateMesh(PstreamBuffers& pBufs) { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::updateMesh : " << name() << endl; - } + DebugPout << "cyclicACMIPolyPatch::updateMesh : " << name() << endl; + cyclicAMIPolyPatch::updateMesh(pBufs); } void Foam::cyclicACMIPolyPatch::clearGeom() { - if (debug) - { - Pout<< "cyclicACMIPolyPatch::clearGeom : " << name() << endl; - } + DebugPout << "cyclicACMIPolyPatch::clearGeom : " << name() << endl; + cyclicAMIPolyPatch::clearGeom(); } @@ -377,16 +348,27 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch const label index, const polyBoundaryMesh& bm, const word& patchType, - const transformType transform + const transformType transform, + const word& defaultAMIMethod ) : - cyclicAMIPolyPatch(name, size, start, index, bm, patchType, transform), + cyclicAMIPolyPatch + ( + name, + size, + start, + index, + bm, + patchType, + transform, + defaultAMIMethod + ), nonOverlapPatchName_(word::null), nonOverlapPatchID_(-1), srcMask_(), tgtMask_() { - AMIRequireMatch_ = false; + AMIPtr_->setRequireMatch(false); // Non-overlapping patch might not be valid yet so cannot determine // associated patchID @@ -399,16 +381,17 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch const dictionary& dict, const label index, const polyBoundaryMesh& bm, - const word& patchType + const word& patchType, + const word& defaultAMIMethod ) : - cyclicAMIPolyPatch(name, dict, index, bm, patchType), - nonOverlapPatchName_(dict.lookup("nonOverlapPatch")), + cyclicAMIPolyPatch(name, dict, index, bm, patchType, defaultAMIMethod), + nonOverlapPatchName_(dict.get<word>("nonOverlapPatch")), nonOverlapPatchID_(-1), srcMask_(), tgtMask_() { - AMIRequireMatch_ = false; + AMIPtr_->setRequireMatch(false); if (nonOverlapPatchName_ == name) { @@ -435,7 +418,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch srcMask_(), tgtMask_() { - AMIRequireMatch_ = false; + AMIPtr_->setRequireMatch(false); // Non-overlapping patch might not be valid yet so cannot determine // associated patchID @@ -459,7 +442,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch srcMask_(), tgtMask_() { - AMIRequireMatch_ = false; + AMIPtr_->setRequireMatch(false); if (nonOverlapPatchName_ == name()) { @@ -489,7 +472,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch srcMask_(), tgtMask_() { - AMIRequireMatch_ = false; + AMIPtr_->setRequireMatch(false); } diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H index c858c14e0485d103ab3930621c137ff830a9d987..51b0f40b5c3d2e92fcd1eb1a7476ceba05ee4e0b 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2018-2019 OpenCFD Ltd. + Copyright (C) 2018-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -41,6 +41,7 @@ SourceFiles #include "cyclicAMIPolyPatch.H" #include "AMIPatchToPatchInterpolation.H" #include "polyBoundaryMesh.H" +#include "partialFaceAreaWeightAMI.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -86,19 +87,10 @@ protected: ) const; //- Reset the AMI interpolator, supply patch points - virtual void resetAMI - ( - const UList<point>& points, - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod = - AMIPatchToPatchInterpolation::imFaceAreaWeight - ) const; + virtual void resetAMI(const UList<point>& points) const; //- Reset the AMI interpolator, use current patch points - virtual void resetAMI - ( - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod = - AMIPatchToPatchInterpolation::imFaceAreaWeight - ) const; + virtual void resetAMI() const; //- Scale patch face areas to maintain physical area virtual void scalePatchFaceAreas(); @@ -145,7 +137,8 @@ public: const label index, const polyBoundaryMesh& bm, const word& patchType, - const transformType transform = UNKNOWN + const transformType transform = UNKNOWN, + const word& defaultAMIMethod = partialFaceAreaWeightAMI::typeName ); //- Construct from dictionary @@ -155,7 +148,8 @@ public: const dictionary& dict, const label index, const polyBoundaryMesh& bm, - const word& patchType + const word& patchType, + const word& defaultAMIMethod = partialFaceAreaWeightAMI::typeName ); //- Construct as copy, resetting the boundary mesh diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C index 04cb5b6bb3f63d06358d959f03f0e17184122a68..85f6170f7185e854219d3992f863303992ac0c42 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C @@ -29,7 +29,9 @@ License #include "cyclicAMIPolyPatch.H" #include "SubField.H" #include "Time.H" -#include "faceAreaIntersect.H" +#include "unitConversion.H" +#include "OFstream.H" +#include "meshTools.H" #include "addToRunTimeSelectionTable.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -59,14 +61,12 @@ Foam::vector Foam::cyclicAMIPolyPatch::findFaceNormalMaxRadius label facei = findMax(magRadSqr); - if (debug) - { - Info<< "findFaceMaxRadius(const pointField&) : patch: " << name() << nl - << " rotFace = " << facei << nl - << " point = " << faceCentres[facei] << nl - << " distance = " << Foam::sqrt(magRadSqr[facei]) - << endl; - } + DebugInFunction + << "Patch: " << name() << nl + << " rotFace = " << facei << nl + << " point = " << faceCentres[facei] << nl + << " distance = " << Foam::sqrt(magRadSqr[facei]) + << endl; return n[facei]; } @@ -288,20 +288,45 @@ void Foam::cyclicAMIPolyPatch::calcTransforms // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * // -void Foam::cyclicAMIPolyPatch::resetAMI -( - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod -) const +const Foam::autoPtr<Foam::searchableSurface>& +Foam::cyclicAMIPolyPatch::surfPtr() const { - resetAMI(boundaryMesh().mesh().points(), AMIMethod); + const word surfType(surfDict_.getOrDefault<word>("type", "none")); + + if (!surfPtr_ && owner() && surfType != "none") + { + word surfName(surfDict_.getOrDefault("name", name())); + + const polyMesh& mesh = boundaryMesh().mesh(); + + surfPtr_ = + searchableSurface::New + ( + surfType, + IOobject + ( + surfName, + mesh.time().constant(), + "triSurface", + mesh, + IOobject::MUST_READ, + IOobject::NO_WRITE + ), + surfDict_ + ); + } + + return surfPtr_; } -void Foam::cyclicAMIPolyPatch::resetAMI -( - const UList<point>& points, - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod -) const +void Foam::cyclicAMIPolyPatch::resetAMI() const +{ + resetAMI(boundaryMesh().mesh().points()); +} + + +void Foam::cyclicAMIPolyPatch::resetAMI(const UList<point>& points) const { DebugInFunction << endl; @@ -310,11 +335,9 @@ void Foam::cyclicAMIPolyPatch::resetAMI return; } - AMIPtr_.clear(); - const cyclicAMIPolyPatch& nbr = neighbPatch(); - pointField srcPoints(points, meshPoints()); - pointField nbrPoints(points, neighbPatch().meshPoints()); + pointField srcPoints(localPoints()); + pointField nbrPoints(nbr.localPoints()); if (debug) { @@ -343,20 +366,12 @@ void Foam::cyclicAMIPolyPatch::resetAMI transformPosition(nbrPoints); primitivePatch nbrPatch0 ( - SubList<face> - ( - nbr.localFaces(), - nbrPatchSize0 - ), + SubList<face>(nbr.localFaces(), nbrPatchSize0), nbrPoints ); primitivePatch patch0 ( - SubList<face> - ( - localFaces(), - patchSize0 - ), + SubList<face>(localFaces(), patchSize0), srcPoints ); @@ -372,23 +387,8 @@ void Foam::cyclicAMIPolyPatch::resetAMI } // Construct/apply AMI interpolation to determine addressing and weights - AMIPtr_.reset - ( - new AMIPatchToPatchInterpolation - ( - patch0, // *this, - nbrPatch0, - surfPtr(), - faceAreaIntersect::tmMesh, - AMIRequireMatch_, - AMIMethod, - AMILowWeightCorrection_, - AMIReverse_ - ) - ); - - // Set the updated flag - updated_ = true; + AMIPtr_->upToDate() = false; + AMIPtr_->calculate(patch0, nbrPatch0, surfPtr()); if (debug) { @@ -424,14 +424,12 @@ void Foam::cyclicAMIPolyPatch::calcTransforms() half1Areas ); - if (debug) - { - Pout<< "calcTransforms() : patch: " << name() << nl - << " forwardT = " << forwardT() << nl - << " reverseT = " << reverseT() << nl - << " separation = " << separation() << nl - << " collocated = " << collocated() << nl << endl; - } + DebugPout + << "calcTransforms() : patch: " << name() << nl + << " forwardT = " << forwardT() << nl + << " reverseT = " << reverseT() << nl + << " separation = " << separation() << nl + << " collocated = " << collocated() << nl << endl; } @@ -439,13 +437,8 @@ void Foam::cyclicAMIPolyPatch::initGeometry(PstreamBuffers& pBufs) { DebugInFunction << endl; - // Only resetting the AMI if not creating additional AMI faces - // - When creating additional faces the AMI is updated externally by the - // dynamic mesh via the setTopology function. - if (!createAMIFaces_ && canResetAMI()) - { - resetAMI(AMIMethod_); - } + // Flag AMI as needing update + AMIPtr_->upToDate() = false; polyPatch::initGeometry(pBufs); @@ -491,7 +484,7 @@ void Foam::cyclicAMIPolyPatch::initMovePoints } else { - resetAMI(p, AMIMethod_); + AMIPtr_->upToDate() = false; } // Early calculation of transforms. See above. @@ -551,9 +544,9 @@ void Foam::cyclicAMIPolyPatch::clearGeom() if (!updatingAMI_) { -//DEBUG("*** CLEARING AMI ***"); - AMIPtr_.clear(); + AMIPtr_->upToDate() = false; } + polyPatch::clearGeom(); } @@ -568,11 +561,11 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch const label index, const polyBoundaryMesh& bm, const word& patchType, - const transformType transform + const transformType transform, + const word& defaultAMIMethod ) : coupledPolyPatch(name, size, start, index, bm, patchType, transform), - updated_(false), nbrPatchName_(word::null), nbrPatchID_(-1), rotationAxis_(Zero), @@ -580,14 +573,11 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch rotationAngleDefined_(false), rotationAngle_(0.0), separationVector_(Zero), - AMIPtr_(nullptr), - AMIMethod_(AMIPatchToPatchInterpolation::imFaceAreaWeight), - AMIReverse_(false), - AMIRequireMatch_(true), - AMILowWeightCorrection_(-1.0), - surfPtr_(nullptr), + AMIPtr_(AMIInterpolation::New(defaultAMIMethod)), surfDict_(fileName("surface")), + surfPtr_(nullptr), createAMIFaces_(false), + moveFaceCentres_(false), updatingAMI_(true), srcFaceIDs_(), tgtFaceIDs_(), @@ -605,12 +595,12 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch const dictionary& dict, const label index, const polyBoundaryMesh& bm, - const word& patchType + const word& patchType, + const word& defaultAMIMethod ) : coupledPolyPatch(name, dict, index, bm, patchType), - updated_(false), - nbrPatchName_(dict.getOrDefault<word>("neighbourPatch", "")), + nbrPatchName_(dict.getOrDefault<word>("neighbourPatch", word::null)), coupleGroup_(dict), nbrPatchID_(-1), rotationAxis_(Zero), @@ -618,27 +608,19 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch rotationAngleDefined_(false), rotationAngle_(0.0), separationVector_(Zero), - AMIPtr_(nullptr), - AMIMethod_ + AMIPtr_ ( - AMIPatchToPatchInterpolation::interpolationMethodNames_ - [ - dict.getOrDefault - ( - "method", - AMIPatchToPatchInterpolation::interpolationMethodNames_ - [ - AMIPatchToPatchInterpolation::imFaceAreaWeight - ] - ) - ] + AMIInterpolation::New + ( + dict.getOrDefault<word>("AMIMethod", defaultAMIMethod), + dict, + dict.getOrDefault("flipNormals", false) + ) ), - AMIReverse_(dict.getOrDefault("flipNormals", false)), - AMIRequireMatch_(true), - AMILowWeightCorrection_(dict.getOrDefault("lowWeightCorrection", -1.0)), - surfPtr_(nullptr), surfDict_(dict.subOrEmptyDict("surface")), - createAMIFaces_(false), + surfPtr_(nullptr), + createAMIFaces_(dict.getOrDefault("createAMIFaces", false)), + moveFaceCentres_(false), updatingAMI_(true), srcFaceIDs_(), tgtFaceIDs_(), @@ -704,18 +686,13 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch // Neighbour patch might not be valid yet so no transformation // calculation possible - // If topology change, recover the sizes of the original patches - label srcSize0 = 0; - if (dict.readIfPresent("srcSize", srcSize0)) - { - srcFaceIDs_.setSize(srcSize0); - createAMIFaces_ = true; - } - label tgtSize0 = 0; - if (dict.readIfPresent("tgtSize", tgtSize0)) + // If topology change, recover the sizes of the original patches and + // read additional controls + if (createAMIFaces_) { - tgtFaceIDs_.setSize(tgtSize0); - createAMIFaces_ = true; + srcFaceIDs_.setSize(dict.get<label>("srcSize")); + tgtFaceIDs_.setSize(dict.get<label>("tgtSize")); + moveFaceCentres_ = dict.getOrDefault("moveFaceCentres", true); } } @@ -727,7 +704,6 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch ) : coupledPolyPatch(pp, bm), - updated_(false), nbrPatchName_(pp.nbrPatchName_), coupleGroup_(pp.coupleGroup_), nbrPatchID_(-1), @@ -736,14 +712,11 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch rotationAngleDefined_(pp.rotationAngleDefined_), rotationAngle_(pp.rotationAngle_), separationVector_(pp.separationVector_), - AMIPtr_(nullptr), - AMIMethod_(pp.AMIMethod_), - AMIReverse_(pp.AMIReverse_), - AMIRequireMatch_(pp.AMIRequireMatch_), - AMILowWeightCorrection_(pp.AMILowWeightCorrection_), - surfPtr_(nullptr), + AMIPtr_(pp.AMIPtr_->clone()), surfDict_(pp.surfDict_), - createAMIFaces_(false), + surfPtr_(nullptr), + createAMIFaces_(pp.createAMIFaces_), + moveFaceCentres_(pp.moveFaceCentres_), updatingAMI_(true), srcFaceIDs_(), tgtFaceIDs_(), @@ -766,7 +739,6 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch ) : coupledPolyPatch(pp, bm, index, newSize, newStart), - updated_(false), nbrPatchName_(nbrPatchName), coupleGroup_(pp.coupleGroup_), nbrPatchID_(-1), @@ -775,14 +747,11 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch rotationAngleDefined_(pp.rotationAngleDefined_), rotationAngle_(pp.rotationAngle_), separationVector_(pp.separationVector_), - AMIPtr_(nullptr), - AMIMethod_(pp.AMIMethod_), - AMIReverse_(pp.AMIReverse_), - AMIRequireMatch_(pp.AMIRequireMatch_), - AMILowWeightCorrection_(pp.AMILowWeightCorrection_), - surfPtr_(nullptr), + AMIPtr_(pp.AMIPtr_->clone()), surfDict_(pp.surfDict_), - createAMIFaces_(false), + surfPtr_(nullptr), + createAMIFaces_(pp.createAMIFaces_), + moveFaceCentres_(pp.moveFaceCentres_), updatingAMI_(true), srcFaceIDs_(), tgtFaceIDs_(), @@ -812,7 +781,6 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch ) : coupledPolyPatch(pp, bm, index, mapAddressing, newStart), - updated_(false), nbrPatchName_(pp.nbrPatchName_), coupleGroup_(pp.coupleGroup_), nbrPatchID_(-1), @@ -821,14 +789,11 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch rotationAngleDefined_(pp.rotationAngleDefined_), rotationAngle_(pp.rotationAngle_), separationVector_(pp.separationVector_), - AMIPtr_(nullptr), - AMIMethod_(pp.AMIMethod_), - AMIReverse_(pp.AMIReverse_), - AMIRequireMatch_(pp.AMIRequireMatch_), - AMILowWeightCorrection_(pp.AMILowWeightCorrection_), - surfPtr_(nullptr), + AMIPtr_(pp.AMIPtr_->clone()), surfDict_(pp.surfDict_), - createAMIFaces_(false), + surfPtr_(nullptr), + createAMIFaces_(pp.createAMIFaces_), + moveFaceCentres_(pp.moveFaceCentres_), updatingAMI_(true), srcFaceIDs_(), tgtFaceIDs_(), @@ -888,38 +853,6 @@ const Foam::cyclicAMIPolyPatch& Foam::cyclicAMIPolyPatch::neighbPatch() const } -const Foam::autoPtr<Foam::searchableSurface>& -Foam::cyclicAMIPolyPatch::surfPtr() const -{ - const word surfType(surfDict_.getOrDefault<word>("type", "none")); - - if (!surfPtr_.valid() && owner() && surfType != "none") - { - word surfName(surfDict_.getOrDefault("name", name())); - - const polyMesh& mesh = boundaryMesh().mesh(); - - surfPtr_ = - searchableSurface::New - ( - surfType, - IOobject - ( - surfName, - mesh.time().constant(), - "triSurface", - mesh, - IOobject::MUST_READ, - IOobject::NO_WRITE - ), - surfDict_ - ); - } - - return surfPtr_; -} - - const Foam::AMIPatchToPatchInterpolation& Foam::cyclicAMIPolyPatch::AMI() const { if (!owner()) @@ -929,9 +862,9 @@ const Foam::AMIPatchToPatchInterpolation& Foam::cyclicAMIPolyPatch::AMI() const << abort(FatalError); } - if (!AMIPtr_.valid()) + if (!AMIPtr_->upToDate()) { - resetAMI(AMIMethod_); + resetAMI(); } return *AMIPtr_; @@ -1206,27 +1139,7 @@ void Foam::cyclicAMIPolyPatch::write(Ostream& os) const } } - if (AMIMethod_ != AMIPatchToPatchInterpolation::imFaceAreaWeight) - { - os.writeEntry - ( - "method", - AMIPatchToPatchInterpolation::interpolationMethodNames_ - [ - AMIMethod_ - ] - ); - } - - if (AMIReverse_) - { - os.writeEntry("flipNormals", AMIReverse_); - } - - if (AMILowWeightCorrection_ > 0) - { - os.writeEntry("lowWeightCorrection", AMILowWeightCorrection_); - } + AMIPtr_->write(os); if (!surfDict_.empty()) { @@ -1235,8 +1148,10 @@ void Foam::cyclicAMIPolyPatch::write(Ostream& os) const if (createAMIFaces_) { + os.writeEntry("createAMIFaces", createAMIFaces_); os.writeEntry("srcSize", srcFaceIDs_.size()); os.writeEntry("tgtSize", tgtFaceIDs_.size()); + os.writeEntry("moveFaceCentres", moveFaceCentres_); } } diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H index 84fab30b99bae4209b920ab3b154c18723092cc4..de1c09d335ed4f80cd6335e50382b939906800eb 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2018-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -41,6 +42,7 @@ SourceFiles #include "AMIPatchToPatchInterpolation.H" #include "polyBoundaryMesh.H" #include "coupleGroupIdentifier.H" +#include "faceAreaWeightAMI.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -74,9 +76,6 @@ protected: // Protected data - //- Flag to indicate that AMI has been updated - mutable bool updated_; - //- Name of other half mutable word nbrPatchName_; @@ -113,24 +112,12 @@ protected: //- AMI interpolation class mutable autoPtr<AMIPatchToPatchInterpolation> AMIPtr_; - //- AMI method - const AMIPatchToPatchInterpolation::interpolationMethod AMIMethod_; - - //- Flag to indicate that slave patch should be reversed for AMI - const bool AMIReverse_; - - //- Flag to indicate that patches should match/overlap - bool AMIRequireMatch_; - - //- Low weight correction threshold for AMI - const scalar AMILowWeightCorrection_; + //- Dictionary used during projection surface construction + const dictionary surfDict_; //- Projection surface mutable autoPtr<searchableSurface> surfPtr_; - //- Dictionary used during projection surface construction - const dictionary surfDict_; - // Change of topology as AMI is updated @@ -138,6 +125,9 @@ protected: // Set by the call to changeTopology mutable bool createAMIFaces_; + //- Move face centres (default = no) + bool moveFaceCentres_; + mutable bool updatingAMI_; labelListList srcFaceIDs_; @@ -163,21 +153,14 @@ protected: virtual void restoreScaledGeometry(); + //- Create and return pointer to the projection surface + const autoPtr<searchableSurface>& surfPtr() const; //- Reset the AMI interpolator, supply patch points - virtual void resetAMI - ( - const UList<point>& points, - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod = - AMIPatchToPatchInterpolation::imFaceAreaWeight - ) const; + virtual void resetAMI(const UList<point>& points) const; //- Reset the AMI interpolator, use current patch points - virtual void resetAMI - ( - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod = - AMIPatchToPatchInterpolation::imFaceAreaWeight - ) const; + virtual void resetAMI() const; //- Recalculate the transformation tensors virtual void calcTransforms(); @@ -221,7 +204,8 @@ public: const label index, const polyBoundaryMesh& bm, const word& patchType, - const transformType transform = UNKNOWN + const transformType transform = UNKNOWN, + const word& defaultAMIMethod = faceAreaWeightAMI::typeName ); //- Construct from dictionary @@ -231,7 +215,8 @@ public: const dictionary& dict, const label index, const polyBoundaryMesh& bm, - const word& patchType + const word& patchType, + const word& defaultAMIMethod = faceAreaWeightAMI::typeName ); //- Construct as copy, resetting the boundary mesh @@ -328,23 +313,12 @@ public: //- Flag to indicate whether the AMI can be reset inline bool canResetAMI() const; - //- Reset the updated flag - inline void setUpdated(bool flag) const; - - //- Return access to the updated flag - inline bool updated() const; - //- Return access to the createAMIFaces flag inline bool createAMIFaces() const; //- Return access to the updated flag inline bool updatingAMI() const; -void clearAMI() const -{ - AMIPtr_.clear(); -} - //- Return true if this patch changes the mesh topology // True when createAMIFaces is true virtual bool changeTopology() const; @@ -371,9 +345,6 @@ void clearAMI() const //- Return a reference to the neighbour patch virtual const cyclicAMIPolyPatch& neighbPatch() const; - //- Return a reference to the projection surface - const autoPtr<searchableSurface>& surfPtr() const; - //- Return a reference to the AMI interpolator const AMIPatchToPatchInterpolation& AMI() const; diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H index 6d55bf1841dce5b9711ab772ee2f7d55e82e2f1e..f4c640d7c95356ef10562c585c1f6ad7cb91a51c 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2013 OpenFOAM Foundation - Copyright (C) 2019 OpenCFD Ltd. + Copyright (C) 2019-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -34,18 +34,6 @@ inline bool Foam::cyclicAMIPolyPatch::canResetAMI() const } -inline void Foam::cyclicAMIPolyPatch::setUpdated(const bool flag) const -{ - updated_ = flag; -} - - -inline bool Foam::cyclicAMIPolyPatch::updated() const -{ - return updated_; -} - - inline bool Foam::cyclicAMIPolyPatch::createAMIFaces() const { return createAMIFaces_; diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C index f7e33910da639eb3b5bb1dabb87360aba8f9fde3..a1c6634e2831b47a34adaefe0889cba00797afdf 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C @@ -2,8 +2,10 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2019 OpenCFD Ltd. + \\ / A nd | www.openfoam.com \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2019-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -53,29 +55,27 @@ void Foam::cyclicAMIPolyPatch::restoreScaledGeometry() vectorField::subField faceAreas = this->faceAreas(); vectorField::subField faceCentres = this->faceCentres(); - if (debug) + DebugInfo + << "Patch:" << name() << " before: sum(mag(faceAreas)):" + << gSum(mag(faceAreas)) << nl + << "Patch:" << name() << " before: sum(mag(faceAreas0)):" + << gSum(mag(faceAreas0_)) << endl; + + faceAreas = faceAreas0_; + if (moveFaceCentres_) { - Info<< "Patch:" << name() << " before: sum(mag(faceAreas)):" - << gSum(mag(faceAreas)) << endl; - Info<< "Patch:" << name() << " before: sum(mag(faceAreas0)):" - << gSum(mag(faceAreas0_)) << endl; + DebugInfo << "Moving face centres" << endl; + faceCentres = faceCentres0_; } - faceAreas = faceAreas0_; - Info<< "WARNING: ==============================" << endl; - Info<< "WARNING: Not updating face cell centres" << endl; - Info<< "WARNING: ==============================" << endl; - //faceCentres = faceCentres0_; faceAreas0_.clear(); faceCentres0_.clear(); - if (debug) - { - Info<< "Patch:" << name() << " after: sum(mag(faceAreas)):" - << gSum(mag(faceAreas)) << endl; - Info<< "Patch:" << name() << " after: sum(mag(faceAreas0)):" - << gSum(mag(faceAreas0_)) << endl; - } + DebugInfo + << "Patch:" << name() << " after: sum(mag(faceAreas)):" + << gSum(mag(faceAreas)) << nl + << "Patch:" << name() << " after: sum(mag(faceAreas0)):" + << gSum(mag(faceAreas0_)) << endl; } @@ -614,9 +614,9 @@ void Foam::cyclicAMIPolyPatch::setAMIFaces() } } - // Update the AMI addressing and weights to reflect the new 1-to-1 + // Reset the AMI addressing and weights to reflect the new 1-to-1 // correspondence - AMIPtr_->update + AMIPtr_->reset ( std::move(srcToTgtMap1), std::move(tgtToSrcMap1), @@ -626,7 +626,9 @@ void Foam::cyclicAMIPolyPatch::setAMIFaces() std::move(newTgtToSrcWeights) ); - AMIPtr_->setAreas(mag(faceAreas0_), mag(nbrFaceAreas0)); + // Need to set areas, e.g. for agglomeration to (re-)normalisation weights + AMIPtr_->srcMagSf() = mag(faceAreas0_); + AMIPtr_->tgtMagSf() = mag(nbrFaceAreas0); if (debug) { @@ -660,7 +662,7 @@ bool Foam::cyclicAMIPolyPatch::setTopology(polyTopoChange& topoChange) { // Calculate the AMI using the new points // Note: mesh still has old points - resetAMI(topoChange.points(), AMIMethod_); + resetAMI(topoChange.points()); removeAMIFaces(topoChange); diff --git a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C index 28c6b98f6d48420a7e1bfe3bd21928ca86eeb78c..d8ade6e0cc2e1538157ee023229f0c64522e44f5 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.C @@ -26,6 +26,7 @@ License \*---------------------------------------------------------------------------*/ #include "cyclicPeriodicAMIPolyPatch.H" +#include "partialFaceAreaWeightAMI.H" #include "addToRunTimeSelectionTable.H" // For debugging @@ -228,10 +229,7 @@ void Foam::cyclicPeriodicAMIPolyPatch::writeOBJ } -void Foam::cyclicPeriodicAMIPolyPatch::resetAMI -( - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod -) const +void Foam::cyclicPeriodicAMIPolyPatch::resetAMI() const { if (owner()) { @@ -326,20 +324,8 @@ void Foam::cyclicPeriodicAMIPolyPatch::resetAMI ); // Construct a new AMI interpolation between the initial patch locations - AMIPtr_.reset - ( - new AMIPatchToPatchInterpolation - ( - thisPatch0, - nbrPatch0, - surfPtr(), - faceAreaIntersect::tmMesh, - false, - AMIPatchToPatchInterpolation::imPartialFaceAreaWeight, - AMILowWeightCorrection_, - AMIReverse_ - ) - ); + AMIPtr_->setRequireMatch(false); + AMIPtr_->calculate(thisPatch0, nbrPatch0, surfPtr()); // Number of geometry replications label iter(0); @@ -358,7 +344,6 @@ void Foam::cyclicPeriodicAMIPolyPatch::resetAMI // Weight sum averages scalar srcSum(gAverage(AMIPtr_->srcWeightsSum())); scalar tgtSum(gAverage(AMIPtr_->tgtWeightsSum())); - // Direction (or rather side of AMI : this or nbr patch) of // geometry replication bool direction = nTransforms_ >= 0; @@ -520,17 +505,13 @@ void Foam::cyclicPeriodicAMIPolyPatch::resetAMI << "The current matchTolerance : " << matchTolerance() << ", sum of owner weights : " << srcSum << ", sum of neighbour weights : " << tgtSum - << "." << nl + << "." << nl << "This is only acceptable during post-processing" << "; not during running. Improve your mesh or increase" << " the 'matchTolerance' setting in the patch specification." << endl; } - // Normalise the weights. Disable printing since weights are - // still areas. - AMIPtr_->normaliseWeights(true, false); - // Print some statistics const label nFace = returnReduce(size(), sumOp<label>()); @@ -577,7 +558,17 @@ Foam::cyclicPeriodicAMIPolyPatch::cyclicPeriodicAMIPolyPatch const transformType transform ) : - cyclicAMIPolyPatch(name, size, start, index, bm, patchType, transform), + cyclicAMIPolyPatch + ( + name, + size, + start, + index, + bm, + patchType, + transform, + partialFaceAreaWeightAMI::typeName + ), periodicPatchName_(word::null), periodicPatchID_(-1), nTransforms_(0), @@ -595,7 +586,15 @@ Foam::cyclicPeriodicAMIPolyPatch::cyclicPeriodicAMIPolyPatch const word& patchType ) : - cyclicAMIPolyPatch(name, dict, index, bm, patchType), + cyclicAMIPolyPatch + ( + name, + dict, + index, + bm, + patchType, + partialFaceAreaWeightAMI::typeName + ), periodicPatchName_(dict.lookup("periodicPatch")), periodicPatchID_(-1), nTransforms_(dict.getOrDefault<label>("nTransforms", 0)), diff --git a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H index 22fbc088455d3b3f30373dcd552e1c769d7e59bf..297c5c3189310d8ef90faf64244c4b0332005075 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicPeriodicAMI/cyclicPeriodicAMIPolyPatch/cyclicPeriodicAMIPolyPatch.H @@ -83,11 +83,7 @@ private: void writeOBJ(const primitivePatch& p, OBJstream& str) const; //- Reset the AMI interpolator - virtual void resetAMI - ( - const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod = - AMIPatchToPatchInterpolation::imFaceAreaWeight - ) const; + virtual void resetAMI() const; public: diff --git a/src/meshTools/Make/files b/src/meshTools/Make/files index c89a1afb64f04ff8dc8dda1d62905eedb4807436..e428fc63dcbdf3ce93c549d8e2c35b97e428872e 100644 --- a/src/meshTools/Make/files +++ b/src/meshTools/Make/files @@ -259,8 +259,13 @@ processorLOD/cellBox/cellBox.C processorLOD/faceBox/faceBox.C AMI=AMIInterpolation -$(AMI)/AMIInterpolation/AMIInterpolationName.C -$(AMI)/AMIInterpolation/AMIPatchToPatchInterpolation.C +$(AMI)/AMIInterpolation/AMIInterpolation.C +$(AMI)/AMIInterpolation/AMIInterpolationNew.C +$(AMI)/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.C +$(AMI)/AMIInterpolation/advancingFrontAMI/advancingFrontAMIParallelOps.C +$(AMI)/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.C +$(AMI)/AMIInterpolation/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C +$(AMI)/AMIInterpolation/nearestFaceAMI/nearestFaceAMI.C $(AMI)/faceAreaIntersect/faceAreaIntersect.C $(AMI)/GAMG/interfaces/cyclicAMIGAMGInterface/cyclicAMIGAMGInterface.C $(AMI)/GAMG/interfaceFields/cyclicAMIGAMGInterfaceField/cyclicAMIGAMGInterfaceField.C diff --git a/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C b/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C index 8dac4e801c5cc156b8fb66c6617ed3241dd70da3..a88173d6bce6befe9fdbbef70ebb264b1f9b46a3 100644 --- a/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C +++ b/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C @@ -45,6 +45,7 @@ License #include "syncTools.H" #include "treeDataCell.H" #include "DynamicField.H" +#include "faceAreaWeightAMI.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -791,15 +792,13 @@ const void Foam::mappedPatchBase::calcAMI() const { - if (AMIPtr_.valid()) + if (AMIPtr_->upToDate()) { - FatalErrorInFunction - << "AMI already calculated" << exit(FatalError); + WarningInFunction + << "AMI already up-to-date" + << endl; } - AMIPtr_.clear(); - - const polyPatch& nbr = samplePolyPatch(); // Transform neighbour patch to local system @@ -829,29 +828,13 @@ void Foam::mappedPatchBase::calcAMI() const } // Construct/apply AMI interpolation to determine addressing and weights - AMIPtr_.reset - ( - new AMIPatchToPatchInterpolation - ( - patch_, - nbrPatch0, - surfPtr(), - faceAreaIntersect::tmMesh, - true, - AMIPatchToPatchInterpolation::imFaceAreaWeight, - -1, - AMIReverse_ - ) - ); + AMIPtr_->calculate(patch_, nbrPatch0, surfPtr()); } // * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * * * // -Foam::mappedPatchBase::mappedPatchBase -( - const polyPatch& pp -) +Foam::mappedPatchBase::mappedPatchBase(const polyPatch& pp) : patch_(pp), sampleRegion_(patch_.boundaryMesh().mesh().name()), @@ -864,8 +847,8 @@ Foam::mappedPatchBase::mappedPatchBase distance_(0), sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(false), + AMIPtr_(new faceAreaWeightAMI(true, AMIReverse_)), surfPtr_(nullptr), surfDict_(fileName("surface")) {} @@ -891,8 +874,8 @@ Foam::mappedPatchBase::mappedPatchBase distance_(0), sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(false), + AMIPtr_(new faceAreaWeightAMI(true, AMIReverse_)), surfPtr_(nullptr), surfDict_(fileName("surface")) {} @@ -918,8 +901,8 @@ Foam::mappedPatchBase::mappedPatchBase distance_(0), sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(false), + AMIPtr_(new faceAreaWeightAMI(true, AMIReverse_)), surfPtr_(nullptr), surfDict_(fileName("surface")) {} @@ -945,8 +928,8 @@ Foam::mappedPatchBase::mappedPatchBase distance_(distance), sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(false), + AMIPtr_(new faceAreaWeightAMI(true, AMIReverse_)), surfPtr_(nullptr), surfDict_(fileName("surface")) {} @@ -969,8 +952,16 @@ Foam::mappedPatchBase::mappedPatchBase distance_(0), sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(dict.getOrDefault("flipNormals", false)), + AMIPtr_ + ( + AMIInterpolation::New + ( + dict.getOrDefault("AMIMethod", faceAreaWeightAMI::typeName), + dict, + AMIReverse_ + ) + ), surfPtr_(nullptr), surfDict_(dict.subOrEmptyDict("surface")) { @@ -1044,8 +1035,16 @@ Foam::mappedPatchBase::mappedPatchBase distance_(0), sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(dict.getOrDefault("flipNormals", false)), + AMIPtr_ + ( + AMIInterpolation::New + ( + dict.lookupOrDefault("AMIMethod", faceAreaWeightAMI::typeName), + dict, + AMIReverse_ + ) + ), surfPtr_(nullptr), surfDict_(dict.subOrEmptyDict("surface")) { @@ -1089,8 +1088,8 @@ Foam::mappedPatchBase::mappedPatchBase distance_(mpb.distance_), sameRegion_(mpb.sameRegion_), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(mpb.AMIReverse_), + AMIPtr_(mpb.AMIPtr_->clone()), surfPtr_(nullptr), surfDict_(mpb.surfDict_) {} @@ -1119,8 +1118,8 @@ Foam::mappedPatchBase::mappedPatchBase distance_(mpb.distance_), sameRegion_(mpb.sameRegion_), mapPtr_(nullptr), - AMIPtr_(nullptr), AMIReverse_(mpb.AMIReverse_), + AMIPtr_(mpb.AMIPtr_->clone()), surfPtr_(nullptr), surfDict_(mpb.surfDict_) {} @@ -1137,8 +1136,8 @@ Foam::mappedPatchBase::~mappedPatchBase() void Foam::mappedPatchBase::clearOut() { mapPtr_.clear(); - AMIPtr_.clear(); surfPtr_.clear(); + AMIPtr_->upToDate() = false; } diff --git a/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.H b/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.H index ac4f6c7fb514960d71a0d9d195482878a8660445..d6b76023dba3cc5da7cf434398d2819e71a740f2 100644 --- a/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.H +++ b/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.H @@ -231,12 +231,12 @@ protected: // AMI interpolator (only for NEARESTPATCHFACEAMI) - //- Pointer to AMI interpolator - mutable autoPtr<AMIPatchToPatchInterpolation> AMIPtr_; - //- Flag to indicate that slave patch should be reversed for AMI const bool AMIReverse_; + //- Pointer to AMI interpolator + mutable autoPtr<AMIPatchToPatchInterpolation> AMIPtr_; + //- Pointer to projection surface employed by AMI interpolator mutable autoPtr<searchableSurface> surfPtr_;