From f30e0ab80887507c903c72399eaa724a9bad7ad0 Mon Sep 17 00:00:00 2001 From: Will Bainbridge <http://cfd.direct> Date: Mon, 6 Nov 2017 14:34:49 +0000 Subject: [PATCH] ENH: reconstructParMesh: Match face point averages on coupled patches In the event that matching centroids across a coupled patch pair fails, we fall back to matching the face point average. The latter can be obtained more reliably on degenerate faces as the calculation does not involve division by the face area. This fallback was already implemented as part of processorPolyPatch. This change also applies it to the faceCoupleInfo class used by reconstructParMesh. --- .../polyMeshAdder/faceCoupleInfo.C | 41 +++++++++++++++++-- .../polyMeshAdder/faceCoupleInfo.H | 12 +++++- .../polyMeshAdder/faceCoupleInfoTemplates.C | 27 +++++++++++- 3 files changed, 75 insertions(+), 5 deletions(-) diff --git a/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.C b/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.C index 3c3edde47d7..79efcbc1be1 100644 --- a/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.C +++ b/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.C @@ -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 | ------------------------------------------------------------------------------- License @@ -1484,12 +1484,13 @@ void Foam::faceCoupleInfo::perfectPointMatch // Faces do not have to be ordered (but all have // to match). Note: Faces will be already ordered if we enter here from // construct from meshes. + matchedAllFaces = matchPoints ( calcFaceCentres<List> ( cutFaces(), - cutPoints_, + cutFaces().points(), 0, cutFaces().size() ), @@ -1501,9 +1502,43 @@ void Foam::faceCoupleInfo::perfectPointMatch slavePatch().size() ), scalarField(slavePatch().size(), absTol), - true, + false, cutToSlaveFaces_ ); + + // If some of the face centres did not match, then try to match the + // point averages instead. There is no division by the face area in + // calculating the point average, so this is more stable when faces + // collapse onto a line or point. + if (!matchedAllFaces) + { + labelList cutToSlaveFacesTemp(cutToSlaveFaces_.size(), -1); + + matchPoints + ( + calcFacePointAverages<List> + ( + cutFaces(), + cutFaces().points(), + 0, + cutFaces().size() + ), + calcFacePointAverages<IndirectList> + ( + slavePatch(), + slavePatch().points(), + 0, + slavePatch().size() + ), + scalarField(slavePatch().size(), absTol), + true, + cutToSlaveFacesTemp + ); + + cutToSlaveFaces_ = max(cutToSlaveFaces_, cutToSlaveFacesTemp); + + matchedAllFaces = min(cutToSlaveFaces_) != -1; + } } diff --git a/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.H b/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.H index e406bb852f8..45810798371 100644 --- a/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.H +++ b/src/dynamicMesh/polyMeshAdder/faceCoupleInfo.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 | ------------------------------------------------------------------------------- License @@ -218,6 +218,16 @@ class faceCoupleInfo const label size ); + //- Calculate face point averages from (subset of) faces. + template<template<class> class FaceList> + static pointField calcFacePointAverages + ( + const FaceList<face>&, + const pointField&, + const label start, + const label size + ); + //- Write edges static void writeOBJ ( diff --git a/src/dynamicMesh/polyMeshAdder/faceCoupleInfoTemplates.C b/src/dynamicMesh/polyMeshAdder/faceCoupleInfoTemplates.C index 3b25a6da516..2158c0c92a5 100644 --- a/src/dynamicMesh/polyMeshAdder/faceCoupleInfoTemplates.C +++ b/src/dynamicMesh/polyMeshAdder/faceCoupleInfoTemplates.C @@ -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 | ------------------------------------------------------------------------------- License @@ -49,4 +49,29 @@ Foam::pointField Foam::faceCoupleInfo::calcFaceCentres } +template<template<class> class FaceList> +Foam::pointField Foam::faceCoupleInfo::calcFacePointAverages +( + const FaceList<face>& faces, + const pointField& points, + const label start, + const label size +) +{ + pointField fpa(size, Zero); + + label facei = start; + + forAll(fpa, i) + { + forAll(faces[facei], j) + { + fpa[i] += points[faces[facei][j]]; + } + fpa[i] /= faces[facei++].size(); + } + return fpa; +} + + // ************************************************************************* // -- GitLab