From 0839258f63b4498f31b4e96e8a4974682cae690d Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 19 Dec 2013 16:05:21 +0000 Subject: [PATCH] STYLE: medialAxisMover: moved routines into meshRefinement --- .../medialAxisMeshMover.C | 80 +++------------- .../medialAxisMeshMover.H | 31 ------ .../medialAxisMeshMoverTemplates.C | 93 ------------------ .../meshRefinement/meshRefinement.C | 62 +++++++++++- .../meshRefinement/meshRefinement.H | 27 ++++++ .../meshRefinement/meshRefinementRefine.C | 36 +++++-- .../meshRefinement/meshRefinementTemplates.C | 96 ++++++++++++++++--- 7 files changed, 206 insertions(+), 219 deletions(-) delete mode 100644 src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMoverTemplates.C diff --git a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C index a6db1ff0b9..0aa74be6c9 100644 --- a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C +++ b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C @@ -51,7 +51,6 @@ namespace Foam // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // - Foam::labelList Foam::medialAxisMeshMover::getFixedValueBCs ( const pointVectorField& fld @@ -118,63 +117,6 @@ Foam::medialAxisMeshMover::getPatch } -void Foam::medialAxisMeshMover::calculateEdgeWeights -( - const PackedBoolList& isMasterEdge, - const labelList& meshEdges, - const labelList& meshPoints, - const edgeList& edges, - scalarField& edgeWeights, - scalarField& invSumWeight -) const -{ - const pointField& pts = mesh().points(); - - // Calculate edgeWeights and inverse sum of edge weights - edgeWeights.setSize(meshEdges.size()); - invSumWeight.setSize(meshPoints.size()); - - forAll(edges, edgeI) - { - const edge& e = edges[edgeI]; - scalar eMag = max - ( - VSMALL, - mag - ( - pts[meshPoints[e[1]]] - - pts[meshPoints[e[0]]] - ) - ); - edgeWeights[edgeI] = 1.0/eMag; - } - - // Sum per point all edge weights - weightedSum - ( - mesh(), - isMasterEdge, - meshEdges, - meshPoints, - edges, - edgeWeights, - scalarField(meshPoints.size(), 1.0), // data - invSumWeight - ); - - // Inplace invert - forAll(invSumWeight, pointI) - { - scalar w = invSumWeight[pointI]; - - if (w > 0.0) - { - invSumWeight[pointI] = 1.0/w; - } - } -} - - void Foam::medialAxisMeshMover::smoothPatchNormals ( const label nSmoothDisp, @@ -193,8 +135,9 @@ void Foam::medialAxisMeshMover::smoothPatchNormals scalarField edgeWeights(meshEdges.size()); scalarField invSumWeight(meshPoints.size()); - calculateEdgeWeights + meshRefinement::calculateEdgeWeights ( + mesh(), isMasterEdge, meshEdges, meshPoints, @@ -207,7 +150,7 @@ void Foam::medialAxisMeshMover::smoothPatchNormals vectorField average; for (label iter = 0; iter < nSmoothDisp; iter++) { - weightedSum + meshRefinement::weightedSum ( mesh(), isMasterEdge, @@ -284,8 +227,9 @@ void Foam::medialAxisMeshMover::smoothNormals scalarField edgeWeights(meshEdges.size()); scalarField invSumWeight(meshPoints.size()); - calculateEdgeWeights + meshRefinement::calculateEdgeWeights ( + mesh(), isMasterEdge, meshEdges, meshPoints, @@ -297,7 +241,7 @@ void Foam::medialAxisMeshMover::smoothNormals vectorField average; for (label iter = 0; iter < nSmoothDisp; iter++) { - weightedSum + meshRefinement::weightedSum ( mesh(), isMasterEdge, @@ -1059,8 +1003,9 @@ void Foam::medialAxisMeshMover::minSmoothField scalarField edgeWeights(meshEdges.size()); scalarField invSumWeight(meshPoints.size()); - calculateEdgeWeights + meshRefinement::calculateEdgeWeights ( + mesh(), isMasterEdge, meshEdges, meshPoints, @@ -1075,7 +1020,7 @@ void Foam::medialAxisMeshMover::minSmoothField for (label iter = 0; iter < nSmoothDisp; iter++) { scalarField average(pp.nPoints()); - weightedSum + meshRefinement::weightedSum ( mesh(), isMasterEdge, @@ -1534,8 +1479,9 @@ void Foam::medialAxisMeshMover::smoothLambdaMuDisplacement // Calculate inverse sum of weights scalarField edgeWeights(meshEdges.size()); scalarField invSumWeight(meshPoints.size()); - calculateEdgeWeights + meshRefinement::calculateEdgeWeights ( + mesh(), isMasterEdge, meshEdges, meshPoints, @@ -1554,7 +1500,7 @@ void Foam::medialAxisMeshMover::smoothLambdaMuDisplacement for (label iter = 0; iter < nSmoothDisp; iter++) { - weightedSum + meshRefinement::weightedSum ( mesh(), isMasterEdge, @@ -1575,7 +1521,7 @@ void Foam::medialAxisMeshMover::smoothLambdaMuDisplacement } } - weightedSum + meshRefinement::weightedSum ( mesh(), isMasterEdge, diff --git a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.H b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.H index b402ecba25..1c270de274 100644 --- a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.H +++ b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.H @@ -105,31 +105,6 @@ class medialAxisMeshMover const labelList& ); - //- Weighted sum (over all subset of mesh points) by - // summing contribution from (master) edges - template - static void weightedSum - ( - const polyMesh& mesh, - const PackedBoolList& isMasterEdge, - const labelList& meshEdges, - const labelList& meshPoints, - const edgeList& edges, - const scalarField& edgeWeights, - const Field& data, - Field& sum - ); - - void calculateEdgeWeights - ( - const PackedBoolList& isMasterEdge, - const labelList& meshEdges, - const labelList& meshPoints, - const edgeList& edges, - scalarField& edgeWeights, - scalarField& invSumWeight - ) const; - // Calculation of medial axis information @@ -314,12 +289,6 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -#ifdef NoRepository -# include "medialAxisMeshMoverTemplates.C" -#endif - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - #endif // ************************************************************************* // diff --git a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMoverTemplates.C b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMoverTemplates.C deleted file mode 100644 index edd053f63f..0000000000 --- a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMoverTemplates.C +++ /dev/null @@ -1,93 +0,0 @@ -/*---------------------------------------------------------------------------*\ - ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox - \\ / O peration | - \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation - \\/ 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 . - -\*---------------------------------------------------------------------------*/ - -#include "medialAxisMeshMover.H" -#include "syncTools.H" - -// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // - -template -void Foam::medialAxisMeshMover::weightedSum -( - const polyMesh& mesh, - const PackedBoolList& isMasterEdge, - const labelList& meshEdges, - const labelList& meshPoints, - const edgeList& edges, - const scalarField& edgeWeights, - const Field& pointData, - Field& sum -) -{ - if - ( - mesh.nEdges() != isMasterEdge.size() - || edges.size() != meshEdges.size() - || edges.size() != edgeWeights.size() - || meshPoints.size() != pointData.size() - ) - { - FatalErrorIn("medialAxisMeshMover::weightedSum(..)") - << "Inconsistent sizes for edge or point data:" - << " meshEdges:" << meshEdges.size() - << " isMasterEdge:" << isMasterEdge.size() - << " edgeWeights:" << edgeWeights.size() - << " edges:" << edges.size() - << " pointData:" << pointData.size() - << " meshPoints:" << meshPoints.size() - << abort(FatalError); - } - - sum.setSize(meshPoints.size()); - sum = pTraits::zero; - - forAll(edges, edgeI) - { - if (isMasterEdge.get(meshEdges[edgeI]) == 1) - { - const edge& e = edges[edgeI]; - - scalar eWeight = edgeWeights[edgeI]; - - label v0 = e[0]; - label v1 = e[1]; - - sum[v0] += eWeight*pointData[v1]; - sum[v1] += eWeight*pointData[v0]; - } - } - - syncTools::syncPointList - ( - mesh, - meshPoints, - sum, - plusEqOp(), - pTraits::zero // null value - ); -} - - -// ************************************************************************* // diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.C b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.C index d93a5ca6d3..8d41eebd25 100644 --- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.C +++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.C @@ -1949,6 +1949,64 @@ void Foam::meshRefinement::checkCoupledFaceZones(const polyMesh& mesh) } +void Foam::meshRefinement::calculateEdgeWeights +( + const polyMesh& mesh, + const PackedBoolList& isMasterEdge, + const labelList& meshEdges, + const labelList& meshPoints, + const edgeList& edges, + scalarField& edgeWeights, + scalarField& invSumWeight +) +{ + const pointField& pts = mesh.points(); + + // Calculate edgeWeights and inverse sum of edge weights + edgeWeights.setSize(meshEdges.size()); + invSumWeight.setSize(meshPoints.size()); + + forAll(edges, edgeI) + { + const edge& e = edges[edgeI]; + scalar eMag = max + ( + VSMALL, + mag + ( + pts[meshPoints[e[1]]] + - pts[meshPoints[e[0]]] + ) + ); + edgeWeights[edgeI] = 1.0/eMag; + } + + // Sum per point all edge weights + weightedSum + ( + mesh, + isMasterEdge, + meshEdges, + meshPoints, + edges, + edgeWeights, + scalarField(meshPoints.size(), 1.0), // data + invSumWeight + ); + + // Inplace invert + forAll(invSumWeight, pointI) + { + scalar w = invSumWeight[pointI]; + + if (w > 0.0) + { + invSumWeight[pointI] = 1.0/w; + } + } +} + + Foam::label Foam::meshRefinement::appendPatch ( fvMesh& mesh, @@ -2167,9 +2225,7 @@ Foam::label Foam::meshRefinement::addMeshedPatch // ); // Store - label sz = meshedPatches_.size(); - meshedPatches_.setSize(sz+1); - meshedPatches_[sz] = name; + meshedPatches_.append(name); return patchI; } diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.H b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.H index f8d6258ba7..57114bd740 100644 --- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.H +++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.H @@ -729,6 +729,33 @@ public: //- Helper function: check that face zones are synced static void checkCoupledFaceZones(const polyMesh&); + //- Helper: calculate edge weights (1/length) + static void calculateEdgeWeights + ( + const polyMesh& mesh, + const PackedBoolList& isMasterEdge, + const labelList& meshEdges, + const labelList& meshPoints, + const edgeList& edges, + scalarField& edgeWeights, + scalarField& invSumWeight + ); + + //- Helper: weighted sum (over all subset of mesh points) by + // summing contribution from (master) edges + template + static void weightedSum + ( + const polyMesh& mesh, + const PackedBoolList& isMasterEdge, + const labelList& meshEdges, + const labelList& meshPoints, + const edgeList& edges, + const scalarField& edgeWeights, + const Field& data, + Field& sum + ); + // Refinement diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementRefine.C b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementRefine.C index a7084e9fc6..701c0dbbc0 100644 --- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementRefine.C +++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementRefine.C @@ -458,10 +458,18 @@ void Foam::meshRefinement::markFeatureCellLevel // Database to pass into trackedParticle::move trackedParticle::trackingData td(startPointCloud, maxFeatureLevel); + // Track all particles to their end position (= starting feature point) // Note that the particle might have started on a different processor // so this will transfer across nicely until we can start tracking proper. scalar maxTrackLen = 2.0*mesh_.bounds().mag(); + + if (debug&meshRefinement::FEATURESEEDS) + { + Pout<< "Tracking " << startPointCloud.size() + << " particles over distance " << maxTrackLen + << " to find the starting cell" << endl; + } startPointCloud.move(td, maxTrackLen); @@ -485,6 +493,11 @@ void Foam::meshRefinement::markFeatureCellLevel IDLList() ); + if (debug&meshRefinement::FEATURESEEDS) + { + Pout<< "Constructing cloud for cell marking" << endl; + } + forAllIter(Cloud, startPointCloud, iter) { const trackedParticle& startTp = iter(); @@ -532,11 +545,14 @@ void Foam::meshRefinement::markFeatureCellLevel while (true) { // Track all particles to their end position. + if (debug&meshRefinement::FEATURESEEDS) + { + Pout<< "Tracking " << cloud.size() + << " particles over distance " << maxTrackLen + << " to mark cells" << endl; + } cloud.move(td, maxTrackLen); - - label nParticles = 0; - // Make particle follow edge. forAllIter(Cloud, cloud, iter) { @@ -578,15 +594,15 @@ void Foam::meshRefinement::markFeatureCellLevel // seeded. Delete particle. cloud.deleteParticle(tp); } - else - { - // Keep particle - nParticles++; - } } - reduce(nParticles, sumOp