diff --git a/src/foam/interpolations/GGIInterpolation/GGIInterpolate.C b/src/foam/interpolations/GGIInterpolation/GGIInterpolate.C index 2e15382a6..b5fa94572 100644 --- a/src/foam/interpolations/GGIInterpolation/GGIInterpolate.C +++ b/src/foam/interpolations/GGIInterpolation/GGIInterpolate.C @@ -97,26 +97,13 @@ void GGIInterpolation::bridge ( const Field& bridgeField, Field& ff, - const labelList& addr, - const labelList& partiallyCoveredAddr, - const scalarField& coveredFractions + const labelList& addr ) { - // Fully uncovered faces forAll (addr, faceI) { ff[addr[faceI]] = bridgeField[addr[faceI]]; } - - // Loop through partially covered faces and correct them. Note the - // operator+= since we assume that the interpolation part is carried out - // before bridging (see e.g. ggiFvPatchField::patchNeighbourField()) using - // weights that do not sum up to 1 - forAll (partiallyCoveredAddr, pcfI) - { - ff[partiallyCoveredAddr[pcfI]] += - coveredFractions[pcfI]*bridgeField[partiallyCoveredAddr[pcfI]]; - } } @@ -127,9 +114,7 @@ void GGIInterpolation::maskedBridge const Field& bridgeField, Field& ff, const labelList& mask, - const labelList& uncoveredFaces, - const labelList& partiallyCoveredAddr, - const scalarField& coveredFractions + const labelList& uncoveredFaces ) { // Note: tricky algorithm @@ -146,7 +131,7 @@ void GGIInterpolation::maskedBridge const label faceI = uncoveredFaces[uncoI]; // Search through the mask - for (; maskAddrI < mask.size(); ++maskAddrI) + for (; maskAddrI < mask.size(); maskAddrI++) { if (faceI == mask[maskAddrI]) { @@ -162,38 +147,7 @@ void GGIInterpolation::maskedBridge // Go one back and check for next uncovered face if (maskAddrI > 0) { - --maskAddrI; - } - - break; - } - } - } - - // Reset maskAddrI - maskAddrI = 0; - - forAll (partiallyCoveredAddr, pcfI) - { - // Pick partially covered face - const label faceI = partiallyCoveredAddr[pcfI]; - - for (; maskAddrI < mask.size(); ++maskAddrI) - { - if (faceI == mask[maskAddrI]) - { - // Found masked partially covered face - ff[maskAddrI] += coveredFractions[pcfI]*bridgeField[maskAddrI]; - - break; - } - else if (mask[maskAddrI] > faceI) - { - // Gone beyond my index: my face is not present in the mask - // Go one back and check for next uncovered face - if (maskAddrI > 0) - { - --maskAddrI; + maskAddrI--; } break; @@ -557,8 +511,7 @@ void GGIInterpolation::bridgeMaster if ( bridgeField.size() != masterPatch_.size() - || ff.size() != masterPatch_.size() - ) + || ff.size() != masterPatch_.size()) { FatalErrorIn ( @@ -574,14 +527,7 @@ void GGIInterpolation::bridgeMaster << abort(FatalError); } - bridge - ( - bridgeField, - ff, - uncoveredMasterFaces(), - partiallyCoveredMasterFaces(), - masterFaceCoveredFractions() - ); + bridge(bridgeField, ff, uncoveredMasterFaces()); } @@ -605,7 +551,7 @@ void GGIInterpolation::maskedBridgeMaster " Field& ff,\n" " const labelList& mask\n" ") const" - ) << "given field does not correspond to patch. Patch (mask) size: " + ) << "given field does not correspond to patch. Patch size: " << masterPatch_.size() << " bridge field size: " << bridgeField.size() << " field size: " << ff.size() @@ -613,15 +559,7 @@ void GGIInterpolation::maskedBridgeMaster << abort(FatalError); } - maskedBridge - ( - bridgeField, - ff, - mask, - uncoveredMasterFaces(), - partiallyCoveredMasterFaces(), - masterFaceCoveredFractions() - ); + maskedBridge(bridgeField, ff, mask, uncoveredMasterFaces()); } @@ -654,14 +592,7 @@ void GGIInterpolation::bridgeSlave << abort(FatalError); } - bridge - ( - bridgeField, - ff, - uncoveredSlaveFaces(), - partiallyCoveredSlaveFaces(), - slaveFaceCoveredFractions() - ); + bridge(bridgeField, ff, uncoveredSlaveFaces()); } @@ -685,7 +616,7 @@ void GGIInterpolation::maskedBridgeSlave " Field& ff\n," " const labelList& mask\n" ") const" - ) << "given field does not correspond to patch. Patch (mask) size: " + ) << "given field does not correspond to patch. Patch size: " << slavePatch_.size() << " bridge field size: " << bridgeField.size() << " field size: " << ff.size() @@ -693,18 +624,12 @@ void GGIInterpolation::maskedBridgeSlave << abort(FatalError); } - maskedBridge - ( - bridgeField, - ff, - mask, - uncoveredSlaveFaces(), - partiallyCoveredSlaveFaces(), - slaveFaceCoveredFractions() - ); + maskedBridge(bridgeField, ff, mask, uncoveredSlaveFaces()); } + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam diff --git a/src/foam/interpolations/GGIInterpolation/GGIInterpolation.C b/src/foam/interpolations/GGIInterpolation/GGIInterpolation.C index 4bdc148ca..49d86233e 100644 --- a/src/foam/interpolations/GGIInterpolation/GGIInterpolation.C +++ b/src/foam/interpolations/GGIInterpolation/GGIInterpolation.C @@ -104,11 +104,7 @@ void GGIInterpolation::clearOut() deleteDemandDrivenData(slaveWeightsPtr_); deleteDemandDrivenData(uncoveredMasterAddrPtr_); - deleteDemandDrivenData(partiallyCoveredMasterAddrPtr_); - deleteDemandDrivenData(masterFaceCoveredFractionsPtr_); deleteDemandDrivenData(uncoveredSlaveAddrPtr_); - deleteDemandDrivenData(partiallyCoveredSlaveAddrPtr_); - deleteDemandDrivenData(slaveFaceCoveredFractionsPtr_); deleteDemandDrivenData(masterPointAddressingPtr_); deleteDemandDrivenData(masterPointWeightsPtr_); @@ -159,11 +155,7 @@ GGIInterpolation::GGIInterpolation slavePointWeightsPtr_(NULL), slavePointDistancePtr_(NULL), uncoveredMasterAddrPtr_(NULL), - partiallyCoveredMasterAddrPtr_(NULL), - masterFaceCoveredFractionsPtr_(NULL), - uncoveredSlaveAddrPtr_(NULL), - partiallyCoveredSlaveAddrPtr_(NULL), - slaveFaceCoveredFractionsPtr_(NULL) + uncoveredSlaveAddrPtr_(NULL) { // Check size of transform. They should be equal to slave patch size // if the transform is not constant @@ -264,32 +256,6 @@ GGIInterpolation::uncoveredMasterFaces() const } -template -const labelList& -GGIInterpolation::partiallyCoveredMasterFaces() const -{ - if (!partiallyCoveredMasterAddrPtr_) - { - calcAddressing(); - } - - return *partiallyCoveredMasterAddrPtr_; -} - - -template -const scalarField& -GGIInterpolation::masterFaceCoveredFractions() const -{ - if (!masterFaceCoveredFractionsPtr_) - { - calcAddressing(); - } - - return *masterFaceCoveredFractionsPtr_; -} - - template const labelList& GGIInterpolation::uncoveredSlaveFaces() const @@ -303,32 +269,6 @@ GGIInterpolation::uncoveredSlaveFaces() const } -template -const labelList& -GGIInterpolation::partiallyCoveredSlaveFaces() const -{ - if (!partiallyCoveredSlaveAddrPtr_) - { - calcAddressing(); - } - - return *partiallyCoveredSlaveAddrPtr_; -} - - -template -const scalarField& -GGIInterpolation::slaveFaceCoveredFractions() const -{ - if (!slaveFaceCoveredFractionsPtr_) - { - calcAddressing(); - } - - return *slaveFaceCoveredFractionsPtr_; -} - - template bool GGIInterpolation::movePoints ( diff --git a/src/foam/interpolations/GGIInterpolation/GGIInterpolationTemplate.H b/src/foam/interpolations/GGIInterpolation/GGIInterpolationTemplate.H index 4dd21035c..3e4b3d305 100644 --- a/src/foam/interpolations/GGIInterpolation/GGIInterpolationTemplate.H +++ b/src/foam/interpolations/GGIInterpolation/GGIInterpolationTemplate.H @@ -218,29 +218,9 @@ class GGIInterpolation //- List of uncovered master faces mutable labelList* uncoveredMasterAddrPtr_; - //- List of partially covered master faces. Contains faces with - // weights between masterNonOverlapFaceTol_ and - // uncoveredFaceAreaTol_. Used in bridging - mutable labelList* partiallyCoveredMasterAddrPtr_; - - //- List of face fractions for partially covered master faces. Range - // anywhere from 1 to uncoveredFaceAreaTol_ (e.g. 0.2 means - // that the face area is covered by 20% on the other side) - mutable scalarField* masterFaceCoveredFractionsPtr_; - //- List of uncovered slave faces mutable labelList* uncoveredSlaveAddrPtr_; - //- List of partially covered slave faces. Contains faces with - // weights between masterNonOverlapFaceTol_ and - // uncoveredFaceAreaTol_. Used in bridging - mutable labelList* partiallyCoveredSlaveAddrPtr_; - - //- List of face fractions for partially covered slave faces. Range - // anywhere from 1 to uncoveredFaceAreaTol_ (e.g. 0.2 means - // that the face area is covered by 20% on the other side) - mutable scalarField* slaveFaceCoveredFractionsPtr_; - // Private static data @@ -253,9 +233,6 @@ class GGIInterpolation //- Facet bound box extension factor static const debug::tolerancesSwitch faceBoundBoxExtendSpanFraction_; - //- Tolerance for uncovered face areas - static const debug::tolerancesSwitch uncoveredFaceAreaTol_; - //- The next 3 attributes are parameters controlling the creation // of an octree search engine for the GGI facets neighbourhood // determination. @@ -440,14 +417,6 @@ class GGIInterpolation const scalar& nonOverlapFaceTol ) const; - //- Calculate partially covered faces - void calcPartiallyCoveredFaces - ( - const scalarListList& patchWeights, - const scalar& nonOverlapFaceTo, - const bool isMaster - ) const; - //- Clear all geometry and addressing void clearOut(); @@ -482,9 +451,7 @@ class GGIInterpolation ( const Field& bridgeField, Field& ff, - const labelList& addr, - const labelList& partiallyCoveredAddr, - const scalarField& coveredFractions + const labelList& addr ); //- Bridge uncovered faces given addressing @@ -495,9 +462,7 @@ class GGIInterpolation const Field& bridgeField, Field& ff, const labelList& mask, - const labelList& uncoveredFaces, - const labelList& partiallyCoveredAddr, - const scalarField& coveredFractions + const labelList& uncoveredFaces ); //- Is a transform required? @@ -568,21 +533,9 @@ public: //- Return reference to the master list of non-overlap faces const labelList& uncoveredMasterFaces() const; - //- Return reference to master list of partially covered faces - const labelList& partiallyCoveredMasterFaces() const; - - //- Return covered fractions of partially covered master faces - const scalarField& masterFaceCoveredFractions() const; - //- Return reference to the slave list of non-overlap faces const labelList& uncoveredSlaveFaces() const; - //- Return reference to slave list of partially covered faces - const labelList& partiallyCoveredSlaveFaces() const; - - //- Return covered fractions of partially covered slave faces - const scalarField& slaveFaceCoveredFractions() const; - //- Return reference to point addressing const List& masterPointAddr() const; @@ -601,7 +554,6 @@ public: //- Return distance to intersection for patch points const scalarField& slavePointDistanceToIntersection() const; - // Interpolation functions //- Interpolate from master to slave @@ -688,6 +640,12 @@ public: const Field& pf ) const; + template + tmp > slaveToMasterPointInterpolate + ( + const Field& pf + ) const; + // Edit diff --git a/src/foam/interpolations/GGIInterpolation/GGIInterpolationWeights.C b/src/foam/interpolations/GGIInterpolation/GGIInterpolationWeights.C index 4b3225498..48845d702 100644 --- a/src/foam/interpolations/GGIInterpolation/GGIInterpolationWeights.C +++ b/src/foam/interpolations/GGIInterpolation/GGIInterpolationWeights.C @@ -68,16 +68,6 @@ GGIInterpolation::featureCosTol_ ); -template -const Foam::debug::tolerancesSwitch -GGIInterpolation::uncoveredFaceAreaTol_ -( - "GGIUncoveredFaceAreaTol", - 0.999, - "Fraction of face area mismatch (sum of weights) to consider a face " - "as uncovered, i.e. not to rescale weights." -); - // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // template @@ -90,11 +80,7 @@ void GGIInterpolation::calcAddressing() const || slaveAddrPtr_ || slaveWeightsPtr_ || uncoveredMasterAddrPtr_ - || partiallyCoveredMasterAddrPtr_ - || masterFaceCoveredFractionsPtr_ || uncoveredSlaveAddrPtr_ - || partiallyCoveredSlaveAddrPtr_ - || slaveFaceCoveredFractionsPtr_ ) { FatalErrorIn @@ -445,7 +431,7 @@ void GGIInterpolation::calcAddressing() const ); // Fix: previously checked for VSMALL. - // HJ, 19/Sep/2016 + // HJ, 19/Sep2/106 if (intersectionArea > intersectionTestArea) { // We compute the GGI weights based on this @@ -569,7 +555,7 @@ void GGIInterpolation::calcAddressing() const const labelList& curMa = ma[mfI]; const scalarList& cursmaW = smaW[mfI]; - // Current master face index + // Current master face indes const label faceMaster = mfI; forAll (curMa, mAI) @@ -615,29 +601,6 @@ void GGIInterpolation::calcAddressing() const findNonOverlappingFaces(saW, slaveNonOverlapFaceTol_) ); - // Calculate master and slave partially covered addressing - - // Note: This function allocates: - // 1. partiallyCoveredMasterAddrPtr_ - // 2. masterFaceCoveredFractionsPtr_ - calcPartiallyCoveredFaces - ( - maW, - masterNonOverlapFaceTol_, - true // This is master - ); - - // Note: this function allocates: - // 1. partiallyCoveredSlaveAddrPtr_ - // 2. slaveFaceCoveredFractionsPtr_ - calcPartiallyCoveredFaces - ( - saW, - slaveNonOverlapFaceTol_, - false // This is not master - ); - - // Rescaling the weighting factors so they will sum up to 1.0 // See the comment for the method ::rescaleWeightingFactors() for // more information. By default, we always rescale. But for some @@ -645,8 +608,7 @@ void GGIInterpolation::calcAddressing() const // then we need the brute values, so no rescaling in that // case. Hence the little flag rescaleGGIWeightingFactors_ - // Not parallelised. HJ, 27/Apr/2016. Rescaling is not performed for - // partially overlapping faces for their correct treatment. VV, 16/Oct/2017. + // Not parallelised. HJ, 27/Apr/2016 if (rescaleGGIWeightingFactors_) { rescaleWeightingFactors(); @@ -694,35 +656,6 @@ void GGIInterpolation::rescaleWeightingFactors() const scalar sumMWC = 0; scalar curMWC = 0; - // Note: do not rescale weighting factors for partially covered faces - if (!partiallyCoveredMasterAddrPtr_ || !partiallyCoveredSlaveAddrPtr_) - { - FatalErrorIn - ( - "void GGIInterpolation::" - "rescaleWeightingFactors() const" - ) << "Master or slave partially covered faces are not calculated." - << abort(FatalError); - } - - const labelList& partiallyCoveredMasterFaces = - *partiallyCoveredMasterAddrPtr_; - const labelList& partiallyCoveredSlaveFaces = - *partiallyCoveredSlaveAddrPtr_; - - // Create a mask for partially covered master/slave faces - boolList masterPCMask(maW.size(), false); - boolList slavePCMask(saW.size(), false); - - forAll (partiallyCoveredMasterFaces, pfmI) - { - masterPCMask[partiallyCoveredMasterFaces[pfmI]] = true; - } - forAll (partiallyCoveredSlaveFaces, pfsI) - { - slavePCMask[partiallyCoveredSlaveFaces[pfsI]] = true; - } - // Rescaling the slave weights if (debug) { @@ -746,7 +679,7 @@ void GGIInterpolation::rescaleWeightingFactors() const { scalar slaveWeightSum = Foam::sum(saW[saWi]); - if (saW[saWi].size() > 0 && !slavePCMask[saWi]) + if (saW[saWi].size() > 0) { saW[saWi] = saW[saWi]/slaveWeightSum; @@ -763,7 +696,7 @@ void GGIInterpolation::rescaleWeightingFactors() const { scalar masterWeightSum = Foam::sum(maW[maWi]); - if (maW[maWi].size() > 0 && !masterPCMask[maWi]) + if (maW[maWi].size() > 0) { maW[maWi] = maW[maWi]/masterWeightSum; @@ -834,119 +767,6 @@ GGIInterpolation::findNonOverlappingFaces return tpatchFaceNonOverlapAddr; } - -template -void GGIInterpolation::calcPartiallyCoveredFaces -( - const scalarListList& patchWeights, - const scalar& nonOverlapFaceTol, - const bool isMaster -) const -{ - // Sanity checks first - if - ( - isMaster - && (partiallyCoveredMasterAddrPtr_ || masterFaceCoveredFractionsPtr_) - ) - { - FatalErrorIn - ( - "void GGIInterpolation::" - "calcPartiallyCoveredFaces() const" - ) << "Already calculated master partially covered faces" - << abort(FatalError); - } - - if - ( - !isMaster - && (partiallyCoveredSlaveAddrPtr_ || slaveFaceCoveredFractionsPtr_) - ) - { - FatalErrorIn - ( - "void GGIInterpolation::" - "calcPartiallyCoveredFaces() const" - ) << "Already calculated slave partially covered faces" - << abort(FatalError); - } - - // Temporary storage - DynamicList patchFacePartialOverlap(patchWeights.size()); - DynamicList uncoveredFaceFractions(patchWeights.size()); - - // Scan the list of patch weights and collect ones inbetween - // nonOverlapFaceTol and uncoveredFaceAreaTol_ - forAll (patchWeights, paWi) - { - const scalar sumWeightsFace = sum(patchWeights[paWi]); - - if - ( - sumWeightsFace > nonOverlapFaceTol - && sumWeightsFace <= uncoveredFaceAreaTol_() - ) - { - // This is considered partially overlapped face, store the index and - // the non-overlapping area (1 - sum of weights) - patchFacePartialOverlap.append(paWi); - uncoveredFaceFractions.append(1.0 - sumWeightsFace); - } - } - - // Transfer the storage - if (isMaster) - { - // Allocate master side - partiallyCoveredMasterAddrPtr_ = - new labelList(patchFacePartialOverlap.xfer()); - masterFaceCoveredFractionsPtr_ = - new scalarField(uncoveredFaceFractions.xfer()); - - if (debug) - { - InfoIn("GGIInterpolation::calcPartiallyCoveredFaces") - << " : Found " << partiallyCoveredMasterAddrPtr_->size() - << " partially overlapping faces for master GGI patch" << endl; - - if (partiallyCoveredMasterAddrPtr_->size()) - { - Info<< "Max coverage: " - << max(*masterFaceCoveredFractionsPtr_) - << ", min coverage: " - << min(*masterFaceCoveredFractionsPtr_) - << endl; - } - } - } - else - { - // Allocate slave side - partiallyCoveredSlaveAddrPtr_ = - new labelList(patchFacePartialOverlap.xfer()); - slaveFaceCoveredFractionsPtr_ = - new scalarField(uncoveredFaceFractions.xfer()); - - if (debug) - { - InfoIn("GGIInterpolation::calcPartiallyCoveredFaces") - << " : Found " << partiallyCoveredSlaveAddrPtr_->size() - << " partially overlapping faces for slave GGI patch" << endl; - - if (partiallyCoveredSlaveAddrPtr_->size()) - { - Info<< "Max coverage: " - << max(*slaveFaceCoveredFractionsPtr_) - << ", min coverage: " - << min(*slaveFaceCoveredFractionsPtr_) - << endl; - } - } - } -} - - template void GGIInterpolation:: calcMasterPointAddressing() const