Commit 1a87cd16 authored by mattijs's avatar mattijs
Browse files

ENH: Add averaging to all iso-surface based sampledSurfaces.

The optional 'average' switch causes use of the average-of-pointvalues
instead of the original volField.
parent fcc442e3
......@@ -284,6 +284,7 @@ Foam::distanceSurface::distanceSurface
distance_(readScalar(dict.lookup("distance"))),
signed_(readBool(dict.lookup("signed"))),
regularise_(dict.lookupOrDefault("regularise", true)),
average_(dict.lookupOrDefault("average", false)),
zoneName_(word::null),
needsUpdate_(true),
isoSurfPtr_(NULL),
......
......@@ -68,6 +68,9 @@ class distanceSurface
//- Whether to coarsen
const Switch regularise_;
//- Whether to recalculate cell values as average of point values
const Switch average_;
//- zone name (if restricted to zones)
word zoneName_;
......
......@@ -62,7 +62,15 @@ Foam::distanceSurface::interpolateField
);
// Sample.
return surface().interpolate(volFld, pointFld());
return surface().interpolate
(
(
average_
? pointAverage(pointFld())()
: volFld
),
pointFld()
);
}
......
......@@ -164,7 +164,10 @@ void Foam::sampledIsoSurface::getIsoFields() const
// point field.
if (average_)
{
storedVolFieldPtr_.reset(average(fvm, *pointFieldPtr_).ptr());
storedVolFieldPtr_.reset
(
pointAverage(*pointFieldPtr_).ptr()
);
volFieldPtr_ = storedVolFieldPtr_.operator->();
}
......@@ -265,7 +268,7 @@ void Foam::sampledIsoSurface::getIsoFields() const
{
storedVolSubFieldPtr_.reset
(
average(subFvm, *pointSubFieldPtr_).ptr()
pointAverage(*pointSubFieldPtr_).ptr()
);
volSubFieldPtr_ = storedVolSubFieldPtr_.operator->();
}
......@@ -286,99 +289,6 @@ void Foam::sampledIsoSurface::getIsoFields() const
}
Foam::tmp<Foam::volScalarField> Foam::sampledIsoSurface::average
(
const fvMesh& mesh,
const pointScalarField& pfld
) const
{
tmp<volScalarField> tcellAvg
(
new volScalarField
(
IOobject
(
"cellAvg",
mesh.time().timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh,
dimensionedScalar("zero", dimless, scalar(0.0))
)
);
volScalarField& cellAvg = tcellAvg();
labelField nPointCells(mesh.nCells(), 0);
{
for (label pointI = 0; pointI < mesh.nPoints(); pointI++)
{
const labelList& pCells = mesh.pointCells(pointI);
forAll(pCells, i)
{
label cellI = pCells[i];
cellAvg[cellI] += pfld[pointI];
nPointCells[cellI]++;
}
}
}
forAll(cellAvg, cellI)
{
cellAvg[cellI] /= nPointCells[cellI];
}
// Give value to calculatedFvPatchFields
cellAvg.correctBoundaryConditions();
return tcellAvg;
}
Foam::tmp<Foam::pointScalarField> Foam::sampledIsoSurface::average
(
const pointMesh& pMesh,
const volScalarField& fld
) const
{
tmp<pointScalarField> tpointAvg
(
new pointScalarField
(
IOobject
(
"pointAvg",
fld.time().timeName(),
fld.db(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
pMesh,
dimensionedScalar("zero", dimless, scalar(0.0))
)
);
pointScalarField& pointAvg = tpointAvg();
for (label pointI = 0; pointI < fld.mesh().nPoints(); pointI++)
{
const labelList& pCells = fld.mesh().pointCells(pointI);
forAll(pCells, i)
{
pointAvg[pointI] += fld[pCells[i]];
}
pointAvg[pointI] /= pCells.size();
}
// Give value to calculatedFvPatchFields
pointAvg.correctBoundaryConditions();
return tpointAvg;
}
bool Foam::sampledIsoSurface::updateGeometry() const
{
const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
......
......@@ -118,18 +118,6 @@ class sampledIsoSurface
//- Get fields needed to recreate iso surface.
void getIsoFields() const;
tmp<volScalarField> average
(
const fvMesh&,
const pointScalarField&
) const;
tmp<pointScalarField> average
(
const pointMesh&,
const volScalarField& fld
) const;
//- Create iso surface (if time has changed)
// Do nothing (and return false) if no update was needed
bool updateGeometry() const;
......
......@@ -71,7 +71,15 @@ Foam::sampledIsoSurface::interpolateField
volPointInterpolation::New(volSubFld.mesh()).interpolate(volSubFld);
// Sample.
return surface().interpolate(volSubFld, tpointSubFld());
return surface().interpolate
(
(
average_
? pointAverage(tpointSubFld())()
: volSubFld
),
tpointSubFld()
);
}
else
{
......@@ -79,7 +87,15 @@ Foam::sampledIsoSurface::interpolateField
volPointInterpolation::New(volFld.mesh()).interpolate(volFld);
// Sample.
return surface().interpolate(volFld, tpointFld());
return surface().interpolate
(
(
average_
? pointAverage(tpointFld())()
: volFld
),
tpointFld()
);
}
}
......
......@@ -254,6 +254,7 @@ Foam::sampledCuttingPlane::sampledCuttingPlane
plane_(dict),
mergeTol_(dict.lookupOrDefault("mergeTol", 1E-6)),
regularise_(dict.lookupOrDefault("regularise", true)),
average_(dict.lookupOrDefault("average", false)),
zoneID_(dict.lookupOrDefault("zone", word::null), mesh.cellZones()),
exposedPatchName_(word::null),
needsUpdate_(true),
......
......@@ -66,6 +66,9 @@ class sampledCuttingPlane
//- Whether to coarsen
const Switch regularise_;
//- Whether to recalculate cell values as average of point values
const Switch average_;
//- zone name/index (if restricted to zones)
mutable cellZoneID zoneID_;
......
......@@ -67,7 +67,15 @@ Foam::sampledCuttingPlane::interpolateField
volPointInterpolation::New(volSubFld.mesh()).interpolate(volSubFld);
// Sample.
return surface().interpolate(volSubFld, tpointSubFld());
return surface().interpolate
(
(
average_
? pointAverage(tpointSubFld())()
: volSubFld
),
tpointSubFld()
);
}
else
{
......@@ -77,7 +85,15 @@ Foam::sampledCuttingPlane::interpolateField
);
// Sample.
return surface().interpolate(volFld, tpointFld());
return surface().interpolate
(
(
average_
? pointAverage(tpointFld())()
: volFld
),
tpointFld()
);
}
}
......
......@@ -312,6 +312,13 @@ public:
//- Project field onto surface
tmp<Field<vector> > project(const Field<tensor>&) const;
//- Interpolate from points to cell centre
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh> > pointAverage
(
const GeometricField<Type, pointPatchField, pointMesh>& pfld
) const;
//- Sample field on surface
virtual tmp<scalarField> sample
(
......
......@@ -155,4 +155,59 @@ Foam::sampledSurface::project
}
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh> >
Foam::sampledSurface::pointAverage
(
const GeometricField<Type, pointPatchField, pointMesh>& pfld
) const
{
const fvMesh& mesh = dynamic_cast<const fvMesh&>(pfld.mesh()());
tmp<GeometricField<Type, fvPatchField, volMesh> > tcellAvg
(
new GeometricField<Type, fvPatchField, volMesh>
(
IOobject
(
"cellAvg",
mesh.time().timeName(),
pfld.db(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh,
dimensioned<Type>("zero", dimless, pTraits<Type>::zero)
)
);
GeometricField<Type, fvPatchField, volMesh>& cellAvg = tcellAvg();
labelField nPointCells(mesh.nCells(), 0);
{
for (label pointI = 0; pointI < mesh.nPoints(); pointI++)
{
const labelList& pCells = mesh.pointCells(pointI);
forAll(pCells, i)
{
label cellI = pCells[i];
cellAvg[cellI] += pfld[pointI];
nPointCells[cellI]++;
}
}
}
forAll(cellAvg, cellI)
{
cellAvg[cellI] /= nPointCells[cellI];
}
// Give value to calculatedFvPatchFields
cellAvg.correctBoundaryConditions();
return tcellAvg;
}
// ************************************************************************* //
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment