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' into cvm

Browse Source
This commit is contained in:
graham
2011-01-24 13:24:07 +00:00
parent b328499d02 4f30c291bd
commit 3171dd78a8
8 changed files with 421 additions and 112 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
applications/utilities/postProcessing/sampling/sample/sampleDict
View File

@ -231,6 +231,7 @@ surfaces
// Sampling on triSurface
type sampledTriSurfaceMesh;
surface integrationPlane.stl;
source boundaryFaces; // sample cells or boundaryFaces
interpolate true;
}
);

6
src/finiteVolume/interpolation/interpolation/interpolationCellPoint/cellPointWeight/cellPointWeight.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -155,8 +155,8 @@ void Foam::cellPointWeight::findTriangle
List<tetIndices> faceTets = polyMeshTetDecomposition::faceTetIndices
(
mesh,
mesh.faceOwner()[faceI],
faceI
faceI,
mesh.faceOwner()[faceI]
);
const scalar faceAreaSqr = magSqr(mesh.faceAreas()[faceI]);

108
src/meshTools/indexedOctree/indexedOctree.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -1623,58 +1623,76 @@ void Foam::indexedOctree<Type>::traverseNode
{
const labelList& indices = contents_[getContent(index)];
if (findAny)
if (indices.size())
{
// Find any intersection
forAll(indices, elemI)
if (findAny)
{
label shapeI = indices[elemI];
// Find any intersection
point pt;
bool hit = shapes_.intersects(shapeI, start, end, pt);
if (hit)
forAll(indices, elemI)
{
// Hit so pt is nearer than nearestPoint.
// Update hit info
hitInfo.setHit();
hitInfo.setIndex(shapeI);
hitInfo.setPoint(pt);
label shapeI = indices[elemI];
point pt;
bool hit = shapes_.intersects(shapeI, start, end, pt);
// Note that intersection of shape might actually be
// in a neighbouring box. For findAny this is not important.
if (hit)
{
// Hit so pt is nearer than nearestPoint.
// Update hit info
hitInfo.setHit();
hitInfo.setIndex(shapeI);
hitInfo.setPoint(pt);
return;
}
}
}
else
{
// Find nearest intersection
const treeBoundBox octantBb(subBbox(nodeI, octant));
point nearestPoint(end);
forAll(indices, elemI)
{
label shapeI = indices[elemI];
point pt;
bool hit = shapes_.intersects
(
shapeI,
start,
nearestPoint,
pt
);
// Note that intersection of shape might actually be
// in a neighbouring box. Since we need to maintain strict
// (findAny=false) ordering skip such an intersection. It
// will be found when we are doing the next box.
if (hit && octantBb.contains(pt))
{
// Hit so pt is nearer than nearestPoint.
nearestPoint = pt;
// Update hit info
hitInfo.setHit();
hitInfo.setIndex(shapeI);
hitInfo.setPoint(pt);
}
}
if (hitInfo.hit())
{
// Found intersected shape.
return;
}
}
}
else
{
// Find nearest intersection.
point nearestPoint(end);
forAll(indices, elemI)
{
label shapeI = indices[elemI];
point pt;
bool hit = shapes_.intersects(shapeI, start, nearestPoint, pt);
if (hit)
{
// Hit so pt is nearer than nearestPoint.
nearestPoint = pt;
// Update hit info
hitInfo.setHit();
hitInfo.setIndex(shapeI);
hitInfo.setPoint(pt);
}
}
if (hitInfo.hit())
{
// Found intersected shape.
return;
}
}
}
// Nothing intersected in this node

