exaFOAM/foam-extend4.1-coherent-io
Archived
1
1
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

Revert GGI change

Browse source
This commit is contained in:
Hrvoje Jasak 2018-01-02 11:24:24 +00:00
parent 5460dcb963
commit d1b82b1a66
4 changed files with 27 additions and 384 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

101
src/foam/interpolations/GGIInterpolation/GGIInterpolate.C
View file

@ -97,26 +97,13 @@ void GGIInterpolation<MasterPatch, SlavePatch>::bridge
(
const Field<Type>& bridgeField,
Field<Type>& 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<MasterPatch, SlavePatch>::maskedBridge
const Field<Type>& bridgeField,
Field<Type>& 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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::bridgeMaster
if
(
bridgeField.size() != masterPatch_.size()
|| ff.size() != masterPatch_.size()
)
|| ff.size() != masterPatch_.size())
{
FatalErrorIn
(
@ -574,14 +527,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::bridgeMaster
<< abort(FatalError);
}
bridge
(
bridgeField,
ff,
uncoveredMasterFaces(),
partiallyCoveredMasterFaces(),
masterFaceCoveredFractions()
);
bridge(bridgeField, ff, uncoveredMasterFaces());
}
@ -605,7 +551,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::maskedBridgeMaster
" Field<Type>& 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<MasterPatch, SlavePatch>::maskedBridgeMaster
<< abort(FatalError);
}
maskedBridge
(
bridgeField,
ff,
mask,
uncoveredMasterFaces(),
partiallyCoveredMasterFaces(),
masterFaceCoveredFractions()
);
maskedBridge(bridgeField, ff, mask, uncoveredMasterFaces());
}
@ -654,14 +592,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::bridgeSlave
<< abort(FatalError);
}
bridge
(
bridgeField,
ff,
uncoveredSlaveFaces(),
partiallyCoveredSlaveFaces(),
slaveFaceCoveredFractions()
);
bridge(bridgeField, ff, uncoveredSlaveFaces());
}
@ -685,7 +616,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::maskedBridgeSlave
" Field<Type>& 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<MasterPatch, SlavePatch>::maskedBridgeSlave
<< abort(FatalError);
}
maskedBridge
(
bridgeField,
ff,
mask,
uncoveredSlaveFaces(),
partiallyCoveredSlaveFaces(),
slaveFaceCoveredFractions()
);
maskedBridge(bridgeField, ff, mask, uncoveredSlaveFaces());
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

62
src/foam/interpolations/GGIInterpolation/GGIInterpolation.C
View file

@ -104,11 +104,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::clearOut()
deleteDemandDrivenData(slaveWeightsPtr_);
deleteDemandDrivenData(uncoveredMasterAddrPtr_);
deleteDemandDrivenData(partiallyCoveredMasterAddrPtr_);
deleteDemandDrivenData(masterFaceCoveredFractionsPtr_);
deleteDemandDrivenData(uncoveredSlaveAddrPtr_);
deleteDemandDrivenData(partiallyCoveredSlaveAddrPtr_);
deleteDemandDrivenData(slaveFaceCoveredFractionsPtr_);
deleteDemandDrivenData(masterPointAddressingPtr_);
deleteDemandDrivenData(masterPointWeightsPtr_);
@ -159,11 +155,7 @@ GGIInterpolation<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::uncoveredMasterFaces() const
}
template<class MasterPatch, class SlavePatch>
const labelList&
GGIInterpolation<MasterPatch, SlavePatch>::partiallyCoveredMasterFaces() const
{
if (!partiallyCoveredMasterAddrPtr_)
{
calcAddressing();
}
return *partiallyCoveredMasterAddrPtr_;
}
template<class MasterPatch, class SlavePatch>
const scalarField&
GGIInterpolation<MasterPatch, SlavePatch>::masterFaceCoveredFractions() const
{
if (!masterFaceCoveredFractionsPtr_)
{
calcAddressing();
}
return *masterFaceCoveredFractionsPtr_;
}
template<class MasterPatch, class SlavePatch>
const labelList&
GGIInterpolation<MasterPatch, SlavePatch>::uncoveredSlaveFaces() const
@ -303,32 +269,6 @@ GGIInterpolation<MasterPatch, SlavePatch>::uncoveredSlaveFaces() const
}
template<class MasterPatch, class SlavePatch>
const labelList&
GGIInterpolation<MasterPatch, SlavePatch>::partiallyCoveredSlaveFaces() const
{
if (!partiallyCoveredSlaveAddrPtr_)
{
calcAddressing();
}
return *partiallyCoveredSlaveAddrPtr_;
}
template<class MasterPatch, class SlavePatch>
const scalarField&
GGIInterpolation<MasterPatch, SlavePatch>::slaveFaceCoveredFractions() const
{
if (!slaveFaceCoveredFractionsPtr_)
{
calcAddressing();
}
return *slaveFaceCoveredFractionsPtr_;
}
template<class MasterPatch, class SlavePatch>
bool GGIInterpolation<MasterPatch, SlavePatch>::movePoints
(

58
src/foam/interpolations/GGIInterpolation/GGIInterpolationTemplate.H
View file

@ -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<Type>& bridgeField,
Field<Type>& 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<Type>& bridgeField,
Field<Type>& 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<labelPair>& 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<Type>& pf
) const;
template<class Type>
tmp<Field<Type> > slaveToMasterPointInterpolate
(
const Field<Type>& pf
) const;
// Edit

