sampling: Removed deprecated isoSurfaceCell

This commit is contained in:
Henry Weller
2018-11-13 12:14:13 +00:00
parent 7b3cc689a6
commit 0b99b5af1a
15 changed files with 25 additions and 3332 deletions

View File

@ -20,9 +20,6 @@ EXE_INC = \
EXE_LIBS = \
$(CGAL_LIBS) \
-lboost_thread \
-lmpfr \
-lgmp \
-lconformalVoronoiMesh \
-ldecompositionMethods /* -L$(FOAM_LIBBIN)/dummy -lscotchDecomp */ \
-ltriSurface \

View File

@ -41,7 +41,7 @@ Description
#include "cellShape.H"
#include "cellModeller.H"
#include "DynamicField.H"
#include "isoSurfaceCell.H"
#include "isoSurface.H"
#include "vtkSurfaceWriter.H"
#include "syncTools.H"
@ -716,9 +716,8 @@ int main(int argc, char *argv[])
}
}
isoSurfaceCell iso
isoSurface iso
(
fvm,
cellDistance,
pointDistance,
0, // distance,

View File

@ -25,8 +25,6 @@ sampledSurface/sampledPatchInternalField/sampledPatchInternalField.C
sampledSurface/sampledPlane/sampledPlane.C
sampledSurface/isoSurface/isoSurface.C
sampledSurface/isoSurface/sampledIsoSurface.C
sampledSurface/isoSurface/isoSurfaceCell.C
sampledSurface/isoSurface/sampledIsoSurfaceCell.C
sampledSurface/distanceSurface/distanceSurface.C
sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C
sampledSurface/sampledSurface/sampledSurface.C

View File

@ -50,7 +50,6 @@ void Foam::distanceSurface::createGeometry()
// Clear any stored topologies
facesPtr_.clear();
isoSurfCellPtr_.clear();
isoSurfPtr_.clear();
// Clear derived data
@ -265,22 +264,6 @@ void Foam::distanceSurface::createGeometry()
//- Direct from cell field and point field.
if (cell_)
{
isoSurfCellPtr_.reset
(
new isoSurfaceCell
(
fvm,
cellDistance,
pointDistance_,
distance_,
regularise_
)
);
}
else
{
isoSurfPtr_.reset
(
new isoSurface
@ -291,7 +274,6 @@ void Foam::distanceSurface::createGeometry()
regularise_
)
);
}
if (debug)
{
@ -330,24 +312,13 @@ Foam::distanceSurface::distanceSurface
),
distance_(readScalar(dict.lookup("distance"))),
signed_(readBool(dict.lookup("signed"))),
cell_(dict.lookupOrDefault("cell", true)),
regularise_(dict.lookupOrDefault("regularise", true)),
average_(dict.lookupOrDefault("average", false)),
zoneKey_(keyType::null),
needsUpdate_(true),
isoSurfCellPtr_(nullptr),
isoSurfPtr_(nullptr),
facesPtr_(nullptr)
{
// dict.readIfPresent("zone", zoneKey_);
//
// if (debug && zoneKey_.size() && mesh.cellZones().findZoneID(zoneKey_) < 0)
// {
// Info<< "cellZone " << zoneKey_
// << " not found - using entire mesh" << endl;
// }
}
{}
Foam::distanceSurface::distanceSurface
@ -384,12 +355,10 @@ Foam::distanceSurface::distanceSurface
),
distance_(distance),
signed_(signedDistance),
cell_(cell),
regularise_(regularise),
average_(average),
zoneKey_(keyType::null),
needsUpdate_(true),
isoSurfCellPtr_(nullptr),
isoSurfPtr_(nullptr),
facesPtr_(nullptr)
{}

View File

@ -27,8 +27,6 @@ Class
Description
A sampledSurface defined by a distance to a surface.
Uses either isoSurfaceCell or isoSurface.
SourceFiles
distanceSurface.C
@ -39,7 +37,6 @@ SourceFiles
#include "sampledSurface.H"
#include "searchableSurface.H"
#include "isoSurfaceCell.H"
#include "isoSurface.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -66,9 +63,6 @@ class distanceSurface
//- Signed distance
const bool signed_;
//- Whether to use isoSurfaceCell or isoSurface
const bool cell_;
//- Whether to coarsen
const Switch regularise_;
@ -88,9 +82,6 @@ class distanceSurface
//- Distance to points
scalarField pointDistance_;
//- Constructed iso surface
autoPtr<isoSurfaceCell> isoSurfCellPtr_;
//- Constructed iso surface
autoPtr<isoSurface> isoSurfPtr_;
@ -190,16 +181,9 @@ public:
}
const triSurface& surface() const
{
if (cell_)
{
return isoSurfCellPtr_();
}
else
{
return isoSurfPtr_();
}
}
//- Sample field on surface
virtual tmp<scalarField> sample

