Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev

This commit is contained in:
andy
2011-12-12 09:47:50 +00:00
115 changed files with 9225 additions and 951 deletions

View File

@ -45,6 +45,7 @@ primitives/Tensor/lists/symmTensorList.C
primitives/Tensor/lists/tensorList.C
primitives/Vector/complexVector/complexVector.C
primitives/Vector/floatVector/floatVector.C
primitives/Vector/labelVector/labelVector.C
primitives/Vector/vector/vector.C
primitives/Vector/lists/vectorList.C
@ -501,7 +502,6 @@ meshes/treeBoundBox/treeBoundBox.C
meshTools = meshes/meshTools
$(meshTools)/matchPoints.C
$(meshTools)/mergePoints.C
fields/UniformDimensionedFields/uniformDimensionedFields.C
fields/cloud/cloud.C

View File

@ -26,6 +26,7 @@ Description
matrix is reduced. The algorithm uses a simple search through the
neighbour list
See http://en.wikipedia.org/wiki/Cuthill-McKee_algorithm
\*---------------------------------------------------------------------------*/

View File

@ -46,6 +46,8 @@ namespace Foam
//- Renumbers the addressing to reduce the band of the matrix.
// The algorithm uses a simple search through the neighbour list
// Returns the order in which the cells need to be visited (i.e. ordered to
// original)
labelList bandCompression(const labelListList& addressing);
} // End namespace Foam

View File