22
src/meshTools/meshSearch/meshSearch.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -101,9 +101,11 @@ Foam::label Foam::meshSearch::findNearestCellTree(const point& location) const
{
const indexedOctree<treeDataCell>& tree = cellTree();
scalar span = tree.bb().mag();
pointIndexHit info = tree.findNearest(location, Foam::sqr(span));
pointIndexHit info = tree.findNearest
(
location,
magSqr(tree.bb().max()-tree.bb().min())
);
if (!info.hit())
{
@ -178,10 +180,12 @@ Foam::label Foam::meshSearch::findNearestFaceTree(const point& location) const
// Search nearest cell centre.
const indexedOctree<treeDataCell>& tree = cellTree();
scalar span = tree.bb().mag();
// Search with decent span
pointIndexHit info = tree.findNearest(location, Foam::sqr(span));
pointIndexHit info = tree.findNearest
(
location,
magSqr(tree.bb().max()-tree.bb().min())
);
if (!info.hit())
{
@ -767,12 +771,10 @@ Foam::label Foam::meshSearch::findNearestBoundaryFace
{
const indexedOctree<treeDataFace>& tree = boundaryTree();
scalar span = tree.bb().mag();
pointIndexHit info = boundaryTree().findNearest
(
location,
Foam::sqr(span)
magSqr(tree.bb().max()-tree.bb().min())
);
if (!info.hit())

284
src/sampling/sampledSurface/sampledTriSurfaceMesh/sampledTriSurfaceMesh.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -28,6 +28,7 @@ License
#include "Tuple2.H"
#include "globalIndex.H"
#include "treeDataCell.H"
#include "treeDataFace.H"
#include "addToRunTimeSelectionTable.H"
@ -43,6 +44,17 @@ namespace Foam
word
);
template<>
const char* NamedEnum<sampledTriSurfaceMesh::samplingSource, 2>::names[] =
{
"cells",
"boundaryFaces"
};
const NamedEnum<sampledTriSurfaceMesh::samplingSource, 2>
sampledTriSurfaceMesh::samplingSourceNames_;
//- Private class for finding nearest
// - global index
// - sqr(distance)
@ -70,7 +82,8 @@ Foam::sampledTriSurfaceMesh::sampledTriSurfaceMesh
(
const word& name,
const polyMesh& mesh,
const word& surfaceName
const word& surfaceName,
const samplingSource sampleSource
)
:
sampledSurface(name, mesh),
@ -78,7 +91,7 @@ Foam::sampledTriSurfaceMesh::sampledTriSurfaceMesh
(
IOobject
(
name,
surfaceName,
mesh.time().constant(), // instance
"triSurface", // local
mesh, // registry
@ -87,9 +100,10 @@ Foam::sampledTriSurfaceMesh::sampledTriSurfaceMesh
false
)
),
sampleSource_(sampleSource),
needsUpdate_(true),
cellLabels_(0),
pointToFace_(0)
sampleElements_(0),
samplePoints_(0)
{}
@ -114,9 +128,10 @@ Foam::sampledTriSurfaceMesh::sampledTriSurfaceMesh
false
)
),
sampleSource_(samplingSourceNames_[dict.lookup("source")]),
needsUpdate_(true),
cellLabels_(0),
pointToFace_(0)
sampleElements_(0),
samplePoints_(0)
{}
@ -144,8 +159,8 @@ bool Foam::sampledTriSurfaceMesh::expire()
sampledSurface::clearGeom();
MeshStorage::clear();
cellLabels_.clear();
pointToFace_.clear();
sampleElements_.clear();
samplePoints_.clear();
needsUpdate_ = true;
return true;
@ -164,43 +179,79 @@ bool Foam::sampledTriSurfaceMesh::update()
// Does approximation by looking at the face centres only
const pointField& fc = surface_.faceCentres();
// Mesh search engine, no triangulation of faces.
meshSearch meshSearcher(mesh(), false);
const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
// Global numbering for cells - only used to uniquely identify local cells.
globalIndex globalCells(mesh().nCells());
List<nearInfo> nearest(fc.size());
// Global numbering for cells/faces - only used to uniquely identify local
// elements
globalIndex globalCells
(
sampleSource_ == cells
? mesh().nCells()
: mesh().nFaces()
);
forAll(nearest, i)
{
nearest[i].first() = GREAT;
nearest[i].second() = labelMax;
}
// Search triangles using nearest on local mesh
forAll(fc, triI)
if (sampleSource_ == cells)
{
pointIndexHit nearInfo = cellTree.findNearest
(
fc[triI],
sqr(GREAT)
);
if (nearInfo.hit())
// Search for nearest cell
const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
forAll(fc, triI)
{
nearest[triI].first() = magSqr(nearInfo.hitPoint()-fc[triI]);
nearest[triI].second() = globalCells.toGlobal(nearInfo.index());
pointIndexHit nearInfo = cellTree.findNearest
(
fc[triI],
sqr(GREAT)
);
if (nearInfo.hit())
{
nearest[triI].first() = magSqr(nearInfo.hitPoint()-fc[triI]);
nearest[triI].second() = globalCells.toGlobal(nearInfo.index());
}
}
}
else
{
// Search for nearest boundaryFace
const indexedOctree<treeDataFace>& bTree = meshSearcher.boundaryTree();
forAll(fc, triI)
{
pointIndexHit nearInfo = bTree.findNearest
(
fc[triI],
sqr(GREAT)
);
if (nearInfo.hit())
{
nearest[triI].first() = magSqr(nearInfo.hitPoint()-fc[triI]);
nearest[triI].second() = globalCells.toGlobal
(
bTree.shapes().faceLabels()[nearInfo.index()]
);
}
}
}
// See which processor has the nearest.
// See which processor has the nearest. Mark and subset
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Pstream::listCombineGather(nearest, nearestEqOp());
Pstream::listCombineScatter(nearest);
boolList include(surface_.size(), false);
cellLabels_.setSize(fc.size());
cellLabels_ = -1;
labelList cellOrFaceLabels(fc.size(), -1);
label nFound = 0;
forAll(nearest, triI)
@ -211,9 +262,10 @@ bool Foam::sampledTriSurfaceMesh::update()
}
else if (globalCells.isLocal(nearest[triI].second()))
{
cellLabels_[triI] = globalCells.toLocal(nearest[triI].second());
include[triI] = true;
cellOrFaceLabels[triI] = globalCells.toLocal
(
nearest[triI].second()
);
nFound++;
}
}
@ -221,7 +273,7 @@ bool Foam::sampledTriSurfaceMesh::update()
if (debug)
{
Pout<< "Local out of faces:" << cellLabels_.size()
Pout<< "Local out of faces:" << cellOrFaceLabels.size()
<< " keeping:" << nFound << endl;
}
@ -243,7 +295,7 @@ bool Foam::sampledTriSurfaceMesh::update()
forAll(s, faceI)
{
if (include[faceI])
if (cellOrFaceLabels[faceI] != -1)
{
faceMap[newFaceI++] = faceI;
@ -262,11 +314,12 @@ bool Foam::sampledTriSurfaceMesh::update()
pointMap.setSize(newPointI);
}
// Subset cellLabels
cellLabels_ = UIndirectList<label>(cellLabels_, faceMap)();
// Store any face per point
pointToFace_.setSize(pointMap.size());
// Subset cellOrFaceLabels
cellOrFaceLabels = UIndirectList<label>(cellOrFaceLabels, faceMap)();
// Store any face per point (without using pointFaces())
labelList pointToFace(pointMap.size());
// Create faces and points for subsetted surface
faceList& faces = this->storedFaces();
@ -284,7 +337,7 @@ bool Foam::sampledTriSurfaceMesh::update()
forAll(newF, fp)
{
pointToFace_[newF[fp]] = i;
pointToFace[newF[fp]] = i;
}
}
@ -296,6 +349,161 @@ bool Foam::sampledTriSurfaceMesh::update()
Pout<< endl;
}
// Collect the samplePoints and sampleElements
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
if (sampledSurface::interpolate())
{
samplePoints_.setSize(pointMap.size());
sampleElements_.setSize(pointMap.size(), -1);
if (sampleSource_ == cells)
{
// samplePoints_ : per surface point a location inside the cell
// sampleElements_ : per surface point the cell
forAll(points(), pointI)
{
const point& pt = points()[pointI];
label cellI = cellOrFaceLabels[pointToFace[pointI]];
sampleElements_[pointI] = cellI;
// Check if point inside cell
if (mesh().pointInCell(pt, sampleElements_[pointI]))
{
samplePoints_[pointI] = pt;
}
else
{
// Find nearest point on faces of cell
const cell& cFaces = mesh().cells()[cellI];
scalar minDistSqr = VGREAT;
forAll(cFaces, i)
{
const face& f = mesh().faces()[cFaces[i]];
pointHit info = f.nearestPoint(pt, mesh().points());
if (info.distance() < minDistSqr)
{
minDistSqr = info.distance();
samplePoints_[pointI] = info.rawPoint();
}
}
}
}
}
else
{
// samplePoints_ : per surface point a location on the boundary
// sampleElements_ : per surface point the boundary face containing
// the location
forAll(points(), pointI)
{
const point& pt = points()[pointI];
label faceI = cellOrFaceLabels[pointToFace[pointI]];
sampleElements_[pointI] = faceI;
samplePoints_[pointI] = mesh().faces()[faceI].nearestPoint
(
pt,
mesh().points()
).rawPoint();
}
}
}
else
{
// if sampleSource_ == cells:
// samplePoints_ : n/a
// sampleElements_ : per surface triangle the cell
// else:
// samplePoints_ : n/a
// sampleElements_ : per surface triangle the boundary face
samplePoints_.clear();
sampleElements_.transfer(cellOrFaceLabels);
}
if (debug)
{
this->clearOut();
OFstream str(mesh().time().path()/"surfaceToMesh.obj");
Info<< "Dumping correspondence from local surface (points or faces)"
<< " to mesh (cells or faces) to " << str.name() << endl;
label vertI = 0;
if (sampledSurface::interpolate())
{
if (sampleSource_ == cells)
{
forAll(samplePoints_, pointI)
{
meshTools::writeOBJ(str, points()[pointI]);
vertI++;
meshTools::writeOBJ(str, samplePoints_[pointI]);
vertI++;
label cellI = sampleElements_[pointI];
meshTools::writeOBJ(str, mesh().cellCentres()[cellI]);
vertI++;
str << "l " << vertI-2 << ' ' << vertI-1 << ' ' << vertI
<< nl;
}
}
else
{
forAll(samplePoints_, pointI)
{
meshTools::writeOBJ(str, points()[pointI]);
vertI++;
meshTools::writeOBJ(str, samplePoints_[pointI]);
vertI++;
label faceI = sampleElements_[pointI];
meshTools::writeOBJ(str, mesh().faceCentres()[faceI]);
vertI++;
str << "l " << vertI-2 << ' ' << vertI-1 << ' ' << vertI
<< nl;
}
}
}
else
{
if (sampleSource_ == cells)
{
forAll(sampleElements_, triI)
{
meshTools::writeOBJ(str, faceCentres()[triI]);
vertI++;
label cellI = sampleElements_[triI];
meshTools::writeOBJ(str, mesh().cellCentres()[cellI]);
vertI++;
str << "l " << vertI-1 << ' ' << vertI << nl;
}
}
else
{
forAll(sampleElements_, triI)
{
meshTools::writeOBJ(str, faceCentres()[triI]);
vertI++;
label faceI = sampleElements_[triI];
meshTools::writeOBJ(str, mesh().faceCentres()[faceI]);
vertI++;
str << "l " << vertI-1 << ' ' << vertI << nl;
}
}
}
}
needsUpdate_ = false;
return true;
}

