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 molecularDynamics

Browse Source
This commit is contained in:
graham
2008-11-14 16:10:41 +00:00
parent 64f3eafd1e a622a9a0ec
commit f22f0af757
42 changed files with 4377 additions and 1073 deletions
Download Patch File Download Diff File Expand all files Collapse all files

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

@ -40,7 +40,7 @@ surfaceFormat vtk;
// 1] vertex values determined from neighbouring cell-centre values
// 2] face values determined using the current face interpolation scheme
// for the field (linear, gamma, etc.)
interpolationScheme cellPointFace;
interpolationScheme cellPoint;
// Fields to sample.
fields
@ -155,21 +155,25 @@ surfaces
// Optional: whether to leave as faces (=default) or triangulate
}
/* not yet (re)implemented --
interpolatedIso
{
// Iso surface for interpolated values only
type isoSurface;
isoField rho;
isoValue 0.5;
interpolate true;
//regularise false; //optional: do not simplify
}
constantIso
{
name iso;
field rho;
value 0.5;
// Iso surface for constant values. Guarantees triangles to not
// cross cells.
type isoSurfaceCell;
isoField rho;
isoValue 0.5;
interpolate false;
//regularise false; //optional: do not simplify
}
someIso
{
type iso;
field rho;
value 0.5;
interpolate true;
}
*/
);

56
src/OpenFOAM/db/runTimeSelection/runTimeSelectionTables.H
View File

@ -28,6 +28,12 @@ Class
Description
Macros to enable the easy declaration of run-time selection tables.
declareRunTimeSelectionTable is used to create a run-time selection table
for a base-class which holds constructor pointers on the table.
declareRunTimeNewSelectionTable is used to create a run-time selection
table for a derived-class which holds "New" pointers on the table.
\*---------------------------------------------------------------------------*/
#include "token.H"
@ -86,6 +92,56 @@ Description
static void destroy##argNames##ConstructorTables()
#define declareRunTimeNewSelectionTable\
(autoPtr,baseType,argNames,argList,parList) \
\
/* Construct from argList function pointer type */ \
typedef autoPtr<baseType> (*argNames##ConstructorPtr)argList; \
\
/* Construct from argList function table type */ \
typedef HashTable<argNames##ConstructorPtr, word, string::hash> \
argNames##ConstructorTable; \
\
/* Construct from argList function pointer table pointer */ \
static argNames##ConstructorTable* argNames##ConstructorTablePtr_; \
\
/* Class to add constructor from argList to table */ \
template<class baseType##Type> \
class add##argNames##ConstructorToTable \
{ \
public: \
\
static autoPtr<baseType> New##baseType argList \
{ \
return autoPtr<baseType>(baseType##Type::New parList.ptr());\
} \
\
add##argNames##ConstructorToTable \
( \
const word& lookup = baseType##Type::typeName \
) \
{ \
construct##argNames##ConstructorTables(); \
argNames##ConstructorTablePtr_->insert \
( \
lookup, \
New##baseType \
); \
} \
\
~add##argNames##ConstructorToTable() \
{ \
destroy##argNames##ConstructorTables(); \
} \
}; \
\
/* Table Constructor called from the table add function */ \
static void construct##argNames##ConstructorTables(); \
\
/* Table destructor called from the table add function destructor */\
static void destroy##argNames##ConstructorTables()
#define defineRunTimeSelectionTableConstructor(baseType,argNames) \
\
/* Table Constructor called from the table add function */ \

11
src/OpenFOAM/dimensionedTypes/dimensionedSymmTensor/dimensionedSymmTensor.C
View File

@ -125,6 +125,17 @@ dimensionedScalar det(const dimensionedSymmTensor& dt)
}
dimensionedSymmTensor cof(const dimensionedSymmTensor& dt)
{
return dimensionedSymmTensor
(
"cof("+dt.name()+')',
dt.dimensions(),
cof(dt.value())
);
}
dimensionedSymmTensor inv(const dimensionedSymmTensor& dt)
{
return dimensionedSymmTensor

1
src/OpenFOAM/dimensionedTypes/dimensionedSymmTensor/dimensionedSymmTensor.H
View File

@ -59,6 +59,7 @@ dimensionedSymmTensor twoSymm(const dimensionedSymmTensor&);
dimensionedSymmTensor dev(const dimensionedSymmTensor&);
dimensionedSymmTensor dev2(const dimensionedSymmTensor&);
dimensionedScalar det(const dimensionedSymmTensor&);
dimensionedSymmTensor cof(const dimensionedSymmTensor&);
dimensionedSymmTensor inv(const dimensionedSymmTensor&);

11
src/OpenFOAM/dimensionedTypes/dimensionedTensor/dimensionedTensor.C
View File

@ -92,6 +92,17 @@ dimensionedScalar det(const dimensionedTensor& dt)
}
dimensionedTensor cof(const dimensionedTensor& dt)
{
return dimensionedTensor
(
"cof("+dt.name()+')',
dt.dimensions(),
cof(dt.value())
);
}
dimensionedTensor inv(const dimensionedTensor& dt)
{
return dimensionedTensor

1
src/OpenFOAM/dimensionedTypes/dimensionedTensor/dimensionedTensor.H
View File

@ -56,6 +56,7 @@ dimensionedScalar tr(const dimensionedTensor&);
dimensionedTensor dev(const dimensionedTensor&);
dimensionedTensor dev2(const dimensionedTensor&);
dimensionedScalar det(const dimensionedTensor&);
dimensionedTensor cof(const dimensionedTensor&);
dimensionedTensor inv(const dimensionedTensor&);
dimensionedSymmTensor symm(const dimensionedTensor&);
dimensionedSymmTensor twoSymm(const dimensionedTensor&);

1
src/OpenFOAM/fields/DimensionedFields/DimensionedSymmTensorField/DimensionedSymmTensorField.C
View File

@ -46,6 +46,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev2, transform)
UNARY_FUNCTION(scalar, symmTensor, det, transform)
UNARY_FUNCTION(symmTensor, symmTensor, cof, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv, inv)

1
src/OpenFOAM/fields/DimensionedFields/DimensionedSymmTensorField/DimensionedSymmTensorField.H
View File

@ -58,6 +58,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev2, transform)
UNARY_FUNCTION(scalar, symmTensor, det, transform)
UNARY_FUNCTION(symmTensor, symmTensor, cof, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv, inv)

1
src/OpenFOAM/fields/DimensionedFields/DimensionedTensorField/DimensionedTensorField.C
View File

@ -45,6 +45,7 @@ UNARY_FUNCTION(tensor, tensor, skew, transform)
UNARY_FUNCTION(tensor, tensor, dev, transform)
UNARY_FUNCTION(tensor, tensor, dev2, transform)
UNARY_FUNCTION(scalar, tensor, det, transform)
UNARY_FUNCTION(tensor, tensor, cof, cof)
UNARY_FUNCTION(tensor, tensor, inv, inv)
UNARY_FUNCTION(vector, tensor, eigenValues, sign)
UNARY_FUNCTION(tensor, tensor, eigenVectors, transform)