View File

@ -37,20 +37,10 @@ Foam::distanceSurface::sampleField
const GeometricField<Type, fvPatchField, volMesh>& vField
) const
{
if (cell_)
{
return tmp<Field<Type>>
(
new Field<Type>(vField, isoSurfCellPtr_().meshCells())
);
}
else
{
return tmp<Field<Type>>
(
new Field<Type>(vField, isoSurfPtr_().meshCells())
);
}
}
@ -72,21 +62,6 @@ Foam::distanceSurface::interpolateField
volPointInterpolation::New(fvm).interpolate(volFld)
);
// Sample.
if (cell_)
{
return isoSurfCellPtr_().interpolate
(
(
average_
? pointAverage(pointFld())()
: volFld
),
pointFld()
);
}
else
{
return isoSurfPtr_().interpolate
(
(
@ -96,7 +71,6 @@ Foam::distanceSurface::interpolateField
),
pointFld()
);
}
}

View File

@ -248,8 +248,6 @@ class isoSurface
) const;
//- Note: cannot use simpler isoSurfaceCell::generateTriPoints since
// the need here to sometimes pass in remote 'snappedPoints'
template<class Type>
void generateTriPoints
(
@ -345,10 +343,6 @@ class isoSurface
public:
//- Declare friendship with isoSurfaceCell to share some functionality
friend class isoSurfaceCell;
//- Runtime type information
TypeName("isoSurface");

File diff suppressed because it is too large Load Diff

View File

@ -1,352 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 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::isoSurfaceCell
Description
A surface formed by the iso value.
After "Polygonising A Scalar Field Using Tetrahedrons", Paul Bourke
(http://paulbourke.net/geometry/polygonise) and
"Regularised Marching Tetrahedra: improved iso-surface extraction",
G.M. Treece, R.W. Prager and A.H. Gee.
See isoSurface. This is a variant which does tetrahedrisation from
triangulation of face to cell centre instead of edge of face to two
neighbouring cell centres. This gives much lower quality triangles
but they are local to a cell.
SourceFiles
isoSurfaceCell.C
\*---------------------------------------------------------------------------*/
#ifndef isoSurfaceCell_H
#define isoSurfaceCell_H
#include "triSurface.H"
#include "labelPair.H"
#include "pointIndexHit.H"
#include "PackedBoolList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class polyMesh;
/*---------------------------------------------------------------------------*\
Class isoSurfaceCell Declaration
\*---------------------------------------------------------------------------*/
class isoSurfaceCell
:
public triSurface
{
// Private data
enum segmentCutType
{
NEARFIRST, // intersection close to e.first()
NEARSECOND, // ,, e.second()
NOTNEAR // no intersection
};
enum cellCutType
{
NOTCUT, // not cut
SPHERE, // all edges to cell centre cut
CUT // normal cut
};
//- Reference to mesh
const polyMesh& mesh_;
const scalarField& cVals_;
const scalarField& pVals_;
//- isoSurfaceCell value
const scalar iso_;
//- When to merge points
const scalar mergeDistance_;
//- Whether cell might be cut
List<cellCutType> cellCutType_;
//- Estimated number of cut cells
label nCutCells_;
//- For every triangle the original cell in mesh
labelList meshCells_;
//- For every unmerged triangle point the point in the triSurface
labelList triPointMergeMap_;
// Private Member Functions
//- Get location of iso value as fraction in between s0,s1
scalar isoFraction
(
const scalar s0,
const scalar s1
) const;
//- Does any edge of triangle cross iso value?
bool isTriCut
(
const triFace& tri,
const scalarField& pointValues
) const;
//- Determine whether cell is cut
cellCutType calcCutType
(
const PackedBoolList&,
const scalarField& cellValues,
const scalarField& pointValues,
const label
) const;
void calcCutTypes
(
const PackedBoolList&,
const scalarField& cellValues,
const scalarField& pointValues
);
//- Return the two common points between two triangles
static labelPair findCommonPoints
(
const labelledTri&,
const labelledTri&
);
//- Calculate centre of surface.
static point calcCentre(const triSurface&);
//- Replace surface (localPoints, localTris) with single point.
// Returns point. Destroys arguments.
pointIndexHit collapseSurface
(
const label celli,
pointField& localPoints,
DynamicList<labelledTri, 64>& localTris
) const;
//- Determine per cc whether all near cuts can be snapped to single
// point.
void calcSnappedCc
(
const PackedBoolList& isTet,
const scalarField& cVals,
const scalarField& pVals,
DynamicList<point>& triPoints,
labelList& snappedCc
) const;
//- Generate triangles for face connected to pointi
void genPointTris
(
const scalarField& cellValues,
const scalarField& pointValues,
const label pointi,
const label facei,
const label celli,
DynamicList<point, 64>& localTriPoints
) const;
//- Generate triangles for tet connected to pointi
void genPointTris
(
const scalarField& pointValues,
const label pointi,
const label facei,
const label celli,
DynamicList<point, 64>& localTriPoints
) const;
//- Determine per point whether all near cuts can be snapped to single
// point.
void calcSnappedPoint
(
const PackedBoolList& isTet,
const scalarField& cVals,
const scalarField& pVals,
DynamicList<point>& triPoints,
labelList& snappedPoint
) const;
//- Generate single point by interpolation or snapping
template<class Type>
Type generatePoint
(
const DynamicList<Type>& snappedPoints,
const scalar s0,
const Type& p0,
const label p0Index,
const scalar s1,
const Type& p1,
const label p1Index
) const;
template<class Type>
void generateTriPoints
(
const DynamicList<Type>& snapped,
const scalar s0,
const Type& p0,
const label p0Index,
const scalar s1,
const Type& p1,
const label p1Index,
const scalar s2,
const Type& p2,
const label p2Index,
const scalar s3,
const Type& p3,
const label p3Index,
DynamicList<Type>& points
) const;
template<class Type>
void generateTriPoints
(
const scalarField& cVals,
const scalarField& pVals,
const Field<Type>& cCoords,
const Field<Type>& pCoords,
const DynamicList<Type>& snappedPoints,
const labelList& snappedCc,
const labelList& snappedPoint,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
) const;
triSurface stitchTriPoints
(
const bool checkDuplicates,
const List<point>& triPoints,
labelList& triPointReverseMap, // unmerged to merged point
labelList& triMap // merged to unmerged triangle
) const;
//- Check single triangle for (topological) validity
static bool validTri(const triSurface&, const label);
//- Determine edge-face addressing
void calcAddressing
(
const triSurface& surf,
List<FixedList<label, 3>>& faceEdges,
labelList& edgeFace0,
labelList& edgeFace1,
Map<labelList>& edgeFacesRest
) const;
//- Is triangle (given by 3 edges) not fully connected?
static bool danglingTriangle
(
const FixedList<label, 3>& fEdges,
const labelList& edgeFace1
);
//- Mark all non-fully connected triangles
static label markDanglingTriangles
(
const List<FixedList<label, 3>>& faceEdges,
const labelList& edgeFace0,
const labelList& edgeFace1,
const Map<labelList>& edgeFacesRest,
boolList& keepTriangles
);
static triSurface subsetMesh
(
const triSurface& s,
const labelList& newToOldFaces,
labelList& oldToNewPoints,
labelList& newToOldPoints
);
public:
//- Runtime type information
TypeName("isoSurfaceCell");
// Constructors
//- Construct from dictionary
isoSurfaceCell
(
const polyMesh& mesh,
const scalarField& cellValues,
const scalarField& pointValues,
const scalar iso,
const bool regularise,
const scalar mergeTol = 1e-6 // fraction of bounding box
);
// Member Functions
//- For every face original cell in mesh
const labelList& meshCells() const
{
return meshCells_;
}
//- Interpolates cCoords,pCoords.
template<class Type>
tmp<Field<Type>> interpolate
(
const Field<Type>& cCoords,
const Field<Type>& pCoords
) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "isoSurfaceCellTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -1,566 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 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 "isoSurfaceCell.H"
#include "polyMesh.H"
#include "tetMatcher.H"
#include "isoSurface.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
Type Foam::isoSurfaceCell::generatePoint
(
const DynamicList<Type>& snappedPoints,
const scalar s0,
const Type& p0,
const label p0Index,
const scalar s1,
const Type& p1,
const label p1Index
) const
{
scalar d = s1-s0;
if (mag(d) > vSmall)
{
scalar s = (iso_-s0)/d;
if (s >= 0.5 && s <= 1 && p1Index != -1)
{
return snappedPoints[p1Index];
}
else if (s >= 0.0 && s <= 0.5 && p0Index != -1)
{
return snappedPoints[p0Index];
}
else
{
return s*p1 + (1.0-s)*p0;
}
}
else
{
scalar s = 0.4999;
return s*p1 + (1.0-s)*p0;
}
}
template<class Type>
void Foam::isoSurfaceCell::generateTriPoints
(
const DynamicList<Type>& snapped,
const scalar s0,
const Type& p0,
const label p0Index,
const scalar s1,
const Type& p1,
const label p1Index,
const scalar s2,
const Type& p2,
const label p2Index,
const scalar s3,
const Type& p3,
const label p3Index,
DynamicList<Type>& pts
) const
{
int triIndex = 0;
if (s0 < iso_)
{
triIndex |= 1;
}
if (s1 < iso_)
{
triIndex |= 2;
}
if (s2 < iso_)
{
triIndex |= 4;
}
if (s3 < iso_)
{
triIndex |= 8;
}
/* Form the vertices of the triangles for each case */
switch (triIndex)
{
case 0x00:
case 0x0F:
break;
case 0x01:
case 0x0E:
{
pts.append
(
generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index)
);
pts.append
(
generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index)
);
pts.append
(
generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index)
);
if (triIndex == 0x0E)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
case 0x02:
case 0x0D:
{
pts.append
(
generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index)
);
pts.append
(
generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index)
);
pts.append
(
generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index)
);
if (triIndex == 0x0D)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
case 0x03:
case 0x0C:
{
Type p0p2 =
generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index);
Type p1p3 =
generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index);
pts.append
(
generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index)
);
pts.append(p1p3);
pts.append(p0p2);
pts.append(p1p3);
pts.append
(
generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index)
);
pts.append(p0p2);
if (triIndex == 0x0C)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-5], pts[sz-4]);
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
case 0x04:
case 0x0B:
{
pts.append
(
generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index)
);
pts.append
(
generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index)
);
pts.append
(
generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index)
);
if (triIndex == 0x0B)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
case 0x05:
case 0x0A:
{
Type p0p1 =
generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
Type p2p3 =
generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
pts.append(p0p1);
pts.append(p2p3);
pts.append
(
generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index)
);
pts.append(p0p1);
pts.append
(
generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index)
);
pts.append(p2p3);
if (triIndex == 0x0A)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-5], pts[sz-4]);
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
case 0x06:
case 0x09:
{
Type p0p1 =
generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
Type p2p3 =
generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
pts.append(p0p1);
pts.append
(
generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index)
);
pts.append(p2p3);
pts.append(p0p1);
pts.append(p2p3);
pts.append
(
generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index)
);
if (triIndex == 0x09)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-5], pts[sz-4]);
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
case 0x08:
case 0x07:
{
pts.append
(
generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index)
);
pts.append
(
generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index)
);
pts.append
(
generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index)
);
if (triIndex == 0x07)
{
// Flip normals
label sz = pts.size();
Swap(pts[sz-2], pts[sz-1]);
}
}
break;
}
}
template<class Type>
void Foam::isoSurfaceCell::generateTriPoints
(
const scalarField& cVals,
const scalarField& pVals,
const Field<Type>& cCoords,
const Field<Type>& pCoords,
const DynamicList<Type>& snappedPoints,
const labelList& snappedCc,
const labelList& snappedPoint,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
) const
{
tetMatcher tet;
label countNotFoundTets = 0;
forAll(mesh_.cells(), celli)
{
if (cellCutType_[celli] != NOTCUT)
{
label oldNPoints = triPoints.size();
const cell& cFaces = mesh_.cells()[celli];
if (tet.isA(mesh_, celli))
{
// For tets don't do cell-centre decomposition, just use the
// tet points and values
const face& f0 = mesh_.faces()[cFaces[0]];
// Get the other point
const face& f1 = mesh_.faces()[cFaces[1]];
label oppositeI = -1;
forAll(f1, fp)
{
oppositeI = f1[fp];
if (findIndex(f0, oppositeI) == -1)
{
break;
}
}
// Start off from positive volume tet to make sure we
// generate outwards pointing tets
if (mesh_.faceOwner()[cFaces[0]] == celli)
{
generateTriPoints
(
snappedPoints,
pVals[f0[1]],
pCoords[f0[1]],
snappedPoint[f0[1]],
pVals[f0[0]],
pCoords[f0[0]],
snappedPoint[f0[0]],
pVals[f0[2]],
pCoords[f0[2]],
snappedPoint[f0[2]],
pVals[oppositeI],
pCoords[oppositeI],
snappedPoint[oppositeI],
triPoints
);
}
else
{
generateTriPoints
(
snappedPoints,
pVals[f0[0]],
pCoords[f0[0]],
snappedPoint[f0[0]],
pVals[f0[1]],
pCoords[f0[1]],
snappedPoint[f0[1]],
pVals[f0[2]],
pCoords[f0[2]],
snappedPoint[f0[2]],
pVals[oppositeI],
pCoords[oppositeI],
snappedPoint[oppositeI],
triPoints
);
}
}
else
{
forAll(cFaces, cFacei)
{
label facei = cFaces[cFacei];
const face& f = mesh_.faces()[facei];
label fp0 = mesh_.tetBasePtIs()[facei];
// Skip undefined tets
if (fp0 < 0)
{
fp0 = 0;
countNotFoundTets++;
}
label fp = f.fcIndex(fp0);
for (label i = 2; i < f.size(); i++)
{
label nextFp = f.fcIndex(fp);
triFace tri(f[fp0], f[fp], f[nextFp]);
// Start off from positive volume tet to make sure we
// generate outwards pointing tets
if (mesh_.faceOwner()[facei] == celli)
{
generateTriPoints
(
snappedPoints,
pVals[tri[1]],
pCoords[tri[1]],
snappedPoint[tri[1]],
pVals[tri[0]],
pCoords[tri[0]],
snappedPoint[tri[0]],
pVals[tri[2]],
pCoords[tri[2]],
snappedPoint[tri[2]],
cVals[celli],
cCoords[celli],
snappedCc[celli],
triPoints
);
}
else
{
generateTriPoints
(
snappedPoints,
pVals[tri[0]],
pCoords[tri[0]],
snappedPoint[tri[0]],
pVals[tri[1]],
pCoords[tri[1]],
snappedPoint[tri[1]],
pVals[tri[2]],
pCoords[tri[2]],
snappedPoint[tri[2]],
cVals[celli],
cCoords[celli],
snappedCc[celli],
triPoints
);
}
fp = nextFp;
}
}
}
// Every three triPoints is a cell
label nCells = (triPoints.size()-oldNPoints)/3;
for (label i = 0; i < nCells; i++)
{
triMeshCells.append(celli);
}
}
}
if (countNotFoundTets > 0)
{
WarningIn("Foam::isoSurfaceCell::generateTriPoints(..)")
<< "Could not find " << countNotFoundTets
<< " tet base points, which may lead to inverted triangles."
<< endl;
}
triPoints.shrink();
triMeshCells.shrink();
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::isoSurfaceCell::interpolate
(
const Field<Type>& cCoords,
const Field<Type>& pCoords
) const
{
DynamicList<Type> triPoints(3*nCutCells_);
DynamicList<label> triMeshCells(nCutCells_);
// Dummy snap data
DynamicList<Type> snappedPoints;
labelList snappedCc(mesh_.nCells(), -1);
labelList snappedPoint(mesh_.nPoints(), -1);
generateTriPoints
(
cVals_,
pVals_,
cCoords,
pCoords,
snappedPoints,
snappedCc,
snappedPoint,
triPoints,
triMeshCells
);
return isoSurface::interpolate
(
points().size(),
triPointMergeMap_,
triPoints
);
}
// ************************************************************************* //

View File

@ -1,368 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 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 "sampledIsoSurfaceCell.H"
#include "dictionary.H"
#include "volFields.H"
#include "volPointInterpolation.H"
#include "addToRunTimeSelectionTable.H"
#include "fvMesh.H"
#include "isoSurfaceCell.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(sampledIsoSurfaceCell, 0);
addNamedToRunTimeSelectionTable
(
sampledSurface,
sampledIsoSurfaceCell,
word,
isoSurfaceCell
);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
bool Foam::sampledIsoSurfaceCell::updateGeometry() const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
// no update needed
if (fvm.time().timeIndex() == prevTimeIndex_)
{
return false;
}
prevTimeIndex_ = fvm.time().timeIndex();
// Clear any stored topo
facesPtr_.clear();
// Clear derived data
sampledSurface::clearGeom();
// Optionally read volScalarField
autoPtr<volScalarField> readFieldPtr_;
// 1. see if field in database
// 2. see if field can be read
const volScalarField* cellFldPtr = nullptr;
if (fvm.foundObject<volScalarField>(isoField_))
{
if (debug)
{
InfoInFunction << "Lookup " << isoField_ << endl;
}
cellFldPtr = &fvm.lookupObject<volScalarField>(isoField_);
}
else
{
// Bit of a hack. Read field and store.
if (debug)
{
InfoInFunction
<< "Reading " << isoField_
<< " from time " <<fvm.time().timeName()
<< endl;
}
readFieldPtr_.reset
(
new volScalarField
(
IOobject
(
isoField_,
fvm.time().timeName(),
fvm,
IOobject::MUST_READ,
IOobject::NO_WRITE,
false
),
fvm
)
);
cellFldPtr = readFieldPtr_.operator->();
}
const volScalarField& cellFld = *cellFldPtr;
tmp<pointScalarField> pointFld
(
volPointInterpolation::New(fvm).interpolate(cellFld)
);
if (average_)
{
//- From point field and interpolated cell.
scalarField cellAvg(fvm.nCells(), scalar(0));
labelField nPointCells(fvm.nCells(), 0);
{
for (label pointi = 0; pointi < fvm.nPoints(); pointi++)
{
const labelList& pCells = fvm.pointCells(pointi);
forAll(pCells, i)
{
label celli = pCells[i];
cellAvg[celli] += pointFld().primitiveField()[pointi];
nPointCells[celli]++;
}
}
}
forAll(cellAvg, celli)
{
cellAvg[celli] /= nPointCells[celli];
}
const isoSurfaceCell iso
(
fvm,
cellAvg,
pointFld().primitiveField(),
isoVal_,
regularise_
);
const_cast<sampledIsoSurfaceCell&>
(
*this
).triSurface::operator=(iso);
meshCells_ = iso.meshCells();
}
else
{
//- Direct from cell field and point field. Gives bad continuity.
const isoSurfaceCell iso
(
fvm,
cellFld.primitiveField(),
pointFld().primitiveField(),
isoVal_,
regularise_
);
const_cast<sampledIsoSurfaceCell&>
(
*this
).triSurface::operator=(iso);
meshCells_ = iso.meshCells();
}
if (debug)
{
Pout<< "sampledIsoSurfaceCell::updateGeometry() : constructed iso:"
<< nl
<< " regularise : " << regularise_ << nl
<< " average : " << average_ << nl
<< " isoField : " << isoField_ << nl
<< " isoValue : " << isoVal_ << nl
<< " points : " << points().size() << nl
<< " tris : " << triSurface::size() << nl
<< " cut cells : " << meshCells_.size() << endl;
}
return true;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::sampledIsoSurfaceCell::sampledIsoSurfaceCell
(
const word& name,
const polyMesh& mesh,
const dictionary& dict
)
:
sampledSurface(name, mesh, dict),
isoField_(dict.lookup("isoField")),
isoVal_(readScalar(dict.lookup("isoValue"))),
regularise_(dict.lookupOrDefault("regularise", true)),
average_(dict.lookupOrDefault("average", true)),
zoneKey_(keyType::null),
facesPtr_(nullptr),
prevTimeIndex_(-1),
meshCells_(0)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::sampledIsoSurfaceCell::~sampledIsoSurfaceCell()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::sampledIsoSurfaceCell::needsUpdate() const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
return fvm.time().timeIndex() != prevTimeIndex_;
}
bool Foam::sampledIsoSurfaceCell::expire()
{
facesPtr_.clear();
// Clear derived data
sampledSurface::clearGeom();
// already marked as expired
if (prevTimeIndex_ == -1)
{
return false;
}
// force update
prevTimeIndex_ = -1;
return true;
}
bool Foam::sampledIsoSurfaceCell::update()
{
return updateGeometry();
}
Foam::tmp<Foam::scalarField>
Foam::sampledIsoSurfaceCell::sample
(
const volScalarField& vField
) const
{
return sampleField(vField);
}
Foam::tmp<Foam::vectorField>
Foam::sampledIsoSurfaceCell::sample
(
const volVectorField& vField
) const
{
return sampleField(vField);
}
Foam::tmp<Foam::sphericalTensorField>
Foam::sampledIsoSurfaceCell::sample
(
const volSphericalTensorField& vField
) const
{
return sampleField(vField);
}
Foam::tmp<Foam::symmTensorField>
Foam::sampledIsoSurfaceCell::sample
(
const volSymmTensorField& vField
) const
{
return sampleField(vField);
}
Foam::tmp<Foam::tensorField>
Foam::sampledIsoSurfaceCell::sample
(
const volTensorField& vField
) const
{
return sampleField(vField);
}
Foam::tmp<Foam::scalarField>
Foam::sampledIsoSurfaceCell::interpolate
(
const interpolation<scalar>& interpolator
) const
{
return interpolateField(interpolator);
}
Foam::tmp<Foam::vectorField>
Foam::sampledIsoSurfaceCell::interpolate
(
const interpolation<vector>& interpolator
) const
{
return interpolateField(interpolator);
}
Foam::tmp<Foam::sphericalTensorField>
Foam::sampledIsoSurfaceCell::interpolate
(
const interpolation<sphericalTensor>& interpolator
) const
{
return interpolateField(interpolator);
}
Foam::tmp<Foam::symmTensorField>
Foam::sampledIsoSurfaceCell::interpolate
(
const interpolation<symmTensor>& interpolator
) const
{
return interpolateField(interpolator);
}
Foam::tmp<Foam::tensorField>
Foam::sampledIsoSurfaceCell::interpolate
(
const interpolation<tensor>& interpolator
) const
{
return interpolateField(interpolator);
}
void Foam::sampledIsoSurfaceCell::print(Ostream& os) const
{
os << "sampledIsoSurfaceCell: " << name() << " :"
<< " field:" << isoField_
<< " value:" << isoVal_;
//<< " faces:" << faces().size() // possibly no geom yet
//<< " points:" << points().size();
}
// ************************************************************************* //

