From b314fb887e71c1ed6eb7e4526edb69352148503a Mon Sep 17 00:00:00 2001 From: andy Date: Wed, 20 Nov 2013 16:22:48 +0000 Subject: [PATCH] ENH: ACMI - code updates to avoid FPEs --- .../cyclicACMIPolyPatch/cyclicACMIPolyPatch.C | 52 ++++++++++++++----- .../cyclicACMIPolyPatch/cyclicACMIPolyPatch.H | 15 ++++++ .../cyclicACMIPolyPatchI.H | 6 +-- .../cyclicACMIPolyPatchTemplates.C | 8 +-- .../cyclicAMIPolyPatch/cyclicAMIPolyPatch.C | 7 +++ .../cyclicAMIPolyPatch/cyclicAMIPolyPatch.H | 2 +- .../cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H | 8 --- 7 files changed, 68 insertions(+), 30 deletions(-) diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C index 30ae0adf37..4c36dc0be4 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C @@ -90,16 +90,16 @@ void Foam::cyclicACMIPolyPatch::resetAMI AMIPatchToPatchInterpolation::imPartialFaceAreaWeight ); - const scalarField& srcWeightSum = AMI().srcWeightsSum(); + srcMask_ = + min(1.0 - tolerance_, max(tolerance_, AMI().srcWeightsSum())); + + tgtMask_ = + min(1.0 - tolerance_, max(tolerance_, AMI().tgtWeightsSum())); - // set patch face areas based on sum of AMI weights per face forAll(Sf, faceI) { - scalar w = srcWeightSum[faceI]; - w = min(1.0 - tolerance_, max(tolerance_, w)); - - Sf[faceI] *= w; - noSf[faceI] *= 1.0 - w; + Sf[faceI] *= srcMask_[faceI]; + noSf[faceI] *= 1.0 - srcMask_[faceI]; } setNeighbourFaceAreas(); @@ -116,8 +116,6 @@ void Foam::cyclicACMIPolyPatch::setNeighbourFaceAreas() const refCast(this->neighbPatch()); const polyPatch& pp = cp.nonOverlapPatch(); - const scalarField& tgtWeightSum = AMI().tgtWeightsSum(); - const vectorField& faceAreas0 = cp.faceAreas0(); vectorField::subField Sf = cp.faceAreas(); @@ -125,11 +123,8 @@ void Foam::cyclicACMIPolyPatch::setNeighbourFaceAreas() const forAll(Sf, faceI) { - scalar w = tgtWeightSum[faceI]; - w = min(1.0 - tolerance_, max(tolerance_, w)); - - Sf[faceI] = w*faceAreas0[faceI]; - noSf[faceI] = (1.0 - w)*faceAreas0[faceI]; + Sf[faceI] = tgtMask_[faceI]*faceAreas0[faceI]; + noSf[faceI] = (1.0 - tgtMask_[faceI])*faceAreas0[faceI]; } } @@ -184,6 +179,18 @@ void Foam::cyclicACMIPolyPatch::clearGeom() } +const Foam::scalarField& Foam::cyclicACMIPolyPatch::srcMask() const +{ + return srcMask_; +} + + +const Foam::scalarField& Foam::cyclicACMIPolyPatch::tgtMask() const +{ + return tgtMask_; +} + + // * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * // Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch @@ -201,6 +208,8 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch faceAreas0_(), nonOverlapPatchName_(word::null), nonOverlapPatchID_(-1), + srcMask_(), + tgtMask_(), updated_(false) { // Non-overlapping patch might not be valid yet so cannot determine @@ -221,6 +230,8 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch faceAreas0_(), nonOverlapPatchName_(dict.lookup("nonOverlapPatch")), nonOverlapPatchID_(-1), + srcMask_(), + tgtMask_(), updated_(false) { if (nonOverlapPatchName_ == name) @@ -256,6 +267,8 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch faceAreas0_(), nonOverlapPatchName_(pp.nonOverlapPatchName_), nonOverlapPatchID_(-1), + srcMask_(), + tgtMask_(), updated_(false) { // Non-overlapping patch might not be valid yet so cannot determine @@ -278,6 +291,8 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch faceAreas0_(), nonOverlapPatchName_(nonOverlapPatchName), nonOverlapPatchID_(-1), + srcMask_(), + tgtMask_(), updated_(false) { if (nonOverlapPatchName_ == name()) @@ -314,6 +329,8 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch faceAreas0_(), nonOverlapPatchName_(pp.nonOverlapPatchName_), nonOverlapPatchID_(-1), + srcMask_(), + tgtMask_(), updated_(false) {} @@ -326,6 +343,13 @@ Foam::cyclicACMIPolyPatch::~cyclicACMIPolyPatch() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +const Foam::cyclicACMIPolyPatch& Foam::cyclicACMIPolyPatch::neighbPatch() const +{ + const polyPatch& pp = this->boundaryMesh()[neighbPatchID()]; + return refCast(pp); +} + + Foam::label Foam::cyclicACMIPolyPatch::nonOverlapPatchID() const { if (nonOverlapPatchID_ == -1) diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H index b092c14f04..b153034d87 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H @@ -66,6 +66,12 @@ private: //- Index of non-overlapping patch mutable label nonOverlapPatchID_; + //- Mask/weighting for source patch + mutable scalarField srcMask_; + + //- Mask/weighting for target patch + mutable scalarField tgtMask_; + //- Flag to indicate that AMI has been updated mutable bool updated_; @@ -111,6 +117,12 @@ protected: //- Clear geometry virtual void clearGeom(); + //- Return the mask/weighting for the source patch + virtual const scalarField& srcMask() const; + + //- Return the mask/weighting for the target patch + virtual const scalarField& tgtMask() const; + public: @@ -245,6 +257,9 @@ public: //- Return access to the original patch face areas inline const vectorField& faceAreas0() const; + //- Return a reference to the neighbour patch + virtual const cyclicACMIPolyPatch& neighbPatch() const; + //- Non-overlapping patch name inline const word& nonOverlapPatchName() const; diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H index a245f779c0..fbb049a725 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H @@ -67,15 +67,15 @@ inline Foam::polyPatch& Foam::cyclicACMIPolyPatch::nonOverlapPatch() } -inline const Foam::scalarField& Foam::cyclicACMIPolyPatch::mask() const +inline const scalarField& Foam::cyclicACMIPolyPatch::mask() const { if (owner()) { - return AMI().srcWeightsSum(); + return srcMask_; } else { - return neighbPatch().AMI().tgtWeightsSum(); + return neighbPatch().tgtMask(); } } diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchTemplates.C b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchTemplates.C index 127085ac44..b64842bf12 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchTemplates.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchTemplates.C @@ -34,7 +34,7 @@ Foam::tmp > Foam::cyclicACMIPolyPatch::interpolate { if (owner()) { - const scalarField& w = AMI().srcWeightsSum(); + const scalarField& w = srcMask_; return w*AMI().interpolateToSource(fldCouple) @@ -42,7 +42,7 @@ Foam::tmp > Foam::cyclicACMIPolyPatch::interpolate } else { - const scalarField& w = neighbPatch().AMI().tgtWeightsSum(); + const scalarField& w = neighbPatch().tgtMask(); return w*neighbPatch().AMI().interpolateToTarget(fldCouple) @@ -73,14 +73,14 @@ void Foam::cyclicACMIPolyPatch::interpolate { if (owner()) { - const scalarField& w = AMI().srcWeightsSum(); + const scalarField& w = srcMask_; AMI().interpolateToSource(fldCouple, cop, result); result = w*result + (1.0 - w)*fldNonOverlap; } else { - const scalarField& w = neighbPatch().AMI().tgtWeightsSum(); + const scalarField& w = neighbPatch().tgtMask(); neighbPatch().AMI().interpolateToTarget(fldCouple, cop, result); result = w*result + (1.0 - w)*fldNonOverlap; diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C index e906689415..52373c438d 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C @@ -642,6 +642,13 @@ bool Foam::cyclicAMIPolyPatch::owner() const } +const Foam::cyclicAMIPolyPatch& Foam::cyclicAMIPolyPatch::neighbPatch() const +{ + const polyPatch& pp = this->boundaryMesh()[neighbPatchID()]; + return refCast(pp); +} + + const Foam::autoPtr& Foam::cyclicAMIPolyPatch::surfPtr() const { diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H index 9c0f5fc316..e3bb1d35d9 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H @@ -284,7 +284,7 @@ public: virtual bool owner() const; //- Return a reference to the neighbour patch - inline const cyclicAMIPolyPatch& neighbPatch() const; + virtual const cyclicAMIPolyPatch& neighbPatch() const; //- Return a reference to the projection surface const autoPtr& surfPtr() const; diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H index 217d189ae3..b136f90ad0 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H +++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H @@ -38,14 +38,6 @@ inline const Foam::word& Foam::cyclicAMIPolyPatch::neighbPatchName() const } -inline const Foam::cyclicAMIPolyPatch& -Foam::cyclicAMIPolyPatch::neighbPatch() const -{ - const polyPatch& pp = this->boundaryMesh()[neighbPatchID()]; - return refCast(pp); -} - - inline const Foam::vector& Foam::cyclicAMIPolyPatch::rotationAxis() const { return rotationAxis_;