From d76923d8511cedb018acb2c1544b056a474486b7 Mon Sep 17 00:00:00 2001 From: Will Bainbridge <http://cfd.direct> Date: Tue, 31 Oct 2017 08:59:21 +0000 Subject: [PATCH] ENH: AMI: Consistency between overlap/normalisation areas The patch magSf calculation has been changed so that it uses the same triangulation as the overlap algorithm. This improves consistency and means that for exactly conforming patches (typically before any mesh motion) the weights do not require normalisation. --- .../AMIInterpolation/AMIInterpolation.C | 54 +++++++++++++---- .../AMIInterpolation/AMIInterpolation.H | 10 +++- .../faceAreaIntersect/faceAreaIntersect.C | 60 ++++++++++++------- .../faceAreaIntersect/faceAreaIntersect.H | 18 +++--- .../faceAreaIntersect/faceAreaIntersectI.H | 25 +------- 5 files changed, 104 insertions(+), 63 deletions(-) diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C index ae79db7a641..48bc0230f79 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C @@ -119,6 +119,48 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::wordTointerpolationMethod } +template<class SourcePatch, class TargetPatch> +template<class Patch> +Foam::tmp<Foam::scalarField> +Foam::AMIInterpolation<SourcePatch, TargetPatch>::patchMagSf +( + const Patch& patch, + const faceAreaIntersect::triangulationMode triMode +) +{ + tmp<scalarField> tResult(new scalarField(patch.size(), Zero)); + scalarField& result = tResult.ref(); + + const pointField& patchPoints = patch.localPoints(); + + faceList patchFaceTris; + + forAll(result, patchFacei) + { + faceAreaIntersect::triangulate + ( + patch.localFaces()[patchFacei], + patchPoints, + triMode, + patchFaceTris + ); + + forAll(patchFaceTris, i) + { + result[patchFacei] += + triPointRef + ( + patchPoints[patchFaceTris[i][0]], + patchPoints[patchFaceTris[i][1]], + patchPoints[patchFaceTris[i][2]] + ).mag(); + } + } + + return tResult; +} + + // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // template<class SourcePatch, class TargetPatch> @@ -884,16 +926,8 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update << endl; // Calculate face areas - srcMagSf_.setSize(srcPatch.size()); - forAll(srcMagSf_, facei) - { - srcMagSf_[facei] = srcPatch[facei].mag(srcPatch.points()); - } - tgtMagSf_.setSize(tgtPatch.size()); - forAll(tgtMagSf_, facei) - { - tgtMagSf_[facei] = tgtPatch[facei].mag(tgtPatch.points()); - } + srcMagSf_ = patchMagSf(srcPatch, triMode_); + tgtMagSf_ = patchMagSf(tgtPatch, triMode_); // Calculate if patches present on multiple processors singlePatchProc_ = calcDistribution(srcPatch, tgtPatch); diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H index 7723262779a..d9a6fc9defc 100644 --- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H +++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd. ------------------------------------------------------------------------------- License @@ -104,6 +104,14 @@ public: const word& method ); + //- Calculate the patch face magnitudes for the given tri-mode + template<class Patch> + static tmp<scalarField> patchMagSf + ( + const Patch& patch, + const faceAreaIntersect::triangulationMode triMode + ); + private: diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C index 4d215af6de3..d52d19cedc8 100644 --- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C +++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -323,41 +323,61 @@ Foam::faceAreaIntersect::faceAreaIntersect // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -Foam::scalar Foam::faceAreaIntersect::calc +void Foam::faceAreaIntersect::triangulate ( - const face& faceA, - const face& faceB, - const vector& n, - const triangulationMode& triMode + const face& f, + const pointField& points, + const triangulationMode& triMode, + faceList& faceTris ) { - // split faces into triangles - DynamicList<face> trisA; - DynamicList<face> trisB; + faceTris.resize(f.nTriangles()); switch (triMode) { case tmFan: { - triangleFan(faceA, trisA); - triangleFan(faceB, trisB); + for (label i = 0; i < f.nTriangles(); ++ i) + { + faceTris[i] = face(3); + faceTris[i][0] = f[0]; + faceTris[i][1] = f[i + 1]; + faceTris[i][2] = f[i + 2]; + } break; } case tmMesh: { - faceA.triangles(pointsA_, trisA); - faceB.triangles(pointsB_, trisB); + const label nFaceTris = f.nTriangles(); - break; - } - default: - { - FatalErrorInFunction - << "Unknown triangulation mode enumeration" - << abort(FatalError); + label nFaceTris1 = 0; + const label nFaceTris2 = f.triangles(points, nFaceTris1, faceTris); + + if (nFaceTris != nFaceTris1 || nFaceTris != nFaceTris2) + { + FatalErrorInFunction + << "The numbers of reported triangles in the face do not " + << "match that generated by the triangulation" + << exit(FatalError); + } } } +} + + +Foam::scalar Foam::faceAreaIntersect::calc +( + const face& faceA, + const face& faceB, + const vector& n, + const triangulationMode& triMode +) +{ + // split faces into triangles + faceList trisA, trisB; + triangulate(faceA, pointsA_, triMode, trisA); + triangulate(faceB, pointsB_, triMode, trisB); // intersect triangles scalar totalArea = 0.0; diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H index 450d343d644..b2fe05b23f6 100644 --- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H +++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -105,13 +105,6 @@ private: FixedList<triPoints, 10>& tris ) const; - //- Decompose face into triangle fan - inline void triangleFan - ( - const face& f, - DynamicList<face>& faces - ) const; - //- Return point of intersection between plane and triangle edge inline point planeIntersection ( @@ -162,6 +155,15 @@ public: //- Fraction of local length scale to use as intersection tolerance inline static scalar& tolerance(); + //- Triangulate a face using the given triangulation mode + static void triangulate + ( + const face& f, + const pointField& points, + const triangulationMode& triMode, + faceList& faceTris + ); + //- Return area of intersection of faceA with faceB scalar calc ( diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H index 76bf413833e..329c75d821f 100644 --- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H +++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | Copyright (C) 2015 OpenCFD Ltd. ------------------------------------------------------------------------------- License @@ -67,29 +67,6 @@ inline Foam::triPoints Foam::faceAreaIntersect::getTriPoints } -inline void Foam::faceAreaIntersect::triangleFan -( - const face& f, - DynamicList<face>& faces -) const -{ - if (f.size() > 2) - { - const label v0 = 0; - - labelList indices(3); - - for (label i = 1; i < f.size() - 1; i++) - { - indices[0] = f[v0]; - indices[1] = f[i]; - indices[2] = f[i + 1]; - faces.append(face(indices)); - } - } -} - - inline Foam::point Foam::faceAreaIntersect::planeIntersection ( const FixedList<scalar, 3>& d, -- GitLab