View File

@ -1,246 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 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::sampledIsoSurfaceCell
Description
A sampledSurface defined by a surface of iso value. Always triangulated.
To be used in sampleSurfaces / functionObjects. Recalculates iso surface
only if time changes.
SourceFiles
sampledIsoSurfaceCell.C
\*---------------------------------------------------------------------------*/
#ifndef sampledIsoSurfaceCell_H
#define sampledIsoSurfaceCell_H
#include "sampledSurface.H"
#include "triSurface.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class sampledIsoSurfaceCell Declaration
\*---------------------------------------------------------------------------*/
class sampledIsoSurfaceCell
:
public sampledSurface,
public triSurface
{
// Private data
//- Field to get isoSurface of
const word isoField_;
//- Iso value
const scalar isoVal_;
//- Whether to coarse
const Switch regularise_;
//- Whether to recalculate cell values as average of point values
const Switch average_;
//- If restricted to zones, name of this zone or a regular expression
keyType zoneKey_;
//- Triangles converted to faceList
mutable autoPtr<faceList> facesPtr_;
// Recreated for every isoSurface
//- Time at last call, also track it surface needs an update
mutable label prevTimeIndex_;
//- For every triangle the original cell in mesh
mutable labelList meshCells_;
// Private Member Functions
//- Create iso surface (if time has changed)
// Do nothing (and return false) if no update was needed
bool updateGeometry() const;
//- Sample field on faces
template<class Type>
tmp<Field<Type>> sampleField
(
const GeometricField<Type, fvPatchField, volMesh>& vField
) const;
template<class Type>
tmp<Field<Type>>
interpolateField(const interpolation<Type>&) const;
public:
//- Runtime type information
TypeName("sampledIsoSurfaceCell");
// Constructors
//- Construct from dictionary
sampledIsoSurfaceCell
(
const word& name,
const polyMesh& mesh,
const dictionary& dict
);
//- Destructor
virtual ~sampledIsoSurfaceCell();
// Member Functions
//- Does the surface need an update?
virtual bool needsUpdate() const;
//- Mark the surface as needing an update.
// May also free up unneeded data.
// Return false if surface was already marked as expired.
virtual bool expire();
//- Update the surface as required.
// Do nothing (and return false) if no update was needed
virtual bool update();
//- Points of surface
virtual const pointField& points() const
{
return triSurface::points();
}
//- Faces of surface
virtual const faceList& faces() const
{
if (facesPtr_.empty())
{
const triSurface& s = *this;
facesPtr_.reset(new faceList(s.size()));
forAll(s, i)
{
facesPtr_()[i] = s[i].triFaceFace();
}
}
return facesPtr_;
}
//- Sample field on surface
virtual tmp<scalarField> sample
(
const volScalarField&
) const;
//- Sample field on surface
virtual tmp<vectorField> sample
(
const volVectorField&
) const;
//- Sample field on surface
virtual tmp<sphericalTensorField> sample
(
const volSphericalTensorField&
) const;
//- Sample field on surface
virtual tmp<symmTensorField> sample
(
const volSymmTensorField&
) const;
//- Sample field on surface
virtual tmp<tensorField> sample
(
const volTensorField&
) const;
//- Interpolate field on surface
virtual tmp<scalarField> interpolate
(
const interpolation<scalar>&
) const;
//- Interpolate field on surface
virtual tmp<vectorField> interpolate
(
const interpolation<vector>&
) const;
//- Interpolate field on surface
virtual tmp<sphericalTensorField> interpolate
(
const interpolation<sphericalTensor>&
) const;
//- Interpolate field on surface
virtual tmp<symmTensorField> interpolate
(
const interpolation<symmTensor>&
) const;
//- Interpolate field on surface
virtual tmp<tensorField> interpolate
(
const interpolation<tensor>&
) const;
//- Write
virtual void print(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "sampledIsoSurfaceCellTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -1,88 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 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 "sampledIsoSurfaceCell.H"
#include "isoSurface.H"
#include "volFieldsFwd.H"
#include "pointFields.H"
#include "volPointInterpolation.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::sampledIsoSurfaceCell::sampleField
(
const GeometricField<Type, fvPatchField, volMesh>& vField
) const
{
// Recreate geometry if time has changed
updateGeometry();
return tmp<Field<Type>>(new Field<Type>(vField, meshCells_));
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::sampledIsoSurfaceCell::interpolateField
(
const interpolation<Type>& interpolator
) const
{
// Recreate geometry if time has changed
updateGeometry();
// One value per point
tmp<Field<Type>> tvalues(new Field<Type>(points().size()));
Field<Type>& values = tvalues.ref();
boolList pointDone(points().size(), false);
forAll(faces(), cutFacei)
{
const face& f = faces()[cutFacei];
forAll(f, faceVertI)
{
label pointi = f[faceVertI];
if (!pointDone[pointi])
{
values[pointi] = interpolator.interpolate
(
points()[pointi],
meshCells_[cutFacei]
);
pointDone[pointi] = true;
}
}
}
return tvalues;
}
// ************************************************************************* //

View File

@ -233,15 +233,6 @@ void Foam::sampledCuttingPlane::createGeometry()
regularise_,
mergeTol_
)
// new isoSurfaceCell
//(
// mesh,
// cellDistance,
// pointDistance_,
// 0.0,
// regularise_,
// mergeTol_
//)
);
if (debug)

View File

@ -38,7 +38,6 @@ SourceFiles
#include "sampledSurface.H"
#include "isoSurface.H"
//#include "isoSurfaceCell.H"
#include "plane.H"
#include "ZoneIDs.H"
#include "fvMeshSubset.H"
@ -90,7 +89,6 @@ class sampledCuttingPlane
scalarField pointDistance_;
//- Constructed iso surface
// autoPtr<isoSurfaceCell> isoSurfPtr_;
autoPtr<isoSurface> isoSurfPtr_;
//- Triangles converted to faceList
@ -173,8 +171,6 @@ public:
return facesPtr_;
}
// const isoSurfaceCell& surface() const
const isoSurface& surface() const
{
return isoSurfPtr_();