@ -23,25 +23,25 @@ License
\*---------------------------------------------------------------------------*/
#include "mergePoints.H"
#include "SortableList.H"
#include "ListOps.H"
#include "point.H"
#include "Field.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::mergePoints
template<class Type>
Foam::label Foam::mergePoints
(
const UList<point>& points,
const UList<Type>& points,
const scalar mergeTol,
const bool verbose,
labelList& pointMap,
List<point>& newPoints,
const point& origin
const Type& origin
)
{
point compareOrigin = origin;
Type compareOrigin = origin;
if (origin == point(VGREAT, VGREAT, VGREAT))
if (origin == Type::max)
{
if (points.size())
{
@ -53,12 +53,9 @@ bool Foam::mergePoints
pointMap.setSize(points.size());
pointMap = -1;
// Storage for merged points
newPoints.setSize(points.size());
if (points.empty())
{
return false;
return points.size();
}
// We're comparing distance squared to origin first.
@ -70,33 +67,56 @@ bool Foam::mergePoints
// x^2+y^2+z^2 + 2*mergeTol*(x+z+y) + mergeTol^2*...
// so the difference will be 2*mergeTol*(x+y+z)
const scalar mergeTolSqr = sqr(mergeTol);
const scalar mergeTolSqr = Foam::sqr(scalar(mergeTol));
// Sort points by magSqr
const pointField d(points - compareOrigin);
SortableList<scalar> sortedMagSqr(magSqr(d));
scalarField sortedTol(points.size());
forAll(sortedMagSqr.indices(), sortI)
const Field<Type> d(points - compareOrigin);
List<scalar> magSqrD(d.size());
forAll(d, pointI)
{
const point& pt = d[sortedMagSqr.indices()[sortI]];
magSqrD[pointI] = magSqr(d[pointI]);
}
labelList order;
sortedOrder(magSqrD, order);
Field<scalar> sortedTol(points.size());
forAll(order, sortI)
{
label pointI = order[sortI];
// Convert to scalar precision
const point pt
(
scalar(d[pointI].x()),
scalar(d[pointI].y()),
scalar(d[pointI].z())
);
sortedTol[sortI] = 2*mergeTol*(mag(pt.x())+mag(pt.y())+mag(pt.z()));
}
bool hasMerged = false;
label newPointI = 0;
// Handle 0th point separately (is always unique)
label pointI = sortedMagSqr.indices()[0];
pointMap[pointI] = newPointI;
newPoints[newPointI++] = points[pointI];
label pointI = order[0];
pointMap[pointI] = newPointI++;
for (label sortI = 1; sortI < sortedMagSqr.size(); sortI++)
for (label sortI = 1; sortI < order.size(); sortI++)
{
// Get original point index
label pointI = sortedMagSqr.indices()[sortI];
label pointI = order[sortI];
const scalar mag2 = magSqrD[order[sortI]];
// Convert to scalar precision
const point pt
(
scalar(points[pointI].x()),
scalar(points[pointI].y()),
scalar(points[pointI].z())
);
// Compare to previous points to find equal one.
label equalPointI = -1;
@ -105,17 +125,19 @@ bool Foam::mergePoints
(
label prevSortI = sortI - 1;
prevSortI >= 0
&& mag
(
sortedMagSqr[prevSortI]
- sortedMagSqr[sortI]
) <= sortedTol[sortI];
&& (mag(magSqrD[order[prevSortI]] - mag2) <= sortedTol[sortI]);
prevSortI--
)
{
label prevPointI = sortedMagSqr.indices()[prevSortI];
label prevPointI = order[prevSortI];
const point prevPt
(
scalar(points[prevPointI].x()),
scalar(points[prevPointI].y()),
scalar(points[prevPointI].z())
);
if (magSqr(points[pointI] - points[prevPointI]) <= mergeTolSqr)
if (magSqr(pt - prevPt) <= mergeTolSqr)
{
// Found match.
equalPointI = prevPointI;
@ -130,8 +152,6 @@ bool Foam::mergePoints
// Same coordinate as equalPointI. Map to same new point.
pointMap[pointI] = pointMap[equalPointI];
hasMerged = true;
if (verbose)
{
Pout<< "Foam::mergePoints : Merging points "
@ -144,14 +164,41 @@ bool Foam::mergePoints
else
{
// Differs. Store new point.
pointMap[pointI] = newPointI;
newPoints[newPointI++] = points[pointI];
pointMap[pointI] = newPointI++;
}
}
newPoints.setSize(newPointI);
return newPointI;
}
return hasMerged;
template<class Type>
bool Foam::mergePoints
(
const UList<Type>& points,
const scalar mergeTol,
const bool verbose,
labelList& pointMap,
List<Type>& newPoints,
const Type& origin = Type::zero
)
{
label nUnique = mergePoints
(
points,
mergeTol,
verbose,
pointMap,
origin
);
newPoints.setSize(nUnique);
forAll(pointMap, pointI)
{
newPoints[pointMap[pointI]] = points[pointI];
}
return (nUnique != points.size());
}

View File

@ -32,8 +32,8 @@ SourceFiles
#ifndef mergePoints_H
#define mergePoints_H
#include "scalarField.H"
#include "pointField.H"
#include "scalar.H"
#include "labelList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -44,23 +44,41 @@ namespace Foam
Function mergePoints Declaration
\*---------------------------------------------------------------------------*/
//- Sort & merge points. All points closer than/equal mergeTol get merged.
// Outputs the new unique points and a map from old to new. Returns
// true if anything merged, false otherwise.
bool mergePoints
//- Sorts and merges points. All points closer than/equal mergeTol get merged.
// Returns the number of unique points and a map from old to new.
template<class Type>
label mergePoints
(
const UList<point>& points,
const UList<Type>& points,
const scalar mergeTol,
const bool verbose,
labelList& pointMap,
List<point>& newPoints,
const point& origin = point::zero
const Type& origin = Type::zero
);
//- Sorts and merges points. Determines new points. Returns true if anything
// merged (though newPoints still sorted even if not merged).
template<class Type>
bool mergePoints
(
const UList<Type>& points,
const scalar mergeTol,
const bool verbose,
labelList& pointMap,
List<Type>& newPoints,
const Type& origin = Type::zero
);
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "mergePoints.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -27,14 +27,22 @@ License
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::globalIndex::globalIndex(const label localSize, const int tag)
Foam::globalIndex::globalIndex
(
const label localSize,
const int tag,
const bool parallel
)
:
offsets_(Pstream::nProcs()+1)
{
labelList localSizes(Pstream::nProcs());
labelList localSizes(Pstream::nProcs(), 0);
localSizes[Pstream::myProcNo()] = localSize;
Pstream::gatherList(localSizes, tag);
Pstream::scatterList(localSizes, tag);
if (parallel)
{
Pstream::gatherList(localSizes, tag);
Pstream::scatterList(localSizes, tag);
}
label offset = 0;
offsets_[0] = 0;

View File

@ -72,7 +72,12 @@ public:
// Constructors
//- Construct from local max size
globalIndex(const label localSize, const int tag = Pstream::msgType());
globalIndex
(
const label localSize,
const int tag = Pstream::msgType(),
const bool parallel = true // use parallel comms
);
//- Construct from Istream
globalIndex(Istream& is);

View File

@ -159,7 +159,7 @@ public:
labelList findIndices
(
const keyType&,
const bool usePatchGroups = false
const bool usePatchGroups = true
) const;
//- Return patch index for the first match, return -1 if not found
@ -184,7 +184,7 @@ public:
(
const UList<wordRe>& patchNames,
const bool warnNotFound = true,
const bool usePatchGroups = false
const bool usePatchGroups = true
) const;
//- Check whether all procs have all patches and in same order. Return

View File

@ -0,0 +1,70 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Description
Vector of floats.
\*---------------------------------------------------------------------------*/
#include "floatVector.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
template<>
const char* const floatVector::typeName = "floatVector";
template<>
const char* floatVector::componentNames[] = {"x", "y", "z"};
template<>
const floatVector floatVector::zero(0, 0, 0);
template<>
const floatVector floatVector::one(1, 1, 1);
template<>
const floatVector floatVector::max
(
floatScalarVGREAT,
floatScalarVGREAT,
floatScalarVGREAT
);
template<>
const floatVector floatVector::min
(
-floatScalarVGREAT,
-floatScalarVGREAT,
-floatScalarVGREAT
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

View File

@ -0,0 +1,64 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Typedef
Foam::floatVector
Description
A float version of vector
SourceFiles
floatVector.C
\*---------------------------------------------------------------------------*/
#ifndef floatVector_H
#define floatVector_H
#include "Vector.H"
#include "contiguous.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
typedef Vector<float> floatVector;
//- Data associated with floatVector type are contiguous
template<>
inline bool contiguous<floatVector>() {return true;}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -36,7 +36,6 @@ License
#include "removeCells.H"
#include "polyModifyFace.H"
#include "polyRemovePoint.H"
#include "mergePoints.H"
#include "mapDistributePolyMesh.H"
#include "surfaceFields.H"
#include "syncTools.H"

View File

@ -1859,24 +1859,22 @@ Foam::Map<Foam::label> Foam::polyMeshAdder::findSharedPoints
);
labelList toMergedPoints;
pointField mergedPoints;
bool hasMerged = Foam::mergePoints
label nUnique = Foam::mergePoints
(
connectedPoints,
mergeDist,
false,
toMergedPoints,
mergedPoints
toMergedPoints
);
if (hasMerged)
if (nUnique < connectedPoints.size())
{
// Invert toMergedPoints
const labelListList mergeSets
(
invertOneToMany
(
mergedPoints.size(),
nUnique,
toMergedPoints
)
);
@ -1919,8 +1917,7 @@ Foam::Map<Foam::label> Foam::polyMeshAdder::findSharedPoints
//- Old: geometric merging. Causes problems for two close shared points.
//labelList sharedToMerged;
//pointField mergedPoints;
//bool hasMerged = Foam::mergePoints
//label nUnique = Foam::mergePoints
//(
// pointField
// (
@ -1929,8 +1926,7 @@ Foam::Map<Foam::label> Foam::polyMeshAdder::findSharedPoints
// ),
// mergeDist,
// false,
// sharedToMerged,
// mergedPoints
// sharedToMerged
//);
//
//// Find out which sets of points get merged and create a map from
@ -1938,7 +1934,7 @@ Foam::Map<Foam::label> Foam::polyMeshAdder::findSharedPoints
//
//Map<label> pointToMaster(10*sharedToMerged.size());
//
//if (hasMerged)
//if (nUnique < sharedPointLabels.size())
//{
// labelListList mergeSets
// (

View File

@ -34,7 +34,6 @@ License
#include "polyModifyFace.H"
#include "polyAddCell.H"
#include "globalIndex.H"
#include "dummyTransform.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

View File

@ -38,7 +38,9 @@ Description
code
#{
operator==(min(10, 0.1*this->db().time().value()));
this->refValue() = min(10, 0.1*this->db().time().value());
this->refGrad() = vector::zero;
this->valueFraction() = 1.0;
#};
//codeInclude
@ -58,13 +60,15 @@ Description
which would have a corresponding entry
\verbatim
rampedMixed
{
code
#{
operator==(min(10, 0.1*this->db().time().value()));
#};
}
rampedMixed
{
code
#{
this->refValue() = min(10, 0.1*this->db().time().value());
this->refGrad() = vector::zero;
this->valueFraction() = 1.0;
#};
}
\endverbatim
SeeAlso

View File

@ -47,6 +47,7 @@ Description
#include "IOmanip.H"
#include "globalIndex.H"
#include "DynamicField.H"
#include "PatchTools.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@ -369,7 +370,16 @@ void Foam::autoLayerDriver::handleFeatureAngle
autoPtr<OFstream> str;
if (debug)
{
str.reset(new OFstream(mesh.time().path()/"featureEdges.obj"));
str.reset
(
new OFstream
(
mesh.time().path()
/ "featureEdges_"
+ meshRefiner_.timeName()
+ ".obj"
)
);
Info<< "Writing feature edges to " << str().name() << endl;
}
@ -1220,49 +1230,19 @@ void Foam::autoLayerDriver::getPatchDisplacement
const vectorField& faceNormals = pp.faceNormals();
const labelListList& pointFaces = pp.pointFaces();
const pointField& localPoints = pp.localPoints();
const labelList& meshPoints = pp.meshPoints();
// Determine pointNormal
// ~~~~~~~~~~~~~~~~~~~~~
pointField pointNormals(pp.nPoints(), vector::zero);
{
labelList nPointFaces(pp.nPoints(), 0);
forAll(faceNormals, faceI)
{
const face& f = pp.localFaces()[faceI];
forAll(f, fp)
{
pointNormals[f[fp]] += faceNormals[faceI];
nPointFaces[f[fp]] ++;
}
}
syncTools::syncPointList
pointField pointNormals
(
PatchTools::pointNormals
(
mesh,
meshPoints,
pointNormals,
plusEqOp<vector>(),
vector::zero // null value
);
syncTools::syncPointList
(
mesh,
meshPoints,
nPointFaces,
plusEqOp<label>(),
label(0) // null value
);
forAll(pointNormals, i)
{
pointNormals[i] /= nPointFaces[i];
}
}
pp,
pp.addressing()
)
);
// Determine local length scale on patch
@ -2266,7 +2246,11 @@ void Foam::autoLayerDriver::addLayers
{
const_cast<Time&>(mesh.time())++;
Info<< "Writing baffled mesh to " << meshRefiner_.timeName() << endl;
mesh.write();
meshRefiner_.write
(
debug,
mesh.time().path()/meshRefiner_.timeName()
);
}
@ -2732,7 +2716,7 @@ void Foam::autoLayerDriver::addLayers
{
dumpDisplacement
(
mesh.time().path()/"layer",
mesh.time().path()/"layer_" + meshRefiner_.timeName(),
pp(),
patchDisp,
extrudeStatus
@ -2744,7 +2728,11 @@ void Foam::autoLayerDriver::addLayers
// See comment in autoSnapDriver why we should not remove meshPhi
// using mesh.clearPout().
mesh.write();
meshRefiner_.write
(
debug,
mesh.time().path()/meshRefiner_.timeName()
);
}

View File

@ -33,6 +33,9 @@ Description
#include "motionSmoother.H"
#include "pointData.H"
#include "PointEdgeWave.H"
#include "OFstream.H"
#include "meshTools.H"
#include "PatchTools.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
@ -334,11 +337,25 @@ bool Foam::autoLayerDriver::isMaxEdge
return false;
}
v0 /= magV0;
v1 /= magV1;
// Test angle.
if ((v0 & v1) < minCos)
//- Detect based on vector to nearest point differing for both endpoints
//v0 /= magV0;
//v1 /= magV1;
//
//// Test angle.
//if ((v0 & v1) < minCos)
//{
// return true;
//}
//else
//{
// return false;
//}
//- Detect based on extrusion vector differing for both endpoints
// the idea is that e.g. a sawtooth wall can still be extruded
// successfully as long as it is done all to the same direction.
if ((pointWallDist[e[0]].v() & pointWallDist[e[1]].v()) < minCos)
{
return true;
}
@ -670,7 +687,6 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
const pointField& points = mesh.points();
const indirectPrimitivePatch& pp = meshMover.patch();
const vectorField& faceNormals = pp.faceNormals();
const labelList& meshPoints = pp.meshPoints();
// Predetermine mesh edges
@ -700,44 +716,15 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
// Determine pointNormal
// ~~~~~~~~~~~~~~~~~~~~~
pointField pointNormals(pp.nPoints(), vector::zero);
{
labelList nPointFaces(pp.nPoints(), 0);
forAll(faceNormals, faceI)
{
const face& f = pp.localFaces()[faceI];
forAll(f, fp)
{
pointNormals[f[fp]] += faceNormals[faceI];
nPointFaces[f[fp]] ++;
}
}
syncTools::syncPointList
pointField pointNormals
(
PatchTools::pointNormals
(
mesh,
meshPoints,
pointNormals,
plusEqOp<vector>(),
vector::zero // null value
);
syncTools::syncPointList
(
mesh,
meshPoints,
nPointFaces,
plusEqOp<label>(),
label(0) // null value
);
forAll(pointNormals, i)
{
pointNormals[i] /= nPointFaces[i];
}
}
pp,
pp.addressing()
)
);
// Smooth patch normal vectors
smoothPatchNormals
@ -1012,6 +999,26 @@ void Foam::autoLayerDriver::shrinkMeshMedialDistance
label numThicknessRatioExclude = 0;
// reduce thickness where thickness/medial axis distance large
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
autoPtr<OFstream> str;
label vertI = 0;
if (debug)
{
str.reset
(
new OFstream
(
mesh.time().path()
/ "thicknessRatioExcludePoints_"
+ meshRefiner_.timeName()
+ ".obj"
)
);
Info<< "Writing points with too large a extrusion distance to "
<< str().name() << endl;
}
forAll(meshPoints, patchPointI)
{
if (extrudeStatus[patchPointI] != NOEXTRUDE)
@ -1025,6 +1032,20 @@ void Foam::autoLayerDriver::shrinkMeshMedialDistance
if (thicknessRatio > maxThicknessToMedialRatio)
{
// Truncate thickness.
if (debug)
{
Pout<< "truncating displacement at "
<< mesh.points()[pointI]
<< " from " << thickness[patchPointI]
<< " to "
<< 0.5
*(
minThickness[patchPointI]
+thickness[patchPointI]
)
<< endl;
}
thickness[patchPointI] =
0.5*(minThickness[patchPointI]+thickness[patchPointI]);
@ -1033,6 +1054,16 @@ void Foam::autoLayerDriver::shrinkMeshMedialDistance
* patchDisp[patchPointI]
/ (mag(patchDisp[patchPointI]) + VSMALL);
numThicknessRatioExclude++;
if (str.valid())
{
const point& pt = mesh.points()[pointI];
meshTools::writeOBJ(str(), pt);
vertI++;
meshTools::writeOBJ(str(), pt+patchDisp[patchPointI]);
vertI++;
str()<< "l " << vertI-1 << ' ' << vertI << nl;
}
}
}
}
@ -1121,6 +1152,31 @@ void Foam::autoLayerDriver::shrinkMeshMedialDistance
*dispVec[pointI];
}
if (debug)
{
const_cast<Time&>(mesh.time())++;
Info<< "Writing wanted-displacement mesh (possibly illegal) to "
<< meshRefiner_.timeName() << endl;
pointField oldPoints(mesh.points());
meshMover.movePoints
(
(
mesh.points()
+ (
meshMover.scale().internalField()
* displacement.internalField()
)
)()
);
meshRefiner_.write
(
debug,
mesh.time().path()/meshRefiner_.timeName()
);
meshMover.movePoints(oldPoints);
}
// Current faces to check. Gets modified in meshMover.scaleMesh
labelList checkFaces(identity(mesh.nFaces()));

View File

@ -63,27 +63,26 @@ Foam::label Foam::autoSnapDriver::getCollocatedPoints
)
{
labelList pointMap;
pointField newPoints;
bool hasMerged = mergePoints
label nUnique = mergePoints
(
points, // points
tol, // mergeTol
false, // verbose
pointMap,
newPoints
pointMap
);
bool hasMerged = (nUnique < points.size());
if (!returnReduce(hasMerged, orOp<bool>()))
{
return 0;
}
// Determine which newPoints are referenced more than once
// Determine which merged points are referenced more than once
label nCollocated = 0;
// Per old point the newPoint. Or -1 (not set yet) or -2 (already seen
// twice)
labelList firstOldPoint(newPoints.size(), -1);
labelList firstOldPoint(nUnique, -1);
forAll(pointMap, oldPointI)
{
label newPointI = pointMap[oldPointI];

View File

@ -31,7 +31,6 @@ License
#include "pointSet.H"
#include "faceSet.H"
#include "indirectPrimitivePatch.H"
#include "OFstream.H"
#include "cellSet.H"
#include "searchableSurfaces.H"
#include "polyMeshGeometry.H"

View File

@ -2096,7 +2096,8 @@ Foam::triSurface Foam::triSurfaceTools::mergePoints
(
newTriangles,
surf.patches(),
newPoints
newPoints,
true //reuse storage
);
}
else

View File

@ -96,6 +96,7 @@ Foam::labelList Foam::decompositionMethod::decompose
mesh,
fineToCoarse,
coarsePoints.size(),
true, // use global cell labels
coarseCellCells
);
@ -158,6 +159,7 @@ void Foam::decompositionMethod::calcCellCells
const polyMesh& mesh,
const labelList& agglom,
const label nCoarse,
const bool parallel,
CompactListList<label>& cellCells
)
{
@ -169,7 +171,7 @@ void Foam::decompositionMethod::calcCellCells
// Create global cell numbers
// ~~~~~~~~~~~~~~~~~~~~~~~~~~
globalIndex globalAgglom(nCoarse);
globalIndex globalAgglom(nCoarse, Pstream::msgType(), parallel);
// Get agglomerate owner on other side of coupled faces
@ -181,7 +183,7 @@ void Foam::decompositionMethod::calcCellCells
{
const polyPatch& pp = patches[patchI];
if (pp.coupled())
if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
{
label faceI = pp.start();
label bFaceI = pp.start() - mesh.nInternalFaces();
@ -203,77 +205,6 @@ void Foam::decompositionMethod::calcCellCells
syncTools::swapBoundaryFaceList(mesh, globalNeighbour);
//// Determine the cellCells (in agglomeration numbering) on coupled faces
//// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//
//labelListList globalCellCells(mesh.nFaces()-mesh.nInternalFaces());
//
//// Current set of face neighbours for the current cell
//labelHashSet cCells;
//
//forAll(patches, patchI)
//{
// const polyPatch& pp = patches[patchI];
//
// if (pp.coupled())
// {
// label faceI = pp.start();
// label bFaceI = pp.start() - mesh.nInternalFaces();
//
// forAll(pp, i)
// {
// label cellI = faceOwner[faceI];
// label globalCellI = globalAgglom.toGlobal(agglom[cellI]);
//
// // First check if agglomerated across coupled patches at all
// // so we don't use memory if not needed
// if (globalNeighbour[bFaceI] == globalCellI)
// {
// cCells.clear();
//
// const cell& cFaces = mesh.cells()[cellI];
//
// forAll(cFaces, i)
// {
// if (mesh.isInternalFace(cFaces[i]))
// {
// label otherCellI = faceOwner[cFaces[i]];
// if (otherCellI == cellI)
// {
// otherCellI = faceNeighbour[cFaces[i]];
// }
//
// cCells.insert
// (
// globalAgglom.toGlobal
// (
// agglom[otherCellI]
// )
// );
// }
// }
// globalCellCells[bFaceI] = cCells.toc();
// }
//
// bFaceI++;
// faceI++;
// }
// }
//}
//
//// Get the cell on the other side of coupled patches
//syncTools::syncBoundaryFaceList
//(
// mesh,
// globalCellCells,
// eqOp<labelList>(),
// dummyTransform()
//);
// Count number of faces (internal + coupled)
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
@ -293,7 +224,7 @@ void Foam::decompositionMethod::calcCellCells
{
const polyPatch& pp = patches[patchI];
if (pp.coupled())
if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
{
label faceI = pp.start();
label bFaceI = pp.start()-mesh.nInternalFaces();
@ -302,18 +233,6 @@ void Foam::decompositionMethod::calcCellCells
{
label own = agglom[faceOwner[faceI]];
//const labelList& cCells = globalCellCells[bFaceI];
//
//forAll(cCells, i)
//{
// label globalNei = cCells[i];
//
// // Allow only processor-local agglomeration
// if (globalAgglom.isLocal(globalNei))
// {
// nFacesPerCell[own]++;
// }
//}
label globalNei = globalNeighbour[bFaceI];
if
(
@ -365,18 +284,6 @@ void Foam::decompositionMethod::calcCellCells
{
label own = agglom[faceOwner[faceI]];
//const labelList& cCells = globalCellCells[bFaceI];
//
//forAll(cCells, i)
//{
// label globalNei = cCells[i];
//
// // Allow only processor-local agglomeration
// if (globalAgglom.isLocal(globalNei))
// {
// m[offsets[own] + nFacesPerCell[own]++] = globalNei;
// }
//}
label globalNei = globalNeighbour[bFaceI];
if
(

View File

@ -57,15 +57,6 @@ protected:
label nProcessors_;
//- Helper: determine (global) cellCells from mesh agglomeration.
static void calcCellCells
(
const polyMesh& mesh,
const labelList& agglom,
const label nCoarse,
CompactListList<label>& cellCells
);
private:
// Private Member Functions
@ -224,6 +215,24 @@ public:
const pointField& cc
);
// Other
//- Helper: determine (local or global) cellCells from mesh
// agglomeration.
// local : connections are in local indices. Coupled across
// cyclics but not processor patches.
// global : connections are in global indices. Coupled across
// cyclics and processor patches.
static void calcCellCells
(
const polyMesh& mesh,
const labelList& agglom,
const label nCoarse,
const bool global,
CompactListList<label>& cellCells
);
};

View File

@ -383,7 +383,7 @@ Foam::labelList Foam::multiLevelDecomp::decompose
)
{
CompactListList<label> cellCells;
calcCellCells(mesh, identity(cc.size()), cc.size(), cellCells);
calcCellCells(mesh, identity(cc.size()), cc.size(), true, cellCells);
labelField finalDecomp(cc.size(), 0);
labelList cellMap(identity(cc.size()));

View File

@ -755,7 +755,14 @@ Foam::labelList Foam::ptscotchDecomp::decompose
CompactListList<label> cellCells;
calcCellCells(mesh, identity(mesh.nCells()), mesh.nCells(), cellCells);
calcCellCells
(
mesh,
identity(mesh.nCells()),
mesh.nCells(),
true,
cellCells
);
// Decompose using default weights
List<int> finalDecomp;
@ -807,7 +814,14 @@ Foam::labelList Foam::ptscotchDecomp::decompose
// adjncy : contains neighbours (= edges in graph)
// xadj(celli) : start of information in adjncy for celli
CompactListList<label> cellCells;
calcCellCells(mesh, agglom, agglomPoints.size(), cellCells);
calcCellCells
(
mesh,
agglom,
agglomPoints.size(),
true,
cellCells
);
// Decompose using weights
List<int> finalDecomp;

View File

@ -589,7 +589,14 @@ Foam::labelList Foam::scotchDecomp::decompose
// Calculate local or global (if Pstream::parRun()) connectivity
CompactListList<label> cellCells;
calcCellCells(mesh, identity(mesh.nCells()), mesh.nCells(), cellCells);
calcCellCells
(
mesh,
identity(mesh.nCells()),
mesh.nCells(),
true,
cellCells
);
// Decompose using default weights
List<int> finalDecomp;
@ -634,7 +641,14 @@ Foam::labelList Foam::scotchDecomp::decompose
// Calculate local or global (if Pstream::parRun()) connectivity
CompactListList<label> cellCells;
calcCellCells(mesh, agglom, agglomPoints.size(), cellCells);
calcCellCells
(
mesh,
agglom,
agglomPoints.size(),
true,
cellCells
);
// Decompose using weights
List<int> finalDecomp;

View File

@ -0,0 +1,112 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "CuthillMcKeeRenumber.H"
#include "addToRunTimeSelectionTable.H"
#include "bandCompression.H"
#include "decompositionMethod.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(CuthillMcKeeRenumber, 0);
addToRunTimeSelectionTable
(
renumberMethod,
CuthillMcKeeRenumber,
dictionary
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::CuthillMcKeeRenumber::CuthillMcKeeRenumber(const dictionary& renumberDict)
:
renumberMethod(renumberDict),
reverse_
(
renumberDict.found(typeName + "Coeffs")
? Switch(renumberDict.subDict(typeName + "Coeffs").lookup("reverse"))
: Switch(false)
)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::CuthillMcKeeRenumber::renumber
(
const polyMesh& mesh,
const pointField& points
)
{
CompactListList<label> cellCells;
decompositionMethod::calcCellCells
(
mesh,
identity(mesh.nCells()),
mesh.nCells(),
false, // local only
cellCells
);
labelList orderedToOld = bandCompression(cellCells());
if (reverse_)
{
for (label i = 0; i < orderedToOld.size()/2; i++)
{
Swap(orderedToOld[i], orderedToOld[orderedToOld.size()-i-1]);
}
}
return invert(orderedToOld.size(), orderedToOld);
}
Foam::labelList Foam::CuthillMcKeeRenumber::renumber
(
const labelListList& cellCells,
const pointField& points
)
{
labelList orderedToOld = bandCompression(cellCells);
if (reverse_)
{
for (label i = 0; i < orderedToOld.size()/2; i++)
{
Swap(orderedToOld[i], orderedToOld[orderedToOld.size()-i-1]);
}
}
return invert(orderedToOld.size(), orderedToOld);
}
// ************************************************************************* //

View File

@ -0,0 +1,118 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Class
Foam::CuthillMcKeeRenumber
Description
Cuthill-McKee renumbering
SourceFiles
CuthillMcKeeRenumber.C
\*---------------------------------------------------------------------------*/
#ifndef CuthillMcKeeRenumber_H
#define CuthillMcKeeRenumber_H
#include "renumberMethod.H"
#include "Switch.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class CuthillMcKeeRenumber Declaration
\*---------------------------------------------------------------------------*/
class CuthillMcKeeRenumber
:
public renumberMethod
{
// Private data
const Switch reverse_;
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
void operator=(const CuthillMcKeeRenumber&);
CuthillMcKeeRenumber(const CuthillMcKeeRenumber&);
public:
//- Runtime type information
TypeName("CuthillMcKee");
// Constructors
//- Construct given the renumber dictionary
CuthillMcKeeRenumber(const dictionary& renumberDict);
//- Destructor
virtual ~CuthillMcKeeRenumber()
{}
// Member Functions
//- Return for every coordinate the wanted processor number.
// We need a polyMesh (to be able to load the file)
virtual labelList renumber(const pointField&)
{
notImplemented("CuthillMcKeeRenumber::renumber(const pointField&)");
return labelList(0);
}
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList renumber
(
const polyMesh& mesh,
const pointField& cc
);
//- Return for every cell the new cell label.
// The connectivity is equal to mesh.cellCells() except
// - the connections are across coupled patches
virtual labelList renumber
(
const labelListList& cellCells,
const pointField& cc
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,7 @@
renumberMethod/renumberMethod.C
manualRenumber/manualRenumber.C
CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
randomRenumber/randomRenumber.C
springRenumber/springRenumber.C
LIB = $(FOAM_LIBBIN)/librenumberMethods

View File

@ -0,0 +1,7 @@
EXE_INC = \
-I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
LIB_LIBS = \
-ldecompositionMethods \
-lmeshTools

View File

@ -0,0 +1,137 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "manualRenumber.H"
#include "addToRunTimeSelectionTable.H"
#include "IFstream.H"
#include "labelIOList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(manualRenumber, 0);
addToRunTimeSelectionTable
(
renumberMethod,
manualRenumber,
dictionary
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::manualRenumber::manualRenumber(const dictionary& renumberDict)
:
renumberMethod(renumberDict),
dataFile_
(
renumberDict.subDict(typeName+"Coeffs").lookup("dataFile")
)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::manualRenumber::renumber
(
const polyMesh& mesh,
const pointField& points
)
{
labelIOList oldToNew
(
IOobject
(
dataFile_,
mesh.facesInstance(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE,
false
)
);
// check if the final renumbering is OK
if (oldToNew.size() != points.size())
{
FatalErrorIn
(
"manualRenumber::renumber(const pointField&, const scalarField&)"
) << "Size of renumber list does not correspond "
<< "to the number of points. Size: "
<< oldToNew.size() << " Number of points: "
<< points.size()
<< ".\n" << "Manual renumbering data read from file "
<< dataFile_ << "." << endl
<< exit(FatalError);
}
// Invert to see if one to one
labelList newToOld(points.size(), -1);
forAll(oldToNew, i)
{
label newI = oldToNew[i];
if (newI < 0 || newI >= oldToNew.size())
{
FatalErrorIn
(
"manualRenumber::renumber(const pointField&"
", const scalarField&)"
) << "Renumbering is not one-to-one. Index "
<< i << " maps onto " << newI
<< ".\n" << "Manual renumbering data read from file "
<< dataFile_ << "." << endl
<< exit(FatalError);
}
if (newToOld[newI] == -1)
{
newToOld[newI] = i;
}
else
{
FatalErrorIn
(
"manualRenumber::renumber(const pointField&"
", const scalarField&)"
) << "Renumbering is not one-to-one. Both index "
<< newToOld[newI]
<< " and " << i << " map onto " << newI
<< ".\n" << "Manual renumbering data read from file "
<< dataFile_ << "." << endl
<< exit(FatalError);
}
}
return oldToNew;
}
// ************************************************************************* //

View File

@ -0,0 +1,125 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Class
Foam::manualRenumber
Description
Renumber given a cell-to-new cell association in a file
SourceFiles
manualRenumber.C
\*---------------------------------------------------------------------------*/
#ifndef manualRenumber_H
#define manualRenumber_H
#include "renumberMethod.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class manualRenumber Declaration
\*---------------------------------------------------------------------------*/
class manualRenumber
:
public renumberMethod
{
// Private data
const fileName dataFile_;
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
void operator=(const manualRenumber&);
manualRenumber(const manualRenumber&);
public:
//- Runtime type information
TypeName("manual");
// Constructors
//- Construct given the renumber dictionary
manualRenumber(const dictionary& renumberDict);
//- Destructor
virtual ~manualRenumber()
{}
// Member Functions
//- Return for every coordinate the wanted processor number.
// We need a polyMesh (to be able to load the file)
virtual labelList renumber(const pointField&)
{
notImplemented("manualRenumber::renumber(const pointField&)");
return labelList(0);
}
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList renumber
(
const polyMesh& mesh,
const pointField& cc
);
//- Return for every cell the new cell label.
// The connectivity is equal to mesh.cellCells() except
// - the connections are across coupled patches
virtual labelList renumber
(
const labelListList& cellCells,
const pointField& cc
)
{
notImplemented
(
"manualRenumber::renumber"
"(const labelListList&, const pointField&)"
);
return labelList(0);
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,96 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "randomRenumber.H"
#include "addToRunTimeSelectionTable.H"
#include "Random.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(randomRenumber, 0);
addToRunTimeSelectionTable
(
renumberMethod,
randomRenumber,
dictionary
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::randomRenumber::randomRenumber(const dictionary& renumberDict)
:
renumberMethod(renumberDict)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::randomRenumber::renumber
(
const pointField& points
)
{
Random rndGen(0);
labelList oldToNew(identity(points.size()));
for (label iter = 0; iter < 10; iter++)
{
forAll(oldToNew, i)
{
label j = rndGen.integer(0, oldToNew.size()-1);
Swap(oldToNew[i], oldToNew[j]);
}
}
return oldToNew;
}
Foam::labelList Foam::randomRenumber::renumber
(
const polyMesh& mesh,
const pointField& points
)
{
return renumber(points);
}
Foam::labelList Foam::randomRenumber::renumber
(
const labelListList& cellCells,
const pointField& points
)
{
return renumber(points);
}
// ************************************************************************* //

View File

@ -0,0 +1,107 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Class
Foam::randomRenumber
Description
Random renumber. Just to see effect of renumbering.
SourceFiles
randomRenumber.C
\*---------------------------------------------------------------------------*/
#ifndef randomRenumber_H
#define randomRenumber_H
#include "renumberMethod.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class randomRenumber Declaration
\*---------------------------------------------------------------------------*/
class randomRenumber
:
public renumberMethod
{
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
void operator=(const randomRenumber&);
randomRenumber(const randomRenumber&);
public:
//- Runtime type information
TypeName("random");
// Constructors
//- Construct given the renumber dictionary
randomRenumber(const dictionary& renumberDict);
//- Destructor
virtual ~randomRenumber()
{}
// Member Functions
//- Return for every coordinate the wanted processor number.
// We need a polyMesh (to be able to load the file)
virtual labelList renumber(const pointField&);
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList renumber
(
const polyMesh& mesh,
const pointField& cc
);
//- Return for every cell the new cell label.
// The connectivity is equal to mesh.cellCells() except
// - the connections are across coupled patches
virtual labelList renumber
(
const labelListList& cellCells,
const pointField& cc
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,131 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
InClass
renumberMethod
\*---------------------------------------------------------------------------*/
#include "renumberMethod.H"
#include "decompositionMethod.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(renumberMethod, 0);
defineRunTimeSelectionTable(renumberMethod, dictionary);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::autoPtr<Foam::renumberMethod> Foam::renumberMethod::New
(
const dictionary& renumberDict
)
{
const word methodType(renumberDict.lookup("method"));
//Info<< "Selecting renumberMethod " << methodType << endl;
dictionaryConstructorTable::iterator cstrIter =
dictionaryConstructorTablePtr_->find(methodType);
if (cstrIter == dictionaryConstructorTablePtr_->end())
{
FatalErrorIn
(
"renumberMethod::New"
"(const dictionary& renumberDict)"
) << "Unknown renumberMethod "
<< methodType << nl << nl
<< "Valid renumberMethods are : " << endl
<< dictionaryConstructorTablePtr_->sortedToc()
<< exit(FatalError);
}
return autoPtr<renumberMethod>(cstrIter()(renumberDict));
}
Foam::labelList Foam::renumberMethod::renumber
(
const polyMesh& mesh,
const pointField& points
)
{
CompactListList<label> cellCells;
decompositionMethod::calcCellCells
(
mesh,
identity(mesh.nCells()),
mesh.nCells(),
false, // local only
cellCells
);
// Renumber based on agglomerated points
return renumber(cellCells(), points);
}
Foam::labelList Foam::renumberMethod::renumber
(
const polyMesh& mesh,
const labelList& fineToCoarse,
const pointField& coarsePoints
)
{
CompactListList<label> coarseCellCells;
decompositionMethod::calcCellCells
(
mesh,
fineToCoarse,
coarsePoints.size(),
false, // local only
coarseCellCells
);
// Renumber based on agglomerated points
labelList coarseDistribution
(
renumber
(
coarseCellCells(),
coarsePoints
)
);
// Rework back into renumbering for original mesh_
labelList fineDistribution(fineToCoarse.size());
forAll(fineDistribution, i)
{
fineDistribution[i] = coarseDistribution[fineToCoarse[i]];
}
return fineDistribution;
}
// ************************************************************************* //

View File

@ -0,0 +1,161 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Class
Foam::renumberMethod
Description
Abstract base class for renumbering
SourceFiles
renumberMethod.C
\*---------------------------------------------------------------------------*/
#ifndef renumberMethod_H
#define renumberMethod_H
#include "polyMesh.H"
#include "pointField.H"
#include "CompactListList.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class renumberMethod Declaration
\*---------------------------------------------------------------------------*/
class renumberMethod
{
protected:
// Protected data
const dictionary& renumberDict_;
private:
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
renumberMethod(const renumberMethod&);
void operator=(const renumberMethod&);
public:
//- Runtime type information
TypeName("renumberMethod");
// Declare run-time constructor selection tables
declareRunTimeSelectionTable
(
autoPtr,
renumberMethod,
dictionary,
(
const dictionary& renumberDict
),
(renumberDict)
);
// Selectors
//- Return a reference to the selected renumbering method
static autoPtr<renumberMethod> New
(
const dictionary& renumberDict
);
// Constructors
//- Construct given the renumber dictionary
renumberMethod(const dictionary& renumberDict)
:
renumberDict_(renumberDict)
{}
//- Destructor
virtual ~renumberMethod()
{}
// Member Functions
//- Return for every cell the new cell label.
// This is only defined for geometric renumberMethods.
virtual labelList renumber(const pointField&)
{
notImplemented
(
"renumberMethod:renumber(const pointField&)"
);
return labelList(0);
}
//- Return for every cell the new cell label. Use the
// mesh connectivity (if needed)
virtual labelList renumber(const polyMesh&, const pointField&);
//- Return for every cell the new cell label. Gets
// passed agglomeration map (from fine to coarse cells) and coarse
// cell
// location. Can be overridden by renumberMethods that provide this
// functionality natively. Coarse cells are local to the processor
// (if in parallel). If you want to have coarse cells spanning
// processors use the globalCellCells instead.
virtual labelList renumber
(
const polyMesh& mesh,
const labelList& cellToRegion,
const pointField& regionPoints
);
//- Return for every cell the new cell label.
// The connectivity is equal to mesh.cellCells() except
// - the connections are across coupled patches
virtual labelList renumber
(
const labelListList& cellCells,
const pointField& cc
) = 0;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,172 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "springRenumber.H"
#include "addToRunTimeSelectionTable.H"
#include "decompositionMethod.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(springRenumber, 0);
addToRunTimeSelectionTable
(
renumberMethod,
springRenumber,
dictionary
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::springRenumber::springRenumber(const dictionary& renumberDict)
:
renumberMethod(renumberDict),
dict_(renumberDict.subDict(typeName+"Coeffs")),
maxCo_(readScalar(dict_.lookup("maxCo"))),
maxIter_(readLabel(dict_.lookup("maxIter"))),
freezeFraction_(readScalar(dict_.lookup("freezeFraction")))
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::springRenumber::renumber
(
const polyMesh& mesh,
const pointField& points
)
{
CompactListList<label> cellCells;
decompositionMethod::calcCellCells
(
mesh,
identity(mesh.nCells()),
mesh.nCells(),
false, // local only
cellCells
);
return renumber(cellCells(), points);
}
Foam::labelList Foam::springRenumber::renumber
(
const labelListList& cellCells,
const pointField& points
)
{
// Look at cell index as a 1D position parameter.
// Move cells to the average 'position' of their neighbour.
scalarField position(cellCells.size());
forAll(position, cellI)
{
position[cellI] = cellI;
}
labelList oldToNew(identity(cellCells.size()));
scalar maxCo = maxCo_ * cellCells.size();
for (label iter = 0; iter < maxIter_; iter++)
{
//Pout<< "Iteration : " << iter << nl
// << "------------"
// << endl;
//Pout<< "Position :" << nl
// << " min : " << min(position) << nl
// << " max : " << max(position) << nl
// << " avg : " << average(position) << nl
// << endl;
// Sum force per cell.
scalarField sumForce(cellCells.size(), 0.0);
forAll(cellCells, oldCellI)
{
const labelList& cCells = cellCells[oldCellI];
label cellI = oldToNew[oldCellI];
forAll(cCells, i)
{
label nbrCellI = oldToNew[cCells[i]];
sumForce[cellI] += (position[nbrCellI]-position[cellI]);
}
}
//Pout<< "Force :" << nl
// << " min : " << min(sumForce) << nl
// << " max : " << max(sumForce) << nl
// << " avgMag : " << average(mag(sumForce)) << nl
// << "DeltaT : " << deltaT << nl
// << endl;
// Limit displacement
scalar deltaT = maxCo / max(mag(sumForce));
Info<< "Iter:" << iter
<< " maxCo:" << maxCo
<< " deltaT:" << deltaT
<< " average force:" << average(mag(sumForce)) << endl;
// Determine displacement.
scalarField displacement = deltaT*sumForce;
//Pout<< "Displacement :" << nl
// << " min : " << min(displacement) << nl
// << " max : " << max(displacement) << nl
// << " avgMag : " << average(mag(displacement)) << nl
// << endl;
// Calculate new position and scale to be within original range
// (0..nCells-1) for ease of postprocessing.
position += displacement;
position -= min(position);
position *= (position.size()-1)/max(position);
// Slowly freeze.
maxCo *= freezeFraction_;
}
//writeOBJ("endPosition.obj", cellCells, position);
// Move cells to new position
labelList shuffle;
sortedOrder(position, shuffle);
// Reorder oldToNew
inplaceReorder(shuffle, oldToNew);
return oldToNew;
}
// ************************************************************************* //

View File

@ -0,0 +1,132 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 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 <http://www.gnu.org/licenses/>.
Class
Foam::springRenumber
Description
Use spring analogy - attract neighbouring cells according to the distance
of their cell indices.
// Maximum jump of cell indices. Is fraction of number of cells
maxCo 0.1;
// Limit the amount of movement; the fraction maxCo gets decreased
// with every iteration.
freezeFraction 0.9;
// Maximum number of iterations
maxIter 1000;
SourceFiles
springRenumber.C
\*---------------------------------------------------------------------------*/
#ifndef springRenumber_H
#define springRenumber_H
#include "renumberMethod.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class springRenumber Declaration
\*---------------------------------------------------------------------------*/
class springRenumber
:
public renumberMethod
{
// Private data
const dictionary& dict_;
const scalar maxCo_;
const label maxIter_;
const scalar freezeFraction_;
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
void operator=(const springRenumber&);
springRenumber(const springRenumber&);
public:
//- Runtime type information
TypeName("spring");
// Constructors
//- Construct given the renumber dictionary
springRenumber(const dictionary& renumberDict);
//- Destructor
virtual ~springRenumber()
{}
// Member Functions
//- Return for every coordinate the wanted processor number.
virtual labelList renumber(const pointField&)
{
notImplemented("springRenumber::renumber(const pointField&)");
return labelList(0);
}
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList renumber
(
const polyMesh& mesh,
const pointField& cc
);
//- Return for every cell the new cell label.
// The connectivity is equal to mesh.cellCells() except
// - the connections are across coupled patches
virtual labelList renumber
(
const labelListList& cellCells,
const pointField& cc
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -33,8 +33,10 @@ License
#include "IFstream.H"
#include "triSurface.H"
#include "STLpoint.H"
#include "floatVector.H"
#include "OSspecific.H"
#include "mergePoints.H"
//#include "memInfo.H"
using namespace Foam;
@ -77,8 +79,8 @@ class STLLexer
label lineNo_;
word startError_;
DynamicList<STLpoint> STLpoints_;
//DynamicList<STLpoint> STLnormals_;
DynamicList<floatVector> STLpoints_;
//DynamicList<floatVector > STLnormals_;
DynamicList<label> STLlabels_;
HashTable<label, word> STLsolidNames_;
@ -103,12 +105,12 @@ public:
return nTriangles_;
}
DynamicList<STLpoint>& STLpoints()
DynamicList<floatVector>& STLpoints()
{
return STLpoints_;
}
//DynamicList<STLpoint>& STLnormals()
//DynamicList<floatVector>& STLnormals()
//{
// return STLnormals_;
//}
@ -202,8 +204,8 @@ endsolid {space}("endsolid"|"ENDSOLID")({some_space}{word})*
// End of read character pointer returned by strtof
// char* endPtr;
STLpoint normal;
STLpoint vertex;
floatVector normal;
floatVector vertex;
label cmpt = 0; // component index used for reading vertex
static const char* stateNames[7] =
@ -387,16 +389,25 @@ bool triSurface::readSTLASCII(const fileName& STLfileName)
<< exit(FatalError);
}
//memInfo memStat;
//memStat.update();
//Pout<< "At start:" << memStat.rss() << endl;
// Create the lexer obtaining the approximate number of vertices in the STL
// from the file size
STLLexer lexer(&STLstream.stdStream(), Foam::fileSize(STLfileName)/400);
while (lexer.lex() != 0)
{}
DynamicList<STLpoint>& STLpoints = lexer.STLpoints();
//memStat.update();
//Pout<< "After lexing:" << memStat.rss() << endl;
DynamicList<floatVector>& STLpoints = lexer.STLpoints();
DynamicList<label>& STLlabels = lexer.STLlabels();
/*
DynamicList<STLpoint>& STLnormals = lexer.STLnormals();
DynamicList<floatVector>& STLnormals = lexer.STLnormals();
if (STLpoints.size() != 3*STLnormals.size())
{
@ -410,35 +421,51 @@ bool triSurface::readSTLASCII(const fileName& STLfileName)
}
*/
pointField rawPoints(STLpoints.size());
labelList pointMap;
label nUniquePoints = mergePoints
(
STLpoints,
100*SMALL, // merge distance
false, // verbose
pointMap
);
forAll(rawPoints, i)
{
rawPoints[i] = STLpoints[i];
}
//memStat.update();
//Pout<< "After merging:" << memStat.rss() << endl;
STLpoints.clear();
pointField& sp = storedPoints();
setSize(lexer.nTriangles());
DynamicList<label>& STLlabels = lexer.STLlabels();
sp.setSize(nUniquePoints);
forAll(STLpoints, pointI)
{
const floatVector& pt = STLpoints[pointI];
sp[pointMap[pointI]] = vector
(
scalar(pt.x()),
scalar(pt.y()),
scalar(pt.z())
);
}
// Assign triangles
label pointi = 0;
label pointI = 0;
forAll(*this, i)
{
operator[](i)[0] = pointi++;
operator[](i)[1] = pointi++;
operator[](i)[2] = pointi++;
operator[](i)[0] = pointMap[pointI++];
operator[](i)[1] = pointMap[pointI++];
operator[](i)[2] = pointMap[pointI++];
operator[](i).region() = STLlabels[i];
}
//memStat.update();
//Pout<< "After assigning:" << memStat.rss() << endl;
STLpoints.clear();
STLlabels.clear();
// Assign coordinates
storedPoints().transfer(rawPoints);
// Stitch all points within SMALL meters.
stitchTriangles();
// Convert solidNames into regionNames
patches_.setSize(lexer.STLsolidNames().size());
@ -457,6 +484,9 @@ bool triSurface::readSTLASCII(const fileName& STLfileName)
// Fill in the missing information in the patches
setDefaultPatches();
//memStat.update();
//Pout<< "After patchifying:" << memStat.rss() << endl;
return true;
}

View File

@ -28,7 +28,8 @@ License
#include "IFstream.H"
#include "OSspecific.H"
#include "gzstream.h"
#include "floatVector.H"
#include "mergePoints.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@ -94,40 +95,55 @@ bool Foam::triSurface::readSTLBINARY(const fileName& STLfileName)
// Everything OK so go ahead and read the triangles.
// Allocate storage for raw points
pointField rawPoints(3*nTris);
// Allocate storage for triangles
List<floatVector> STLpoints(3*nTris);
setSize(nTris);
label rawPointI = 0;
label pointI = 0;
// Read the triangles
forAll(*this, i)
for (label i = 0; i < nTris; i++)
{
// Read an STL triangle
STLtriangle stlTri(STLfile);
// Set the rawPoints to the vertices of the STL triangle
// and set the point labels of the labelledTri
rawPoints[rawPointI] = stlTri.a();
operator[](i)[0] = rawPointI++;
rawPoints[rawPointI] = stlTri.b();
operator[](i)[1] = rawPointI++;
rawPoints[rawPointI] = stlTri.c();
operator[](i)[2] = rawPointI++;
// Set the STLpoints to the vertices of the STL triangle
STLpoints[pointI++] = stlTri.a();
STLpoints[pointI++] = stlTri.b();
STLpoints[pointI++] = stlTri.c();
operator[](i).region() = stlTri.region();
}
//STLfile.close();
// Stitch points
labelList pointMap;
label nUniquePoints = mergePoints
(
STLpoints,
10*SMALL, // merge distance
false, // verbose
pointMap // old to new
);
// Assign coordinates
storedPoints().transfer(rawPoints);
pointField& sp = storedPoints();
// Stitch all points within SMALL meters.
stitchTriangles();
sp.setSize(nUniquePoints);
forAll(STLpoints, pointI)
{
const floatVector& pt = STLpoints[pointI];
sp[pointMap[pointI]] = vector
(
scalar(pt.x()),
scalar(pt.y()),
scalar(pt.z())
);
}
// Assign triangles
pointI = 0;
forAll(*this, i)
{
operator[](i)[0] = pointMap[pointI++];
operator[](i)[1] = pointMap[pointI++];
operator[](i)[2] = pointMap[pointI++];
}
return true;
}