190
src/foam/interpolations/GGIInterpolation/GGIInterpolationWeights.C
View file

@ -68,16 +68,6 @@ GGIInterpolation<MasterPatch, SlavePatch>::featureCosTol_
);
template<class MasterPatch, class SlavePatch>
const Foam::debug::tolerancesSwitch
GGIInterpolation<MasterPatch, SlavePatch>::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<class MasterPatch, class SlavePatch>
@ -90,11 +80,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::calcAddressing() const
|| slaveAddrPtr_
|| slaveWeightsPtr_
|| uncoveredMasterAddrPtr_
|| partiallyCoveredMasterAddrPtr_
|| masterFaceCoveredFractionsPtr_
|| uncoveredSlaveAddrPtr_
|| partiallyCoveredSlaveAddrPtr_
|| slaveFaceCoveredFractionsPtr_
)
{
FatalErrorIn
@ -445,7 +431,7 @@ void GGIInterpolation<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::rescaleWeightingFactors() const
scalar sumMWC = 0;
scalar curMWC = 0;
// Note: do not rescale weighting factors for partially covered faces
if (!partiallyCoveredMasterAddrPtr_ || !partiallyCoveredSlaveAddrPtr_)
{
FatalErrorIn
(
"void GGIInterpolation<MasterPatch, SlavePatch>::"
"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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::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<MasterPatch, SlavePatch>::findNonOverlappingFaces
return tpatchFaceNonOverlapAddr;
}
template<class MasterPatch, class SlavePatch>
void GGIInterpolation<MasterPatch, SlavePatch>::calcPartiallyCoveredFaces
(
const scalarListList& patchWeights,
const scalar& nonOverlapFaceTol,
const bool isMaster
) const
{
// Sanity checks first
if
(
isMaster
&& (partiallyCoveredMasterAddrPtr_ || masterFaceCoveredFractionsPtr_)
)
{
FatalErrorIn
(
"void GGIInterpolation<MasterPatch, SlavePatch>::"
"calcPartiallyCoveredFaces() const"
) << "Already calculated master partially covered faces"
<< abort(FatalError);
}
if
(
!isMaster
&& (partiallyCoveredSlaveAddrPtr_ || slaveFaceCoveredFractionsPtr_)
)
{
FatalErrorIn
(
"void GGIInterpolation<MasterPatch, SlavePatch>::"
"calcPartiallyCoveredFaces() const"
) << "Already calculated slave partially covered faces"
<< abort(FatalError);
}
// Temporary storage
DynamicList<label, 64> patchFacePartialOverlap(patchWeights.size());
DynamicList<scalar, 64> 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<class FromPatch, class ToPatch>
void GGIInterpolation<FromPatch, ToPatch>::
calcMasterPointAddressing() const