1
src/OpenFOAM/fields/DimensionedFields/DimensionedTensorField/DimensionedTensorField.H
View File

@ -58,6 +58,7 @@ UNARY_FUNCTION(tensor, tensor, skew, transform)
UNARY_FUNCTION(tensor, tensor, dev, transform)
UNARY_FUNCTION(tensor, tensor, dev2, transform)
UNARY_FUNCTION(scalar, tensor, det, transform)
UNARY_FUNCTION(tensor, tensor, cof, cof)
UNARY_FUNCTION(tensor, tensor, inv, inv)
UNARY_FUNCTION(vector, tensor, eigenValues, sign)
UNARY_FUNCTION(tensor, tensor, eigenVectors, transform)

1
src/OpenFOAM/fields/FieldFields/symmTensorFieldField/symmTensorFieldField.C
View File

@ -48,6 +48,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm)
UNARY_FUNCTION(symmTensor, symmTensor, dev)
UNARY_FUNCTION(symmTensor, symmTensor, dev2)
UNARY_FUNCTION(scalar, symmTensor, det)
UNARY_FUNCTION(symmTensor, symmTensor, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv)

1
src/OpenFOAM/fields/FieldFields/symmTensorFieldField/symmTensorFieldField.H
View File

@ -58,6 +58,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm)
UNARY_FUNCTION(symmTensor, symmTensor, dev)
UNARY_FUNCTION(symmTensor, symmTensor, dev2)
UNARY_FUNCTION(scalar, symmTensor, det)
UNARY_FUNCTION(symmTensor, symmTensor, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv)

1
src/OpenFOAM/fields/FieldFields/tensorFieldField/tensorFieldField.C
View File

@ -47,6 +47,7 @@ UNARY_FUNCTION(tensor, tensor, skew)
UNARY_FUNCTION(tensor, tensor, dev)
UNARY_FUNCTION(tensor, tensor, dev2)
UNARY_FUNCTION(scalar, tensor, det)
UNARY_FUNCTION(tensor, tensor, cof)
UNARY_FUNCTION(tensor, tensor, inv)
UNARY_FUNCTION(vector, tensor, eigenValues)
UNARY_FUNCTION(tensor, tensor, eigenVectors)

1
src/OpenFOAM/fields/FieldFields/tensorFieldField/tensorFieldField.H
View File

@ -58,6 +58,7 @@ UNARY_FUNCTION(tensor, tensor, skew)
UNARY_FUNCTION(tensor, tensor, dev)
UNARY_FUNCTION(tensor, tensor, dev2)
UNARY_FUNCTION(scalar, tensor, det)
UNARY_FUNCTION(tensor, tensor, cof)
UNARY_FUNCTION(tensor, tensor, inv)
UNARY_FUNCTION(vector, tensor, eigenValues)
UNARY_FUNCTION(tensor, tensor, eigenVectors)

1
src/OpenFOAM/fields/Fields/symmTensorField/symmTensorField.C
View File

@ -46,6 +46,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm)
UNARY_FUNCTION(symmTensor, symmTensor, dev)
UNARY_FUNCTION(symmTensor, symmTensor, dev2)
UNARY_FUNCTION(scalar, symmTensor, det)
UNARY_FUNCTION(symmTensor, symmTensor, cof)
void inv(Field<symmTensor>& tf, const UList<symmTensor>& tf1)
{

1
src/OpenFOAM/fields/Fields/symmTensorField/symmTensorField.H
View File

@ -62,6 +62,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm)
UNARY_FUNCTION(symmTensor, symmTensor, dev)
UNARY_FUNCTION(symmTensor, symmTensor, dev2)
UNARY_FUNCTION(scalar, symmTensor, det)
UNARY_FUNCTION(symmTensor, symmTensor, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv)

1
src/OpenFOAM/fields/Fields/tensorField/tensorField.C
View File

@ -45,6 +45,7 @@ UNARY_FUNCTION(tensor, tensor, skew)
UNARY_FUNCTION(tensor, tensor, dev)
UNARY_FUNCTION(tensor, tensor, dev2)
UNARY_FUNCTION(scalar, tensor, det)
UNARY_FUNCTION(tensor, tensor, cof)
void inv(Field<tensor>& tf, const UList<tensor>& tf1)
{

1
src/OpenFOAM/fields/Fields/tensorField/tensorField.H
View File

@ -62,6 +62,7 @@ UNARY_FUNCTION(tensor, tensor, skew)
UNARY_FUNCTION(tensor, tensor, dev)
UNARY_FUNCTION(tensor, tensor, dev2)
UNARY_FUNCTION(scalar, tensor, det)
UNARY_FUNCTION(tensor, tensor, cof)
UNARY_FUNCTION(tensor, tensor, inv)
UNARY_FUNCTION(vector, tensor, eigenValues)
UNARY_FUNCTION(tensor, tensor, eigenVectors)

1
src/OpenFOAM/fields/GeometricFields/GeometricSymmTensorField/GeometricSymmTensorField.C
View File

@ -46,6 +46,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev2, transform)
UNARY_FUNCTION(scalar, symmTensor, det, transform)
UNARY_FUNCTION(symmTensor, symmTensor, cof, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv, inv)

1
src/OpenFOAM/fields/GeometricFields/GeometricSymmTensorField/GeometricSymmTensorField.H
View File

@ -58,6 +58,7 @@ UNARY_FUNCTION(symmTensor, symmTensor, twoSymm, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev, transform)
UNARY_FUNCTION(symmTensor, symmTensor, dev2, transform)
UNARY_FUNCTION(scalar, symmTensor, det, transform)
UNARY_FUNCTION(symmTensor, symmTensor, cof, cof)
UNARY_FUNCTION(symmTensor, symmTensor, inv, inv)

1
src/OpenFOAM/fields/GeometricFields/GeometricTensorField/GeometricTensorField.C
View File

@ -45,6 +45,7 @@ UNARY_FUNCTION(tensor, tensor, skew, transform)
UNARY_FUNCTION(tensor, tensor, dev, transform)
UNARY_FUNCTION(tensor, tensor, dev2, transform)
UNARY_FUNCTION(scalar, tensor, det, transform)
UNARY_FUNCTION(tensor, tensor, cof, cof)
UNARY_FUNCTION(tensor, tensor, inv, inv)
UNARY_FUNCTION(vector, tensor, eigenValues, sign)
UNARY_FUNCTION(tensor, tensor, eigenVectors, transform)

1
src/OpenFOAM/fields/GeometricFields/GeometricTensorField/GeometricTensorField.H
View File