29
src/sampling/sampledSurface/sampledTriSurfaceMesh/sampledTriSurfaceMesh.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -61,22 +61,38 @@ class sampledTriSurfaceMesh
public sampledSurface,
public MeshedSurface<face>
{
public:
//- Types of communications
enum samplingSource
{
cells,
boundaryFaces
};
private:
//- Private typedefs for convenience
typedef MeshedSurface<face> MeshStorage;
// Private data
static const NamedEnum<samplingSource, 2> samplingSourceNames_;
//- Surface to sample on
const triSurfaceMesh surface_;
//- Whether to sample internal cell values or boundary values
const samplingSource sampleSource_;
//- Track if the surface needs an update
mutable bool needsUpdate_;
//- From local surface triangle to mesh cell.
labelList cellLabels_;
//- From local surface triangle to mesh cell/face.
labelList sampleElements_;
//- From local surface back to surface_
labelList pointToFace_;
//- Local points to sample per point
pointField samplePoints_;
// Private Member Functions
@ -106,7 +122,8 @@ public:
(
const word& name,
const polyMesh& mesh,
const word& surfaceName
const word& surfaceName,
const samplingSource sampleSource
);
//- Construct from dictionary

76
src/sampling/sampledSurface/sampledTriSurfaceMesh/sampledTriSurfaceMeshTemplates.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -35,12 +35,48 @@ Foam::sampledTriSurfaceMesh::sampleField
) const
{
// One value per face
tmp<Field<Type> > tvalues(new Field<Type>(cellLabels_.size()));
tmp<Field<Type> > tvalues(new Field<Type>(sampleElements_.size()));
Field<Type>& values = tvalues();
forAll(cellLabels_, triI)
if (sampleSource_ == cells)
{
values[triI] = vField[cellLabels_[triI]];
// Sample cells
forAll(sampleElements_, triI)
{
values[triI] = vField[sampleElements_[triI]];
}
}
else
{
// Sample boundary faces
const polyBoundaryMesh& pbm = mesh().boundaryMesh();
label nBnd = mesh().nFaces()-mesh().nInternalFaces();
// Create flat boundary field
Field<Type> bVals(nBnd, pTraits<Type>::zero);
forAll(vField.boundaryField(), patchI)
{
label bFaceI = pbm[patchI].start() - mesh().nInternalFaces();
SubList<Type>
(
bVals,
vField.boundaryField()[patchI].size(),
bFaceI
).assign(vField.boundaryField()[patchI]);
}
// Sample in flat boundary field
forAll(sampleElements_, triI)
{
label faceI = sampleElements_[triI];
values[triI] = bVals[faceI-mesh().nInternalFaces()];
}
}
return tvalues;
@ -55,15 +91,37 @@ Foam::sampledTriSurfaceMesh::interpolateField
) const
{
// One value per vertex
tmp<Field<Type> > tvalues(new Field<Type>(pointToFace_.size()));
tmp<Field<Type> > tvalues(new Field<Type>(sampleElements_.size()));
Field<Type>& values = tvalues();
forAll(pointToFace_, pointI)
if (sampleSource_ == cells)
{
label triI = pointToFace_[pointI];
label cellI = cellLabels_[triI];
// Sample cells.
values[pointI] = interpolator.interpolate(points()[pointI], cellI);
forAll(sampleElements_, pointI)
{
values[pointI] = interpolator.interpolate
(
samplePoints_[pointI],
sampleElements_[pointI]
);
}
}
else
{
// Sample boundary faces.
forAll(samplePoints_, pointI)
{
label faceI = sampleElements_[pointI];
values[pointI] = interpolator.interpolate
(
samplePoints_[pointI],
mesh().faceOwner()[faceI],
faceI
);
}
}
return tvalues;

7
src/triSurface/triSurface/stitchTriangles.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -151,6 +151,11 @@ bool Foam::triSurface::stitchTriangles
}
}
}
else
{
// Can happen for e.g. single triangle or cloud of unconnected triangles
storedPoints() = rawPoints;
}
return hasMerged;
}