diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C b/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C
index 80232d4498e081d7e3098538278b994bf3426bc5..e5818171d53ac6f85da1ceb6399a5e40ab79df99 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C
@@ -45,6 +45,7 @@ void Foam::polyMesh::updateMesh(const mapPolyMesh& mpm)
<< "Updating addressing and (optional) pointMesh/pointFields" << endl;
// Update boundaryMesh (note that patches themselves already ok)
+// boundary_.updateMesh(mpm);
boundary_.updateMesh();
// Update zones
diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C
index 3e915b3e7e6938eb96cbbd47ae49c16b44365b25..6db7dc22fcb8cb648992b24ff9cbfe31edaeb668 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C
@@ -63,6 +63,7 @@ void Foam::polyPatch::movePoints(PstreamBuffers&, const pointField& p)
primitivePatch::movePoints(p);
}
+
void Foam::polyPatch::updateMesh(PstreamBuffers&)
{
primitivePatch::clearGeom();
diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H
index 08bbf0c732aa34cc498a86359be0b314378086c6..f7f01abb091986e668a1cf66b95ff9fa121ad05b 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H
@@ -55,6 +55,7 @@ namespace Foam
// Forward declarations
class polyBoundaryMesh;
class polyPatch;
+class polyTopoChange;
class PstreamBuffers;
Ostream& operator<<(Ostream&, const polyPatch&);
@@ -419,6 +420,19 @@ public:
labelList& rotation
) const;
+ //- For dynamic mesh cases - return true if this patch will change the
+ //- topology
+ virtual bool changeTopology() const
+ {
+ return false;
+ }
+
+ //- Collect topology changes in a polyTopoChange object
+ virtual bool setTopology(polyTopoChange&)
+ {
+ return false;
+ }
+
// Member operators
diff --git a/src/dynamicFvMesh/Make/files b/src/dynamicFvMesh/Make/files
index b8432c07c905108747f9b94eb4cffefd26d1e0bd..794af6f3d0dd9754f6371c669da5e3c10da2b07d 100644
--- a/src/dynamicFvMesh/Make/files
+++ b/src/dynamicFvMesh/Make/files
@@ -10,4 +10,9 @@ dynamicMotionSolverListFvMesh/dynamicMotionSolverListFvMesh.C
simplifiedDynamicFvMesh/simplifiedDynamicFvMeshes.C
simplifiedDynamicFvMesh/simplifiedDynamicFvMesh.C
+
+
+dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.C
+
+
LIB = $(FOAM_LIBBIN)/libdynamicFvMesh
diff --git a/src/dynamicFvMesh/dynamicMotionSolverFvMesh/dynamicMotionSolverFvMesh.C b/src/dynamicFvMesh/dynamicMotionSolverFvMesh/dynamicMotionSolverFvMesh.C
index ba1d636969a1aeb010e0e0712c4a0a6ffc5cad1d..5b7a90cddc85ae66857042fc5c54bb5347372998 100644
--- a/src/dynamicFvMesh/dynamicMotionSolverFvMesh/dynamicMotionSolverFvMesh.C
+++ b/src/dynamicFvMesh/dynamicMotionSolverFvMesh/dynamicMotionSolverFvMesh.C
@@ -92,6 +92,9 @@ const Foam::motionSolver& Foam::dynamicMotionSolverFvMesh::motion() const
bool Foam::dynamicMotionSolverFvMesh::update()
{
+ // Scan through AMI patches and update
+
+
fvMesh::movePoints(motionPtr_->newPoints());
volVectorField* Uptr = getObjectPtr<volVectorField>("U");
diff --git a/src/dynamicFvMesh/dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.C b/src/dynamicFvMesh/dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.C
new file mode 100644
index 0000000000000000000000000000000000000000..5de3713a53f959125d81181f24c78d6a7596c78d
--- /dev/null
+++ b/src/dynamicFvMesh/dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.C
@@ -0,0 +1,216 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2019-2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "dynamicMotionSolverFvMeshAMI.H"
+#include "addToRunTimeSelectionTable.H"
+#include "motionSolver.H"
+#include "volFields.H"
+#include "surfaceFields.H"
+#include "cyclicAMIPolyPatch.H"
+#include "polyTopoChange.H"
+#include "MeshObject.H"
+#include "lduMesh.H"
+
+#include "processorFvPatch.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(dynamicMotionSolverFvMeshAMI, 0);
+ addToRunTimeSelectionTable
+ (
+ dynamicFvMesh,
+ dynamicMotionSolverFvMeshAMI,
+ IOobject
+ );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::dynamicMotionSolverFvMeshAMI::dynamicMotionSolverFvMeshAMI
+(
+ const IOobject& io
+)
+:
+ dynamicFvMesh(io),
+ motionPtr_(motionSolver::New(*this))
+{}
+
+
+Foam::dynamicMotionSolverFvMeshAMI::dynamicMotionSolverFvMeshAMI
+(
+ const IOobject& io,
+ pointField&& points,
+ faceList&& faces,
+ labelList&& allOwner,
+ labelList&& allNeighbour,
+ const bool syncPar
+)
+:
+ dynamicFvMesh
+ (
+ io,
+ std::move(points),
+ std::move(faces),
+ std::move(allOwner),
+ std::move(allNeighbour),
+ syncPar
+ ),
+ motionPtr_(motionSolver::New(*this))
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+const Foam::motionSolver& Foam::dynamicMotionSolverFvMeshAMI::motion() const
+{
+ return *motionPtr_;
+}
+
+
+bool Foam::dynamicMotionSolverFvMeshAMI::update()
+{
+ // Mesh not moved/changed yet
+ moving(false);
+ topoChanging(false);
+
+ if (debug)
+ {
+ for (const fvPatch& fvp : boundary())
+ {
+ if (!isA<processorFvPatch>(fvp))
+ {
+ Info<< "1 --- patch:" << fvp.patch().name()
+ << " area:" << gSum(fvp.magSf()) << endl;
+ }
+ }
+ }
+
+ pointField newPoints(motionPtr_->curPoints());
+
+ polyBoundaryMesh& pbm = const_cast<polyBoundaryMesh&>(boundaryMesh());
+
+ // Scan all patches and see if we want to apply a mesh topology update
+ bool changeRequired = false;
+ for (polyPatch& pp : pbm)
+ {
+ DebugInfo
+ << "pre-topology change: patch " << pp.name()
+ << " size:" << returnReduce(pp.size(), sumOp<label>())
+ << " mag(faceAreas):" << gSum(mag(pp.faceAreas())) << endl;
+
+ //changeRequired = pp.changeTopology(newPoints) || changeRequired;
+ changeRequired = pp.changeTopology() || changeRequired;
+ }
+
+ reduce(changeRequired, orOp<bool>());
+
+ if (changeRequired)
+ {
+ polyTopoChange polyTopo(*this);
+
+ // Set new point positions in polyTopo object
+ polyTopo.movePoints(newPoints);
+
+ // Accumulate the patch-based mesh changes on the current mesh
+ // Note:
+ // - updates the AMIs using the new points
+ // - creates a topo change object that removes old added faces and
+ // adds the new faces
+ for (polyPatch& pp : pbm)
+ {
+ pp.setTopology(polyTopo);
+ }
+
+ // Update geometry
+ // Note
+ // - changeMesh leads to polyMesh::resetPrimitives which will also
+ // trigger polyBoundaryMesh::updateMesh (init and update) and
+ // ::calcGeometry (with topoChanging = false)
+ // - BUT: mesh still corresponds to original (non-extended mesh) so
+ // we want to bypass these calls...
+ // - after changes topoChanging = true
+ autoPtr<mapPolyMesh> map =
+ polyTopo.changeMesh
+ (
+ *this,
+ true // We will be calling movePoints after this update
+ );
+
+ // Apply topology change - update fv geometry and map fields
+ // - polyMesh::updateMesh
+ // - fires initUpdateMesh and updateMesh in AMI BCs - called before
+ // mapFields
+ // - AMI addressing must be up-to-date - used by, e.g. FaceCellWave
+ // - will trigger (again) polyBoundaryMesh::updateMesh (init and update)
+ updateMesh(map());
+
+ // Move points and update derived properties
+ // Note:
+ // - resets face areas based on raw point locations!
+ // - polyBoundaryMesh::updateMesh (init and update)
+ // Note:
+ // - processorPolyPatches will trigger calculation of faceCentres
+ // (and therefore cell volumes), so need to update faceAreas in
+ // initMovePoints since proc patches will be evaluated later than
+ // AMI patches
+ if (map().hasMotionPoints())
+ {
+ movePoints(map().preMotionPoints());
+ }
+ }
+ else
+ {
+ fvMesh::movePoints(newPoints);
+ }
+
+ volVectorField* Uptr = getObjectPtr<volVectorField>("U");
+
+ if (Uptr)
+ {
+ Uptr->correctBoundaryConditions();
+ }
+
+ if (debug)
+ {
+ for (const fvPatch& fvp : boundary())
+ {
+ if (!isA<processorFvPatch>(fvp))
+ {
+ Info<< "2 --- patch:" << fvp.patch().name()
+ << " area:" << gSum(fvp.magSf()) << endl;
+ }
+ }
+ }
+
+ return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/dynamicFvMesh/dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.H b/src/dynamicFvMesh/dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.H
new file mode 100644
index 0000000000000000000000000000000000000000..e21f5623ed7df6416380459ef8b6470b6ebb42f0
--- /dev/null
+++ b/src/dynamicFvMesh/dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.H
@@ -0,0 +1,121 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2011-2017 OpenFOAM Foundation
+ Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::dynamicMotionSolverFvMeshAMI
+
+Description
+ The dynamicMotionSolverFvMeshAMI
+
+SourceFiles
+ dynamicMotionSolverFvMeshAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef dynamicMotionSolverFvMeshAMI_H
+#define dynamicMotionSolverFvMeshAMI_H
+
+#include "dynamicFvMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class motionSolver;
+
+/*---------------------------------------------------------------------------*\
+ Class dynamicMotionSolverFvMeshAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+class dynamicMotionSolverFvMeshAMI
+:
+ public dynamicFvMesh
+{
+ // Private data
+
+ autoPtr<motionSolver> motionPtr_;
+
+
+ // Private Member Functions
+
+ //- No copy construct
+ dynamicMotionSolverFvMeshAMI
+ (
+ const dynamicMotionSolverFvMeshAMI&
+ ) = delete;
+
+ //- No copy assignment
+ void operator=(const dynamicMotionSolverFvMeshAMI&) = delete;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("dynamicMotionSolverFvMeshAMI");
+
+
+ // Constructors
+
+ //- Construct from IOobject
+ dynamicMotionSolverFvMeshAMI(const IOobject& io);
+
+ //- Construct from components without boundary.
+ // Boundary is added using addFvPatches() member function
+ dynamicMotionSolverFvMeshAMI
+ (
+ const IOobject& io,
+ pointField&& points,
+ faceList&& faces,
+ labelList&& allOwner,
+ labelList&& allNeighbour,
+ const bool syncPar = true
+ );
+
+
+ //- Destructor
+ virtual ~dynamicMotionSolverFvMeshAMI() = default;
+
+
+ // Member Functions
+
+ //- Return the motionSolver
+ const motionSolver& motion() const;
+
+ //- Update the mesh for both mesh motion and topology change
+ virtual bool update();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/Make/options b/src/finiteVolume/Make/options
index 851f94b6b3910365ebdd5ad09a66903da978f118..ea11ab0d45a46d24a226464a14238a684ff6278a 100644
--- a/src/finiteVolume/Make/options
+++ b/src/finiteVolume/Make/options
@@ -1,7 +1,8 @@
EXE_INC = \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude \
- -I$(LIB_SRC)/meshTools/lnInclude
+ -I$(LIB_SRC)/meshTools/lnInclude \
+ -I$(LIB_SRC)/dynamicMesh/lnInclude
LIB_LIBS = \
-lOpenFOAM \
diff --git a/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.H b/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.H
index eae2dd49dad5656f8de4b2d6690f9dc8dad0f3b1..db093b3940174646b4663dd5a93f52491f220d18 100644
--- a/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.H
+++ b/src/finiteVolume/fields/fvPatchFields/constraint/cyclicACMI/cyclicACMIFvPatchField.H
@@ -173,7 +173,7 @@ public:
//- Return true if this patch field fixes a value
// Needed to check if a level has to be specified while solving
- // Poissons equations
+ // Poisson equations
virtual bool fixesValue() const
{
const scalarField& mask =
diff --git a/src/finiteVolume/fvMesh/fvMesh.C b/src/finiteVolume/fvMesh/fvMesh.C
index 9aff27408a42320329b01ab30f83d41d4349f1d5..e02ae63dccb3177dc8a6766dc2b9f2187453afe9 100644
--- a/src/finiteVolume/fvMesh/fvMesh.C
+++ b/src/finiteVolume/fvMesh/fvMesh.C
@@ -726,6 +726,8 @@ void Foam::fvMesh::mapFields(const mapPolyMesh& meshMap)
Foam::tmp<Foam::scalarField> Foam::fvMesh::movePoints(const pointField& p)
{
+ DebugInFunction << endl;
+
// Grab old time volumes if the time has been incremented
// This will update V0, V00
if (curTimeIndex_ < time().timeIndex())
@@ -733,6 +735,14 @@ Foam::tmp<Foam::scalarField> Foam::fvMesh::movePoints(const pointField& p)
storeOldVol(V());
}
+
+ // Move the polyMesh and set the mesh motion fluxes to the swept-volumes
+
+ scalar rDeltaT = 1.0/time().deltaTValue();
+
+ tmp<scalarField> tsweptVols = polyMesh::movePoints(p);
+ scalarField& sweptVols = tsweptVols.ref();
+
if (!phiPtr_)
{
// Create mesh motion flux
@@ -761,14 +771,6 @@ Foam::tmp<Foam::scalarField> Foam::fvMesh::movePoints(const pointField& p)
}
surfaceScalarField& phi = *phiPtr_;
-
- // Move the polyMesh and set the mesh motion fluxes to the swept-volumes
-
- scalar rDeltaT = 1.0/time().deltaTValue();
-
- tmp<scalarField> tsweptVols = polyMesh::movePoints(p);
- scalarField& sweptVols = tsweptVols.ref();
-
phi.primitiveFieldRef() =
scalarField::subField(sweptVols, nInternalFaces());
phi.primitiveFieldRef() *= rDeltaT;
@@ -776,7 +778,6 @@ Foam::tmp<Foam::scalarField> Foam::fvMesh::movePoints(const pointField& p)
const fvPatchList& patches = boundary();
surfaceScalarField::Boundary& phibf = phi.boundaryFieldRef();
-
forAll(patches, patchi)
{
phibf[patchi] = patches[patchi].patchSlice(sweptVols);
@@ -805,6 +806,8 @@ Foam::tmp<Foam::scalarField> Foam::fvMesh::movePoints(const pointField& p)
void Foam::fvMesh::updateMesh(const mapPolyMesh& mpm)
{
+ DebugInFunction << endl;
+
// Update polyMesh. This needs to keep volume existent!
polyMesh::updateMesh(mpm);
diff --git a/src/finiteVolume/fvMesh/fvMesh.H b/src/finiteVolume/fvMesh/fvMesh.H
index 53e57c7f683590c72fa7f0cd51eda9a89fbb46c9..0e35ef9f9d2a4da19d29e39bfedc28a33c41f9fb 100644
--- a/src/finiteVolume/fvMesh/fvMesh.H
+++ b/src/finiteVolume/fvMesh/fvMesh.H
@@ -91,6 +91,8 @@ class fvMesh
public fvSolution,
public data
{
+protected:
+
// Private data
//- Boundary mesh
diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C
index c4427d54849224973829e37955bd6144034c7c95..23a12796ad1d9a2c682e398edeb62f53fae10fda 100644
--- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C
+++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C
@@ -6,6 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2013-2016 OpenFOAM Foundation
+ Copyright (C) 2019 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -26,9 +27,10 @@ License
\*---------------------------------------------------------------------------*/
#include "cyclicACMIFvPatch.H"
-#include "addToRunTimeSelectionTable.H"
#include "fvMesh.H"
#include "transform.H"
+#include "surfaceFields.H"
+#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -41,45 +43,13 @@ namespace Foam
// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
-void Foam::cyclicACMIFvPatch::updateAreas() const
+void Foam::cyclicACMIFvPatch::resetPatchAreas(const fvPatch& fvp) const
{
- if (cyclicACMIPolyPatch_.updated())
- {
- if (debug)
- {
- Pout<< "cyclicACMIFvPatch::updateAreas() : updating fv areas for "
- << name() << " and " << this->nonOverlapPatch().name()
- << endl;
- }
+ const_cast<vectorField&>(fvp.Sf()) = fvp.patch().faceAreas();
+ const_cast<scalarField&>(fvp.magSf()) = mag(fvp.patch().faceAreas());
- // owner couple
- const_cast<vectorField&>(Sf()) = patch().faceAreas();
- const_cast<scalarField&>(magSf()) = mag(patch().faceAreas());
-
- // owner non-overlapping
- const fvPatch& nonOverlapPatch = this->nonOverlapPatch();
- const_cast<vectorField&>(nonOverlapPatch.Sf()) =
- nonOverlapPatch.patch().faceAreas();
- const_cast<scalarField&>(nonOverlapPatch.magSf()) =
- mag(nonOverlapPatch.patch().faceAreas());
-
- // neighbour couple
- const cyclicACMIFvPatch& nbrACMI = neighbPatch();
- const_cast<vectorField&>(nbrACMI.Sf()) =
- nbrACMI.patch().faceAreas();
- const_cast<scalarField&>(nbrACMI.magSf()) =
- mag(nbrACMI.patch().faceAreas());
-
- // neighbour non-overlapping
- const fvPatch& nbrNonOverlapPatch = nbrACMI.nonOverlapPatch();
- const_cast<vectorField&>(nbrNonOverlapPatch.Sf()) =
- nbrNonOverlapPatch.patch().faceAreas();
- const_cast<scalarField&>(nbrNonOverlapPatch.magSf()) =
- mag(nbrNonOverlapPatch.patch().faceAreas());
-
- // set the updated flag
- cyclicACMIPolyPatch_.setUpdated(false);
- }
+ DebugPout
+ << fvp.patch().name() << " area:" << sum(fvp.magSf()) << endl;
}
@@ -100,13 +70,13 @@ void Foam::cyclicACMIFvPatch::makeWeights(scalarField& w) const
)
);
- scalar tol = cyclicACMIPolyPatch::tolerance();
+ const scalar tol = cyclicACMIPolyPatch::tolerance();
forAll(deltas, facei)
{
- scalar di = deltas[facei];
- scalar dni = nbrDeltas[facei];
+ scalar di = mag(deltas[facei]);
+ scalar dni = mag(nbrDeltas[facei]);
if (dni < tol)
{
@@ -147,7 +117,6 @@ Foam::tmp<Foam::vectorField> Foam::cyclicACMIFvPatch::delta() const
vectorField nbrPatchD(interpolate(nbrPatch.coupledFvPatch::delta()));
-
tmp<vectorField> tpdv(new vectorField(patchD.size()));
vectorField& pdv = tpdv.ref();
@@ -201,4 +170,94 @@ Foam::tmp<Foam::labelField> Foam::cyclicACMIFvPatch::internalFieldTransfer
}
+void Foam::cyclicACMIFvPatch::movePoints()
+{
+ if (!cyclicACMIPolyPatch_.owner())
+ {
+ return;
+ }
+
+ // Set the patch face areas to be consistent with the changes made at the
+ // polyPatch level
+
+ const fvPatch& nonOverlapPatch = this->nonOverlapPatch();
+ const cyclicACMIFvPatch& nbrACMI = neighbPatch();
+ const fvPatch& nbrNonOverlapPatch = nbrACMI.nonOverlapPatch();
+
+ resetPatchAreas(*this);
+ resetPatchAreas(nonOverlapPatch);
+ resetPatchAreas(nbrACMI);
+ resetPatchAreas(nbrNonOverlapPatch);
+
+ // Scale the mesh flux
+
+ const labelListList& newSrcAddr = AMI().srcAddress();
+ const labelListList& newTgtAddr = AMI().tgtAddress();
+
+ const fvMesh& mesh = boundaryMesh().mesh();
+ surfaceScalarField& meshPhi = const_cast<fvMesh&>(mesh).setPhi();
+ surfaceScalarField::Boundary& meshPhiBf = meshPhi.boundaryFieldRef();
+
+ // Note: phip and phiNonOverlapp will be different sizes if new faces
+ // have been added
+ scalarField& phip = meshPhiBf[cyclicACMIPolyPatch_.index()];
+ scalarField& phiNonOverlapp =
+ meshPhiBf[nonOverlapPatch.patch().index()];
+
+ forAll(phip, facei)
+ {
+ if (newSrcAddr[facei].empty())
+ {
+ // AMI patch with no connection to other coupled faces
+ phip[facei] = 0.0;
+ }
+ else
+ {
+ // Scale the mesh flux according to the area fraction
+ const face& fAMI = cyclicACMIPolyPatch_.localFaces()[facei];
+
+ // Note: using raw point locations to calculate the geometric
+ // area - faces areas are currently scaled (decoupled from
+ // mesh points)
+ const scalar geomArea = fAMI.mag(cyclicACMIPolyPatch_.localPoints());
+ phip[facei] *= magSf()[facei]/geomArea;
+ }
+ }
+
+ forAll(phiNonOverlapp, facei)
+ {
+ const scalar w = 1.0 - cyclicACMIPolyPatch_.srcMask()[facei];
+ phiNonOverlapp[facei] *= w;
+ }
+
+ scalarField& nbrPhip = meshPhiBf[nbrACMI.patch().index()];
+ scalarField& nbrPhiNonOverlapp =
+ meshPhiBf[nbrNonOverlapPatch.patch().index()];
+
+ forAll(nbrPhip, facei)
+ {
+ if (newTgtAddr[facei].empty())
+ {
+ nbrPhip[facei] = 0.0;
+ }
+ else
+ {
+ const face& fAMI = nbrACMI.patch().localFaces()[facei];
+
+ // Note: using raw point locations to calculate the geometric
+ // area - faces areas are currently scaled (decoupled from
+ // mesh points)
+ const scalar geomArea = fAMI.mag(nbrACMI.patch().localPoints());
+ nbrPhip[facei] *= nbrACMI.magSf()[facei]/geomArea;
+ }
+ }
+
+ forAll(nbrPhiNonOverlapp, facei)
+ {
+ const scalar w = 1.0 - cyclicACMIPolyPatch_.tgtMask()[facei];
+ nbrPhiNonOverlapp[facei] *= w;
+ }
+}
+
+
// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.H b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.H
index 037df53fe3797d7dffdef2af16cafdfc6a325603..c3aed02df707219e3826b30a8ee8c8ae7b140df4 100644
--- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.H
+++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.H
@@ -65,12 +65,15 @@ protected:
// Protected Member functions
- //- Update the patch areas after AMI update
- void updateAreas() const;
+ //- Helper function to reset the FV patch areas from the primitive patch
+ void resetPatchAreas(const fvPatch& fvp) const;
//- Make patch weighting factors
void makeWeights(scalarField&) const;
+ //- Correct patches after moving points
+ virtual void movePoints();
+
public:
@@ -134,12 +137,7 @@ public:
//- Return a reference to the AMI interpolator
virtual const AMIPatchToPatchInterpolation& AMI() const
{
- const AMIPatchToPatchInterpolation& AMI =
- cyclicACMIPolyPatch_.AMI();
-
- updateAreas();
-
- return AMI;
+ return cyclicACMIPolyPatch_.AMI();
}
//- Are the cyclic planes parallel
@@ -169,8 +167,8 @@ public:
}
//- Return true if this patch is coupled. This is equivalent
- // to the coupledPolyPatch::coupled() if parallel running or
- // both sides present, false otherwise
+ //- to the coupledPolyPatch::coupled() if parallel running or
+ //- both sides present, false otherwise
virtual bool coupled() const;
//- Return delta (P to N) vectors across coupled patch
@@ -183,8 +181,6 @@ public:
const Field<Type>& fld
) const
{
- updateAreas();
-
return
cyclicACMIPolyPatch_.cyclicAMIPolyPatch::interpolate
(
@@ -192,7 +188,7 @@ public:
);
}
- //- Interpolate (make sure to have uptodate areas)
+ //- Interpolate (make sure to have up-to-date areas)
template<class Type>
tmp<Field<Type>> interpolate
(
@@ -206,7 +202,7 @@ public:
// Interface transfer functions
//- Return the values of the given internal data adjacent to
- // the interface as a field
+ //- the interface as a field
virtual tmp<labelField> interfaceInternalField
(
const labelUList& internalData
diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C
index 356c61db10f7e9ed761916ed54bb4cbbb417b6bc..a9c03e340c24d43d29e28aea662dc733344e10f2 100644
--- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C
+++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C
@@ -29,6 +29,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "fvMesh.H"
#include "transform.H"
+#include "surfaceFields.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -82,8 +83,9 @@ void Foam::cyclicAMIFvPatch::makeWeights(scalarField& w) const
forAll(deltas, facei)
{
- scalar di = deltas[facei];
- scalar dni = nbrDeltas[facei];
+ // Note use of mag
+ scalar di = mag(deltas[facei]);
+ scalar dni = mag(nbrDeltas[facei]);
w[facei] = dni/(di + dni);
}
@@ -174,4 +176,75 @@ Foam::tmp<Foam::labelField> Foam::cyclicAMIFvPatch::internalFieldTransfer
}
+void Foam::cyclicAMIFvPatch::movePoints()
+{
+ if (!owner() || !cyclicAMIPolyPatch_.createAMIFaces())
+ {
+ // Only manipulating patch face areas and mesh motion flux if the AMI
+ // creates additional faces
+ return;
+ }
+
+ // Update face data based on values set by the AMI manipulations
+ const_cast<vectorField&>(Sf()) = cyclicAMIPolyPatch_.faceAreas();
+ const_cast<vectorField&>(Cf()) = cyclicAMIPolyPatch_.faceCentres();
+ const_cast<scalarField&>(magSf()) = mag(Sf());
+
+ const cyclicAMIFvPatch& nbr = neighbPatch();
+ const_cast<vectorField&>(nbr.Sf()) = nbr.cyclicAMIPatch().faceAreas();
+ const_cast<vectorField&>(nbr.Cf()) = nbr.cyclicAMIPatch().faceCentres();
+ const_cast<scalarField&>(nbr.magSf()) = mag(nbr.Sf());
+
+
+ // Set consitent mesh motion flux
+ // TODO: currently maps src mesh flux to tgt - update to
+ // src = src + mapped(tgt) and tgt = tgt + mapped(src)?
+
+ const fvMesh& mesh = boundaryMesh().mesh();
+ surfaceScalarField& meshPhi = const_cast<fvMesh&>(mesh).setPhi();
+ surfaceScalarField::Boundary& meshPhiBf = meshPhi.boundaryFieldRef();
+
+ if (cyclicAMIPolyPatch_.owner())
+ {
+ scalarField& phip = meshPhiBf[patch().index()];
+ forAll(phip, facei)
+ {
+ const face& f = cyclicAMIPolyPatch_.localFaces()[facei];
+
+ // Note: using raw point locations to calculate the geometric
+ // area - faces areas are currently scaled by the AMI weights
+ // (decoupled from mesh points)
+ const scalar geomArea = f.mag(cyclicAMIPolyPatch_.localPoints());
+
+ const scalar scaledArea = magSf()[facei];
+ phip[facei] *= scaledArea/geomArea;
+ }
+
+ scalarField srcMeshPhi(phip);
+ if (Pstream::parRun())
+ {
+ AMI().srcMap().distribute(srcMeshPhi);
+ }
+
+ const labelListList& tgtToSrcAddr = AMI().tgtAddress();
+ scalarField& nbrPhip = meshPhiBf[nbr.index()];
+
+ forAll(tgtToSrcAddr, tgtFacei)
+ {
+ // Note: now have 1-to-1 mapping so tgtToSrcAddr[tgtFacei] is size 1
+ const label srcFacei = tgtToSrcAddr[tgtFacei][0];
+ nbrPhip[tgtFacei] = -srcMeshPhi[srcFacei];
+ }
+
+ DebugInfo
+ << "patch:" << patch().name()
+ << " sum(area):" << gSum(magSf())
+ << " min(mag(faceAreas):" << gMin(magSf())
+ << " sum(meshPhi):" << gSum(phip) << nl
+ << " sum(nbrMeshPhi):" << gSum(nbrPhip) << nl
+ << endl;
+ }
+}
+
+
// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H
index 2b643895982dbf9bedab765a25894d67dc8a79e4..4245dd4be83c232b14262ca5840853d8c5b3ac28 100644
--- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H
+++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.H
@@ -68,6 +68,9 @@ protected:
//- Make patch weighting factors
void makeWeights(scalarField&) const;
+ //- Correct patches after moving points
+ virtual void movePoints();
+
public:
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolation/surfaceInterpolation.C b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolation/surfaceInterpolation.C
index 2b31e5c66bded5a82ef8aef4f05d36121a1e54ff..c716e1456a060ee1a870bf5af17121018cd7a778 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolation/surfaceInterpolation.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolation/surfaceInterpolation.C
@@ -172,7 +172,6 @@ void Foam::surfaceInterpolation::makeWeights() const
// ... and reference to the internal field of the weighting factors
scalarField& w = weights.primitiveFieldRef();
-
forAll(owner, facei)
{
// Note: mag in the dot-product.
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C
index b7c297829447d505c46765b0afa9119056d508fa..090a3812a2f6105a93362af93ecb9a0bcf7ebc97 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C
@@ -33,6 +33,8 @@ License
#include "flipOp.H"
#include "profiling.H"
+#define DEBUGAMI(msg){Pout<< "[" << __FILE__ << ":" << __LINE__ << "]: " << #msg << "=" << msg << endl;}
+
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@@ -49,10 +51,12 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolationMethodNames_
{ interpolationMethod::imPartialFaceAreaWeight, "partialFaceAreaWeightAMI" }
});
+
template<class SourcePatch, class TargetPatch>
bool Foam::AMIInterpolation<SourcePatch, TargetPatch>::cacheIntersections_ =
false;
+
template<class SourcePatch, class TargetPatch>
template<class Patch>
Foam::tmp<Foam::scalarField>
@@ -264,7 +268,6 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::agglomerate
}
}
-
// Agglomerate weights and indices
if (targetMapPtr.valid())
{
@@ -298,7 +301,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::agglomerate
// the slots are equal to face indices.
// A mapDistribute has:
// - a subMap : these are face indices
- // - a constructMap : these are from 'transferred-date' to slots
+ // - a constructMap : these are from 'transferred-data' to slots
labelListList tgtSubMap(Pstream::nProcs());
@@ -616,10 +619,12 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
+ srcCentroids_(),
tgtMagSf_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
+ tgtCentroids_(),
triMode_(triMode),
srcMapPtr_(nullptr),
tgtMapPtr_(nullptr)
@@ -649,10 +654,12 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
+ srcCentroids_(),
tgtMagSf_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
+ tgtCentroids_(),
triMode_(triMode),
srcMapPtr_(nullptr),
tgtMapPtr_(nullptr)
@@ -683,10 +690,12 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
+ srcCentroids_(),
tgtMagSf_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
+ tgtCentroids_(),
triMode_(triMode),
srcMapPtr_(nullptr),
tgtMapPtr_(nullptr)
@@ -717,10 +726,12 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
+ srcCentroids_(),
tgtMagSf_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
+ tgtCentroids_(),
triMode_(triMode),
srcMapPtr_(nullptr),
tgtMapPtr_(nullptr)
@@ -834,13 +845,6 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
}
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-template<class SourcePatch, class TargetPatch>
-Foam::AMIInterpolation<SourcePatch, TargetPatch>::~AMIInterpolation()
-{}
-
-
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@@ -932,6 +936,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
(
srcAddress_,
srcWeights_,
+ srcCentroids_,
tgtAddress_,
tgtWeights_
);
@@ -1009,9 +1014,16 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
AMIPtr->normaliseWeights(true, *this);
// Cache maps and reset addresses
- List<Map<label>> cMap;
- srcMapPtr_.reset(new mapDistribute(globalSrcFaces, tgtAddress_, cMap));
- tgtMapPtr_.reset(new mapDistribute(globalTgtFaces, srcAddress_, cMap));
+ List<Map<label>> cMapSrc;
+ srcMapPtr_.reset
+ (
+ new mapDistribute(globalSrcFaces, tgtAddress_, cMapSrc)
+ );
+ List<Map<label>> cMapTgt;
+ tgtMapPtr_.reset
+ (
+ new mapDistribute(globalTgtFaces, srcAddress_, cMapTgt)
+ );
}
else
{
@@ -1033,6 +1045,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
(
srcAddress_,
srcWeights_,
+ srcCentroids_,
tgtAddress_,
tgtWeights_
);
@@ -1056,6 +1069,42 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
}
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
+(
+ autoPtr<mapDistribute>&& srcToTgtMap,
+ autoPtr<mapDistribute>&& tgtToSrcMap,
+ labelListList&& srcAddress,
+ scalarListList&& srcWeights,
+ labelListList&& tgtAddress,
+ scalarListList&& tgtWeights
+)
+{
+ DebugInFunction<< endl;
+
+ srcAddress_.transfer(srcAddress);
+ srcWeights_.transfer(srcWeights);
+ tgtAddress_.transfer(tgtAddress);
+ tgtWeights_.transfer(tgtWeights);
+
+ // Reset the sums of the weights
+ srcWeightsSum_.setSize(srcWeights_.size());
+ forAll(srcWeights_, facei)
+ {
+ srcWeightsSum_[facei] = sum(srcWeights_[facei]);
+ }
+
+ tgtWeightsSum_.setSize(tgtWeights_.size());
+ forAll(tgtWeights_, facei)
+ {
+ tgtWeightsSum_[facei] = sum(tgtWeights_[facei]);
+ }
+
+ srcMapPtr_ = srcToTgtMap;
+ tgtMapPtr_ = tgtToSrcMap;
+}
+
+
template<class SourcePatch, class TargetPatch>
void Foam::AMIInterpolation<SourcePatch, TargetPatch>::append
(
@@ -1397,7 +1446,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToSource
FatalErrorInFunction
<< "Employing default values when sum of weights falls below "
<< lowWeightCorrection_
- << " but supplied default field size is not equal to target "
+ << " but supplied default field size is not equal to source "
<< "patch size" << nl
<< " default values = " << defaultValues.size() << nl
<< " source patch = " << srcAddress_.size() << nl
@@ -1693,6 +1742,73 @@ const
}
+template<class SourcePatch, class TargetPatch>
+bool Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkSymmetricWeights
+(
+ const bool log
+) const
+{
+ if (Pstream::parRun() && (singlePatchProc_ == -1))
+ {
+ Log << "Checks only valid for serial running (currently)" << endl;
+
+ return true;
+ }
+
+ bool symmetricSrc = true;
+
+ Log << " Checking for missing src face in tgt lists" << nl;
+
+ forAll(srcAddress_, srcFacei)
+ {
+ const labelList& tgtIds = srcAddress_[srcFacei];
+ for (const label tgtFacei : tgtIds)
+ {
+ if (!tgtAddress_[tgtFacei].found(srcFacei))
+ {
+ symmetricSrc = false;
+
+ Log << " srcFacei:" << srcFacei
+ << " not found in tgtToSrc list for tgtFacei:"
+ << tgtFacei << nl;
+ }
+ }
+ }
+
+ if (symmetricSrc)
+ {
+ Log << " - symmetric" << endl;
+ }
+
+ bool symmetricTgt = true;
+
+ Log << " Checking for missing tgt face in src lists" << nl;
+
+ forAll(tgtAddress_, tgtFacei)
+ {
+ const labelList& srcIds = tgtAddress_[tgtFacei];
+ for (const label srcFacei : srcIds)
+ {
+ if (!srcAddress_[srcFacei].found(tgtFacei))
+ {
+ symmetricTgt = false;
+
+ Log << " tgtFacei:" << tgtFacei
+ << " not found in srcToTgt list for srcFacei:"
+ << srcFacei << nl;
+ }
+ }
+ }
+
+ if (symmetricTgt)
+ {
+ Log << " - symmetric" << endl;
+ }
+
+ return symmetricSrc && symmetricTgt;
+}
+
+
template<class SourcePatch, class TargetPatch>
void Foam::AMIInterpolation<SourcePatch, TargetPatch>::writeFaceConnectivity
(
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H
index fcf418ff6bf2b67087a6c2cfbcfed0ca3c3fd8fb..2e57a3c7dd3da042b11e146eacdacbe8ae2f7f9d 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H
@@ -61,6 +61,7 @@ SourceFiles
#include "globalIndex.H"
#include "ops.H"
#include "Enum.H"
+#include "pointList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -146,6 +147,9 @@ private:
//- Sum of weights of target faces per source face
scalarField srcWeightsSum_;
+ //- Centroid of target faces per source face
+ pointListList srcCentroids_;
+
// Target patch
@@ -161,6 +165,9 @@ private:
//- Sum of weights of source faces per target face
scalarField tgtWeightsSum_;
+ //- Centroid of source faces per target face
+ pointListList tgtCentroids_;
+
//- Face triangulation mode
const faceAreaIntersect::triangulationMode triMode_;
@@ -347,7 +354,7 @@ public:
//- Destructor
- ~AMIInterpolation();
+ ~AMIInterpolation() = default;
// Typedef to SourcePatch type this AMIInterpolation is instantiated on
typedef SourcePatch sourcePatchType;
@@ -399,6 +406,12 @@ public:
//- patch weights (i.e. the sum before normalisation)
inline scalarField& srcWeightsSum();
+ //- Return const access to source patch face centroids
+ inline const pointListList& srcCentroids() const;
+
+ //- Return access to source patch face centroids
+ inline pointListList& srcCentroids();
+
//- Source map pointer - valid only if singlePatchProc = -1
//- This gets source data into a form to be consumed by
//- tgtAddress, tgtWeights
@@ -448,6 +461,22 @@ public:
const TargetPatch& tgtPatch
);
+ void update
+ (
+ autoPtr<mapDistribute>&& srcToTgtMap,
+ autoPtr<mapDistribute>&& tgtToSrcMap,
+ labelListList&& srcAddress,
+ scalarListList&& srcWeights,
+ labelListList&& tgtAddress,
+ scalarListList&& tgtWeights
+ );
+
+ void setAreas(const scalarList& srcMagSf, const scalarList& tgtMagSf)
+ {
+ srcMagSf_ = srcMagSf;
+ tgtMagSf_ = tgtMagSf;
+ }
+
//- Append additional addressing and weights
void append
(
@@ -582,6 +611,10 @@ public:
// Checks
+ //- Check if src addresses are present in tgt addresses and
+ //- viceversa
+ bool checkSymmetricWeights(const bool log) const;
+
//- Write face connectivity as OBJ file
void writeFaceConnectivity
(
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H
index 9cd5e51ff1dce34480e328ddd6cd7532e31c5fff..027a504521d7bf285941b0e4fe8a414ac7f9a72b 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H
@@ -115,6 +115,22 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeightsSum()
}
+template<class SourcePatch, class TargetPatch>
+inline const Foam::pointListList&
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcCentroids() const
+{
+ return srcCentroids_;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+inline Foam::pointListList&
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcCentroids()
+{
+ return srcCentroids_;
+}
+
+
template<class SourcePatch, class TargetPatch>
inline const Foam::mapDistribute&
Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcMap() const
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
index 3b1f06d405cf74a363d3ae5066bb522fe4500b7b..5a8f11ec74b6d925f1fa0ec55ec074b517a27ceb 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
@@ -31,6 +31,8 @@ License
#include "mapDistribute.H"
#include "unitConversion.H"
+#include "findNearestMaskedOp.H"
+
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@@ -231,17 +233,25 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::resetTree()
template<class SourcePatch, class TargetPatch>
Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::findTargetFace
(
- const label srcFacei
+ const label srcFacei,
+ const UList<label>& excludeFaces,
+ const label srcFacePti
) const
{
label targetFacei = -1;
const pointField& srcPts = srcPatch_.points();
const face& srcFace = srcPatch_[srcFacei];
- const point srcPt = srcFace.centre(srcPts);
- const scalar srcFaceArea = srcMagSf_[srcFacei];
- pointIndexHit sample = treePtr_->findNearest(srcPt, 10.0*srcFaceArea);
+ findNearestMaskedOp<TargetPatch> fnOp(*treePtr_, excludeFaces);
+
+ boundBox bb(srcPts, srcFace, false);
+
+ const point srcPt =
+ srcFacePti == -1 ? bb.centre() : srcPts[srcFace[srcFacePti]];
+
+ const pointIndexHit sample =
+ treePtr_->findNearest(srcPt, magSqr(bb.max() - bb.centre()), fnOp);
if (sample.hit())
{
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H
index bf569aad560a8e392ed97e9277573397b55bead0..e019c384cea9d5997738eeaf75ce2b69046d8199 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H
@@ -45,6 +45,7 @@ SourceFiles
#include "treeDataPrimitivePatch.H"
#include "treeBoundBoxList.H"
#include "runTimeSelectionTables.H"
+#include "pointList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -101,7 +102,7 @@ protected:
List<scalar> tgtMagSf_;
//- Labels of faces that are not overlapped by any target faces
- //- (should be empty for correct functioning)
+ //- (should be empty for correct functioning for fully covered AMIs)
labelList srcNonOverlap_;
//- Octree used to find face seeds
@@ -145,8 +146,12 @@ protected:
//- Reset the octree for the target patch face search
void resetTree();
- //- Find face on target patch that overlaps source face
- label findTargetFace(const label srcFacei) const;
+ label findTargetFace
+ (
+ const label srcFacei,
+ const UList<label>& excludeFaces = UList<label>::null(),
+ const label srcFacePti = -1
+ ) const;
//- Add faces neighbouring facei to the ID list
void appendNbrFaces
@@ -251,6 +256,7 @@ public:
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei = -1,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C
index 894b33fe0f63d5f0292d48e3a46495ca7ced0d30..14989fdfbed58aa035df160d765e266136da3f37 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C
@@ -205,13 +205,6 @@ Foam::directAMI<SourcePatch, TargetPatch>::directAMI
{}
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-template<class SourcePatch, class TargetPatch>
-Foam::directAMI<SourcePatch, TargetPatch>::~directAMI()
-{}
-
-
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@@ -219,6 +212,7 @@ void Foam::directAMI<SourcePatch, TargetPatch>::calculate
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H
index a451cad7da2c5e5c3bb4a15ea42b31eab12245be..a07aa4f617b08122b838da3f1722e1112ef9417f 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H
@@ -118,7 +118,7 @@ public:
//- Destructor
- virtual ~directAMI();
+ virtual ~directAMI() = default;
// Member Functions
@@ -130,6 +130,7 @@ public:
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei = -1,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C
index 3cbe9cd65f8cedcbb26730f9722640790cef792d..f770dfbb58e99f7fef51ee95ab427f15d29aaf7b 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C
@@ -28,6 +28,7 @@ License
#include "faceAreaWeightAMI.H"
#include "profiling.H"
+#include "OBJstream.H"
// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
@@ -36,6 +37,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing
(
List<DynamicList<label>>& srcAddr,
List<DynamicList<scalar>>& srcWght,
+ List<DynamicList<point>>& srcCtr,
List<DynamicList<label>>& tgtAddr,
List<DynamicList<scalar>>& tgtWght,
label srcFacei,
@@ -60,7 +62,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing
seedFaces[srcFacei] = tgtFacei;
// list to keep track of whether src face can be mapped
- boolList mapFlag(nFacesRemaining, true);
+ bitSet mapFlag(nFacesRemaining, true);
// reset starting seed
label startSeedi = 0;
@@ -68,6 +70,9 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing
DynamicList<label> nonOverlapFaces;
do
{
+ nbrFaces.clear();
+ visitedFaces.clear();
+
// Do advancing front starting from srcFacei,tgtFacei
bool faceProcessed = processSourceFace
(
@@ -79,11 +84,12 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcAddressing
srcAddr,
srcWght,
+ srcCtr,
tgtAddr,
tgtWght
);
- mapFlag[srcFacei] = false;
+ mapFlag.unset(srcFacei);
nFacesRemaining--;
@@ -122,9 +128,10 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace
// list of faces currently visited for srcFacei to avoid multiple hits
DynamicList<label>& visitedFaces,
- // temporary storage for addressing and weights
+ // temporary storage for addressing, weights and centroid
List<DynamicList<label>>& srcAddr,
List<DynamicList<scalar>>& srcWght,
+ List<DynamicList<point>>& srcCtr,
List<DynamicList<label>>& tgtAddr,
List<DynamicList<scalar>>& tgtWght
)
@@ -136,9 +143,6 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace
return false;
}
- nbrFaces.clear();
- visitedFaces.clear();
-
// append initial target face and neighbours
nbrFaces.append(tgtStartFacei);
this->appendNbrFaces
@@ -158,16 +162,24 @@ bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace
// process new target face
label tgtFacei = nbrFaces.remove();
visitedFaces.append(tgtFacei);
- scalar area = interArea(srcFacei, tgtFacei);
+
+ scalar interArea = 0;
+ vector interCentroid(Zero);
+ calcInterArea(srcFacei, tgtFacei, interArea, interCentroid);
// store when intersection fractional area > tolerance
- if (area/this->srcMagSf_[srcFacei] > faceAreaIntersect::tolerance())
+ if
+ (
+ interArea/this->srcMagSf_[srcFacei]
+ > faceAreaIntersect::tolerance()
+ )
{
srcAddr[srcFacei].append(tgtFacei);
- srcWght[srcFacei].append(area);
+ srcWght[srcFacei].append(interArea);
+ srcCtr[srcFacei].append(interCentroid);
tgtAddr[tgtFacei].append(srcFacei);
- tgtWght[tgtFacei].append(area);
+ tgtWght[tgtFacei].append(interArea);
this->appendNbrFaces
(
@@ -199,10 +211,10 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
label& startSeedi,
label& srcFacei,
label& tgtFacei,
- const boolList& mapFlag,
+ const bitSet& mapFlag,
labelList& seedFaces,
const DynamicList<label>& visitedFaces,
- bool errorOnNotFound
+ const bool errorOnNotFound
) const
{
addProfiling(ami, "faceAreaWeightAMI::setNextFaces");
@@ -216,15 +228,14 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
bool valuesSet = false;
for (label faceS: srcNbrFaces)
{
- if (mapFlag[faceS] && seedFaces[faceS] == -1)
+ if (mapFlag.test(faceS) && seedFaces[faceS] == -1)
{
for (label faceT : visitedFaces)
{
const scalar threshold =
this->srcMagSf_[faceS]*faceAreaIntersect::tolerance();
- // store when intersection fractional area > tolerance
- // Check that faces have enough overlap for robust walking
+ // Store when intersection fractional area > threshold
if (overlaps(faceS, faceT, threshold))
{
seedFaces[faceS] = faceT;
@@ -248,25 +259,31 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
else
{
// try to use existing seed
+ label facei = startSeedi;
+ if (!mapFlag.test(startSeedi))
+ {
+ facei = mapFlag.find_next(facei);
+ }
+ const label startSeedi0 = facei;
+
bool foundNextSeed = false;
- for (label facei = startSeedi; facei < mapFlag.size(); ++facei)
+ while (facei != -1)
{
- if (mapFlag[facei])
+ if (!foundNextSeed)
{
- if (!foundNextSeed)
- {
- startSeedi = facei;
- foundNextSeed = true;
- }
+ startSeedi = facei;
+ foundNextSeed = true;
+ }
- if (seedFaces[facei] != -1)
- {
- srcFacei = facei;
- tgtFacei = seedFaces[facei];
+ if (seedFaces[facei] != -1)
+ {
+ srcFacei = facei;
+ tgtFacei = seedFaces[facei];
- return;
- }
+ return;
}
+
+ facei = mapFlag.find_next(facei);
}
// perform new search to find match
@@ -276,25 +293,18 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
<< "target face" << endl;
}
- foundNextSeed = false;
- for (label facei = startSeedi; facei < mapFlag.size(); ++facei)
+ facei = startSeedi0;
+ while (facei != -1)
{
- if (mapFlag[facei])
- {
- if (!foundNextSeed)
- {
- startSeedi = facei + 1;
- foundNextSeed = true;
- }
-
- srcFacei = facei;
- tgtFacei = this->findTargetFace(srcFacei);
+ srcFacei = facei;
+ tgtFacei = this->findTargetFace(srcFacei, visitedFaces);
- if (tgtFacei >= 0)
- {
- return;
- }
+ if (tgtFacei >= 0)
+ {
+ return;
}
+
+ facei = mapFlag.find_next(facei);
}
if (errorOnNotFound)
@@ -307,16 +317,16 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
template<class SourcePatch, class TargetPatch>
-Foam::scalar Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::interArea
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calcInterArea
(
const label srcFacei,
- const label tgtFacei
+ const label tgtFacei,
+ scalar& area,
+ vector& centroid
) const
{
addProfiling(ami, "faceAreaWeightAMI::interArea");
- scalar area = 0;
-
const pointField& srcPoints = this->srcPatch_.points();
const pointField& tgtPoints = this->tgtPatch_.points();
@@ -331,7 +341,7 @@ Foam::scalar Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::interArea
|| (this->tgtMagSf_[tgtFacei] < ROOTVSMALL)
)
{
- return area;
+ return;
}
// create intersection object
@@ -359,7 +369,7 @@ Foam::scalar Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::interArea
if (magN > ROOTVSMALL)
{
- area = inter.calc(src, tgt, n/magN);
+ inter.calc(src, tgt, n/magN, area, centroid);
}
else
{
@@ -371,13 +381,24 @@ Foam::scalar Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::interArea
<< endl;
}
+ if
+ (
+ AMIInterpolation<SourcePatch, TargetPatch>::cacheIntersections_
+ && debug
+ )
+ {
+ static OBJstream tris("intersectionTris.obj");
+ const auto& triPts = inter.triangles();
+ for (const auto& tp : triPts)
+ {
+ tris.write(triPointRef(tp[0], tp[1], tp[2]), false);
+ }
+ }
if ((debug > 1) && (area > 0))
{
this->writeIntersectionOBJ(area, src, tgt, srcPoints, tgtPoints);
}
-
- return area;
}
@@ -452,96 +473,71 @@ restartUncoveredSourceFace
(
List<DynamicList<label>>& srcAddr,
List<DynamicList<scalar>>& srcWght,
+ List<DynamicList<point>>& srcCtr,
List<DynamicList<label>>& tgtAddr,
List<DynamicList<scalar>>& tgtWght
)
{
addProfiling(ami, "faceAreaWeightAMI::restartUncoveredSourceFace");
- // Collect all src faces with a low weight
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+ // Note: exclude faces in srcNonOverlap_ for ACMI?
- labelHashSet lowWeightFaces(100);
- forAll(srcWght, srcFacei)
- {
- scalar s = sum(srcWght[srcFacei]);
- scalar t = s/this->srcMagSf_[srcFacei];
+ label nBelowMinWeight = 0;
+ const scalar minWeight = 0.95;
- if (t < 0.5)
- {
- lowWeightFaces.insert(srcFacei);
- }
- }
+ // list of tgt face neighbour faces
+ DynamicList<label> nbrFaces(10);
- if (debug)
- {
- Pout<< "faceAreaWeightAMI: restarting search on "
- << lowWeightFaces.size() << " faces since sum of weights < 0.5"
- << endl;
- }
+ // list of faces currently visited for srcFacei to avoid multiple hits
+ DynamicList<label> visitedFaces(10);
- if (lowWeightFaces.size() > 0)
+ forAll(srcWght, srcFacei)
{
- // Erase all the lowWeight source faces from the target
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+ const scalar s = sum(srcWght[srcFacei]);
+ const scalar t = s/this->srcMagSf_[srcFacei];
- DynamicList<label> okSrcFaces(10);
- DynamicList<scalar> okSrcWeights(10);
- forAll(tgtAddr, tgtFacei)
+ if (t < minWeight)
{
- okSrcFaces.clear();
- okSrcWeights.clear();
- DynamicList<label>& srcFaces = tgtAddr[tgtFacei];
- DynamicList<scalar>& srcWeights = tgtWght[tgtFacei];
- forAll(srcFaces, i)
- {
- if (!lowWeightFaces.found(srcFaces[i]))
- {
- okSrcFaces.append(srcFaces[i]);
- okSrcWeights.append(srcWeights[i]);
- }
- }
- if (okSrcFaces.size() < srcFaces.size())
- {
- srcFaces.transfer(okSrcFaces);
- srcWeights.transfer(okSrcWeights);
- }
- }
-
-
+ ++nBelowMinWeight;
- // Restart search from best hit
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+ const face& f = this->srcPatch_[srcFacei];
- // list of tgt face neighbour faces
- DynamicList<label> nbrFaces(10);
-
- // list of faces currently visited for srcFacei to avoid multiple hits
- DynamicList<label> visitedFaces(10);
-
- for (const label srcFacei : lowWeightFaces)
- {
- const label tgtFacei = this->findTargetFace(srcFacei);
- if (tgtFacei != -1)
+ forAll(f, fpi)
{
- //bool faceProcessed =
- processSourceFace
- (
- srcFacei,
- tgtFacei,
-
- nbrFaces,
- visitedFaces,
+ const label tgtFacei =
+ this->findTargetFace(srcFacei, srcAddr[srcFacei], fpi);
- srcAddr,
- srcWght,
- tgtAddr,
- tgtWght
- );
- // ? Check faceProcessed to see if restarting has worked.
+ if (tgtFacei != -1)
+ {
+ nbrFaces.clear();
+ visitedFaces = srcAddr[srcFacei];
+
+ (void)processSourceFace
+ (
+ srcFacei,
+ tgtFacei,
+
+ nbrFaces,
+ visitedFaces,
+
+ srcAddr,
+ srcWght,
+ srcCtr,
+ tgtAddr,
+ tgtWght
+ );
+ }
}
}
}
+
+ if (debug && nBelowMinWeight)
+ {
+ WarningInFunction
+ << "Restarted search on " << nBelowMinWeight
+ << " faces since sum of weights < " << minWeight
+ << endl;
+ }
}
@@ -572,14 +568,43 @@ Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::faceAreaWeightAMI
{
this->triangulatePatch(srcPatch, srcTris_, this->srcMagSf_);
this->triangulatePatch(tgtPatch, tgtTris_, this->tgtMagSf_);
-}
+ if (debug)
+ {
+ static label nAMI = 0;
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+ // Write out triangulated surfaces as OBJ files
+ OBJstream srcTriObj("srcTris_" + Foam::name(nAMI) + ".obj");
+ const pointField& srcPts = srcPatch.points();
+ forAll(srcTris_, facei)
+ {
+ const DynamicList<face>& faces = srcTris_[facei];
+ for (const face& f : faces)
+ {
+ srcTriObj.write
+ (
+ triPointRef(srcPts[f[0]], srcPts[f[1]], srcPts[f[2]])
+ );
+ }
+ }
-template<class SourcePatch, class TargetPatch>
-Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::~faceAreaWeightAMI()
-{}
+ OBJstream tgtTriObj("tgtTris_" + Foam::name(nAMI) + ".obj");
+ const pointField& tgtPts = tgtPatch.points();
+ forAll(tgtTris_, facei)
+ {
+ const DynamicList<face>& faces = tgtTris_[facei];
+ for (const face& f : faces)
+ {
+ tgtTriObj.write
+ (
+ triPointRef(tgtPts[f[0]], tgtPts[f[1]], tgtPts[f[2]])
+ );
+ }
+ }
+
+ ++nAMI;
+ }
+}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@@ -589,6 +614,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei,
@@ -608,6 +634,9 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
tgtFacei
);
+ srcCentroids.setSize(srcAddress.size());
+
+
if (!ok)
{
return;
@@ -616,6 +645,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
// temporary storage for addressing and weights
List<DynamicList<label>> srcAddr(this->srcPatch_.size());
List<DynamicList<scalar>> srcWght(srcAddr.size());
+ List<DynamicList<point>> srcCtr(srcAddr.size());
List<DynamicList<label>> tgtAddr(this->tgtPatch_.size());
List<DynamicList<scalar>> tgtWght(tgtAddr.size());
@@ -623,6 +653,7 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
(
srcAddr,
srcWght,
+ srcCtr,
tgtAddr,
tgtWght,
srcFacei,
@@ -644,17 +675,18 @@ void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
(
srcAddr,
srcWght,
+ srcCtr,
tgtAddr,
tgtWght
);
}
-
- // transfer data to persistent storage
+ // Transfer data to persistent storage
forAll(srcAddr, i)
{
srcAddress[i].transfer(srcAddr[i]);
srcWeights[i].transfer(srcWght[i]);
+ srcCentroids[i].transfer(srcCtr[i]);
}
forAll(tgtAddr, i)
{
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H
index 0c6dec88d6ce98d25f95c8efa330720eea418278..10622d727a13580a49d044f5d0264c91242d45e9 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H
@@ -80,11 +80,13 @@ protected:
// Marching front
- //- Calculate addressing and weights using temporary storage
+ //- Calculate addressing, weights and centroids using temporary
+ //- storage
virtual void calcAddressing
(
List<DynamicList<label>>& srcAddress,
List<DynamicList<scalar>>& srcWeights,
+ List<DynamicList<point>>& srcCentroids,
List<DynamicList<label>>& tgtAddress,
List<DynamicList<scalar>>& tgtWeights,
label srcFacei,
@@ -100,6 +102,7 @@ protected:
DynamicList<label>& visitedFaces,
List<DynamicList<label>>& srcAddr,
List<DynamicList<scalar>>& srcWght,
+ List<DynamicList<point>>& srcCtr,
List<DynamicList<label>>& tgtAddr,
List<DynamicList<scalar>>& tgtWght
);
@@ -109,6 +112,7 @@ protected:
(
List<DynamicList<label>>& srcAddr,
List<DynamicList<scalar>>& srcWght,
+ List<DynamicList<point>>& srcCtr,
List<DynamicList<label>>& tgtAddr,
List<DynamicList<scalar>>& tgtWght
);
@@ -119,20 +123,22 @@ protected:
label& startSeedi,
label& srcFacei,
label& tgtFacei,
- const boolList& mapFlag,
+ const bitSet& mapFlag,
labelList& seedFaces,
const DynamicList<label>& visitedFaces,
- bool errorOnNotFound = true
+ const bool errorOnNotFound = true
) const;
// Evaluation
//- Area of intersection between source and target faces
- virtual scalar interArea
+ virtual void calcInterArea
(
const label srcFacei,
- const label tgtFacei
+ const label tgtFacei,
+ scalar& area,
+ vector& centroid
) const;
//- Return true if faces overlap
@@ -165,18 +171,19 @@ public:
//- Destructor
- virtual ~faceAreaWeightAMI();
+ virtual ~faceAreaWeightAMI() = default;
// Member Functions
// Manipulation
- //- Update addressing and weights
+ //- Update addressing, weights and centroids
virtual void calculate
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei = -1,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C
index 5bfaa347ca68747b5ae725a9fe2c8f63b7761a61..a159fae7a0ba328d906ee66fce64556855ae0472 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C
@@ -190,13 +190,6 @@ Foam::mapNearestAMI<SourcePatch, TargetPatch>::mapNearestAMI
{}
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-template<class SourcePatch, class TargetPatch>
-Foam::mapNearestAMI<SourcePatch, TargetPatch>::~mapNearestAMI()
-{}
-
-
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@@ -204,6 +197,7 @@ void Foam::mapNearestAMI<SourcePatch, TargetPatch>::calculate
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H
index b59e0682d3bcc5f21b8d301e561c64371e9b7b25..f347fcc7676455e875d77ee1a0ff82d3cc8cef3f 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H
@@ -123,7 +123,7 @@ public:
//- Destructor
- virtual ~mapNearestAMI();
+ virtual ~mapNearestAMI() = default;
// Member Functions
@@ -135,6 +135,7 @@ public:
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei = -1,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C
index 78275e0111c64cf05b958dcfeee5f80f35aea498..5692a656baaf96b738d1c99b6bf2235694aea5e6 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.C
@@ -35,7 +35,7 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
label& startSeedi,
label& srcFacei,
label& tgtFacei,
- const boolList& mapFlag,
+ const bitSet& mapFlag,
labelList& seedFaces,
const DynamicList<label>& visitedFaces,
const bool errorOnNotFound
@@ -78,14 +78,6 @@ partialFaceAreaWeightAMI
{}
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-template<class SourcePatch, class TargetPatch>
-Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::
-~partialFaceAreaWeightAMI()
-{}
-
-
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@@ -100,6 +92,7 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei,
@@ -122,9 +115,12 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
return;
}
+ srcCentroids.setSize(srcAddress.size());
+
// temporary storage for addressing and weights
List<DynamicList<label>> srcAddr(this->srcPatch_.size());
List<DynamicList<scalar>> srcWght(srcAddr.size());
+ List<DynamicList<point>> srcCtr(srcAddr.size());
List<DynamicList<label>> tgtAddr(this->tgtPatch_.size());
List<DynamicList<scalar>> tgtWght(tgtAddr.size());
@@ -132,6 +128,7 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
(
srcAddr,
srcWght,
+ srcCtr,
tgtAddr,
tgtWght,
srcFacei,
@@ -143,6 +140,7 @@ void Foam::partialFaceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
{
srcAddress[i].transfer(srcAddr[i]);
srcWeights[i].transfer(srcWght[i]);
+ srcCentroids[i].transfer(srcCtr[i]);
}
forAll(tgtAddr, i)
{
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H
index becb1c69fb4f3e1cf08f68db6e40c8850c54f5a1..2d0255d6dcb4ca659736e23f0c6756dfda3e1891 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/partialFaceAreaWeightAMI/partialFaceAreaWeightAMI.H
@@ -65,19 +65,22 @@ private:
//- No copy assignment
void operator=(const partialFaceAreaWeightAMI&) = delete;
- // Marching front
- //- Set the source and target seed faces
- virtual void setNextFaces
- (
- label& startSeedi,
- label& srcFacei,
- label& tgtFacei,
- const boolList& mapFlag,
- labelList& seedFaces,
- const DynamicList<label>& visitedFaces,
- bool errorOnNotFound = true
- ) const;
+protected:
+
+ // Protected Member Functions
+
+ //- Set the source and target seed faces
+ virtual void setNextFaces
+ (
+ label& startSeedi,
+ label& srcFacei,
+ label& tgtFacei,
+ const bitSet& mapFlag,
+ labelList& seedFaces,
+ const DynamicList<label>& visitedFaces,
+ const bool errorOnNotFound = true
+ ) const;
public:
@@ -100,7 +103,7 @@ public:
//- Destructor
- virtual ~partialFaceAreaWeightAMI();
+ virtual ~partialFaceAreaWeightAMI() = default;
// Member Functions
@@ -118,6 +121,7 @@ public:
(
labelListList& srcAddress,
scalarListList& srcWeights,
+ pointListList& srcCentroids,
labelListList& tgtAddress,
scalarListList& tgtWeights,
label srcFacei = -1,
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/findNearestMaskedOp.H b/src/meshTools/AMIInterpolation/AMIInterpolation/findNearestMaskedOp.H
new file mode 100644
index 0000000000000000000000000000000000000000..4bc2507cc998b85706a7b372031a32402704034f
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/findNearestMaskedOp.H
@@ -0,0 +1,62 @@
+#include "indexedOctree.H"
+#include "treeDataPrimitivePatch.H"
+
+namespace Foam
+{
+
+template<class PatchType>
+class findNearestMaskedOp
+{
+ const indexedOctree<treeDataPrimitivePatch<PatchType>>& tree_;
+ const labelUList& excludeIndices_;
+
+public:
+
+ findNearestMaskedOp
+ (
+ const indexedOctree<treeDataPrimitivePatch<PatchType>>& tree,
+ const labelUList& excludeIndices
+ )
+ :
+ tree_(tree),
+ excludeIndices_(excludeIndices)
+ {}
+
+ void operator()
+ (
+ const labelUList& indices,
+ const point& sample,
+
+ scalar& nearestDistSqr,
+ label& minIndex,
+ point& nearestPoint
+ ) const
+ {
+ const treeDataPrimitivePatch<PatchType>& shape = tree_.shapes();
+ const PatchType& patch = shape.patch();
+
+ const pointField& points = patch.points();
+
+ forAll(indices, i)
+ {
+ const label index = indices[i];
+
+ if (!excludeIndices_.found(index))
+ {
+ const typename PatchType::FaceType& f = patch[index];
+
+ pointHit nearHit = f.nearestPoint(sample, points);
+ scalar distSqr = sqr(nearHit.distance());
+
+ if (distSqr < nearestDistSqr)
+ {
+ nearestDistSqr = distSqr;
+ minIndex = index;
+ nearestPoint = nearHit.rawPoint();
+ }
+ }
+ }
+ }
+};
+
+} // End namespace Foam
diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C
index 97428ec0e3bf66b9772c21e025ec297a8d54455d..8976e57b586dc365a0b6dd4366d35a19fc567a82 100644
--- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C
+++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C
@@ -219,13 +219,15 @@ void Foam::faceAreaIntersect::triSliceWithPlane
}
-Foam::scalar Foam::faceAreaIntersect::triangleIntersect
+void Foam::faceAreaIntersect::triangleIntersect
(
const triPoints& src,
const point& tgt0,
const point& tgt1,
const point& tgt2,
- const vector& n
+ const vector& n,
+ scalar& area,
+ vector& centroid
) const
{
// Work storage
@@ -241,7 +243,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
const scalar srcArea(triArea(src));
if (srcArea < ROOTVSMALL)
{
- return 0.0;
+ return;
}
// Typical length scale
@@ -255,7 +257,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
scalar s = mag(tgt1 - tgt0);
if (s < ROOTVSMALL)
{
- return 0.0;
+ return;
}
// Note: outer product with n pre-scaled with edge length. This is
@@ -268,7 +270,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
// Triangle either zero edge length (s == 0) or
// perpendicular to face normal n. In either case zero
// overlap area
- return 0.0;
+ return;
}
else
{
@@ -279,7 +281,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
if (nWorkTris1 == 0)
{
- return 0.0;
+ return;
}
// Edge 1
@@ -290,7 +292,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
scalar s = mag(tgt2 - tgt1);
if (s < ROOTVSMALL)
{
- return 0.0;
+ return;
}
const vector n1((tgt1 - tgt2)^(-s*n));
const scalar magSqrN1(magSqr(n1));
@@ -300,7 +302,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
// Triangle either zero edge length (s == 0) or
// perpendicular to face normal n. In either case zero
// overlap area
- return 0.0;
+ return;
}
else
{
@@ -315,7 +317,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
if (nWorkTris2 == 0)
{
- return 0.0;
+ return;
}
}
}
@@ -328,7 +330,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
scalar s = mag(tgt2 - tgt0);
if (s < ROOTVSMALL)
{
- return 0.0;
+ return;
}
const vector n2((tgt2 - tgt0)^(-s*n));
const scalar magSqrN2(magSqr(n2));
@@ -338,7 +340,7 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
// Triangle either zero edge length (s == 0) or
// perpendicular to face normal n. In either case zero
// overlap area
- return 0.0;
+ return;
}
else
{
@@ -353,23 +355,25 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
if (nWorkTris1 == 0)
{
- return 0.0;
+ return;
}
else
{
// Calculate area of sub-triangles
- scalar area = 0.0;
for (label i = 0; i < nWorkTris1; ++i)
{
- area += triArea(workTris1[i]);
+ // Area of intersection
+ const scalar currArea = triArea(workTris1[i]);
+ area += currArea;
+
+ // Area-weighted centroid of intersection
+ centroid += currArea*triCentroid(workTris1[i]);
if (cacheTriangulation_)
{
triangles_.append(workTris1[i]);
}
}
-
- return area;
}
}
}
@@ -443,12 +447,13 @@ void Foam::faceAreaIntersect::triangulate
}
-
-Foam::scalar Foam::faceAreaIntersect::calc
+void Foam::faceAreaIntersect::calc
(
const face& faceA,
const face& faceB,
- const vector& n
+ const vector& n,
+ scalar& area,
+ vector& centroid
) const
{
if (cacheTriangulation_)
@@ -456,8 +461,10 @@ Foam::scalar Foam::faceAreaIntersect::calc
triangles_.clear();
}
+ area = 0.0;
+ centroid = vector::zero;
+
// Intersect triangles
- scalar totalArea = 0.0;
for (const face& triA : trisA_)
{
triPoints tpA = getTriPoints(pointsA_, triA, false);
@@ -466,32 +473,38 @@ Foam::scalar Foam::faceAreaIntersect::calc
{
if (reverseB_)
{
- totalArea +=
- triangleIntersect
- (
- tpA,
- pointsB_[triB[0]],
- pointsB_[triB[1]],
- pointsB_[triB[2]],
- n
- );
+ triangleIntersect
+ (
+ tpA,
+ pointsB_[triB[0]],
+ pointsB_[triB[1]],
+ pointsB_[triB[2]],
+ n,
+ area,
+ centroid
+ );
}
else
{
- totalArea +=
- triangleIntersect
- (
- tpA,
- pointsB_[triB[2]],
- pointsB_[triB[1]],
- pointsB_[triB[0]],
- n
- );
+ triangleIntersect
+ (
+ tpA,
+ pointsB_[triB[2]],
+ pointsB_[triB[1]],
+ pointsB_[triB[0]],
+ n,
+ area,
+ centroid
+ );
}
}
}
- return totalArea;
+ // Area weighed centroid
+ if (area > 0)
+ {
+ centroid /= area;
+ }
}
@@ -503,8 +516,10 @@ bool Foam::faceAreaIntersect::overlaps
const scalar threshold
) const
{
+ scalar area = 0.0;
+ vector centroid(Zero);
+
// Intersect triangles
- scalar totalArea = 0.0;
for (const face& triA : trisA_)
{
const triPoints tpA = getTriPoints(pointsA_, triA, false);
@@ -513,30 +528,32 @@ bool Foam::faceAreaIntersect::overlaps
{
if (reverseB_)
{
- totalArea +=
- triangleIntersect
- (
- tpA,
- pointsB_[triB[0]],
- pointsB_[triB[1]],
- pointsB_[triB[2]],
- n
- );
+ triangleIntersect
+ (
+ tpA,
+ pointsB_[triB[0]],
+ pointsB_[triB[1]],
+ pointsB_[triB[2]],
+ n,
+ area,
+ centroid
+ );
}
else
{
- totalArea +=
- triangleIntersect
- (
- tpA,
- pointsB_[triB[2]],
- pointsB_[triB[1]],
- pointsB_[triB[0]],
- n
- );
+ triangleIntersect
+ (
+ tpA,
+ pointsB_[triB[2]],
+ pointsB_[triB[1]],
+ pointsB_[triB[0]],
+ n,
+ area,
+ centroid
+ );
}
- if (totalArea > threshold)
+ if (area > threshold)
{
return true;
}
diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H
index f0ba9b3e92510c4ad9a11f1b8ed8116982b6b817..916d40e36f7cf3d02504b8ac97eedc30f5260b12 100644
--- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H
+++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H
@@ -97,6 +97,7 @@ private:
// Static data members
+ //- Tolerance
static scalar tol;
@@ -132,6 +133,9 @@ private:
//- Return triangle area
inline scalar triArea(const triPoints& t) const;
+ //- Return triangle centre
+ inline vector triCentroid(const triPoints& t) const;
+
//- Slice triangle with plane and generate new cut sub-triangles
void triSliceWithPlane
(
@@ -143,13 +147,15 @@ private:
) const;
//- Return area of intersection of triangles src and tgt
- scalar triangleIntersect
+ void triangleIntersect
(
const triPoints& src,
const point& tgt0,
const point& tgt1,
const point& tgt2,
- const vector& n
+ const vector& n,
+ scalar& area,
+ vector& centroid
) const;
@@ -199,12 +205,15 @@ public:
DynamicList<face>& faces
);
- //- Return area of intersection of faceA with faceB
- scalar calc
+ //- Return area of intersection of faceA with faceB and effective
+ //- face centre
+ void calc
(
const face& faceA,
const face& faceB,
- const vector& n
+ const vector& n,
+ scalar& area,
+ vector& centroid
) const;
//- Return area of intersection of faceA with faceB
diff --git a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H
index 12113a594d6cc71da2ad12ad5c0a85b7417dc46b..836b87d76200a29f5ec519dfca87d7559d1a50ad 100644
--- a/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H
+++ b/src/meshTools/AMIInterpolation/faceAreaIntersect/faceAreaIntersectI.H
@@ -113,6 +113,15 @@ inline Foam::scalar Foam::faceAreaIntersect::triArea(const triPoints& t) const
}
+inline Foam::vector Foam::faceAreaIntersect::triCentroid
+(
+ const triPoints& t
+) const
+{
+ return (t[0] + t[1] + t[2])/3;
+}
+
+
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
Foam::scalar& Foam::faceAreaIntersect::tolerance()
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C
index 3cd26ef2b82f8069e73f5bb0215a4719d4968ad4..cbba50fb341b641aeee43b308a25d338107b6c38 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C
+++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2013-2016 OpenFOAM Foundation
- Copyright (C) 2017-2018 OpenCFD Ltd.
+ Copyright (C) 2017-2019 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -41,195 +41,225 @@ namespace Foam
addToRunTimeSelectionTable(polyPatch, cyclicACMIPolyPatch, dictionary);
}
-const Foam::scalar Foam::cyclicACMIPolyPatch::tolerance_ = 1e-10;
-
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
+void Foam::cyclicACMIPolyPatch::reportCoverage
+(
+ const word& name,
+ const scalarField& weightSum
+) const
+{
+ label nUncovered = 0;
+ label nCovered = 0;
+ for (const scalar sum : weightSum)
+ {
+ if (sum < tolerance_)
+ {
+ ++nUncovered;
+ }
+ else if (sum > scalar(1) - tolerance_)
+ {
+ ++nCovered;
+ }
+ }
+ reduce(nUncovered, sumOp<label>());
+ reduce(nCovered, sumOp<label>());
+ label nTotal = returnReduce(weightSum.size(), sumOp<label>());
+
+ Info<< "ACMI: Patch " << name << " uncovered/blended/covered = "
+ << nUncovered << ", " << nTotal-nUncovered-nCovered
+ << ", " << nCovered << endl;
+}
+
+
+void Foam::cyclicACMIPolyPatch::resetAMI
+(
+ const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod
+) const
+{
+ resetAMI(boundaryMesh().mesh().points(), AMIMethod);
+}
+
+
void Foam::cyclicACMIPolyPatch::resetAMI
(
- const AMIPatchToPatchInterpolation::interpolationMethod&
+ const UList<point>& points,
+ const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod
) const
{
- if (owner())
+ if (!owner())
{
- const polyPatch& nonOverlapPatch = this->nonOverlapPatch();
+ return;
+ }
+
+ const polyPatch& nonOverlapPatch = this->nonOverlapPatch();
+ if (debug)
+ {
+ Pout<< "cyclicACMIPolyPatch::resetAMI : recalculating weights"
+ << " for " << name() << " and " << nonOverlapPatch.name()
+ << endl;
+ }
+
+
+WarningInFunction<< "DEACTIVATED clearGeom()" << endl;
+// if (boundaryMesh().mesh().hasCellCentres())
+ if (0 && boundaryMesh().mesh().hasCellCentres())
+ {
if (debug)
{
- Pout<< "cyclicACMIPolyPatch::resetAMI : recalculating weights"
+ Pout<< "cyclicACMIPolyPatch::resetAMI : clearing cellCentres"
<< " for " << name() << " and " << nonOverlapPatch.name()
<< endl;
}
- if (boundaryMesh().mesh().hasCellCentres())
- {
- if (debug)
- {
- Pout<< "cyclicACMIPolyPatch::resetAMI : clearing cellCentres"
- << " for " << name() << " and " << nonOverlapPatch.name()
- << endl;
- }
-
- //WarningInFunction
- // << "The mesh already has cellCentres calculated when"
- // << " resetting ACMI " << name() << "." << endl
- // << "This is a problem since ACMI adapts the face areas"
- // << " (to close cells) so this has" << endl
- // << "to be done before cell centre calculation." << endl
- // << "This can happen if e.g. the cyclicACMI is after"
- // << " any processor patches in the boundary." << endl;
- const_cast<polyMesh&>
- (
- boundaryMesh().mesh()
- ).primitiveMesh::clearGeom();
- }
+ //WarningInFunction
+ // << "The mesh already has cellCentres calculated when"
+ // << " resetting ACMI " << name() << "." << endl
+ // << "This is a problem since ACMI adapts the face areas"
+ // << " (to close cells) so this has" << endl
+ // << "to be done before cell centre calculation." << endl
+ // << "This can happen if e.g. the cyclicACMI is after"
+ // << " any processor patches in the boundary." << endl;
+ const_cast<polyMesh&>
+ (
+ boundaryMesh().mesh()
+ ).primitiveMesh::clearGeom();
+ }
- // Trigger re-building of faceAreas
- (void)boundaryMesh().mesh().faceAreas();
+ // Trigger re-building of faceAreas
+ (void)boundaryMesh().mesh().faceAreas();
- // Calculate the AMI using partial face-area-weighted. This leaves
- // the weights as fractions of local areas (sum(weights) = 1 means
- // face is fully covered)
- cyclicAMIPolyPatch::resetAMI
- (
- AMIPatchToPatchInterpolation::imPartialFaceAreaWeight
- );
+ // Calculate the AMI using partial face-area-weighted. This leaves
+ // the weights as fractions of local areas (sum(weights) = 1 means
+ // face is fully covered)
+ cyclicAMIPolyPatch::resetAMI
+ (
+ points,
+ AMIPatchToPatchInterpolation::imPartialFaceAreaWeight
+ );
- AMIPatchToPatchInterpolation& AMI =
- const_cast<AMIPatchToPatchInterpolation&>(this->AMI());
+ const AMIPatchToPatchInterpolation& AMI = this->AMI();
- // Output some stats. AMIInterpolation will have already output the
- // average weights ("sum(weights) min:1 max:1 average:1")
- {
- const scalarField& wghtsSum = AMI.srcWeightsSum();
+ // Output some statistics
+ reportCoverage("source", AMI.srcWeightsSum());
+ reportCoverage("target", AMI.tgtWeightsSum());
- label nUncovered = 0;
- label nCovered = 0;
- forAll(wghtsSum, facei)
- {
- scalar sum = wghtsSum[facei];
- if (sum < tolerance_)
- {
- nUncovered++;
- }
- else if (sum > scalar(1)-tolerance_)
- {
- nCovered++;
- }
- }
- reduce(nUncovered, sumOp<label>());
- reduce(nCovered, sumOp<label>());
- label nTotal = returnReduce(wghtsSum.size(), sumOp<label>());
+ // Set the mask fields
+ // Note:
+ // - assumes that the non-overlap patches are decomposed using the same
+ // decomposition as the coupled patches (per side)
+ srcMask_ = min(scalar(1), max(scalar(0), AMI.srcWeightsSum()));
+ tgtMask_ = min(scalar(1), max(scalar(0), AMI.tgtWeightsSum()));
- Info<< "ACMI: Patch source uncovered/blended/covered = "
- << nUncovered << ", " << nTotal-nUncovered-nCovered
- << ", " << nCovered << endl;
+ if (debug)
+ {
+ Pout<< "resetAMI" << endl;
+ {
+ const cyclicACMIPolyPatch& patch = *this;
+ Pout<< "patch:" << patch.name() << " size:" << patch.size()
+ << " non-overlap patch: " << patch.nonOverlapPatch().name()
+ << " size:" << patch.nonOverlapPatch().size()
+ << " mask size:" << patch.srcMask().size() << endl;
}
{
- const scalarField& wghtsSum = AMI.tgtWeightsSum();
+ const cyclicACMIPolyPatch& patch = this->neighbPatch();
+ Pout<< "patch:" << patch.name() << " size:" << patch.size()
+ << " non-overlap patch: " << patch.nonOverlapPatch().name()
+ << " size:" << patch.nonOverlapPatch().size()
+ << " mask size:" << patch.neighbPatch().tgtMask().size()
+ << endl;
+ }
+ }
+}
- label nUncovered = 0;
- label nCovered = 0;
- forAll(wghtsSum, facei)
- {
- scalar sum = wghtsSum[facei];
- if (sum < tolerance_)
- {
- nUncovered++;
- }
- else if (sum > scalar(1)-tolerance_)
- {
- nCovered++;
- }
- }
- reduce(nUncovered, sumOp<label>());
- reduce(nCovered, sumOp<label>());
- label nTotal = returnReduce(wghtsSum.size(), sumOp<label>());
- Info<< "ACMI: Patch target uncovered/blended/covered = "
- << nUncovered << ", " << nTotal-nUncovered-nCovered
- << ", " << nCovered << endl;
- }
+void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas()
+{
+ if (!owner() || !canResetAMI())
+ {
+ return;
+ }
- srcMask_ =
- min(scalar(1) - tolerance_, max(tolerance_, AMI.srcWeightsSum()));
+ scalePatchFaceAreas(*this);
+ scalePatchFaceAreas(this->neighbPatch());
+}
- tgtMask_ =
- min(scalar(1) - tolerance_, max(tolerance_, AMI.tgtWeightsSum()));
+void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas
+(
+ const cyclicACMIPolyPatch& acmipp
+)
+{
+ // Primitive patch face areas have been cleared/reset based on the raw
+ // points - need to reset to avoid double-accounting of face areas
- // Adapt owner side areas. Note that in uncoupled situations (e.g.
- // decomposePar) srcMask, tgtMask can be zero size.
- if (srcMask_.size())
- {
- vectorField::subField Sf = faceAreas();
- vectorField::subField noSf = nonOverlapPatch.faceAreas();
+ DebugPout
+ << "rescaling non-overlap patch areas" << endl;
- forAll(Sf, facei)
- {
- Sf[facei] *= srcMask_[facei];
- noSf[facei] *= 1.0 - srcMask_[facei];
- }
- }
- // Adapt slave side areas
- if (tgtMask_.size())
- {
- const cyclicACMIPolyPatch& cp =
- refCast<const cyclicACMIPolyPatch>(this->neighbPatch());
- const polyPatch& pp = cp.nonOverlapPatch();
+ const scalar maxTol = scalar(1) - tolerance_;
+ const scalarField& mask = acmipp.mask();
- vectorField::subField Sf = cp.faceAreas();
- vectorField::subField noSf = pp.faceAreas();
+ const polyPatch& nonOverlapPatch = acmipp.nonOverlapPatch();
+ vectorField::subField noSf = nonOverlapPatch.faceAreas();
- forAll(Sf, facei)
- {
- Sf[facei] *= tgtMask_[facei];
- noSf[facei] *= 1.0 - tgtMask_[facei];
- }
+ if (mask.size() != noSf.size())
+ {
+ WarningInFunction
+ << "Inconsistent sizes for patch: " << acmipp.name()
+ << " - not manipulating patches" << nl
+ << " - size: " << size() << nl
+ << " - non-overlap patch size: " << noSf.size() << nl
+ << " - mask size: " << mask.size() << nl
+ << "This is OK for decomposition but should be considered fatal "
+ << "at run-time" << endl;
+
+ return;
+ }
+
+ forAll(noSf, facei)
+ {
+ const scalar w = min(maxTol, max(tolerance_, mask[facei]));
+ noSf[facei] *= scalar(1) - w;
+ }
+
+ if (!createAMIFaces_)
+ {
+ // Note: for topological update (createAMIFaces_ = true)
+ // AMI coupled patch face areas are updated as part of the topological
+ // updates, e.g. by the calls to cyclicAMIPolyPatch's setTopology and
+ // initMovePoints
+ DebugPout
+ << "scaling coupled patch areas" << endl;
+
+ // Scale the coupled patch face areas
+ vectorField::subField Sf = acmipp.faceAreas();
+
+ forAll(Sf, facei)
+ {
+ Sf[facei] *= max(tolerance_, mask[facei]);
}
// Re-normalise the weights since the effect of overlap is already
- // accounted for in the area.
+ // accounted for in the area
+ auto& weights = const_cast<scalarListList&>(acmipp.weights());
+ auto& weightsSum = const_cast<scalarField&>(acmipp.weightsSum());
+ forAll(weights, i)
{
- scalarListList& srcWeights = AMI.srcWeights();
- scalarField& srcWeightsSum = AMI.srcWeightsSum();
- forAll(srcWeights, i)
+ scalarList& wghts = weights[i];
+ if (wghts.size())
{
- scalarList& wghts = srcWeights[i];
- if (wghts.size())
- {
- scalar& sum = srcWeightsSum[i];
+ scalar& sum = weightsSum[i];
- forAll(wghts, j)
- {
- wghts[j] /= sum;
- }
- sum = 1.0;
- }
- }
- }
- {
- scalarListList& tgtWeights = AMI.tgtWeights();
- scalarField& tgtWeightsSum = AMI.tgtWeightsSum();
- forAll(tgtWeights, i)
- {
- scalarList& wghts = tgtWeights[i];
- if (wghts.size())
+ forAll(wghts, j)
{
- scalar& sum = tgtWeightsSum[i];
- forAll(wghts, j)
- {
- wghts[j] /= sum;
- }
- sum = 1.0;
+ wghts[j] /= sum;
}
+ sum = 1.0;
}
}
-
- // Set the updated flag
- updated_ = true;
}
}
@@ -244,9 +274,9 @@ void Foam::cyclicACMIPolyPatch::initGeometry(PstreamBuffers& pBufs)
// Note: calculates transformation and triggers face centre calculation
cyclicAMIPolyPatch::initGeometry(pBufs);
- // Initialise the AMI early to make sure we adapt the face areas before the
- // cell centre calculation gets triggered.
- resetAMI();
+ // On start-up there are no topological updates so scale the face areas
+ // - Note: resetAMI called by cyclicAMIPolyPatch::initGeometry
+ scalePatchFaceAreas();
}
@@ -272,10 +302,10 @@ void Foam::cyclicACMIPolyPatch::initMovePoints
}
// Note: calculates transformation and triggers face centre calculation
+ // - Note: resetAMI called by cyclicAMIPolyPatch::initMovePoints
cyclicAMIPolyPatch::initMovePoints(pBufs, p);
- // Initialise the AMI early. See initGeometry.
- resetAMI();
+ scalePatchFaceAreas();
}
@@ -289,6 +319,8 @@ void Foam::cyclicACMIPolyPatch::movePoints
{
Pout<< "cyclicACMIPolyPatch::movePoints : " << name() << endl;
}
+
+ // When topology is changing, this will scale the duplicate AMI faces
cyclicAMIPolyPatch::movePoints(pBufs, p);
}
@@ -352,8 +384,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch
nonOverlapPatchName_(word::null),
nonOverlapPatchID_(-1),
srcMask_(),
- tgtMask_(),
- updated_(false)
+ tgtMask_()
{
AMIRequireMatch_ = false;
@@ -375,8 +406,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch
nonOverlapPatchName_(dict.lookup("nonOverlapPatch")),
nonOverlapPatchID_(-1),
srcMask_(),
- tgtMask_(),
- updated_(false)
+ tgtMask_()
{
AMIRequireMatch_ = false;
@@ -403,8 +433,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch
nonOverlapPatchName_(pp.nonOverlapPatchName_),
nonOverlapPatchID_(-1),
srcMask_(),
- tgtMask_(),
- updated_(false)
+ tgtMask_()
{
AMIRequireMatch_ = false;
@@ -428,8 +457,7 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch
nonOverlapPatchName_(nonOverlapPatchName),
nonOverlapPatchID_(-1),
srcMask_(),
- tgtMask_(),
- updated_(false)
+ tgtMask_()
{
AMIRequireMatch_ = false;
@@ -459,19 +487,12 @@ Foam::cyclicACMIPolyPatch::cyclicACMIPolyPatch
nonOverlapPatchName_(pp.nonOverlapPatchName_),
nonOverlapPatchID_(-1),
srcMask_(),
- tgtMask_(),
- updated_(false)
+ tgtMask_()
{
AMIRequireMatch_ = false;
}
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-Foam::cyclicACMIPolyPatch::~cyclicACMIPolyPatch()
-{}
-
-
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
const Foam::cyclicACMIPolyPatch& Foam::cyclicACMIPolyPatch::neighbPatch() const
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H
index 5f2bd9dfad989c9a68b8e1d394ce0263c165fcca..c858c14e0485d103ab3930621c137ff830a9d987 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H
+++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.H
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2013-2016 OpenFOAM Foundation
- Copyright (C) 2018 OpenCFD Ltd.
+ Copyright (C) 2018-2019 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -60,9 +60,6 @@ private:
// Private data
- //- Fraction of face area below which face is considered disconnected
- static const scalar tolerance_;
-
//- Name of non-overlapping patch
const word nonOverlapPatchName_;
@@ -75,21 +72,40 @@ private:
//- Mask/weighting for target patch
mutable scalarField tgtMask_;
- //- Flag to indicate that AMI has been updated
- mutable bool updated_;
-
protected:
// Protected Member Functions
- //- Reset the AMI interpolator
+ //- Report coverage statics, e.g. number of uncovered/blended/covered
+ //- faces
+ void reportCoverage
+ (
+ const word& name,
+ const scalarField& weightSum
+ ) const;
+
+ //- Reset the AMI interpolator, supply patch points
+ virtual void resetAMI
+ (
+ const UList<point>& points,
+ const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod =
+ AMIPatchToPatchInterpolation::imFaceAreaWeight
+ ) const;
+
+ //- Reset the AMI interpolator, use current patch points
virtual void resetAMI
(
const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod =
AMIPatchToPatchInterpolation::imFaceAreaWeight
) const;
+ //- Scale patch face areas to maintain physical area
+ virtual void scalePatchFaceAreas();
+
+ //- Scale patch face areas to maintain physical area
+ virtual void scalePatchFaceAreas(const cyclicACMIPolyPatch& acmipp);
+
//- Initialise the calculation of the patch geometry
virtual void initGeometry(PstreamBuffers&);
@@ -111,12 +127,6 @@ 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:
@@ -235,19 +245,13 @@ public:
//- Destructor
- virtual ~cyclicACMIPolyPatch();
+ virtual ~cyclicACMIPolyPatch() = default;
// Member Functions
// Access
- //- Reset the updated flag
- inline void setUpdated(bool flag) const;
-
- //- Return access to the updated flag
- inline bool updated() const;
-
//- Return a reference to the neighbour patch
virtual const cyclicACMIPolyPatch& neighbPatch() const;
@@ -263,9 +267,15 @@ public:
//- Return a reference to the non-overlapping patch
inline polyPatch& nonOverlapPatch();
- //- Mask field where 1 = overlap, 0 = no-overlap
+ //- Mask field where 1 = overlap(coupled), 0 = no-overlap
inline const scalarField& mask() const;
+ //- 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;
+
//- Overlap tolerance
inline static scalar tolerance();
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H
index 42b41196673b4f316e4bcb331b14a4a93b3ea646..50f4022ed63e2ce5204ebda9f9b33b1f653218fa 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H
+++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatchI.H
@@ -27,18 +27,6 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-inline void Foam::cyclicACMIPolyPatch::setUpdated(const bool flag) const
-{
- updated_ = flag;
-}
-
-
-inline bool Foam::cyclicACMIPolyPatch::updated() const
-{
- return updated_;
-}
-
-
inline const Foam::word& Foam::cyclicACMIPolyPatch::nonOverlapPatchName() const
{
return nonOverlapPatchName_;
@@ -69,10 +57,8 @@ inline const Foam::scalarField& Foam::cyclicACMIPolyPatch::mask() const
{
return srcMask_;
}
- else
- {
- return neighbPatch().tgtMask();
- }
+
+ return neighbPatch().tgtMask();
}
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
index 51f14ccb5ed67a11988c894ce82fe7259b128cdd..04cb5b6bb3f63d06358d959f03f0e17184122a68 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
+++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
@@ -27,13 +27,10 @@ License
\*---------------------------------------------------------------------------*/
#include "cyclicAMIPolyPatch.H"
-#include "transformField.H"
#include "SubField.H"
-#include "polyMesh.H"
#include "Time.H"
-#include "addToRunTimeSelectionTable.H"
#include "faceAreaIntersect.H"
-#include "ops.H"
+#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -45,6 +42,7 @@ namespace Foam
addToRunTimeSelectionTable(polyPatch, cyclicAMIPolyPatch, dictionary);
}
+const Foam::scalar Foam::cyclicAMIPolyPatch::tolerance_ = 1e-10;
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
@@ -295,77 +293,114 @@ void Foam::cyclicAMIPolyPatch::resetAMI
const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod
) const
{
- if (owner())
+ resetAMI(boundaryMesh().mesh().points(), AMIMethod);
+}
+
+
+void Foam::cyclicAMIPolyPatch::resetAMI
+(
+ const UList<point>& points,
+ const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod
+) const
+{
+ DebugInFunction << endl;
+
+ if (!owner())
{
- AMIPtr_.clear();
+ return;
+ }
- const polyPatch& nbr = neighbPatch();
- pointField nbrPoints
- (
- neighbPatch().boundaryMesh().mesh().points(),
- neighbPatch().meshPoints()
- );
+ AMIPtr_.clear();
+
+ const cyclicAMIPolyPatch& nbr = neighbPatch();
+ pointField srcPoints(points, meshPoints());
+ pointField nbrPoints(points, neighbPatch().meshPoints());
+
+ if (debug)
+ {
+ const Time& t = boundaryMesh().mesh().time();
+ OFstream os(t.path()/name() + "_neighbourPatch-org.obj");
+ meshTools::writeOBJ(os, neighbPatch().localFaces(), nbrPoints);
+ }
- if (debug)
+ label patchSize0 = size();
+ label nbrPatchSize0 = nbr.size();
+
+ if (createAMIFaces_)
+ {
+ // AMI is created based on the original patch faces (non-extended patch)
+ if (srcFaceIDs_.size())
+ {
+ patchSize0 = srcFaceIDs_.size();
+ }
+ if (tgtFaceIDs_.size())
{
- const Time& t = boundaryMesh().mesh().time();
- OFstream os(t.path()/name() + "_neighbourPatch-org.obj");
- meshTools::writeOBJ(os, neighbPatch().localFaces(), nbrPoints);
+ nbrPatchSize0 = tgtFaceIDs_.size();
}
+ }
- // Transform neighbour patch to local system
- transformPosition(nbrPoints);
- primitivePatch nbrPatch0
+ // Transform neighbour patch to local system
+ transformPosition(nbrPoints);
+ primitivePatch nbrPatch0
+ (
+ SubList<face>
(
- SubList<face>
- (
- nbr.localFaces(),
- nbr.size()
- ),
- nbrPoints
- );
+ nbr.localFaces(),
+ nbrPatchSize0
+ ),
+ nbrPoints
+ );
+ primitivePatch patch0
+ (
+ SubList<face>
+ (
+ localFaces(),
+ patchSize0
+ ),
+ srcPoints
+ );
- if (debug)
- {
- const Time& t = boundaryMesh().mesh().time();
- OFstream osN(t.path()/name() + "_neighbourPatch-trans.obj");
- meshTools::writeOBJ(osN, nbrPatch0.localFaces(), nbrPoints);
- OFstream osO(t.path()/name() + "_ownerPatch.obj");
- meshTools::writeOBJ(osO, this->localFaces(), localPoints());
- }
+ if (debug)
+ {
+ const Time& t = boundaryMesh().mesh().time();
+ OFstream osN(t.path()/name() + "_neighbourPatch-trans.obj");
+ meshTools::writeOBJ(osN, nbrPatch0.localFaces(), nbrPoints);
+
+ OFstream osO(t.path()/name() + "_ownerPatch.obj");
+ meshTools::writeOBJ(osO, this->localFaces(), localPoints());
+ }
- // Construct/apply AMI interpolation to determine addressing and weights
- AMIPtr_.reset
+ // Construct/apply AMI interpolation to determine addressing and weights
+ AMIPtr_.reset
+ (
+ new AMIPatchToPatchInterpolation
(
- new AMIPatchToPatchInterpolation
- (
- *this,
- nbrPatch0,
- surfPtr(),
- faceAreaIntersect::tmMesh,
- AMIRequireMatch_,
- AMIMethod,
- AMILowWeightCorrection_,
- AMIReverse_
- )
- );
+ patch0, // *this,
+ nbrPatch0,
+ surfPtr(),
+ faceAreaIntersect::tmMesh,
+ AMIRequireMatch_,
+ AMIMethod,
+ AMILowWeightCorrection_,
+ AMIReverse_
+ )
+ );
- if (debug)
- {
- Pout<< "cyclicAMIPolyPatch : " << name()
- << " constructed AMI with " << nl
- << " " << "srcAddress:" << AMIPtr_().srcAddress().size()
- << nl
- << " " << "tgAddress :" << AMIPtr_().tgtAddress().size()
- << nl << endl;
- }
+ // Set the updated flag
+ updated_ = true;
+
+ if (debug)
+ {
+ AMIPtr_->checkSymmetricWeights(true);
}
}
void Foam::cyclicAMIPolyPatch::calcTransforms()
{
+ DebugInFunction << endl;
+
const cyclicAMIPolyPatch& half0 = *this;
vectorField half0Areas(half0.size());
forAll(half0, facei)
@@ -402,22 +437,29 @@ void Foam::cyclicAMIPolyPatch::calcTransforms()
void Foam::cyclicAMIPolyPatch::initGeometry(PstreamBuffers& pBufs)
{
- // The AMI is no longer valid. Leave it up to demand-driven calculation
- AMIPtr_.clear();
+ DebugInFunction << endl;
+
+ // Only resetting the AMI if not creating additional AMI faces
+ // - When creating additional faces the AMI is updated externally by the
+ // dynamic mesh via the setTopology function.
+ if (!createAMIFaces_ && canResetAMI())
+ {
+ resetAMI(AMIMethod_);
+ }
polyPatch::initGeometry(pBufs);
// Early calculation of transforms so e.g. cyclicACMI can use them.
// Note: also triggers primitiveMesh face centre. Note that cell
// centres should -not- be calculated
- // since e.g. cyclicACMI override face areas
+ // since e.g. cyclicACMI overrides face areas
calcTransforms();
}
void Foam::cyclicAMIPolyPatch::calcGeometry(PstreamBuffers& pBufs)
{
- // All geometry done inside initGeometry
+ DebugInFunction << endl;
}
@@ -427,14 +469,31 @@ void Foam::cyclicAMIPolyPatch::initMovePoints
const pointField& p
)
{
- // The AMI is no longer valid. Leave it up to demand-driven calculation
- AMIPtr_.clear();
-
- polyPatch::initMovePoints(pBufs, p);
+ DebugInFunction << endl;
// See below. Clear out any local geometry
primitivePatch::movePoints(p);
+ // Note: processorPolyPatch::initMovePoints calls
+ // processorPolyPatch::initGeometry which will trigger calculation of
+ // patch faceCentres() and cell volumes...
+
+ if (createAMIFaces_)
+ {
+ // Note: AMI should have been updated in setTopology
+
+ // faceAreas() and faceCentres() have been reset and will be
+ // recalculated on-demand using the mesh points and no longer
+ // correspond to the scaled areas!
+ restoreScaledGeometry();
+
+ // deltas need to be recalculated to use new face centres!
+ }
+ else
+ {
+ resetAMI(p, AMIMethod_);
+ }
+
// Early calculation of transforms. See above.
calcTransforms();
}
@@ -446,30 +505,55 @@ void Foam::cyclicAMIPolyPatch::movePoints
const pointField& p
)
{
- polyPatch::movePoints(pBufs, p);
-
- // All transformation tensors already done in initMovePoints
+ DebugInFunction << endl;
+
+// Note: not calling movePoints since this will undo our manipulations!
+// polyPatch::movePoints(pBufs, p);
+/*
+ polyPatch::movePoints
+ -> primitivePatch::movePoints
+ -> primitivePatch::clearGeom:
+ deleteDemandDrivenData(localPointsPtr_);
+ deleteDemandDrivenData(faceCentresPtr_);
+ deleteDemandDrivenData(faceAreasPtr_);
+ deleteDemandDrivenData(magFaceAreasPtr_);
+ deleteDemandDrivenData(faceNormalsPtr_);
+ deleteDemandDrivenData(pointNormalsPtr_);
+*/
}
void Foam::cyclicAMIPolyPatch::initUpdateMesh(PstreamBuffers& pBufs)
{
- // The AMI is no longer valid. Leave it up to demand-driven calculation
- AMIPtr_.clear();
+ DebugInFunction << endl;
polyPatch::initUpdateMesh(pBufs);
+
+ if (createAMIFaces_ && boundaryMesh().mesh().topoChanging() && owner())
+ {
+ setAMIFaces();
+ }
}
void Foam::cyclicAMIPolyPatch::updateMesh(PstreamBuffers& pBufs)
{
+ DebugInFunction << endl;
+
+ // Note: this clears out cellCentres(), faceCentres() and faceAreas()
polyPatch::updateMesh(pBufs);
}
void Foam::cyclicAMIPolyPatch::clearGeom()
{
- AMIPtr_.clear();
+ DebugInFunction << endl;
+
+ if (!updatingAMI_)
+ {
+//DEBUG("*** CLEARING AMI ***");
+ AMIPtr_.clear();
+ }
polyPatch::clearGeom();
}
@@ -488,6 +572,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
)
:
coupledPolyPatch(name, size, start, index, bm, patchType, transform),
+ updated_(false),
nbrPatchName_(word::null),
nbrPatchID_(-1),
rotationAxis_(Zero),
@@ -501,7 +586,13 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
AMIRequireMatch_(true),
AMILowWeightCorrection_(-1.0),
surfPtr_(nullptr),
- surfDict_(fileName("surface"))
+ surfDict_(fileName("surface")),
+ createAMIFaces_(false),
+ updatingAMI_(true),
+ srcFaceIDs_(),
+ tgtFaceIDs_(),
+ faceAreas0_(),
+ faceCentres0_()
{
// Neighbour patch might not be valid yet so no transformation
// calculation possible
@@ -518,6 +609,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
)
:
coupledPolyPatch(name, dict, index, bm, patchType),
+ updated_(false),
nbrPatchName_(dict.getOrDefault<word>("neighbourPatch", "")),
coupleGroup_(dict),
nbrPatchID_(-1),
@@ -545,7 +637,13 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
AMIRequireMatch_(true),
AMILowWeightCorrection_(dict.getOrDefault("lowWeightCorrection", -1.0)),
surfPtr_(nullptr),
- surfDict_(dict.subOrEmptyDict("surface"))
+ surfDict_(dict.subOrEmptyDict("surface")),
+ createAMIFaces_(false),
+ updatingAMI_(true),
+ srcFaceIDs_(),
+ tgtFaceIDs_(),
+ faceAreas0_(),
+ faceCentres0_()
{
if (nbrPatchName_ == word::null && !coupleGroup_.valid())
{
@@ -605,6 +703,20 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
// Neighbour patch might not be valid yet so no transformation
// calculation possible
+
+ // If topology change, recover the sizes of the original patches
+ label srcSize0 = 0;
+ if (dict.readIfPresent("srcSize", srcSize0))
+ {
+ srcFaceIDs_.setSize(srcSize0);
+ createAMIFaces_ = true;
+ }
+ label tgtSize0 = 0;
+ if (dict.readIfPresent("tgtSize", tgtSize0))
+ {
+ tgtFaceIDs_.setSize(tgtSize0);
+ createAMIFaces_ = true;
+ }
}
@@ -615,6 +727,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
)
:
coupledPolyPatch(pp, bm),
+ updated_(false),
nbrPatchName_(pp.nbrPatchName_),
coupleGroup_(pp.coupleGroup_),
nbrPatchID_(-1),
@@ -629,7 +742,13 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
AMIRequireMatch_(pp.AMIRequireMatch_),
AMILowWeightCorrection_(pp.AMILowWeightCorrection_),
surfPtr_(nullptr),
- surfDict_(pp.surfDict_)
+ surfDict_(pp.surfDict_),
+ createAMIFaces_(false),
+ updatingAMI_(true),
+ srcFaceIDs_(),
+ tgtFaceIDs_(),
+ faceAreas0_(),
+ faceCentres0_()
{
// Neighbour patch might not be valid yet so no transformation
// calculation possible
@@ -647,6 +766,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
)
:
coupledPolyPatch(pp, bm, index, newSize, newStart),
+ updated_(false),
nbrPatchName_(nbrPatchName),
coupleGroup_(pp.coupleGroup_),
nbrPatchID_(-1),
@@ -661,7 +781,13 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
AMIRequireMatch_(pp.AMIRequireMatch_),
AMILowWeightCorrection_(pp.AMILowWeightCorrection_),
surfPtr_(nullptr),
- surfDict_(pp.surfDict_)
+ surfDict_(pp.surfDict_),
+ createAMIFaces_(false),
+ updatingAMI_(true),
+ srcFaceIDs_(),
+ tgtFaceIDs_(),
+ faceAreas0_(),
+ faceCentres0_()
{
if (nbrPatchName_ == name())
{
@@ -686,6 +812,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
)
:
coupledPolyPatch(pp, bm, index, mapAddressing, newStart),
+ updated_(false),
nbrPatchName_(pp.nbrPatchName_),
coupleGroup_(pp.coupleGroup_),
nbrPatchID_(-1),
@@ -700,13 +827,13 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
AMIRequireMatch_(pp.AMIRequireMatch_),
AMILowWeightCorrection_(pp.AMILowWeightCorrection_),
surfPtr_(nullptr),
- surfDict_(pp.surfDict_)
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-Foam::cyclicAMIPolyPatch::~cyclicAMIPolyPatch()
+ surfDict_(pp.surfDict_),
+ createAMIFaces_(false),
+ updatingAMI_(true),
+ srcFaceIDs_(),
+ tgtFaceIDs_(),
+ faceAreas0_(),
+ faceCentres0_()
{}
@@ -1105,6 +1232,12 @@ void Foam::cyclicAMIPolyPatch::write(Ostream& os) const
{
surfDict_.writeEntry(surfDict_.dictName(), os);
}
+
+ if (createAMIFaces_)
+ {
+ os.writeEntry("srcSize", srcFaceIDs_.size());
+ os.writeEntry("tgtSize", tgtFaceIDs_.size());
+ }
}
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
index 25259ffacdfc8dff92aa2ae82363658b16f99f5e..84fab30b99bae4209b920ab3b154c18723092cc4 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
+++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
@@ -74,6 +74,9 @@ protected:
// Protected data
+ //- Flag to indicate that AMI has been updated
+ mutable bool updated_;
+
//- Name of other half
mutable word nbrPatchName_;
@@ -129,9 +132,47 @@ protected:
const dictionary surfDict_;
+ // Change of topology as AMI is updated
+
+ //- Flag to indicate that new AMI faces will created
+ // Set by the call to changeTopology
+ mutable bool createAMIFaces_;
+
+ mutable bool updatingAMI_;
+
+ labelListList srcFaceIDs_;
+
+ labelListList tgtFaceIDs_;
+
+ //- Temporary storage for AMI face areas
+ mutable vectorField faceAreas0_;
+
+ //- Temporary storage for AMI face centres
+ mutable vectorField faceCentres0_;
+
+
// Protected Member Functions
- //- Reset the AMI interpolator
+ // Topology change
+
+ virtual bool removeAMIFaces(polyTopoChange& topoChange);
+
+ virtual bool addAMIFaces(polyTopoChange& topoChange);
+
+ virtual void setAMIFaces();
+
+ virtual void restoreScaledGeometry();
+
+
+ //- Reset the AMI interpolator, supply patch points
+ virtual void resetAMI
+ (
+ const UList<point>& points,
+ const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod =
+ AMIPatchToPatchInterpolation::imFaceAreaWeight
+ ) const;
+
+ //- Reset the AMI interpolator, use current patch points
virtual void resetAMI
(
const AMIPatchToPatchInterpolation::interpolationMethod& AMIMethod =
@@ -274,15 +315,45 @@ public:
//- Destructor
- virtual ~cyclicAMIPolyPatch();
+ virtual ~cyclicAMIPolyPatch() = default;
// Member Functions
// Access
+ //- Tolerance used e.g. for area calculations/limits
+ static const scalar tolerance_;
+
+ //- Flag to indicate whether the AMI can be reset
+ inline bool canResetAMI() const;
+
+ //- Reset the updated flag
+ inline void setUpdated(bool flag) const;
+
+ //- Return access to the updated flag
+ inline bool updated() const;
+
+ //- Return access to the createAMIFaces flag
+ inline bool createAMIFaces() const;
+
+ //- Return access to the updated flag
+ inline bool updatingAMI() const;
+
+void clearAMI() const
+{
+ AMIPtr_.clear();
+}
+
+ //- Return true if this patch changes the mesh topology
+ // True when createAMIFaces is true
+ virtual bool changeTopology() const;
+
+ //- Set topology changes in the polyTopoChange object
+ virtual bool setTopology(polyTopoChange& topoChange);
+
//- Is patch 'coupled'. Note that on AMI the geometry is not
- // coupled but the fields are!
+ //- coupled but the fields are!
virtual bool coupled() const
{
return false;
@@ -306,9 +377,23 @@ public:
//- Return a reference to the AMI interpolator
const AMIPatchToPatchInterpolation& AMI() const;
+ //- Helper function to return the weights
+ inline const scalarListList& weights() const;
+
+ //- Helper function to return the weights sum
+ inline const scalarField& weightsSum() const;
+
//- Return true if applying the low weight correction
bool applyLowWeightCorrection() const;
+ //- Return access to the initial face areas
+ // Used for topology change
+ inline vectorField& faceAreas0() const;
+
+ //- Return access to the initial face centres
+ // Used for topology change
+ inline vectorField& faceCentres0() const;
+
// Transformations
@@ -388,7 +473,7 @@ public:
);
//- Initialize ordering for primitivePatch. Does not
- // refer to *this (except for name() and type() etc.)
+ //- refer to *this (except for name() and type() etc.)
virtual void initOrder
(
PstreamBuffers&,
@@ -409,7 +494,7 @@ public:
) const;
//- Return face index on neighbour patch which shares point p
- // following trajectory vector n
+ //- following trajectory vector n
label pointFace
(
const label facei,
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H
index ccc36fc3007d3fdc7f7768c85c61edc63250ba42..6d55bf1841dce5b9711ab772ee2f7d55e82e2f1e 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H
+++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H
@@ -6,6 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2013 OpenFOAM Foundation
+ Copyright (C) 2019 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -27,6 +28,36 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+inline bool Foam::cyclicAMIPolyPatch::canResetAMI() const
+{
+ return Pstream::parRun() || !boundaryMesh().mesh().time().processorCase();
+}
+
+
+inline void Foam::cyclicAMIPolyPatch::setUpdated(const bool flag) const
+{
+ updated_ = flag;
+}
+
+
+inline bool Foam::cyclicAMIPolyPatch::updated() const
+{
+ return updated_;
+}
+
+
+inline bool Foam::cyclicAMIPolyPatch::createAMIFaces() const
+{
+ return createAMIFaces_;
+}
+
+
+inline bool Foam::cyclicAMIPolyPatch::updatingAMI() const
+{
+ return updatingAMI_;
+}
+
+
inline const Foam::word& Foam::cyclicAMIPolyPatch::neighbPatchName() const
{
if (nbrPatchName_.empty())
@@ -39,6 +70,38 @@ inline const Foam::word& Foam::cyclicAMIPolyPatch::neighbPatchName() const
return nbrPatchName_;
}
+inline const Foam::scalarListList& Foam::cyclicAMIPolyPatch::weights() const
+{
+ if (owner())
+ {
+ return AMI().srcWeights();
+ }
+
+ return neighbPatch().AMI().tgtWeights();
+}
+
+
+inline const Foam::scalarField& Foam::cyclicAMIPolyPatch::weightsSum() const
+{
+ if (owner())
+ {
+ return AMI().srcWeightsSum();
+ }
+
+ return neighbPatch().AMI().tgtWeightsSum();
+}
+
+
+inline Foam::vectorField& Foam::cyclicAMIPolyPatch::faceAreas0() const
+{
+ return faceAreas0_;
+}
+
+
+inline Foam::vectorField& Foam::cyclicAMIPolyPatch::faceCentres0() const
+{
+ return faceCentres0_;
+}
inline const Foam::vector& Foam::cyclicAMIPolyPatch::rotationAxis() const
{
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C
new file mode 100644
index 0000000000000000000000000000000000000000..f7e33910da639eb3b5bb1dabb87360aba8f9fde3
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C
@@ -0,0 +1,676 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2019 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cyclicAMIPolyPatch.H"
+#include "SubField.H"
+#include "vectorList.H"
+#include "polyTopoChange.H"
+
+// * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
+
+void Foam::cyclicAMIPolyPatch::restoreScaledGeometry()
+{
+ DebugInFunction << endl;
+
+ // Note: only used for topology update (createAMIFaces_ flag)
+ if (!createAMIFaces_)
+ {
+ FatalErrorInFunction
+ << "Attempted to perform topology update when createAMIFaces_ "
+ << "flag is set to false"
+ << abort(FatalError);
+ }
+
+ if (boundaryMesh().mesh().hasCellVolumes())
+ {
+ WarningInFunction
+ << "Mesh already has volumes set!"
+ << endl;
+ }
+
+ vectorField::subField faceAreas = this->faceAreas();
+ vectorField::subField faceCentres = this->faceCentres();
+
+ if (debug)
+ {
+ Info<< "Patch:" << name() << " before: sum(mag(faceAreas)):"
+ << gSum(mag(faceAreas)) << endl;
+ Info<< "Patch:" << name() << " before: sum(mag(faceAreas0)):"
+ << gSum(mag(faceAreas0_)) << endl;
+ }
+
+ faceAreas = faceAreas0_;
+ Info<< "WARNING: ==============================" << endl;
+ Info<< "WARNING: Not updating face cell centres" << endl;
+ Info<< "WARNING: ==============================" << endl;
+ //faceCentres = faceCentres0_;
+ faceAreas0_.clear();
+ faceCentres0_.clear();
+
+ if (debug)
+ {
+ Info<< "Patch:" << name() << " after: sum(mag(faceAreas)):"
+ << gSum(mag(faceAreas)) << endl;
+ Info<< "Patch:" << name() << " after: sum(mag(faceAreas0)):"
+ << gSum(mag(faceAreas0_)) << endl;
+ }
+}
+
+
+bool Foam::cyclicAMIPolyPatch::removeAMIFaces(polyTopoChange& topoChange)
+{
+ DebugInFunction << endl;
+
+ // Note: only used for topology update (createAMIFaces_ flag)
+ if (!createAMIFaces_)
+ {
+ FatalErrorInFunction
+ << "Attempted to perform topology update when createAMIFaces_ "
+ << "flag is set to false"
+ << abort(FatalError);
+ }
+
+ if (!owner())
+ {
+ return false;
+ }
+
+ bool changeRequired = false;
+
+ // Remove any faces that we inserted to create the 1-to-1 match...
+
+ const cyclicAMIPolyPatch& nbr = neighbPatch();
+
+ const label newSrcFaceStart = srcFaceIDs_.size();
+
+ if (newSrcFaceStart != 0)
+ {
+ for (label facei = newSrcFaceStart; facei < size(); ++facei)
+ {
+ changeRequired = true;
+ label meshFacei = start() + facei;
+ topoChange.removeFace(meshFacei, -1);
+ }
+ }
+
+ const label newTgtFaceStart = tgtFaceIDs_.size();
+
+ if (newTgtFaceStart != 0)
+ {
+ for (label facei = newTgtFaceStart; facei < nbr.size(); ++facei)
+ {
+ changeRequired = true;
+ label meshFacei = nbr.start() + facei;
+ topoChange.removeFace(meshFacei, -1);
+ }
+ }
+
+ srcFaceIDs_.clear();
+ tgtFaceIDs_.clear();
+
+ return changeRequired;
+}
+
+
+bool Foam::cyclicAMIPolyPatch::addAMIFaces(polyTopoChange& topoChange)
+{
+ DebugInFunction << endl;
+
+ // Note: only used for topology update (createAMIFaces_ flag = true)
+ if (!createAMIFaces_)
+ {
+ FatalErrorInFunction
+ << "Attempted to perform topology update when createAMIFaces_ "
+ << "flag is set to false"
+ << abort(FatalError);
+ }
+
+ bool changedFaces = false;
+ const cyclicAMIPolyPatch& nbr = neighbPatch();
+
+ polyMesh& mesh = const_cast<polyMesh&>(boundaryMesh().mesh());
+ const faceZoneMesh& faceZones = mesh.faceZones();
+
+ // First face address and weight are used to manipulate the
+ // original face - all other addresses and weights are used to
+ // create additional faces
+ const labelListList& srcToTgtAddr = AMI().srcAddress();
+ const labelListList& tgtToSrcAddr = AMI().tgtAddress();
+
+ const label nSrcFace = srcToTgtAddr.size();
+ const label nTgtFace = tgtToSrcAddr.size();
+
+ srcFaceIDs_.setSize(nSrcFace);
+ tgtFaceIDs_.setSize(nTgtFace);
+
+ label nNewSrcFaces = 0;
+ forAll(srcToTgtAddr, srcFacei)
+ {
+ const labelList& tgtAddr = srcToTgtAddr[srcFacei];
+
+ // No tgt faces linked to srcFacei (ACMI)
+ if (tgtAddr.empty()) continue;
+
+ srcFaceIDs_[srcFacei].setSize(tgtAddr.size());
+ srcFaceIDs_[srcFacei][0] = srcFacei;
+
+ const label meshFacei = start() + srcFacei;
+ for (label addri = 1; addri < tgtAddr.size(); ++addri)
+ {
+ changedFaces = true;
+
+ // Note: new faces reuse originating face points
+ // - but areas are scaled by the weights (later)
+
+ // New source face for each target face address
+ srcFaceIDs_[srcFacei][addri] = nNewSrcFaces + nSrcFace;
+ ++nNewSrcFaces;
+ (void)topoChange.addFace
+ (
+ mesh.faces()[meshFacei], // modified face
+ mesh.faceOwner()[meshFacei], // owner
+ -1, // neighbour
+ -1, // master point
+ -1, // master edge
+ meshFacei, // master face
+ false, // face flip
+ index(), // patch for face
+ faceZones.whichZone(meshFacei), // zone for original face
+ false // face flip in zone
+ );
+ }
+ }
+
+ label nNewTgtFaces = 0;
+ forAll(tgtToSrcAddr, tgtFacei)
+ {
+ const labelList& srcAddr = tgtToSrcAddr[tgtFacei];
+
+ // No src faces linked to tgtFacei (ACMI)
+ if (srcAddr.empty()) continue;
+
+ tgtFaceIDs_[tgtFacei].setSize(srcAddr.size());
+ tgtFaceIDs_[tgtFacei][0] = tgtFacei;
+
+ const label meshFacei = nbr.start() + tgtFacei;
+ for (label addri = 1; addri < srcAddr.size(); ++addri)
+ {
+ changedFaces = true;
+
+ // Note: new faces reuse originating face points
+ // - but areas are scaled by the weights (later)
+
+ // New target face for each source face address
+ tgtFaceIDs_[tgtFacei][addri] = nNewTgtFaces + nTgtFace;
+ ++nNewTgtFaces;
+
+ (void)topoChange.addFace
+ (
+ mesh.faces()[meshFacei], // modified face
+ mesh.faceOwner()[meshFacei], // owner
+ -1, // neighbour
+ -1, // master point
+ -1, // master edge
+ meshFacei, // master face
+ false, // face flip
+ nbr.index(), // patch for face
+ faceZones.whichZone(meshFacei), // zone for original face
+ false // face flip in zone
+ );
+ }
+ }
+
+ Info<< "AMI: Patch " << name() << " additional faces: "
+ << returnReduce(nNewSrcFaces, sumOp<label>()) << nl
+ << "AMI: Patch " << nbr.name() << " additional faces: "
+ << returnReduce(nNewTgtFaces, sumOp<label>())
+ << endl;
+
+ if (debug)
+ {
+ Pout<< "New faces - " << name() << ": " << nNewSrcFaces
+ << " " << nbr.name() << ": " << nNewTgtFaces << endl;
+ }
+
+ return returnReduce(changedFaces, orOp<bool>());
+}
+
+
+void Foam::cyclicAMIPolyPatch::setAMIFaces()
+{
+ // Create new mesh faces so that there is a 1-to-1 correspondence
+ // between faces on each side of the AMI
+
+ // Note: only used for topology update (createAMIFaces_ flag)
+ if (!createAMIFaces_)
+ {
+ FatalErrorInFunction
+ << "Attempted to perform topology update when createAMIFaces_ "
+ << "flag is set to false"
+ << abort(FatalError);
+ }
+
+
+ DebugInFunction << endl;
+
+ if (!owner())
+ {
+ return;
+ }
+
+ const cyclicAMIPolyPatch& nbr = neighbPatch();
+
+ vectorField& nbrFaceAreas0 = nbr.faceAreas0();
+ vectorField& nbrFaceCentres0 = nbr.faceCentres0();
+
+ // Scale the new face areas and set the centroids
+ // Note:
+ // - storing local copies so that they can be re-applied after the call to
+ // movePoints that will reset any changes to the areas and centroids
+ //
+ // - For AMI, src and tgt patches should be the same
+ // - For ACMI they are likely to be different!
+ faceAreas0_ = faceAreas();
+ faceCentres0_ = faceCentres();
+ nbrFaceAreas0 = nbr.faceAreas();
+ nbrFaceCentres0 = nbr.faceCentres();
+
+ // Original AMI info (based on the mesh state when the AMI was evaluated)
+ const labelListList& srcToTgtAddr0 = AMIPtr_->srcAddress();
+ const labelListList& tgtToSrcAddr0 = AMIPtr_->tgtAddress();
+ const pointListList& srcCtr0 = AMIPtr_->srcCentroids();
+ const scalarListList& srcToTgtWght0 = AMIPtr_->srcWeights();
+
+ // New addressing on new mesh (extended by polyTopoChange)
+ labelListList srcToTgtAddr1(size(), labelList());
+ labelListList tgtToSrcAddr1(nbr.size(), labelList());
+
+ // Need to calc new parallel maps (mesh has changed since AMI was computed)
+ autoPtr<mapDistribute> srcToTgtMap1;
+ autoPtr<mapDistribute> tgtToSrcMap1;
+
+ if (AMIPtr_->singlePatchProc() == -1)
+ {
+ // Parallel running
+
+ // Global index based on old patch sizes (when AMI was computed)
+ globalIndex globalSrcFaces0(srcToTgtAddr0.size());
+ globalIndex globalTgtFaces0(tgtToSrcAddr0.size());
+
+ // Global index based on new patch sizes
+ globalIndex globalSrcFaces1(size());
+ globalIndex globalTgtFaces1(nbr.size());
+
+
+ // Gather source side info
+ // =======================
+
+ // Note: using new global index for addressing, and distributed using
+ // the old AMI map
+ labelListList newTgtGlobalFaces(tgtFaceIDs_);
+ forAll(newTgtGlobalFaces, tgtFacei)
+ {
+ globalTgtFaces1.inplaceToGlobal(newTgtGlobalFaces[tgtFacei]);
+ }
+ AMIPtr_->tgtMap().distribute(newTgtGlobalFaces);
+
+ // Now have new tgt face indices for each src face
+
+ labelList globalSrcFaceIDs(identity(srcToTgtAddr0.size()));
+ globalSrcFaces0.inplaceToGlobal(globalSrcFaceIDs);
+ AMIPtr_->srcMap().distribute(globalSrcFaceIDs);
+ // globalSrcFaceIDs now has remote data for each srcFacei0 known to the
+ // tgt patch
+
+ List<List<point>> globalSrcCtrs0(srcCtr0);
+ AMIPtr_->srcMap().distribute(globalSrcCtrs0);
+
+ labelList globalTgtFaceIDs(identity(tgtToSrcAddr0.size()));
+ globalTgtFaces0.inplaceToGlobal(globalTgtFaceIDs);
+ AMIPtr_->tgtMap().distribute(globalTgtFaceIDs);
+ // globalTgtFaceIDs now has remote data for each tgtFacei0 known to the
+ // src patch
+
+ // For debug - send tgt face centres and compare against mapped src
+ // face centres
+ //List<List<point>> globalTgtCtrs0(tgtCtr0);
+ //AMIPtr_->tgtMap().distribute(globalTgtCtrs0);
+
+ labelListList globalTgtToSrcAddr(tgtToSrcAddr0);
+ forAll(tgtToSrcAddr0, tgtFacei0)
+ {
+ forAll(tgtToSrcAddr0[tgtFacei0], addri)
+ {
+ const label globalSrcFacei =
+ globalSrcFaceIDs[tgtToSrcAddr0[tgtFacei0][addri]];
+ globalTgtToSrcAddr[tgtFacei0][addri] = globalSrcFacei;
+ }
+ }
+ AMIPtr_->tgtMap().distribute(globalTgtToSrcAddr);
+
+ labelListList globalSrcToTgtAddr(srcToTgtAddr0);
+ forAll(srcToTgtAddr0, srcFacei0)
+ {
+ forAll(srcToTgtAddr0[srcFacei0], addri)
+ {
+ const label globalTgtFacei =
+ globalTgtFaceIDs[srcToTgtAddr0[srcFacei0][addri]];
+ globalSrcToTgtAddr[srcFacei0][addri] = globalTgtFacei;
+ }
+ }
+ AMIPtr_->srcMap().distribute(globalSrcToTgtAddr);
+
+ label nError = 0;
+ forAll(srcToTgtAddr0, srcFacei0)
+ {
+ const labelList& newSrcFaces = srcFaceIDs_[srcFacei0];
+
+ forAll(newSrcFaces, i)
+ {
+ const label srcFacei1 = newSrcFaces[i];
+
+ // What index did srcFacei0 appear in tgtToSrc0 list?
+ // - if first index, all ok
+ // - else tgt face has been moved to according to tgtFaceIDs_
+ const label tgtFacei0 = srcToTgtAddr0[srcFacei0][i];
+ const label addri =
+ globalTgtToSrcAddr[tgtFacei0].find
+ (
+ globalSrcFaceIDs[srcFacei0]
+ );
+
+ if (addri == -1)
+ {
+ ++nError;
+ continue;
+
+ if (debug)
+ {
+ Pout<< "Unable to find global source face "
+ << globalSrcFaceIDs[srcFacei0]
+ << " in globalTgtToSrcAddr[" << tgtFacei0 << "]: "
+ << globalTgtToSrcAddr[tgtFacei0]
+ << endl;
+ }
+ }
+
+ const label tgtFacei1 = newTgtGlobalFaces[tgtFacei0][addri];
+
+ // Sanity check to see that we've picked the correct face
+ // point tgtCtr0(globalTgtCtrs0[tgtFacei0][addri]);
+ // Pout<< "srcCtr:" << srcCtr0[srcFacei0][i]
+ // << " tgtCtr:" << tgtCtr0 << endl;
+
+ srcToTgtAddr1[srcFacei1] = labelList(1, tgtFacei1);
+ faceAreas0_[srcFacei1] *= srcToTgtWght0[srcFacei0][i];
+ faceCentres0_[srcFacei1] = srcCtr0[srcFacei0][i];
+ }
+ }
+
+ if (nError)
+ {
+ FatalErrorInFunction
+ << "Unable to find " << nError << " global source faces"
+ << abort(FatalError);
+ }
+
+
+ // Gather Target side info
+ // =======================
+
+ labelListList newSrcGlobalFaces(srcFaceIDs_);
+ forAll(newSrcGlobalFaces, srcFacei)
+ {
+ globalSrcFaces1.inplaceToGlobal(newSrcGlobalFaces[srcFacei]);
+ }
+
+ AMIPtr_->srcMap().distribute(newSrcGlobalFaces);
+
+ // Now have new src face indices for each tgt face
+ forAll(tgtToSrcAddr0, tgtFacei0)
+ {
+ const labelList& newTgtFaces = tgtFaceIDs_[tgtFacei0];
+ forAll(newTgtFaces, i)
+ {
+ const label srcFacei0 = tgtToSrcAddr0[tgtFacei0][i];
+
+ const label addri =
+ globalSrcToTgtAddr[srcFacei0].find
+ (
+ globalTgtFaceIDs[tgtFacei0]
+ );
+
+ if (addri == -1)
+ {
+ ++nError;
+ continue;
+
+ if (debug)
+ {
+ Pout<< "Unable to find global target face "
+ << globalTgtFaceIDs[tgtFacei0]
+ << " in globalSrcToTgtAddr[" << srcFacei0 << "]: "
+ << globalSrcToTgtAddr[srcFacei0]
+ << endl;
+ }
+ }
+
+ const label srcFacei1 = newSrcGlobalFaces[srcFacei0][addri];
+
+ // Sanity check to see that we've picked the correct face
+ point srcCtr0(globalSrcCtrs0[srcFacei0][addri]);
+ reverseTransformPosition(srcCtr0, srcFacei0);
+
+ const label tgtFacei1 = newTgtFaces[i];
+ tgtToSrcAddr1[tgtFacei1] = labelList(1, srcFacei1);
+ nbrFaceCentres0[tgtFacei1] = srcCtr0;
+ }
+ }
+
+ if (nError)
+ {
+ FatalErrorInFunction
+ << "Unable to find " << nError << " global target faces"
+ << abort(FatalError);
+ }
+
+ // Update the maps
+ {
+ List<Map<label>> cMap;
+ srcToTgtMap1.reset
+ (
+ new mapDistribute(globalSrcFaces1, tgtToSrcAddr1, cMap)
+ );
+ }
+ {
+ List<Map<label>> cMap;
+ tgtToSrcMap1.reset
+ (
+ new mapDistribute(globalTgtFaces1, srcToTgtAddr1, cMap)
+ );
+ }
+
+ // Reset tgt patch areas using the new map
+ vectorList newSrcGlobalFaceAreas(faceAreas0_);
+
+ srcToTgtMap1->distribute(newSrcGlobalFaceAreas);
+ forAll(nbrFaceAreas0, tgtFacei)
+ {
+ if (!tgtToSrcAddr1[tgtFacei].empty())
+ {
+ const label srcFacei = tgtToSrcAddr1[tgtFacei][0];
+ nbrFaceAreas0[tgtFacei] = -newSrcGlobalFaceAreas[srcFacei];
+ }
+ }
+ }
+ else
+ {
+ label nError = 0;
+ forAll(srcToTgtAddr0, srcFacei0)
+ {
+ const labelList& srcFaceTgtAddr0 = srcToTgtAddr0[srcFacei0];
+ const scalarList& srcFaceTgtWght0 = srcToTgtWght0[srcFacei0];
+ const pointList& srcFaceTgtCtr0 = srcCtr0[srcFacei0];
+ forAll(srcFaceTgtAddr0, addri)
+ {
+ const label srcFacei1 = srcFaceIDs_[srcFacei0][addri];
+
+ // Find which slot srcFacei0 appears in tgt->src addressing
+ const label tgtFacei0 = srcFaceTgtAddr0[addri];
+ const label tgtAddri0 =
+ tgtToSrcAddr0[tgtFacei0].find(srcFacei0);
+
+ if (tgtAddri0 == -1)
+ {
+ ++nError;
+ continue;
+
+ if (debug)
+ {
+ Pout<< "Unable to find source face " << srcFacei0
+ << " in tgtToSrcAddr0[" << tgtFacei0 << "]: "
+ << tgtToSrcAddr0[tgtFacei0]
+ << endl;
+ }
+ }
+
+ const label tgtFacei1 = tgtFaceIDs_[tgtFacei0][tgtAddri0];
+
+ faceAreas0_[srcFacei1] *= srcFaceTgtWght0[addri];
+ nbrFaceAreas0[tgtFacei1] = -faceAreas0_[srcFacei1];
+
+ point pt(srcFaceTgtCtr0[addri]);
+ faceCentres0_[srcFacei1] = pt;
+ reverseTransformPosition(pt, srcFacei0);
+ nbrFaceCentres0[tgtFacei1] = pt;
+
+ // SANITY CHECK
+ // Info<< "srcPt:" << srcFaceCentres[srcFacei1]
+ // << " tgtPt:" << tgtFaceCentres[tgtFacei1] << endl;
+
+ srcToTgtAddr1[srcFacei1] = labelList(1, tgtFacei1);
+ tgtToSrcAddr1[tgtFacei1] = labelList(1, srcFacei1);
+ }
+ }
+
+ if (nError)
+ {
+ FatalErrorInFunction
+ << "Unable to find " << nError
+ << " source faces in tgtToSrcAddr0"
+ << abort(FatalError);
+ }
+ }
+
+ scalarListList newSrcToTgtWeights(srcToTgtAddr1.size());
+ forAll(srcToTgtAddr1, facei)
+ {
+ if (srcToTgtAddr1[facei].size())
+ {
+ newSrcToTgtWeights[facei] = scalarList(1, scalar(1));
+ }
+ else
+ {
+ // No connection - effect of face removed by setting area to a
+ // a small value
+ faceAreas0_[facei] *= tolerance_;
+ }
+ }
+
+ scalarListList newTgtToSrcWeights(tgtToSrcAddr1.size());
+ forAll(tgtToSrcAddr1, facei)
+ {
+ if (tgtToSrcAddr1[facei].size())
+ {
+ newTgtToSrcWeights[facei] = scalarList(1, scalar(1));
+ }
+ else
+ {
+ // No connection - effect of face removed by setting area to a
+ // a small value
+ nbrFaceAreas0[facei] *= tolerance_;
+ }
+ }
+
+ // Update the AMI addressing and weights to reflect the new 1-to-1
+ // correspondence
+ AMIPtr_->update
+ (
+ std::move(srcToTgtMap1),
+ std::move(tgtToSrcMap1),
+ std::move(srcToTgtAddr1),
+ std::move(newSrcToTgtWeights),
+ std::move(tgtToSrcAddr1),
+ std::move(newTgtToSrcWeights)
+ );
+
+ AMIPtr_->setAreas(mag(faceAreas0_), mag(nbrFaceAreas0));
+
+ if (debug)
+ {
+ Pout<< "cyclicAMIPolyPatch : " << name()
+ << " constructed AMI with " << nl
+ << " " << "srcAddress:" << AMIPtr_().srcAddress().size()
+ << nl
+ << " " << "tgAddress :" << AMIPtr_().tgtAddress().size()
+ << nl << endl;
+ }
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+bool Foam::cyclicAMIPolyPatch::changeTopology() const
+{
+ DebugInFunction << endl;
+
+ createAMIFaces_ = true;
+
+ return true;
+}
+
+
+bool Foam::cyclicAMIPolyPatch::setTopology(polyTopoChange& topoChange)
+{
+ DebugInFunction << endl;
+
+ if (createAMIFaces_ && owner())
+ {
+ // Calculate the AMI using the new points
+ // Note: mesh still has old points
+ resetAMI(topoChange.points(), AMIMethod_);
+
+ removeAMIFaces(topoChange);
+
+ addAMIFaces(topoChange);
+
+ return true;
+ }
+
+ return false;
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/Make/files b/src/meshTools/Make/files
index 22b12ee488cfd3e4ba2d8395e20ddb4ff3c4d0cd..c89a1afb64f04ff8dc8dda1d62905eedb4807436 100644
--- a/src/meshTools/Make/files
+++ b/src/meshTools/Make/files
@@ -271,6 +271,7 @@ AMICycPatches=$(AMI)/patches/cyclicAMI
$(AMICycPatches)/cyclicAMILduInterfaceField/cyclicAMILduInterface.C
$(AMICycPatches)/cyclicAMILduInterfaceField/cyclicAMILduInterfaceField.C
$(AMICycPatches)/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
+$(AMICycPatches)/cyclicAMIPolyPatch/cyclicAMIPolyPatchTopologyChange.C
$(AMICycPatches)/cyclicAMIPointPatch/cyclicAMIPointPatch.C
$(AMICycPatches)/cyclicAMIPointPatchField/cyclicAMIPointPatchFields.C
diff --git a/src/meshTools/Make/options b/src/meshTools/Make/options
index 9ee5884e5908ccdb6b1ba3acf323b8c524ec996a..5b4f9ec84181988556c25393c12e10eb3e3b4531 100644
--- a/src/meshTools/Make/options
+++ b/src/meshTools/Make/options
@@ -1,6 +1,7 @@
EXE_INC = \
-I$(LIB_SRC)/fileFormats/lnInclude \
- -I$(LIB_SRC)/surfMesh/lnInclude
+ -I$(LIB_SRC)/surfMesh/lnInclude \
+ -I$(LIB_SRC)/dynamicMesh/lnInclude
LIB_LIBS = \
-lfileFormats \
diff --git a/src/meshTools/polyTopoChange/polyTopoChange.H b/src/meshTools/polyTopoChange/polyTopoChange.H
index d861a2d7da259af1f9cfc462ecf6663ca794c0ef..613f6d398f3bafb16e051d0ca40b2d44c46da63d 100644
--- a/src/meshTools/polyTopoChange/polyTopoChange.H
+++ b/src/meshTools/polyTopoChange/polyTopoChange.H
@@ -94,7 +94,7 @@ class IOobject;
template<class T, class Container> class CompactListList;
/*---------------------------------------------------------------------------*\
- Class polyTopoChange Declaration
+ Class polyTopoChange Declaration
\*---------------------------------------------------------------------------*/
class polyTopoChange
@@ -121,7 +121,7 @@ class polyTopoChange
DynamicList<label> pointMap_;
//- For all original and added points contains new point label.
- // (used to map return value of addPoint to new mesh point)
+ //- (used to map return value of addPoint to new mesh point)
DynamicList<label> reversePointMap_;
//- Zone of point
@@ -324,10 +324,10 @@ class polyTopoChange
);
//- Remove all unused/removed points/faces/cells and update
- // face ordering (always), cell ordering (bandcompression,
- // orderCells=true),
- // point ordering (sorted into internal and boundary points,
- // orderPoints=true)
+ //- face ordering (always), cell ordering (bandcompression,
+ //- orderCells=true),
+ //- point ordering (sorted into internal and boundary points,
+ //- orderPoints=true)
void compact
(
const bool orderCells,
@@ -433,7 +433,7 @@ public:
// Constructors
//- Construct without mesh. Either specify nPatches or use
- // setNumPatches before trying to make a mesh (makeMesh, changeMesh)
+ //- setNumPatches before trying to make a mesh (makeMesh, changeMesh)
polyTopoChange(const label nPatches, const bool strict = true);
//- Construct from mesh. Adds all points/face/cells from mesh.
@@ -471,15 +471,15 @@ public:
}
//- Is point removed?
- // Considered removed if point is GREAT.
+ //- Considered removed if point is GREAT.
inline bool pointRemoved(const label pointi) const;
//- Is face removed?
- // Considered removed if face is empty
+ //- Considered removed if face is empty
inline bool faceRemoved(const label facei) const;
//- Is cell removed?
- // Considered removed if the cellMap is -2
+ //- Considered removed if the cellMap is -2
inline bool cellRemoved(const label celli) const;
@@ -489,7 +489,7 @@ public:
void clear();
//- Add all points/faces/cells of mesh. Additional offset for patch
- // or zone ids.
+ //- or zone ids.
void addMesh
(
const polyMesh& mesh,
@@ -500,7 +500,7 @@ public:
);
//- Explicitly pre-size the dynamic storage for expected mesh
- // size for if construct-without-mesh
+ //- size for if construct-without-mesh
void setCapacity
(
const label nPoints,
@@ -590,7 +590,7 @@ public:
void removeCell(const label celli, const label mergeCelli);
//- Explicitly set the number of patches if construct-without-mesh
- // used.
+ //- used.
inline void setNumPatches(const label nPatches);
// Other
@@ -628,7 +628,6 @@ public:
const bool orderCells = false,
const bool orderPoints = false
);
-
};