@ -58,6 +58,7 @@ UNARY_FUNCTION(tensor, tensor, skew, transform)
UNARY_FUNCTION(tensor, tensor, dev, transform)
UNARY_FUNCTION(tensor, tensor, dev2, transform)
UNARY_FUNCTION(scalar, tensor, det, transform)
UNARY_FUNCTION(tensor, tensor, cof, cof)
UNARY_FUNCTION(tensor, tensor, inv, inv)
UNARY_FUNCTION(vector, tensor, eigenValues, sign)
UNARY_FUNCTION(tensor, tensor, eigenVectors, transform)

15
src/OpenFOAM/primitives/SymmTensor/SymmTensorI.H
View File

@ -33,13 +33,11 @@ namespace Foam
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
//- Construct null
template <class Cmpt>
inline SymmTensor<Cmpt>::SymmTensor()
{}
//- Construct given VectorSpace
template <class Cmpt>
inline SymmTensor<Cmpt>::SymmTensor
(
@ -50,7 +48,6 @@ inline SymmTensor<Cmpt>::SymmTensor
{}
//- Construct given SphericalTensor
template <class Cmpt>
inline SymmTensor<Cmpt>::SymmTensor(const SphericalTensor<Cmpt>& st)
{
@ -60,7 +57,6 @@ inline SymmTensor<Cmpt>::SymmTensor(const SphericalTensor<Cmpt>& st)
}
//- Construct from components
template <class Cmpt>
inline SymmTensor<Cmpt>::SymmTensor
(
@ -75,7 +71,6 @@ inline SymmTensor<Cmpt>::SymmTensor
}
//- Construct from Istream
template <class Cmpt>
inline SymmTensor<Cmpt>::SymmTensor(Istream& is)
:
@ -159,7 +154,6 @@ inline Cmpt& SymmTensor<Cmpt>::zz()
}
//- Return symmetric tensor transpose
template <class Cmpt>
inline const SymmTensor<Cmpt>& SymmTensor<Cmpt>::T() const
{
@ -323,22 +317,19 @@ inline Cmpt det(const SymmTensor<Cmpt>& st)
}
//- Return the cofactor tensor of a symmetric tensor
//- Return the cofactor symmetric tensor of a symmetric tensor
template <class Cmpt>
inline SymmTensor<Cmpt> cofactors(const SymmTensor<Cmpt>& st)
inline SymmTensor<Cmpt> cof(const SymmTensor<Cmpt>& st)
{
return SymmTensor<Cmpt>
(
st.yy()*st.zz() - st.yz()*st.yz(),
st.xz()*st.yz() - st.xy()*st.zz(),
st.xy()*st.yz() - st.yy()*st.xz(),
st.xy()*st.yz() - st.xz()*st.yy(),
st.xz()*st.yz() - st.xy()*st.zz(),
st.xx()*st.zz() - st.xz()*st.xz(),
st.xy()*st.xz() - st.xx()*st.yz(),
st.xy()*st.yz() - st.xz()*st.yy(),
st.xy()*st.xz() - st.xx()*st.yz(),
st.xx()*st.yy() - st.xy()*st.xy()
);
}

2
src/OpenFOAM/primitives/Tensor/TensorI.H
View File

@ -467,7 +467,7 @@ inline Cmpt det(const Tensor<Cmpt>& t)
//- Return the cofactor tensor of a tensor
template <class Cmpt>
inline Tensor<Cmpt> cofactors(const Tensor<Cmpt>& t)
inline Tensor<Cmpt> cof(const Tensor<Cmpt>& t)
{
return Tensor<Cmpt>
(

22
src/OpenFOAM/primitives/Tensor2D/Tensor2DI.H
View File

@ -31,13 +31,11 @@ namespace Foam
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
//- Construct null
template <class Cmpt>
inline Tensor2D<Cmpt>::Tensor2D()
{}
//- Construct given VectorSpace
template <class Cmpt>
inline Tensor2D<Cmpt>::Tensor2D(const VectorSpace<Tensor2D<Cmpt>, Cmpt, 4>& vs)
:
@ -45,7 +43,6 @@ inline Tensor2D<Cmpt>::Tensor2D(const VectorSpace<Tensor2D<Cmpt>, Cmpt, 4>& vs)
{}
//- Construct given SphericalTensor2D
template <class Cmpt>
inline Tensor2D<Cmpt>::Tensor2D(const SphericalTensor2D<Cmpt>& st)
{
@ -54,7 +51,6 @@ inline Tensor2D<Cmpt>::Tensor2D(const SphericalTensor2D<Cmpt>& st)
}
//- Construct from components
template <class Cmpt>
inline Tensor2D<Cmpt>::Tensor2D
(
@ -67,7 +63,6 @@ inline Tensor2D<Cmpt>::Tensor2D
}
//- Construct from Istream
template <class Cmpt>
inline Tensor2D<Cmpt>::Tensor2D(Istream& is)
:
@ -153,7 +148,6 @@ inline Cmpt& Tensor2D<Cmpt>::yy()
}
//- Return tensor transpose
template <class Cmpt>
inline Tensor2D<Cmpt> Tensor2D<Cmpt>::T() const
{
@ -303,17 +297,27 @@ inline Cmpt det(const Tensor2D<Cmpt>& t)
return(t.xx()*t.yy() - t.xy()*t.yx());
}
//- Return the inverse of a tensor give the determinant
//- Return the cofactor tensor of a tensor
template <class Cmpt>
inline Tensor2D<Cmpt> inv(const Tensor2D<Cmpt>& t, const Cmpt dett)
inline Tensor2D<Cmpt> cof(const Tensor2D<Cmpt>& t)
{
return Tensor2D<Cmpt>
(
t.yy(), -t.xy(),
-t.yx(), t.xx()
)/dett;
);
}
//- Return the inverse of a tensor given the determinant
template <class Cmpt>
inline Tensor2D<Cmpt> inv(const Tensor2D<Cmpt>& t, const Cmpt dett)
{
return cof(t)/dett;
}
//- Return the inverse of a tensor
template <class Cmpt>
inline Tensor2D<Cmpt> inv(const Tensor2D<Cmpt>& t)

8
src/Pstream/mpi/IPread.C
View File

