zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

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

Browse Source
This commit is contained in:
william
2013-12-20 10:32:38 +00:00
parent 005ff02085 01d4c05328
commit b0218cfcb8
29 changed files with 722 additions and 902 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
applications/utilities/surface/surfaceLambdaMuSmooth/surfaceLambdaMuSmooth.C
View File

@ -33,6 +33,9 @@ Description
Provide an edgeMesh file containing points that are not to be moved during
smoothing in order to preserve features.
lambda/mu smoothing: G. Taubin, IBM Research report Rc-19923 (02/01/95)
"A signal processing approach to fair surface design"
\*---------------------------------------------------------------------------*/
#include "argList.H"

10
src/OpenFOAM/db/IOstreams/Fstreams/OFstream.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -104,10 +104,10 @@ Foam::OFstream::OFstream
{
if (debug)
{
Info<< "IFstream::IFstream(const fileName&,"
"streamFormat format=ASCII,"
"versionNumber version=currentVersion) : "
"could not open file for input\n"
Info<< "OFstream::OFstream(const fileName&,"
"streamFormat, versionNumber, compressionType) : "
"could not open file " << pathname
<< "for input\n"
"in stream " << info() << Foam::endl;
}

15
src/OpenFOAM/db/IOstreams/Pstreams/PstreamBuffers.C
View File

@ -145,4 +145,19 @@ void Foam::PstreamBuffers::finishedSends(labelListList& sizes, const bool block)
}
void Foam::PstreamBuffers::clear()
{
forAll(sendBuf_, i)
{
sendBuf_[i].clear();
}
forAll(recvBuf_, i)
{
recvBuf_[i].clear();
}
recvBufPos_ = 0;
finishedSendsCalled_ = false;
}
// ************************************************************************* //

5
src/OpenFOAM/db/IOstreams/Pstreams/PstreamBuffers.H
View File

@ -112,8 +112,6 @@ class PstreamBuffers
bool finishedSendsCalled_;
// Private Member Functions
public:
// Static data
@ -155,6 +153,9 @@ public:
// non-blocking.
void finishedSends(labelListList& sizes, const bool block = true);
//- Clear storage and reset
void clear();
};

2
src/edgeMesh/Make/files
View File

@ -35,8 +35,6 @@ eem = $(em)/extendedEdgeMesh
$(eem)/extendedEdgeMesh.C
$(eem)/extendedEdgeMeshNew.C
$(eem)/extendedEdgeMeshFormats/obj/OBJextendedEdgeFormat.C
$(eem)/extendedEdgeMeshFormats/obj/OBJextendedEdgeFormatRunTime.C
$(eem)/extendedEdgeMeshFormats/extendedEdgeMeshFormat/extendedEdgeMeshFormat.C
$(eem)/extendedEdgeMeshFormats/extendedEdgeMeshFormat/extendedEdgeMeshFormatRunTime.C

63
src/edgeMesh/extendedEdgeMesh/extendedEdgeMesh.C
View File

@ -30,6 +30,7 @@ License
#include "Time.H"
#include "OBJstream.H"
#include "DynamicField.H"
#include "edgeMeshFormatsCore.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@ -104,6 +105,68 @@ Foam::label Foam::extendedEdgeMesh::nPointTypes = 4;
Foam::label Foam::extendedEdgeMesh::nEdgeTypes = 5;
Foam::wordHashSet Foam::extendedEdgeMesh::readTypes()
{
return wordHashSet(*fileExtensionConstructorTablePtr_);
}
Foam::wordHashSet Foam::extendedEdgeMesh::writeTypes()
{
return wordHashSet(*writefileExtensionMemberFunctionTablePtr_);
}
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
bool Foam::extendedEdgeMesh::canReadType
(
const word& ext,
const bool verbose
)
{
return edgeMeshFormatsCore::checkSupport
(
readTypes(),
ext,
verbose,
"reading"
);
}
bool Foam::extendedEdgeMesh::canWriteType
(
const word& ext,
const bool verbose
)
{
return edgeMeshFormatsCore::checkSupport
(
writeTypes(),
ext,
verbose,
"writing"
);
}
bool Foam::extendedEdgeMesh::canRead
(
const fileName& name,
const bool verbose
)
{
word ext = name.ext();
if (ext == "gz")
{
ext = name.lessExt().ext();
}
return canReadType(ext, verbose);
}
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
Foam::extendedEdgeMesh::pointStatus

12
src/edgeMesh/extendedEdgeMesh/extendedEdgeMesh.H
View File

@ -225,6 +225,18 @@ public:
//- Number of possible feature edge types (i.e. number of slices)
static label nEdgeTypes;
//- Can we read this file format?
static bool canRead(const fileName&, const bool verbose=false);
//- Can we read this file format?
static bool canReadType(const word& ext, const bool verbose=false);
//- Can we write this file format type?
static bool canWriteType(const word& ext, const bool verbose=false);
static wordHashSet readTypes();
static wordHashSet writeTypes();
// Constructors

66
src/edgeMesh/extendedEdgeMesh/extendedEdgeMeshFormats/obj/OBJextendedEdgeFormat.C
View File