@ -193,15 +193,13 @@ void Foam::IPstream::waitRequests()
{
if (IPstream_outstandingRequests_.size() > 0)
{
List<MPI_Status> status(IPstream_outstandingRequests_.size());
if
(
MPI_Waitall
(
IPstream_outstandingRequests_.size(),
IPstream_outstandingRequests_.begin(),
status.begin()
MPI_STATUSES_IGNORE
)
)
{
@ -231,9 +229,7 @@ bool Foam::IPstream::finishedRequest(const label i)
}
int flag;
MPI_Status status;
MPI_Test(&IPstream_outstandingRequests_[i], &flag, &status);
MPI_Test(&IPstream_outstandingRequests_[i], &flag, MPI_STATUS_IGNORE);
return flag != 0;
}

8
src/Pstream/mpi/OPwrite.C
View File

@ -131,15 +131,13 @@ void Foam::OPstream::waitRequests()
{
if (OPstream_outstandingRequests_.size() > 0)
{
List<MPI_Status> status(OPstream_outstandingRequests_.size());
if
(
MPI_Waitall
(
OPstream_outstandingRequests_.size(),
OPstream_outstandingRequests_.begin(),
status.begin()
MPI_STATUSES_IGNORE
)
)
{
@ -169,9 +167,7 @@ bool Foam::OPstream::finishedRequest(const label i)
}
int flag;
MPI_Status status;
MPI_Test(&OPstream_outstandingRequests_[i], &flag, &status);
MPI_Test(&OPstream_outstandingRequests_[i], &flag, MPI_STATUS_IGNORE);
return flag != 0;
}

12
src/Pstream/mpi/Pstream.C
View File

@ -157,8 +157,6 @@ void Foam::reduce(scalar& Value, const sumOp<scalar>& bop)
if (Pstream::nProcs() <= Pstream::nProcsSimpleSum)
{
MPI_Status status;
if (Pstream::master())
{
for
@ -180,7 +178,7 @@ void Foam::reduce(scalar& Value, const sumOp<scalar>& bop)
Pstream::procID(slave),
Pstream::msgType(),
MPI_COMM_WORLD,
&status
MPI_STATUS_IGNORE
)
)
{
@ -260,7 +258,7 @@ void Foam::reduce(scalar& Value, const sumOp<scalar>& bop)
Pstream::procID(Pstream::masterNo()),
Pstream::msgType(),
MPI_COMM_WORLD,
&status
MPI_STATUS_IGNORE
)
)
{
@ -279,8 +277,6 @@ void Foam::reduce(scalar& Value, const sumOp<scalar>& bop)
Value = sum;
/*
MPI_Status status;
int myProcNo = Pstream::myProcNo();
int nProcs = Pstream::nProcs();
@ -314,7 +310,7 @@ void Foam::reduce(scalar& Value, const sumOp<scalar>& bop)
Pstream::procID(childProcId),
Pstream::msgType(),
MPI_COMM_WORLD,
&status
MPI_STATUS_IGNORE
)
)
{
@ -370,7 +366,7 @@ void Foam::reduce(scalar& Value, const sumOp<scalar>& bop)
Pstream::procID(parentId),
Pstream::msgType(),
MPI_COMM_WORLD,
&status
MPI_STATUS_IGNORE
)
)
{

21
src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C
View File

@ -2875,6 +2875,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
}
// 2. Extend to 2:1. I don't understand yet why this is not done
// 2. Extend to 2:1. For non-cube cells the scalar distance does not work
// so make sure it at least provides 2:1.
PackedList<1> refineCell(mesh_.nCells(), 0);
forAll(allCellInfo, cellI)
{
@ -2887,6 +2889,25 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
}
faceConsistentRefinement(true, refineCell);
while (true)
{
label nChanged = faceConsistentRefinement(true, refineCell);
reduce(nChanged, sumOp<label>());
if (debug)
{
Pout<< "hexRef8::consistentSlowRefinement2 : Changed " << nChanged
<< " refinement levels due to 2:1 conflicts."
<< endl;
}
if (nChanged == 0)
{
break;
}
}
// 3. Convert back to labelList.
label nRefined = 0;

2
src/sampling/Make/files
View File

@ -24,6 +24,8 @@ sampledSurface/patch/sampledPatch.C
sampledSurface/plane/sampledPlane.C
sampledSurface/isoSurface/isoSurface.C
sampledSurface/isoSurface/sampledIsoSurface.C
sampledSurface/isoSurface/isoSurfaceCell.C
sampledSurface/isoSurface/sampledIsoSurfaceCell.C
sampledSurface/distanceSurface/distanceSurface.C
sampledSurface/sampledSurface/sampledSurface.C
sampledSurface/sampledSurfaces/sampledSurfaces.C

3
src/sampling/sampledSurface/distanceSurface/distanceSurface.C
View File

@ -30,6 +30,7 @@ License
#include "volPointInterpolation.H"
#include "addToRunTimeSelectionTable.H"
#include "fvMesh.H"
#include "isoSurfaceCell.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@ -110,7 +111,7 @@ void Foam::distanceSurface::createGeometry() const
}
//- Direct from cell field and point field.
const isoSurface iso
const isoSurfaceCell iso
(
mesh(),
cellDistance,

1304
src/sampling/sampledSurface/isoSurface/isoSurface.C
View File

File diff suppressed because it is too large Load Diff

181
src/sampling/sampledSurface/isoSurface/isoSurface.H
View File

@ -27,9 +27,7 @@ Class
Description
A surface formed by the iso value.
After "Polygonising A Scalar Field Using Tetrahedrons", Paul Bourke
and
"Regularised Marching Tetrahedra: improved iso-surface extraction",
After "Regularised Marching Tetrahedra: improved iso-surface extraction",
G.M. Treece, R.W. Prager and A.H. Gee.
Note:
@ -39,7 +37,9 @@ Description
- does tets without using cell centres/cell values. Not tested.
- regularisation can create duplicate triangles/non-manifold surfaces.
Current handling of those is bit ad-hoc for now and not perfect.
Uses geometric merge.
- uses geometric merge with fraction of bounding box as distance.
- triangles can be between two cell centres so constant sampling
does not make sense.
- does not do 2D correctly, creates non-flat iso surface.
SourceFiles
@ -54,13 +54,14 @@ SourceFiles
#include "labelPair.H"
#include "pointIndexHit.H"
#include "PackedList.H"
#include "volFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class polyMesh;
class fvMesh;
/*---------------------------------------------------------------------------*\
Class isoSurface Declaration
@ -88,7 +89,7 @@ class isoSurface
//- Reference to mesh
const polyMesh& mesh_;
const fvMesh& mesh_;
//- Isosurface value
const scalar iso_;
@ -97,6 +98,15 @@ class isoSurface
const scalar mergeDistance_;
//- Whether face might be cut
List<cellCutType> faceCutType_;
//- 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_;
@ -113,16 +123,12 @@ class isoSurface
const scalar s1
) const;
//List<triFace> triangulate(const face& f) const;
//- Determine whether cell is cut
cellCutType calcCutType
//- Set faceCutType,cellCutType.
void calcCutTypes
(
const PackedList<1>&,
const scalarField& cellValues,
const scalarField& pointValues,
const label
) const;
const volScalarField& cVals,
const scalarField& pVals
);
static labelPair findCommonPoints
(
@ -138,79 +144,115 @@ class isoSurface
DynamicList<labelledTri, 64>& localTris
);
void getNeighbour
(
const labelList& boundaryRegion,
const volScalarField& cVals,
const label cellI,
const label faceI,
scalar& nbrValue,
point& nbrPoint
) const;
//- Determine per cc whether all near cuts can be snapped to single
// point.
void calcSnappedCc
(
const PackedList<1>& isTet,
const scalarField& cVals,
const labelList& boundaryRegion,
const volScalarField& cVals,
const scalarField& pVals,
const List<cellCutType>& cellCutType,
DynamicList<point>& triPoints,
DynamicList<point>& snappedPoints,
labelList& snappedCc
) const;
//- Generate triangles for face connected to pointI
void genPointTris
(
const scalarField& cellValues,
const scalarField& pointValues,
const label pointI,
const face& f,
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 cellI,
DynamicList<point, 64>& localTriPoints
) const;
//- Determine per point whether all near cuts can be snapped to single
// point.
void calcSnappedPoint
(
const PackedList<1>& isBoundaryPoint,
const PackedList<1>& isTet,
const scalarField& cVals,
const labelList& boundaryRegion,
const volScalarField& cVals,
const scalarField& pVals,
const List<cellCutType>& cellCutType,
DynamicList<point>& triPoints,
DynamicList<point>& snappedPoints,
labelList& snappedPoint
) const;
//- Generate single point by interpolation or snapping
point generatePoint
template<class Type>
Type generatePoint
(
const DynamicList<point>& snappedPoints,
const DynamicList<Type>& snappedPoints,
const scalar s0,
const point& p0,
const Type& p0,
const label p0Index,
const scalar s1,
const point& p1,
const Type& p1,
const label p1Index
) const;
template<class Type>
void generateTriPoints
(
const DynamicList<point>& snapped,
const DynamicList<Type>& snapped,
const scalar s0,
const point& p0,
const Type& p0,
const label p0Index,
const scalar s1,
const point& p1,
const Type& p1,
const label p1Index,
const scalar s2,
const point& p2,
const Type& p2,
const label p2Index,
const scalar s3,
const point& p3,
const Type& p3,
const label p3Index,
DynamicList<point>& points
DynamicList<Type>& points
) const;
template<class Type>
label generateTriPoints
(
const volScalarField& cVals,
const scalarField& pVals,
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
const Field<Type>& pCoords,
const DynamicList<Type>& snappedPoints,
const labelList& snappedCc,
const labelList& snappedPoint,
const label faceI,
const scalar neiVal,
const Type& neiPt,
const label neiSnap,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
) const;
template<class Type>
void generateTriPoints
(
const volScalarField& cVals,
const scalarField& pVals,
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
const Field<Type>& pCoords,
const DynamicList<Type>& snappedPoints,
const labelList& snappedCc,
const labelList& snappedPoint,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
) const;
triSurface stitchTriPoints
@ -276,12 +318,10 @@ class isoSurface
(
const triSurface& s,
const labelList& newToOldFaces,
labelList& oldToNewPoints,
labelList& newToOldPoints
);
//- Combine all triangles inside a cell into a minimal triangulation
void combineCellTriangles();
public:
//- Runtime type information
@ -290,12 +330,13 @@ public:
// Constructors
//- Construct from dictionary
//- Construct from cell values and point values. Uses boundaryField
// for boundary values. Requires on coupled patchfields to be set
// to the opposite cell value.
isoSurface
(
const polyMesh& mesh,
const scalarField& cellValues,
const scalarField& pointValues,
const volScalarField& cellIsoVals,
const scalarField& pointIsoVals,
const scalar iso,
const bool regularise,
const scalar mergeTol = 1E-6 // fraction of bounding box
@ -311,10 +352,22 @@ public:
}
//- For every unmerged triangle point the point in the triSurface
const labelList triPointMergeMap() const
const labelList& triPointMergeMap() const
{
return triPointMergeMap_;
}
//- Interpolates cCoords,pCoords. Takes the original fields
// used to create the iso surface.
template <class Type>
tmp<Field<Type> > interpolate
(
const volScalarField& cellIsoVals,
const scalarField& pointIsoVals,
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
const Field<Type>& pCoords
) const;
};
@ -324,6 +377,12 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "isoSurfaceTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

2000
src/sampling/sampledSurface/isoSurface/isoSurfaceCell.C Normal file
View File

File diff suppressed because it is too large Load Diff

325
src/sampling/sampledSurface/isoSurface/isoSurfaceCell.H Normal file
View File

@ -0,0 +1,325 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::isoSurfaceCell
Description
A surface formed by the iso value.
After "Polygonising A Scalar Field Using Tetrahedrons", Paul Bourke
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 "PackedList.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_;
//- isoSurfaceCell value
const scalar iso_;
//- When to merge points
const scalar mergeDistance_;
//- 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 inbetween s0,s1
scalar isoFraction
(
const scalar s0,
const scalar s1
) const;
//List<triFace> triangulate(const face& f) const;
//- Determine whether cell is cut
cellCutType calcCutType
(
const PackedList<1>&,
const scalarField& cellValues,
const scalarField& pointValues,
const label
) const;
static labelPair findCommonPoints
(
const labelledTri&,
const labelledTri&
);
static point calcCentre(const triSurface&);
static pointIndexHit collapseSurface
(
pointField& localPoints,
DynamicList<labelledTri, 64>& localTris
);
//- Determine per cc whether all near cuts can be snapped to single
// point.
void calcSnappedCc
(
const PackedList<1>& isTet,
const scalarField& cVals,
const scalarField& pVals,
const List<cellCutType>& cellCutType,
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 face& f,
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 cellI,
DynamicList<point, 64>& localTriPoints
) const;
//- Determine per point whether all near cuts can be snapped to single
// point.
void calcSnappedPoint
(
const PackedList<1>& isBoundaryPoint,
const PackedList<1>& isTet,
const scalarField& cVals,
const scalarField& pVals,
const List<cellCutType>& cellCutType,
DynamicList<point>& triPoints,
labelList& snappedPoint
) const;
//- Generate single point by interpolation or snapping
point generatePoint
(
const DynamicList<point>& snappedPoints,
const scalar s0,
const point& p0,
const label p0Index,
const scalar s1,
const point& p1,
const label p1Index
) const;
void generateTriPoints
(
const DynamicList<point>& snapped,
const scalar s0,
const point& p0,
const label p0Index,
const scalar s1,
const point& p1,
const label p1Index,
const scalar s2,
const point& p2,
const label p2Index,
const scalar s3,
const point& p3,
const label p3Index,
DynamicList<point>& points
) 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;
//- Determine orientation
static void walkOrientation
(
const triSurface& surf,
const List<FixedList<label, 3> >& faceEdges,
const labelList& edgeFace0,
const labelList& edgeFace1,
const label seedTriI,
labelList& flipState
);
//- Orient surface
static void orientSurface
(
triSurface&,
const List<FixedList<label, 3> >& faceEdges,
const labelList& edgeFace0,
const labelList& edgeFace1,
const Map<labelList>& edgeFacesRest
);
//- 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
);
//- Combine all triangles inside a cell into a minimal triangulation
void combineCellTriangles();
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_;
}
//- For every unmerged triangle point the point in the triSurface
const labelList triPointMergeMap() const
{
return triPointMergeMap_;
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

485
src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C Normal file
View File

@ -0,0 +1,485 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "isoSurface.H"
#include "polyMesh.H"
#include "syncTools.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
Type Foam::isoSurface::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::isoSurface::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
{
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 0x0E:
case 0x01:
points.append(generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index));
points.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
break;
case 0x0D:
case 0x02:
points.append(generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index));
points.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
break;
case 0x0C:
case 0x03:
{
Type tp1 = generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index);
Type tp2 = generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index);
points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
points.append(tp1);
points.append(tp2);
points.append(tp2);
points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
points.append(tp1);
}
break;
case 0x0B:
case 0x04:
{
points.append(generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index));
points.append(generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index));
points.append(generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index));
}
break;
case 0x0A:
case 0x05:
{
Type tp0 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
Type tp1 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
points.append(tp0);
points.append(tp1);
points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
points.append(tp0);
points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
points.append(tp1);
}
break;
case 0x09:
case 0x06:
{
Type tp0 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
Type tp1 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
points.append(tp0);
points.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
points.append(tp1);
points.append(tp0);
points.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
points.append(tp1);
}
break;
case 0x07:
case 0x08:
points.append(generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index));
points.append(generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index));
points.append(generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index));
break;
}
}
template<class Type>
Foam::label Foam::isoSurface::generateTriPoints
(
const volScalarField& cVals,
const scalarField& pVals,
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
const Field<Type>& pCoords,
const DynamicList<Type>& snappedPoints,
const labelList& snappedCc,
const labelList& snappedPoint,
const label faceI,
const scalar neiVal,
const Type& neiPt,
const label neiSnap,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
) const
{
label own = mesh_.faceOwner()[faceI];
label oldNPoints = triPoints.size();
const face& f = mesh_.faces()[faceI];
forAll(f, fp)
{
label pointI = f[fp];
label nextPointI = f[f.fcIndex(fp)];
generateTriPoints
(
snappedPoints,
pVals[pointI],
pCoords[pointI],
snappedPoint[pointI],
pVals[nextPointI],
pCoords[nextPointI],
snappedPoint[nextPointI],
cVals[own],
cCoords[own],
snappedCc[own],
neiVal,
neiPt,
neiSnap,
triPoints
);
}
// Every three triPoints is a triangle
label nTris = (triPoints.size()-oldNPoints)/3;
for (label i = 0; i < nTris; i++)
{
triMeshCells.append(own);
}
return nTris;
}
template<class Type>
void Foam::isoSurface::generateTriPoints
(
const volScalarField& cVals,
const scalarField& pVals,
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
const Field<Type>& pCoords,
const DynamicList<Type>& snappedPoints,
const labelList& snappedCc,
const labelList& snappedPoint,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
) const
{
const polyBoundaryMesh& patches = mesh_.boundaryMesh();
const labelList& own = mesh_.faceOwner();
const labelList& nei = mesh_.faceNeighbour();
// Determine neighbouring snap status
labelList neiSnappedCc(mesh_.nFaces()-mesh_.nInternalFaces(), -1);
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
if (pp.coupled())
{
label faceI = pp.start();
forAll(pp, i)
{
neiSnappedCc[faceI-mesh_.nInternalFaces()] =
snappedCc[own[faceI]];
faceI++;
}
}
}
syncTools::swapBoundaryFaceList(mesh_, neiSnappedCc, false);
// Generate triangle points
triPoints.clear();
triMeshCells.clear();
for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
{
if (faceCutType_[faceI] != NOTCUT)
{
generateTriPoints
(
cVals,
pVals,
cCoords,
pCoords,
snappedPoints,
snappedCc,
snappedPoint,
faceI,
cVals[nei[faceI]],
cCoords[nei[faceI]],
snappedCc[nei[faceI]],
triPoints,
triMeshCells
);
}
}
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
if (pp.coupled())
{
if (refCast<const processorPolyPatch>(pp).owner())
{
label faceI = pp.start();
forAll(pp, i)
{
if (faceCutType_[faceI] != NOTCUT)
{
generateTriPoints
(
cVals,
pVals,
cCoords,
pCoords,
snappedPoints,
snappedCc,
snappedPoint,
faceI,
cVals.boundaryField()[patchI][i],
cCoords.boundaryField()[patchI][i],
neiSnappedCc[faceI-mesh_.nInternalFaces()],
triPoints,
triMeshCells
);
}
faceI++;
}
}
}
else
{
label faceI = pp.start();
forAll(pp, i)
{
if (faceCutType_[faceI] != NOTCUT)
{
generateTriPoints
(
cVals,
pVals,
cCoords,
pCoords,
snappedPoints,
snappedCc,
snappedPoint,
faceI,
cVals.boundaryField()[patchI][i],
cCoords.boundaryField()[patchI][i],
-1, // fc not snapped
triPoints,
triMeshCells
);
}
faceI++;
}
}
}
triPoints.shrink();
triMeshCells.shrink();
}
//template <class Type>
//Foam::tmp<Foam::Field<Type> >
//Foam::isoSurface::sample(const Field<Type>& vField) const
//{
// return tmp<Field<Type> >(new Field<Type>(vField, meshCells()));
//}
//
//
template <class Type>
Foam::tmp<Foam::Field<Type> >
Foam::isoSurface::interpolate
(
const volScalarField& cVals,
const scalarField& pVals,
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
const Field<Type>& pCoords
) const
{
DynamicList<Type> triPoints(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
);
// One value per point
tmp<Field<Type> > tvalues(new Field<Type>(points().size()));
Field<Type>& values = tvalues();
forAll(triPoints, i)
{
label mergedPointI = triPointMergeMap_[i];
if (mergedPointI >= 0)
{
values[mergedPointI] = triPoints[i];
}
}
return tvalues;
}
// ************************************************************************* //

218
src/sampling/sampledSurface/isoSurface/sampledIsoSurface.C
View File

@ -41,32 +41,22 @@ namespace Foam
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::sampledIsoSurface::createGeometry() const
void Foam::sampledIsoSurface::getIsoFields() const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
if (fvm.time().timeIndex() != storedTimeIndex_)
{
storedTimeIndex_ = fvm.time().timeIndex();
// Get volField
// ~~~~~~~~~~~~
// Clear any stored topo
facesPtr_.clear();
// Optionally read volScalarField
autoPtr<volScalarField> readFieldPtr_;
// 1. see if field in database
// 2. see if field can be read
const volScalarField* cellFldPtr = NULL;
if (fvm.foundObject<volScalarField>(isoField_))
{
if (debug)
{
Info<< "sampledIsoSurface::createGeometry() : lookup "
Info<< "sampledIsoSurface::getIsoField() : lookup "
<< isoField_ << endl;
}
cellFldPtr = &fvm.lookupObject<volScalarField>(isoField_);
storedVolFieldPtr_.clear();
volFieldPtr_ = &fvm.lookupObject<volScalarField>(isoField_);
}
else
{
@ -74,12 +64,12 @@ void Foam::sampledIsoSurface::createGeometry() const
if (debug)
{
Info<< "sampledIsoSurface::createGeometry() : reading "
Info<< "sampledIsoSurface::getIsoField() : reading "
<< isoField_ << " from time " <<fvm.time().timeName()
<< endl;
}
readFieldPtr_.reset
storedVolFieldPtr_.reset
(
new volScalarField
(
@ -95,70 +85,148 @@ void Foam::sampledIsoSurface::createGeometry() const
fvm
)
);
cellFldPtr = readFieldPtr_.operator->();
volFieldPtr_ = storedVolFieldPtr_.operator->();
}
const volScalarField& cellFld = *cellFldPtr;
tmp<pointScalarField> pointFld
// Get pointField
// ~~~~~~~~~~~~~~
word pointFldName = "volPointInterpolate(" + isoField_ + ')';
if (fvm.foundObject<pointScalarField>(pointFldName))
{
if (debug)
{
Info<< "sampledIsoSurface::getIsoField() : lookup "
<< pointFldName << endl;
}
storedPointFieldPtr_.clear();
pointFieldPtr_ = &fvm.lookupObject<pointScalarField>(pointFldName);
}
else
{
// Not in registry. Interpolate.
if (debug)
{
Info<< "sampledIsoSurface::getIsoField() : interpolating "
<< pointFldName << endl;
}
storedPointFieldPtr_.reset
(
volPointInterpolation::New(fvm).interpolate(cellFld)
volPointInterpolation::New(fvm).interpolate(*volFieldPtr_).ptr()
);
pointFieldPtr_ = storedPointFieldPtr_.operator->();
}
//- Direct from cell field and point field. Gives bad continuity.
const isoSurface iso
if (debug)
{
Info<< "sampledIsoSurface::getIsoField() : obtained volField "
<< volFieldPtr_->name() << " min:" << min(*volFieldPtr_).value()
<< " max:" << max(*volFieldPtr_).value() << endl;
Info<< "sampledIsoSurface::getIsoField() : obtained pointField "
<< pointFieldPtr_->name()
<< " min:" << gMin(pointFieldPtr_->internalField())
<< " max:" << gMax(pointFieldPtr_->internalField()) << endl;
}
}
void Foam::sampledIsoSurface::createGeometry() const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
if (fvm.time().timeIndex() != storedTimeIndex_)
{
storedTimeIndex_ = fvm.time().timeIndex();
getIsoFields();
// Clear any stored topo
surfPtr_.clear();
facesPtr_.clear();
if (average_)
{
//- From point field and interpolated cell.
volScalarField cellAvg
(
IOobject
(
"cellAvg",
fvm.time().timeName(),
fvm.time(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
fvm,
cellFld.internalField(),
pointFld().internalField(),
dimensionedScalar("zero", dimless, scalar(0.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] += (*pointFieldPtr_)[pointI];
nPointCells[cellI]++;
}
}
}
forAll(cellAvg, cellI)
{
cellAvg[cellI] /= nPointCells[cellI];
}
// Give value to calculatedFvPatchFields
cellAvg.correctBoundaryConditions();
surfPtr_.reset
(
new isoSurface
(
cellAvg,
*pointFieldPtr_,
isoVal_,
regularise_
)
);
}
else
{
//- Direct from cell field and point field.
surfPtr_.reset
(
new isoSurface
(
*volFieldPtr_,
*pointFieldPtr_,
isoVal_,
regularise_
)
);
}
////- From point field and interpolated cell.
//scalarField cellAvg(fvm.nCells(), scalar(0.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().internalField()[pointI];
// nPointCells[cellI]++;
// }
// }
//}
//forAll(cellAvg, cellI)
//{
// cellAvg[cellI] /= nPointCells[cellI];
//}
//
//const isoSurface iso
//(
// fvm,
// cellAvg,
// pointFld().internalField(),
// isoVal_,
// regularise_
//);
const_cast<sampledIsoSurface&>(*this).triSurface::operator=(iso);
meshCells_ = iso.meshCells();
if (debug)
{
Pout<< "sampledIsoSurface::createGeometry() : 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;
<< " tris : " << surface().size() << nl
<< " cut cells : " << surface().meshCells().size()
<< endl;
}
}
}
@ -177,11 +245,26 @@ Foam::sampledIsoSurface::sampledIsoSurface
isoField_(dict.lookup("isoField")),
isoVal_(readScalar(dict.lookup("isoValue"))),
regularise_(dict.lookupOrDefault("regularise", true)),
average_(dict.lookupOrDefault("average", false)),
zoneName_(word::null),
surfPtr_(NULL),
facesPtr_(NULL),
storedTimeIndex_(-1),
meshCells_(0)
storedVolFieldPtr_(NULL),
volFieldPtr_(NULL),
storedPointFieldPtr_(NULL),
pointFieldPtr_(NULL)
{
if (!sampledSurface::interpolate())
{
FatalErrorIn
(
"sampledIsoSurface::sampledIsoSurface"
"(const word&, const polyMesh&, const dictionary&)"
) << "Non-interpolated iso surface not supported since triangles"
<< " span across cells." << exit(FatalError);
}
// label zoneId = -1;
// if (dict.readIfPresent("zone", zoneName_))
// {
@ -209,6 +292,7 @@ void Foam::sampledIsoSurface::correct(const bool meshChanged)
// Only change of mesh changes plane - zone restriction gets lost
if (meshChanged)
{
surfPtr_.clear();
facesPtr_.clear();
}
}
@ -317,9 +401,9 @@ void Foam::sampledIsoSurface::print(Ostream& os) const
{
os << "sampledIsoSurface: " << name() << " :"
<< " field:" << isoField_
<< " value:" << isoVal_
<< " faces:" << faces().size()
<< " points:" << points().size();
<< " value:" << isoVal_;
//<< " faces:" << faces().size() // note: possibly no geom yet
//<< " points:" << points().size();
}

37
src/sampling/sampledSurface/isoSurface/sampledIsoSurface.H
View File

@ -39,7 +39,7 @@ SourceFiles
#define sampledIsoSurface_H
#include "sampledSurface.H"
#include "triSurface.H"
#include "isoSurface.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -52,8 +52,7 @@ namespace Foam
class sampledIsoSurface
:
public sampledSurface,
public triSurface
public sampledSurface
{
// Private data
@ -66,9 +65,14 @@ class sampledIsoSurface
//- Whether to coarse
const Switch regularise_;
//- Whether to recalculate cell values as average of point values
const Switch average_;
//- zone name (if restricted to zones)
word zoneName_;
mutable autoPtr<isoSurface> surfPtr_;
//- triangles converted to faceList
mutable autoPtr<faceList> facesPtr_;
@ -78,12 +82,21 @@ class sampledIsoSurface
//- Time at last call
mutable label storedTimeIndex_;
//- For every triangle the original cell in mesh
mutable labelList meshCells_;
//- Cached volfield
mutable autoPtr<volScalarField> storedVolFieldPtr_;
mutable const volScalarField* volFieldPtr_;
//- Cached pointfield
mutable autoPtr<pointScalarField> storedPointFieldPtr_;
mutable const pointScalarField* pointFieldPtr_;
// Private Member Functions
//- Get fields needed to recreate iso surface.
void getIsoFields() const;
//- Create iso surface (if time has changed)
void createGeometry() const;
@ -127,7 +140,7 @@ public:
//- Points of surface
virtual const pointField& points() const
{
return triSurface::points();
return surface().points();
}
//- Faces of surface
@ -135,7 +148,7 @@ public:
{
if (!facesPtr_.valid())
{
const triSurface& s = *this;
const triSurface& s = surface();
facesPtr_.reset(new faceList(s.size()));
@ -147,11 +160,19 @@ public:
return facesPtr_;
}
//- Correct for mesh movement and/or field changes
virtual void correct(const bool meshChanged);
const isoSurface& surface() const
{
return surfPtr_();
}
//- Lookup or read isoField. Sets volFieldPtr_ and pointFieldPtr_.
void getIsoField();
//- sample field on surface
virtual tmp<scalarField> sample
(

344
src/sampling/sampledSurface/isoSurface/sampledIsoSurfaceCell.C Normal file
View File

@ -0,0 +1,344 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#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 * * * * * * * * * * * //
void Foam::sampledIsoSurfaceCell::createGeometry() const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
if (fvm.time().timeIndex() != storedTimeIndex_)
{
storedTimeIndex_ = fvm.time().timeIndex();
// Clear any stored topo
facesPtr_.clear();
// Optionally read volScalarField
autoPtr<volScalarField> readFieldPtr_;
// 1. see if field in database
// 2. see if field can be read
const volScalarField* cellFldPtr = NULL;
if (fvm.foundObject<volScalarField>(isoField_))
{
if (debug)
{
Info<< "sampledIsoSurfaceCell::createGeometry() : lookup "
<< isoField_ << endl;
}
cellFldPtr = &fvm.lookupObject<volScalarField>(isoField_);
}
else
{
// Bit of a hack. Read field and store.
if (debug)
{
Info<< "sampledIsoSurfaceCell::createGeometry() : 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.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().internalField()[pointI];
nPointCells[cellI]++;
}
}
}
forAll(cellAvg, cellI)
{
cellAvg[cellI] /= nPointCells[cellI];
}
const isoSurfaceCell iso
(
fvm,
cellAvg,
pointFld().internalField(),
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.internalField(),
pointFld().internalField(),
isoVal_,
regularise_
);
const_cast<sampledIsoSurfaceCell&>
(
*this
).triSurface::operator=(iso);
meshCells_ = iso.meshCells();
}
if (debug)
{
Pout<< "sampledIsoSurfaceCell::createGeometry() : 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;
}
}
}
// * * * * * * * * * * * * * * * * 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)),
zoneName_(word::null),
facesPtr_(NULL),
storedTimeIndex_(-1),
meshCells_(0)
{
// label zoneId = -1;
// if (dict.readIfPresent("zone", zoneName_))
// {
// zoneId = mesh.cellZones().findZoneID(zoneName_);
// if (debug && zoneId < 0)
// {
// Info<< "cellZone \"" << zoneName_
// << "\" not found - using entire mesh"
// << endl;
// }
// }
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::sampledIsoSurfaceCell::~sampledIsoSurfaceCell()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::sampledIsoSurfaceCell::correct(const bool meshChanged)
{
// Only change of mesh changes plane - zone restriction gets lost
if (meshChanged)
{
facesPtr_.clear();
}
}
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();
}
// ************************************************************************* //

238
src/sampling/sampledSurface/isoSurface/sampledIsoSurfaceCell.H Normal file
View File

@ -0,0 +1,238 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
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_;
//- zone name (if restricted to zones)
word zoneName_;
//- triangles converted to faceList
mutable autoPtr<faceList> facesPtr_;
// Recreated for every isoSurface
//- Time at last call
mutable label storedTimeIndex_;
//- For every triangle the original cell in mesh
mutable labelList meshCells_;
// Private Member Functions
//- Create iso surface (if time has changed)
void createGeometry() 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
//- Points of surface
virtual const pointField& points() const
{
return triSurface::points();
}
//- Faces of surface
virtual const faceList& faces() const
{
if (!facesPtr_.valid())
{
const triSurface& s = *this;
facesPtr_.reset(new faceList(s.size()));
forAll(s, i)
{
facesPtr_()[i] = s[i].triFaceFace();
}
}
return facesPtr_;
}
//- Correct for mesh movement and/or field changes
virtual void correct(const bool meshChanged);
//- 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
// ************************************************************************* //

89
src/sampling/sampledSurface/isoSurface/sampledIsoSurfaceCellTemplates.C Normal file
View File

@ -0,0 +1,89 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#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
createGeometry();
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
createGeometry();
// One value per point
tmp<Field<Type> > tvalues(new Field<Type>(points().size()));
Field<Type>& values = tvalues();
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;
}
// ************************************************************************* //

46
src/sampling/sampledSurface/isoSurface/sampledIsoSurfaceTemplates.C
View File

@ -25,7 +25,6 @@ License
\*---------------------------------------------------------------------------*/
#include "sampledIsoSurface.H"
#include "isoSurface.H"
#include "volFieldsFwd.H"
#include "pointFields.H"
#include "volPointInterpolation.H"
@ -41,8 +40,7 @@ Foam::sampledIsoSurface::sampleField
{
// Recreate geometry if time has changed
createGeometry();
return tmp<Field<Type> >(new Field<Type>(vField, meshCells_));
return tmp<Field<Type> >(new Field<Type>(vField, surface().meshCells()));
}
@ -53,36 +51,28 @@ Foam::sampledIsoSurface::interpolateField
const interpolation<Type>& interpolator
) const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
// Get fields to sample. Assume volPointInterpolation!
const GeometricField<Type, fvPatchField, volMesh>& volFld =
interpolator.psi();
tmp<GeometricField<Type, pointPatchField, pointMesh> > pointFld
(
volPointInterpolation::New(fvm).interpolate(volFld)
);
// Recreate geometry if time has changed
createGeometry();
// One value per point
tmp<Field<Type> > tvalues(new Field<Type>(points().size()));
Field<Type>& values = tvalues();
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
// Sample.
return surface().interpolate
(
points()[pointI],
meshCells_[cutFaceI]
*volFieldPtr_,
*pointFieldPtr_,
volFld,
pointFld()
);
pointDone[pointI] = true;
}
}
}
return tvalues;
}