@ -1,66 +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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "OBJextendedEdgeFormat.H"
#include "OBJedgeFormat.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::fileFormats::OBJextendedEdgeFormat::OBJextendedEdgeFormat
(
const fileName& filename
)
{
read(filename);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::fileFormats::OBJextendedEdgeFormat::read(const fileName& filename)
{
edgeMesh em(filename);
clear();
// Note: should transfer here instead
storedPoints() = em.points();
storedEdges() = em.edges();
return true;
}
void Foam::fileFormats::OBJextendedEdgeFormat::write
(
const fileName& filename,
const extendedEdgeMesh& mesh
)
{
OBJedgeFormat::write(filename, mesh);
}
// ************************************************************************* //

119
src/edgeMesh/extendedEdgeMesh/extendedEdgeMeshFormats/obj/OBJextendedEdgeFormat.H
View File

@ -1,119 +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 <http://www.gnu.org/licenses/>.
Class
Foam::fileFormats::OBJextendedEdgeFormat
Description
Provide a means of reading/writing Alias/Wavefront OBJ format.
Does not handle negative vertex indices.
SourceFiles
OBJextendedEdgeFormat.C
\*---------------------------------------------------------------------------*/
#ifndef OBJextendedEdgeFormat_H
#define OBJextendedEdgeFormat_H
#include "extendedEdgeMesh.H"
//#include "IFstream.H"
//#include "Ostream.H"
#include "OFstream.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace fileFormats
{
/*---------------------------------------------------------------------------*\
Class OBJextendedEdgeFormat Declaration
\*---------------------------------------------------------------------------*/
class OBJextendedEdgeFormat
:
public extendedEdgeMesh
{
// Private Member Functions
//- Disallow default bitwise copy construct
OBJextendedEdgeFormat(const OBJextendedEdgeFormat&);
//- Disallow default bitwise assignment
void operator=(const OBJextendedEdgeFormat&);
public:
// Constructors
//- Construct from file name
OBJextendedEdgeFormat(const fileName&);
// // Selectors
//
// //- Read file and return surface
// static autoPtr<extendedEdgeMesh> New(const fileName& name)
// {
// return autoPtr<extendedEdgeMesh>
// (
// new OBJextendedEdgeFormat(name)
// );
// }
//
//- Destructor
virtual ~OBJextendedEdgeFormat()
{}
// Member Functions
//- Write surface mesh components by proxy
static void write(const fileName&, const extendedEdgeMesh&);
//- Read from file
virtual bool read(const fileName&);
//- Write object file
virtual void write(const fileName& name) const
{
write(name, *this);
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace fileFormats
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

60
src/edgeMesh/extendedEdgeMesh/extendedEdgeMeshFormats/obj/OBJextendedEdgeFormatRunTime.C
View File

@ -1,60 +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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "OBJextendedEdgeFormat.H"
#include "addToRunTimeSelectionTable.H"
#include "addToMemberFunctionSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace fileFormats
{
// read extendedEdgeMesh
addNamedToRunTimeSelectionTable
(
extendedEdgeMesh,
OBJextendedEdgeFormat,
fileExtension,
obj
);
//// write extendedEdgeMesh
//addNamedToMemberFunctionSelectionTable
//(
// extendedEdgeMesh,
// OBJextendedEdgeFormat,
// write,
// fileExtension,
// obj
//);
}
}
// ************************************************************************* //

28
src/lagrangian/basic/Cloud/Cloud.C
View File

@ -219,9 +219,7 @@ void Foam::Cloud<ParticleType>::move(TrackData& td, const scalar trackTime)
// Reset nTrackingRescues
nTrackingRescues_ = 0;
// While there are particles to transfer
while (true)
{
// List of lists of particles to be transfered for all of the
// neighbour processors
List<IDLList<ParticleType> > particleTransferLists
@ -236,6 +234,19 @@ void Foam::Cloud<ParticleType>::move(TrackData& td, const scalar trackTime)
neighbourProcs.size()
);
// Allocate transfer buffers
PstreamBuffers pBufs(Pstream::nonBlocking);
// While there are particles to transfer
while (true)
{
particleTransferLists = IDLList<ParticleType>();
forAll(patchIndexTransferLists, i)
{
patchIndexTransferLists[i].clear();
}
// Loop over all particles
forAllIter(typename Cloud<ParticleType>, *this, pIter)
{
@ -288,8 +299,9 @@ void Foam::Cloud<ParticleType>::move(TrackData& td, const scalar trackTime)
break;
}
// Allocate transfer buffers
PstreamBuffers pBufs(Pstream::nonBlocking);
// Clear transfer buffers
pBufs.clear();
// Stream into send buffers
forAll(particleTransferLists, i)
@ -308,12 +320,12 @@ void Foam::Cloud<ParticleType>::move(TrackData& td, const scalar trackTime)
}
}
// Set up transfers when in non-blocking mode. Returns sizes (in bytes)
// to be sent/received.
labelListList allNTrans(Pstream::nProcs());
// Start sending. Sets number of bytes transferred
labelListList allNTrans(Pstream::nProcs());
pBufs.finishedSends(allNTrans);
bool transfered = false;
forAll(allNTrans, i)

5
src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriver.H
View File

@ -310,14 +310,13 @@ class autoSnapDriver
//- Find point on nearest feature edge (within searchDist).
// Return point and feature
// and store feature-edge to mesh-point and vice versa
pointIndexHit findNearFeatureEdge
Tuple2<label, pointIndexHit> findNearFeatureEdge
(
const indirectPrimitivePatch& pp,
const scalarField& snapDist,
const label pointI,
const point& estimatedPt,
label& featI,
List<List<DynamicList<point> > >&,
List<List<DynamicList<pointConstraint> > >&,
vectorField&,
@ -330,7 +329,7 @@ class autoSnapDriver
// If another mesh point already referring to this feature
// point and further away, reset that one to a near feature
// edge (using findNearFeatureEdge above)
labelPair findNearFeaturePoint
Tuple2<label, pointIndexHit> findNearFeaturePoint
(
const indirectPrimitivePatch& pp,
const scalarField& snapDist,

178
src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C
View File

@ -170,7 +170,7 @@ void Foam::autoSnapDriver::smoothAndConstrain
if (isMasterEdge[meshEdges[edgeI]])
{
label nbrPointI = edges[pEdges[i]].otherVertex(pointI);
label nbrPointI = edges[edgeI].otherVertex(pointI);
if (constraints[nbrPointI].first() >= nConstraints)
{
dispSum[pointI] += disp[nbrPointI];
@ -213,80 +213,6 @@ void Foam::autoSnapDriver::smoothAndConstrain
}
}
}
//XXXXXX
//TODO: make proper parallel so coupled edges don't have double influence
//void Foam::autoSnapDriver::smoothAndConstrain2
//(
// const bool applyConstraints,
// const indirectPrimitivePatch& pp,
// const List<pointConstraint>& constraints,
// vectorField& disp
//) const
//{
// const fvMesh& mesh = meshRefiner_.mesh();
//
// for (label avgIter = 0; avgIter < 20; avgIter++)
// {
// vectorField dispSum(pp.nPoints(), vector::zero);
// labelList dispCount(pp.nPoints(), 0);
//
// const labelListList& pointEdges = pp.pointEdges();
// const edgeList& edges = pp.edges();
//
// forAll(pointEdges, pointI)
// {
// const labelList& pEdges = pointEdges[pointI];
//
// forAll(pEdges, i)
// {
// label nbrPointI = edges[pEdges[i]].otherVertex(pointI);
// dispSum[pointI] += disp[nbrPointI];
// dispCount[pointI]++;
// }
// }
//
// syncTools::syncPointList
// (
// mesh,
// pp.meshPoints(),
// dispSum,
// plusEqOp<point>(),
// vector::zero,
// mapDistribute::transform()
// );
// syncTools::syncPointList
// (
// mesh,
// pp.meshPoints(),
// dispCount,
// plusEqOp<label>(),
// 0,
// mapDistribute::transform()
// );
//
// // Constraints
// forAll(constraints, pointI)
// {
// if (dispCount[pointI] > 0)// && constraints[pointI].first() <= 1)
// {
// // Mix my displacement with neighbours' displacement
// disp[pointI] =
// 0.5
// *(disp[pointI] + dispSum[pointI]/dispCount[pointI]);
//
// if (applyConstraints)
// {
// disp[pointI] = transform
// (
// constraints[pointI].constraintTransformation(),
// disp[pointI]
// );
// }
// }
// }
// }
//}
//XXXXXX
void Foam::autoSnapDriver::calcNearestFace
@ -716,7 +642,7 @@ void Foam::autoSnapDriver::calcNearestFacePointProperties
pNormals = List<point>(pNormals, visitOrder);
pDisp = List<point>(pDisp, visitOrder);
pFc = List<point>(pFc, visitOrder);
pFid = UIndirectList<label>(pFid, visitOrder);
pFid = UIndirectList<label>(pFid, visitOrder)();
}
}
@ -1412,14 +1338,14 @@ Foam::labelPair Foam::autoSnapDriver::findDiagonalAttraction
}
Foam::pointIndexHit Foam::autoSnapDriver::findNearFeatureEdge
Foam::Tuple2<Foam::label, Foam::pointIndexHit>
Foam::autoSnapDriver::findNearFeatureEdge
(
const indirectPrimitivePatch& pp,
const scalarField& snapDist,
const label pointI,
const point& estimatedPt,
label& featI,
List<List<DynamicList<point> > >& edgeAttractors,
List<List<DynamicList<pointConstraint> > >& edgeConstraints,
vectorField& patchAttraction,
@ -1430,16 +1356,18 @@ Foam::pointIndexHit Foam::autoSnapDriver::findNearFeatureEdge
labelList nearEdgeFeat;
List<pointIndexHit> nearEdgeInfo;
vectorField nearNormal;
features.findNearestEdge
(
pointField(1, estimatedPt),
scalarField(1, sqr(snapDist[pointI])),
nearEdgeFeat,
nearEdgeInfo
nearEdgeInfo,
nearNormal
);
const pointIndexHit& nearInfo = nearEdgeInfo[0];
featI = nearEdgeFeat[0];
label featI = nearEdgeFeat[0];
if (nearInfo.hit())
{
@ -1449,12 +1377,7 @@ Foam::pointIndexHit Foam::autoSnapDriver::findNearFeatureEdge
(
nearInfo.hitPoint()
);
pointConstraint c;
const edge e = features[featI].edges()[nearInfo.index()];
vector eVec = e.vec(features[featI].points());
eVec /= mag(eVec)+VSMALL;
c.first() = 2;
c.second() = eVec;
pointConstraint c(Tuple2<label, vector>(2, nearNormal[0]));
edgeConstraints[featI][nearInfo.index()].append(c);
// Store for later use
@ -1462,11 +1385,12 @@ Foam::pointIndexHit Foam::autoSnapDriver::findNearFeatureEdge
nearInfo.hitPoint()-pp.localPoints()[pointI];
patchConstraints[pointI] = c;
}
return nearInfo;
return Tuple2<label, pointIndexHit>(featI, nearInfo);
}
Foam::labelPair Foam::autoSnapDriver::findNearFeaturePoint
Foam::Tuple2<Foam::label, Foam::pointIndexHit>
Foam::autoSnapDriver::findNearFeaturePoint
(
const indirectPrimitivePatch& pp,
const scalarField& snapDist,
@ -1487,26 +1411,23 @@ Foam::labelPair Foam::autoSnapDriver::findNearFeaturePoint
const refinementFeatures& features = meshRefiner_.features();
labelList nearFeat;
labelList nearIndex;
List<pointIndexHit> nearInfo;
features.findNearestPoint
(
pointField(1, estimatedPt),
scalarField(1, sqr(snapDist[pointI])),
nearFeat,
nearIndex
nearInfo
);
label featI = nearFeat[0];
label featPointI = -1;
if (featI != -1)
{
const point& pt = pp.localPoints()[pointI];
const treeDataPoint& shapes =
features.pointTrees()[featI].shapes();
featPointI = shapes.pointLabels()[nearIndex[0]];
const point& featPt = shapes.points()[featPointI];
label featPointI = nearInfo[0].index();
const point& featPt = nearInfo[0].hitPoint();
scalar distSqr = magSqr(featPt-pt);
// Check if already attracted
@ -1537,7 +1458,6 @@ Foam::labelPair Foam::autoSnapDriver::findNearFeaturePoint
patchAttraction[oldPointI] = vector::zero;
patchConstraints[oldPointI] = pointConstraint();
label edgeFeatI;
findNearFeatureEdge
(
pp,
@ -1545,7 +1465,6 @@ Foam::labelPair Foam::autoSnapDriver::findNearFeaturePoint
oldPointI,
pp.localPoints()[oldPointI],
edgeFeatI,
edgeAttractors,
edgeConstraints,
patchAttraction,
@ -1566,7 +1485,7 @@ Foam::labelPair Foam::autoSnapDriver::findNearFeaturePoint
}
}
return labelPair(featI, featPointI);
return Tuple2<label, pointIndexHit>(featI, nearInfo[0]);
}
@ -1637,7 +1556,6 @@ void Foam::autoSnapDriver::determineFeatures
<< featurePointStr().name() << endl;
}
const refinementFeatures& features = meshRefiner_.features();
forAll(pp.localPoints(), pointI)
{
@ -1664,6 +1582,7 @@ void Foam::autoSnapDriver::determineFeatures
constraint
);
if
(
(constraint.first() > patchConstraints[pointI].first())
@ -1696,15 +1615,14 @@ void Foam::autoSnapDriver::determineFeatures
// Behave like when having two surface normals so
// attract to nearest feature edge (with a guess for
// the multipatch point as starting point)
label featI = -1;
pointIndexHit nearInfo = findNearFeatureEdge
Tuple2<label, pointIndexHit> nearInfo =
findNearFeatureEdge
(
pp,
snapDist,
pointI,
multiPatchPt.hitPoint(), //estimatedPt
featI,
edgeAttractors,
edgeConstraints,
@ -1712,14 +1630,15 @@ void Foam::autoSnapDriver::determineFeatures
patchConstraints
);
if (nearInfo.hit())
const pointIndexHit& info = nearInfo.second();
if (info.hit())
{
// Dump
if (featureEdgeStr.valid())
{
featureEdgeStr().write
(
linePointRef(pt, nearInfo.hitPoint())
linePointRef(pt, info.hitPoint())
);
}
}
@ -1729,7 +1648,7 @@ void Foam::autoSnapDriver::determineFeatures
{
missedEdgeStr().write
(
linePointRef(pt, nearInfo.missPoint())
linePointRef(pt, info.missPoint())
);
}
}
@ -1746,15 +1665,13 @@ void Foam::autoSnapDriver::determineFeatures
// Determine nearest point on feature edge. Store constraint
// (calculated from feature edge, alternative would be to
// use constraint calculated from both surfaceNormals)
label featI = -1;
pointIndexHit nearInfo = findNearFeatureEdge
Tuple2<label, pointIndexHit> nearInfo = findNearFeatureEdge
(
pp,
snapDist,
pointI,
estimatedPt,
featI,
edgeAttractors,
edgeConstraints,
@ -1762,14 +1679,15 @@ void Foam::autoSnapDriver::determineFeatures
patchConstraints
);
if (nearInfo.hit())
// Dump to obj
const pointIndexHit& info = nearInfo.second();
if (info.hit())
{
// Dump
if (featureEdgeStr.valid())
{
featureEdgeStr().write
(
linePointRef(pt, nearInfo.hitPoint())
linePointRef(pt, info.hitPoint())
);
}
}
@ -1779,7 +1697,7 @@ void Foam::autoSnapDriver::determineFeatures
{
missedEdgeStr().write
(
linePointRef(pt, nearInfo.missPoint())
linePointRef(pt, info.missPoint())
);
}
}
@ -1789,7 +1707,7 @@ void Foam::autoSnapDriver::determineFeatures
// Mark point on the nearest feature point.
const point estimatedPt(pt + patchAttraction[pointI]);
labelPair nearInfo = findNearFeaturePoint
Tuple2<label, pointIndexHit> nearInfo = findNearFeaturePoint
(
pp,
snapDist,
@ -1807,18 +1725,10 @@ void Foam::autoSnapDriver::determineFeatures
patchConstraints
);
if (nearInfo.first() != -1)
const pointIndexHit& info = nearInfo.second();
if (info.hit() && featurePointStr.valid())
{
// Dump
if (featurePointStr.valid())
{
const treeDataPoint& shapes =
features.pointTrees()[nearInfo.first()].shapes();
const point& featPt =
shapes.points()[nearInfo.second()];
featurePointStr().write(linePointRef(pt, featPt));
}
featurePointStr().write(linePointRef(pt, info.hitPoint()));
}
}
}
@ -2029,22 +1939,20 @@ void Foam::autoSnapDriver::featureAttractionUsingFeatureEdges
if (pointStatus[pointI] == 0) // baffle edge
{
label featI;
const pointIndexHit nearInfo = findNearFeatureEdge
Tuple2<label, pointIndexHit> nearInfo = findNearFeatureEdge
(
pp,
snapDist,
pointI,
pt,
featI,
edgeAttractors,
edgeConstraints,
rawPatchAttraction,
rawPatchConstraints
);
if (!nearInfo.hit())
if (!nearInfo.second().hit())
{
//Pout<< "*** Failed to find close edge to point " << pt
// << endl;
@ -2053,23 +1961,21 @@ void Foam::autoSnapDriver::featureAttractionUsingFeatureEdges
else if (pointStatus[pointI] == 1) // baffle point
{
labelList nearFeat;
labelList nearIndex;
List<pointIndexHit> nearInfo;
features.findNearestPoint
(
pointField(1, pt),
scalarField(1, sqr(snapDist[pointI])),
nearFeat,
nearIndex
nearInfo
);
label featI = nearFeat[0];
if (featI != -1)
{
const treeDataPoint& shapes =
features.pointTrees()[featI].shapes();
label featPointI = shapes.pointLabels()[nearIndex[0]];
const point& featPt = shapes.points()[featPointI];
label featPointI = nearInfo[0].index();
const point& featPt = nearInfo[0].hitPoint();
scalar distSqr = magSqr(featPt-pt);
// Check if already attracted
@ -2099,7 +2005,6 @@ void Foam::autoSnapDriver::featureAttractionUsingFeatureEdges
// The current point is closer so wins. Reset
// the old point to attract to nearest edge
// instead.
label edgeFeatI;
findNearFeatureEdge
(
pp,
@ -2107,7 +2012,6 @@ void Foam::autoSnapDriver::featureAttractionUsingFeatureEdges
oldPointI,
pp.localPoints()[oldPointI],
edgeFeatI,
edgeAttractors,
edgeConstraints,
rawPatchAttraction,
@ -2130,7 +2034,6 @@ void Foam::autoSnapDriver::featureAttractionUsingFeatureEdges
// << " for baffle-feature-point " << pt
// << endl;
label featI;
findNearFeatureEdge
(
pp,
@ -2138,7 +2041,6 @@ void Foam::autoSnapDriver::featureAttractionUsingFeatureEdges
pointI,
pt, // starting point
featI,
edgeAttractors,
edgeConstraints,
rawPatchAttraction,

80
src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C
View File

@ -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,

31
src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.H
View File

@ -105,31 +105,6 @@ class medialAxisMeshMover
const labelList&
);
//- Weighted sum (over all subset of mesh points) by
// summing contribution from (master) edges
template<class Type>
static void weightedSum
(
const polyMesh& mesh,
const PackedBoolList& isMasterEdge,
const labelList& meshEdges,
const labelList& meshPoints,
const edgeList& edges,
const scalarField& edgeWeights,
const Field<Type>& data,
Field<Type>& 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
// ************************************************************************* //

93
src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMoverTemplates.C
View File

@ -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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "medialAxisMeshMover.H"
#include "syncTools.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
void Foam::medialAxisMeshMover::weightedSum
(
const polyMesh& mesh,
const PackedBoolList& isMasterEdge,
const labelList& meshEdges,
const labelList& meshPoints,
const edgeList& edges,
const scalarField& edgeWeights,
const Field<Type>& pointData,
Field<Type>& 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<Type>::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<Type>(),
pTraits<Type>::zero // null value
);
}
// ************************************************************************* //

62
src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.C
View File

@ -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;
}

27
src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.H
View File

@ -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<class Type>
static void weightedSum
(
const polyMesh& mesh,
const PackedBoolList& isMasterEdge,
const labelList& meshEdges,
const labelList& meshPoints,
const edgeList& edges,
const scalarField& edgeWeights,
const Field<Type>& data,
Field<Type>& sum
);
// Refinement

36
src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementRefine.C
View File

@ -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<trackedParticle>()
);
if (debug&meshRefinement::FEATURESEEDS)
{
Pout<< "Constructing cloud for cell marking" << endl;
}
forAllIter(Cloud<trackedParticle>, 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<trackedParticle>, cloud, iter)
{
@ -578,15 +594,15 @@ void Foam::meshRefinement::markFeatureCellLevel
// seeded. Delete particle.
cloud.deleteParticle(tp);
}
else
{
// Keep particle
nParticles++;
}
}
reduce(nParticles, sumOp<label>());
if (nParticles == 0)
if (debug&meshRefinement::FEATURESEEDS)
{
Pout<< "Remaining particles " << cloud.size() << endl;
}
if (returnReduce(cloud.size(), sumOp<label>()) == 0)
{
break;
}

96
src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementTemplates.C
View File

@ -26,16 +26,12 @@ License
#include "meshRefinement.H"
#include "fvMesh.H"
#include "globalIndex.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
#include "syncTools.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
// Add a T entry
template<class T> void meshRefinement::updateList
template<class T> void Foam::meshRefinement::updateList
(
const labelList& newToOld,
const T& nullValue,
@ -59,7 +55,7 @@ template<class T> void meshRefinement::updateList
template<class T>
T meshRefinement::gAverage
T Foam::meshRefinement::gAverage
(
const polyMesh& mesh,
const PackedBoolList& isMasterElem,
@ -109,7 +105,7 @@ T meshRefinement::gAverage
template<class T>
T meshRefinement::gAverage
T Foam::meshRefinement::gAverage
(
const polyMesh& mesh,
const PackedBoolList& isMasterElem,
@ -162,7 +158,7 @@ T meshRefinement::gAverage
// Compare two lists over all boundary faces
template<class T>
void meshRefinement::testSyncBoundaryFaceList
void Foam::meshRefinement::testSyncBoundaryFaceList
(
const scalar tol,
const string& msg,
@ -222,7 +218,7 @@ void meshRefinement::testSyncBoundaryFaceList
// Print list sorted by coordinates. Used for comparing non-parallel v.s.
// parallel operation
template<class T>
void meshRefinement::collectAndPrint
void Foam::meshRefinement::collectAndPrint
(
const UList<point>& points,
const UList<T>& data
@ -268,7 +264,11 @@ void meshRefinement::collectAndPrint
//template<class T, class Mesh>
template<class GeoField>
void meshRefinement::addPatchFields(fvMesh& mesh, const word& patchFieldType)
void Foam::meshRefinement::addPatchFields
(
fvMesh& mesh,
const word& patchFieldType
)
{
HashTable<GeoField*> flds
(
@ -298,7 +298,11 @@ void meshRefinement::addPatchFields(fvMesh& mesh, const word& patchFieldType)
// Reorder patch field
template<class GeoField>
void meshRefinement::reorderPatchFields(fvMesh& mesh, const labelList& oldToNew)
void Foam::meshRefinement::reorderPatchFields
(
fvMesh& mesh,
const labelList& oldToNew
)
{
HashTable<GeoField*> flds
(
@ -316,7 +320,11 @@ void meshRefinement::reorderPatchFields(fvMesh& mesh, const labelList& oldToNew)
template<class Enum>
int meshRefinement::readFlags(const Enum& namedEnum, const wordList& words)
int Foam::meshRefinement::readFlags
(
const Enum& namedEnum,
const wordList& words
)
{
int flags = 0;
@ -330,8 +338,66 @@ int meshRefinement::readFlags(const Enum& namedEnum, const wordList& words)
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type>
void Foam::meshRefinement::weightedSum
(
const polyMesh& mesh,
const PackedBoolList& isMasterEdge,
const labelList& meshEdges,
const labelList& meshPoints,
const edgeList& edges,
const scalarField& edgeWeights,
const Field<Type>& pointData,
Field<Type>& 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<Type>::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<Type>(),
pTraits<Type>::zero // null value
);
}
} // End namespace Foam
// ************************************************************************* //

375
src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.C
View File

@ -26,6 +26,7 @@ License
#include "refinementFeatures.H"
#include "Time.H"
#include "Tuple2.H"
#include "DynamicField.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
@ -41,7 +42,40 @@ void Foam::refinementFeatures::read
fileName featFileName(dict.lookup("file"));
// Try reading extendedEdgeMesh first
IOobject extFeatObj
(
featFileName, // name
io.time().constant(), // instance
"extendedFeatureEdgeMesh", // local
io.time(), // registry
IOobject::MUST_READ,
IOobject::NO_WRITE,
false
);
const fileName fName(extFeatObj.filePath());
if (!fName.empty() && extendedEdgeMesh::canRead(fName))
{
autoPtr<extendedEdgeMesh> eMeshPtr = extendedEdgeMesh::New
(
fName
);
Info<< "Read extendedFeatureEdgeMesh " << extFeatObj.name()
<< nl << incrIndent;
eMeshPtr().writeStats(Info);
Info<< decrIndent << endl;
set(featI, new extendedFeatureEdgeMesh(extFeatObj, eMeshPtr()));
}
else
{
// Try reading edgeMesh
IOobject featObj
(
featFileName, // name
@ -53,21 +87,106 @@ void Foam::refinementFeatures::read
false
);
autoPtr<edgeMesh> eMeshPtr = edgeMesh::New(featObj.filePath());
const fileName fName(featObj.filePath());
set
if (fName.empty())
{
FatalIOErrorIn
(
featI,
new featureEdgeMesh
"refinementFeatures::read"
"(const objectRegistry&"
", const PtrList<dictionary>&)",
dict
) << "Could not open " << featObj.objectPath()
<< exit(FatalIOError);
}
// Read as edgeMesh
autoPtr<edgeMesh> eMeshPtr = edgeMesh::New(fName);
const edgeMesh& eMesh = eMeshPtr();
Info<< "Read edgeMesh " << featObj.name() << nl
<< incrIndent;
eMesh.writeStats(Info);
Info<< decrIndent << endl;
// Analyse for feature points. These are all classified as mixed
// points for lack of anything better
const labelListList& pointEdges = eMesh.pointEdges();
labelList oldToNew(eMesh.points().size(), -1);
DynamicField<point> newPoints(eMesh.points().size());
forAll(pointEdges, pointI)
{
if (pointEdges[pointI].size() > 2)
{
oldToNew[pointI] = newPoints.size();
newPoints.append(eMesh.points()[pointI]);
}
//else if (pointEdges[pointI].size() == 2)
//MEJ: do something based on a feature angle?
}
label nFeatures = newPoints.size();
forAll(oldToNew, pointI)
{
if (oldToNew[pointI] == -1)
{
oldToNew[pointI] = newPoints.size();
newPoints.append(eMesh.points()[pointI]);
}
}
const edgeList& edges = eMesh.edges();
edgeList newEdges(edges.size());
forAll(edges, edgeI)
{
const edge& e = edges[edgeI];
newEdges[edgeI] = edge
(
featObj,
eMeshPtr->points(),
eMeshPtr->edges()
)
oldToNew[e[0]],
oldToNew[e[1]]
);
}
const featureEdgeMesh& eMesh = operator[](featI);
// Construct an extendedEdgeMesh with
// - all points on more than 2 edges : mixed feature points
// - all edges : external edges
extendedEdgeMesh eeMesh
(
newPoints, // pts
newEdges, // eds
0, // (point) concaveStart
0, // (point) mixedStart
nFeatures, // (point) nonFeatureStart
edges.size(), // (edge) internalStart
edges.size(), // (edge) flatStart
edges.size(), // (edge) openStart
edges.size(), // (edge) multipleStart
vectorField(0), // normals
List<extendedEdgeMesh::sideVolumeType>(0),// normalVolumeTypes
vectorField(0), // edgeDirections
labelListList(0), // normalDirections
labelListList(0), // edgeNormals
labelListList(0), // featurePointNormals
labelListList(0), // featurePointEdges
labelList(0) // regionEdges
);
Info<< "Constructed extendedFeatureEdgeMesh " << featObj.name()
<< nl << incrIndent;
eeMesh.writeStats(Info);
Info<< decrIndent << endl;
set(featI, new extendedFeatureEdgeMesh(featObj, eeMesh));
}
const edgeMesh& eMesh = operator[](featI);
//eMesh.mergePoints(meshRefiner_.mergeDistance());
@ -145,7 +264,7 @@ void Foam::refinementFeatures::buildTrees
const labelList& featurePoints
)
{
const featureEdgeMesh& eMesh = operator[](featI);
const edgeMesh& eMesh = operator[](featI);
const pointField& points = eMesh.points();
const edgeList& edges = eMesh.edges();
@ -277,7 +396,7 @@ Foam::refinementFeatures::refinementFeatures
const PtrList<dictionary>& featDicts
)
:
PtrList<featureEdgeMesh>(featDicts.size()),
PtrList<extendedFeatureEdgeMesh>(featDicts.size()),
distances_(featDicts.size()),
levels_(featDicts.size()),
edgeTrees_(featDicts.size()),
@ -289,7 +408,7 @@ Foam::refinementFeatures::refinementFeatures
// Search engines
forAll(*this, i)
{
const featureEdgeMesh& eMesh = operator[](i);
const extendedEdgeMesh& eMesh = operator[](i);
const labelListList& pointEdges = eMesh.pointEdges();
DynamicList<label> featurePoints;
@ -303,81 +422,82 @@ Foam::refinementFeatures::refinementFeatures
Info<< "Detected " << featurePoints.size()
<< " featurePoints out of " << pointEdges.size()
<< " points on feature " << eMesh.name() << endl;
<< " points on feature " << operator[](i).name()
<< endl;
buildTrees(i, featurePoints);
}
}
Foam::refinementFeatures::refinementFeatures
(
const objectRegistry& io,
const PtrList<dictionary>& featDicts,
const scalar minCos
)
:
PtrList<featureEdgeMesh>(featDicts.size()),
distances_(featDicts.size()),
levels_(featDicts.size()),
edgeTrees_(featDicts.size()),
pointTrees_(featDicts.size())
{
// Read features
read(io, featDicts);
// Search engines
forAll(*this, i)
{
const featureEdgeMesh& eMesh = operator[](i);
const pointField& points = eMesh.points();
const edgeList& edges = eMesh.edges();
const labelListList& pointEdges = eMesh.pointEdges();
DynamicList<label> featurePoints;
forAll(pointEdges, pointI)
{
const labelList& pEdges = pointEdges[pointI];
if (pEdges.size() > 2)
{
featurePoints.append(pointI);
}
else if (pEdges.size() == 2)
{
// Check the angle
const edge& e0 = edges[pEdges[0]];
const edge& e1 = edges[pEdges[1]];
const point& p = points[pointI];
const point& p0 = points[e0.otherVertex(pointI)];
const point& p1 = points[e1.otherVertex(pointI)];
vector v0 = p-p0;
scalar v0Mag = mag(v0);
vector v1 = p1-p;
scalar v1Mag = mag(v1);
if
(
v0Mag > SMALL
&& v1Mag > SMALL
&& ((v0/v0Mag & v1/v1Mag) < minCos)
)
{
featurePoints.append(pointI);
}
}
}
Info<< "Detected " << featurePoints.size()
<< " featurePoints out of " << points.size()
<< " points on feature " << eMesh.name()
<< " when using feature cos " << minCos << endl;
buildTrees(i, featurePoints);
}
}
//Foam::refinementFeatures::refinementFeatures
//(
// const objectRegistry& io,
// const PtrList<dictionary>& featDicts,
// const scalar minCos
//)
//:
// PtrList<extendedFeatureEdgeMesh>(featDicts.size()),
// distances_(featDicts.size()),
// levels_(featDicts.size()),
// edgeTrees_(featDicts.size()),
// pointTrees_(featDicts.size())
//{
// // Read features
// read(io, featDicts);
//
// // Search engines
// forAll(*this, i)
// {
// const edgeMesh& eMesh = operator[](i);
// const pointField& points = eMesh.points();
// const edgeList& edges = eMesh.edges();
// const labelListList& pointEdges = eMesh.pointEdges();
//
// DynamicList<label> featurePoints;
// forAll(pointEdges, pointI)
// {
// const labelList& pEdges = pointEdges[pointI];
// if (pEdges.size() > 2)
// {
// featurePoints.append(pointI);
// }
// else if (pEdges.size() == 2)
// {
// // Check the angle
// const edge& e0 = edges[pEdges[0]];
// const edge& e1 = edges[pEdges[1]];
//
// const point& p = points[pointI];
// const point& p0 = points[e0.otherVertex(pointI)];
// const point& p1 = points[e1.otherVertex(pointI)];
//
// vector v0 = p-p0;
// scalar v0Mag = mag(v0);
//
// vector v1 = p1-p;
// scalar v1Mag = mag(v1);
//
// if
// (
// v0Mag > SMALL
// && v1Mag > SMALL
// && ((v0/v0Mag & v1/v1Mag) < minCos)
// )
// {
// featurePoints.append(pointI);
// }
// }
// }
//
// Info<< "Detected " << featurePoints.size()
// << " featurePoints out of " << points.size()
// << " points on feature " << i //eMesh.name()
// << " when using feature cos " << minCos << endl;
//
// buildTrees(i, featurePoints);
// }
//}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@ -387,12 +507,16 @@ void Foam::refinementFeatures::findNearestEdge
const pointField& samples,
const scalarField& nearestDistSqr,
labelList& nearFeature,
List<pointIndexHit>& nearInfo
List<pointIndexHit>& nearInfo,
vectorField& nearNormal
) const
{
nearFeature.setSize(samples.size());
nearFeature = -1;
nearInfo.setSize(samples.size());
nearInfo = pointIndexHit();
nearNormal.setSize(samples.size());
nearNormal = vector::zero;
forAll(edgeTrees_, featI)
{
@ -418,8 +542,18 @@ void Foam::refinementFeatures::findNearestEdge
if (info.hit())
{
nearInfo[sampleI] = info;
nearFeature[sampleI] = featI;
nearInfo[sampleI] = pointIndexHit
(
info.hit(),
info.hitPoint(),
tree.shapes().edgeLabels()[info.index()]
);
const treeDataEdge& td = tree.shapes();
const edge& e = td.edges()[nearInfo[sampleI].index()];
nearNormal[sampleI] = e.vec(td.points());
nearNormal[sampleI] /= mag(nearNormal[sampleI])+VSMALL;
}
}
}
@ -427,18 +561,71 @@ void Foam::refinementFeatures::findNearestEdge
}
//void Foam::refinementFeatures::findNearestPoint
//(
// const pointField& samples,
// const scalarField& nearestDistSqr,
// labelList& nearFeature,
// labelList& nearIndex
//) const
//{
// nearFeature.setSize(samples.size());
// nearFeature = -1;
// nearIndex.setSize(samples.size());
// nearIndex = -1;
//
// forAll(pointTrees_, featI)
// {
// const indexedOctree<treeDataPoint>& tree = pointTrees_[featI];
//
// if (tree.shapes().pointLabels().size() > 0)
// {
// forAll(samples, sampleI)
// {
// const point& sample = samples[sampleI];
//
// scalar distSqr;
// if (nearFeature[sampleI] != -1)
// {
// label nearFeatI = nearFeature[sampleI];
// const indexedOctree<treeDataPoint>& nearTree =
// pointTrees_[nearFeatI];
// label featPointI =
// nearTree.shapes().pointLabels()[nearIndex[sampleI]];
// const point& featPt =
// operator[](nearFeatI).points()[featPointI];
// distSqr = magSqr(featPt-sample);
// }
// else
// {
// distSqr = nearestDistSqr[sampleI];
// }
//
// pointIndexHit info = tree.findNearest(sample, distSqr);
//
// if (info.hit())
// {
// nearFeature[sampleI] = featI;
// nearIndex[sampleI] = info.index();
// }
// }
// }
// }
//}
void Foam::refinementFeatures::findNearestPoint
(
const pointField& samples,
const scalarField& nearestDistSqr,
labelList& nearFeature,
labelList& nearIndex
List<pointIndexHit>& nearInfo
) const
{
nearFeature.setSize(samples.size());
nearFeature = -1;
nearIndex.setSize(samples.size());
nearIndex = -1;
nearInfo.setSize(samples.size());
nearInfo = pointIndexHit();
forAll(pointTrees_, featI)
{
@ -453,14 +640,7 @@ void Foam::refinementFeatures::findNearestPoint
scalar distSqr;
if (nearFeature[sampleI] != -1)
{
label nearFeatI = nearFeature[sampleI];
const indexedOctree<treeDataPoint>& nearTree =
pointTrees_[nearFeatI];
label featPointI =
nearTree.shapes().pointLabels()[nearIndex[sampleI]];
const point& featPt =
operator[](nearFeatI).points()[featPointI];
distSqr = magSqr(featPt-sample);
distSqr = magSqr(nearInfo[sampleI].hitPoint()-sample);
}
else
{
@ -472,7 +652,12 @@ void Foam::refinementFeatures::findNearestPoint
if (info.hit())
{
nearFeature[sampleI] = featI;
nearIndex[sampleI] = info.index();
nearInfo[sampleI] = pointIndexHit
(
info.hit(),
info.hitPoint(),
tree.shapes().pointLabels()[info.index()]
);
}
}
}

61
src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -35,7 +35,7 @@ SourceFiles
#ifndef refinementFeatures_H
#define refinementFeatures_H
#include "featureEdgeMesh.H"
#include "extendedFeatureEdgeMesh.H"
#include "indexedOctree.H"
#include "treeDataEdge.H"
#include "treeDataPoint.H"
@ -51,7 +51,7 @@ namespace Foam
class refinementFeatures
:
public PtrList<featureEdgeMesh>
public PtrList<extendedFeatureEdgeMesh>
{
private:
@ -86,6 +86,20 @@ private:
labelList& maxLevel
) const;
protected:
const PtrList<indexedOctree<treeDataEdge> >& edgeTrees() const
{
return edgeTrees_;
}
const PtrList<indexedOctree<treeDataPoint> >& pointTrees() const
{
return pointTrees_;
}
public:
// Constructors
@ -97,15 +111,6 @@ public:
const PtrList<dictionary>& featDicts
);
//- Construct from description and do geometric analysis to determine
// feature points
refinementFeatures
(
const objectRegistry& io,
const PtrList<dictionary>& featDicts,
const scalar minCos
);
// Member Functions
@ -123,42 +128,38 @@ public:
return distances_;
}
const PtrList<indexedOctree<treeDataEdge> >& edgeTrees() const
{
return edgeTrees_;
}
const PtrList<indexedOctree<treeDataPoint> >& pointTrees() const
{
return pointTrees_;
}
// Query
//- Highest distance of all features
scalar maxDistance() const;
//- Find nearest point on nearest feature edge
//- Find nearest point on nearest feature edge. Sets
// - nearFeature: index of feature mesh
// - nearInfo : location on feature edge and edge index
// (note: not feature edge index but index into
// edges() directly)
// - nearNormal : local feature edge normal
void findNearestEdge
(
const pointField& samples,
const scalarField& nearestDistSqr,
labelList& nearFeature,
List<pointIndexHit>& nearInfo
List<pointIndexHit>& nearInfo,
vectorField& nearNormal
) const;
//- Find nearest feature point. Is an index into feature points
// which itself is an index into the edgeMesh points.
// So the point index is
// pointTrees()[nearFeature].shapes().pointLabels()[nearIndex]
// Wip.
//- Find nearest feature point. Sets
// - nearFeature: index of feature mesh
// - nearInfo : location on feature point and point index.
// (note: not index into shapes().pointLabels() but
// index into points() directly)
void findNearestPoint
(
const pointField& samples,
const scalarField& nearestDistSqr,
labelList& nearFeature,
labelList& nearIndex
List<pointIndexHit>& nearInfo
) const;
//- Find shell level higher than ptLevel

3
tutorials/incompressible/icoFoam/Allclean
View File

@ -30,4 +30,7 @@ do
removeCase $caseName
done
(cd elbow && ./Allclean)
# ----------------------------------------------------------------- end-of-file

18
tutorials/incompressible/icoFoam/Allrun
View File

@ -18,12 +18,6 @@ runMapFieldsConsistent()
mapFieldsNew $1 -case $2 -sourceTime latestTime -consistent > $2/log.mapFields 2>&1
}
runFluentMeshToFoam()
{
echo "fluentMeshToFoam: converting mesh $2"
fluentMeshToFoam $2 -case $1 > $1/log.fluentMeshToFoam 2>&1
}
copySolutionDirs()
{
echo "Copying $2/0* directory to $1"
@ -95,16 +89,8 @@ do
( cd $caseName && runApplication `getApplication` )
done
# elbow case for testing Fluent-FOAM conversion tools
runFluentMeshToFoam elbow elbow/elbow.msh
(
cd elbow || exit
runApplication `getApplication`
runApplication foamMeshToFluent
runApplication foamDataToFluent
)
(cd elbow && ./Allrun)
# ----------------------------------------------------------------- end-of-file

11
tutorials/incompressible/icoFoam/elbow/Allclean Executable file
View File

@ -0,0 +1,11 @@
#!/bin/sh
cd ${0%/*} || exit 1 # run from this directory
# Source tutorial clean functions
. $WM_PROJECT_DIR/bin/tools/CleanFunctions
rm -f constant/polyMesh/boundary
rm -rf fluentInterface
cleanCase
# ----------------------------------------------------------------- end-of-file

15
tutorials/incompressible/icoFoam/elbow/Allrun Executable file
View File

@ -0,0 +1,15 @@
#!/bin/sh
cd ${0%/*} || exit 1 # run from this directory
# Source tutorial run functions
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# Get application directory
application=`getApplication`
runApplication fluentMeshToFoam elbow.msh
runApplication "$application"
runApplication foamMeshToFluent
runApplication foamDataToFluent
# ----------------------------------------------------------------- end-of-file

58
tutorials/incompressible/icoFoam/elbow/constant/polyMesh/boundary
View File

@ -1,58 +0,0 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: dev |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class polyBoundaryMesh;
location "constant/polyMesh";
object boundary;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
6
(
wall-4
{
type wall;
nFaces 100;
startFace 1300;
}
velocity-inlet-5
{
type patch;
nFaces 8;
startFace 1400;
}
velocity-inlet-6
{
type patch;
nFaces 4;
startFace 1408;
}
pressure-outlet-7
{
type patch;
nFaces 8;
startFace 1412;
}
wall-8
{
type wall;
nFaces 34;
startFace 1420;
}
frontAndBackPlanes
{
type empty;
nFaces 1836;
startFace 1454;
}
)
// ************************************************************************* //

35
tutorials/incompressible/simpleFoam/pitzDailyExptInlet/0/epsilon
View File

@ -26,41 +26,6 @@ boundaryField
type timeVaryingMappedFixedValue;
setAverage 0;
offset 0;
value nonuniform List<scalar>
30
(
9813.84
8665.24
1866.31
755.118
205.654
76.6694
28.4518
16.0868
15.9867
11.0187
7.95753
5.26064
3.44136
2.55317
2.27183
2.33608
2.9115
3.59492
3.0497
2.716
2.9325
3.88456
6.91821
14.9754
37.5461
217.022
2043.58
4864.22
6244
6334.7
)
;
}
outlet
{

35
tutorials/incompressible/simpleFoam/pitzDailyExptInlet/0/k
View File

@ -26,41 +26,6 @@ boundaryField
type timeVaryingMappedFixedValue;
setAverage 0;
offset 0;
value nonuniform List<scalar>
30
(
2.95219
2.54219
0.725449
0.486465
0.353566
0.240375
0.172984
0.147052
0.146827
0.135658
0.12147
0.0942189
0.0833465
0.0828453
0.0955983
0.0920838
0.0967682
0.0990811
0.100866
0.101556
0.0967155
0.0841739
0.0904567
0.130411
0.194046
0.219327
0.975528
2.22578
3.12421
2.28104
)
;
}
outlet
{