diff --git a/applications/utilities/postProcessing/sampling/sample/sampleDict b/applications/utilities/postProcessing/sampling/sample/sampleDict index 1616bbb662..b605f29e57 100644 --- a/applications/utilities/postProcessing/sampling/sample/sampleDict +++ b/applications/utilities/postProcessing/sampling/sample/sampleDict @@ -177,6 +177,25 @@ surfaces interpolate false; regularise false; // do not simplify } + + triangleCut + { + // Cutingplaneusing iso surface + type cuttingPlane; + planeType pointAndNormal; + pointAndNormalDict + { + basePoint (0.4 0 0.4); + normalVector (1 0.2 0.2); + } + interpolate true; + + //zone ABC; // Optional: zone only + //exposedPatchName fixedWalls; // Optional: zone only + + // regularise false; // Optional: do not simplify + } + ); diff --git a/src/sampling/Make/files b/src/sampling/Make/files index 43bc910a3c..004f81d435 100644 --- a/src/sampling/Make/files +++ b/src/sampling/Make/files @@ -29,6 +29,7 @@ sampledSurface/isoSurface/sampledIsoSurface.C sampledSurface/isoSurface/isoSurfaceCell.C sampledSurface/isoSurface/sampledIsoSurfaceCell.C sampledSurface/distanceSurface/distanceSurface.C +sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C sampledSurface/sampledSurface/sampledSurface.C sampledSurface/sampledSurfaces/sampledSurfaces.C sampledSurface/sampledSurfacesFunctionObject/sampledSurfacesFunctionObject.C diff --git a/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C b/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C new file mode 100644 index 0000000000..85dbeca595 --- /dev/null +++ b/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C @@ -0,0 +1,420 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 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 "sampledCuttingPlane.H" +#include "dictionary.H" +#include "volFields.H" +#include "volPointInterpolation.H" +#include "addToRunTimeSelectionTable.H" +#include "fvMesh.H" +#include "isoSurface.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(sampledCuttingPlane, 0); + addNamedToRunTimeSelectionTable(sampledSurface, sampledCuttingPlane, word, cuttingPlane); +} + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::sampledCuttingPlane::createGeometry() +{ + if (debug) + { + Pout<< "sampledCuttingPlane::createGeometry :updating geometry." + << endl; + } + + // Clear any stored topologies + facesPtr_.clear(); + isoSurfPtr_.ptr(); + pointDistance_.clear(); + cellDistancePtr_.clear(); + + + // Get any subMesh + if (zoneID_.index() != -1 && !subMeshPtr_.valid()) + { + const polyBoundaryMesh& patches = mesh().boundaryMesh(); + + // Patch to put exposed internal faces into + label exposedPatchI = patches.findPatchID(exposedPatchName_); + + if (debug) + { + Info<< "Allocating subset of size " + << mesh().cellZones()[zoneID_.index()].size() + << " with exposed faces into patch " + << patches[exposedPatchI].name() << endl; + } + + subMeshPtr_.reset + ( + new fvMeshSubset(static_cast(mesh())) + ); + subMeshPtr_().setLargeCellSubset + ( + labelHashSet(mesh().cellZones()[zoneID_.index()]), + exposedPatchI + ); + } + + + // Select either the submesh or the underlying mesh + const fvMesh& fvm = + ( + subMeshPtr_.valid() + ? subMeshPtr_().subMesh() + : static_cast(mesh()) + ); + + + // Distance to cell centres + // ~~~~~~~~~~~~~~~~~~~~~~~~ + + cellDistancePtr_.reset + ( + new volScalarField + ( + IOobject + ( + "cellDistance", + fvm.time().timeName(), + fvm.time(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + fvm, + dimensionedScalar("zero", dimLength, 0) + ) + ); + volScalarField& cellDistance = cellDistancePtr_(); + + // Internal field + { + const pointField& cc = fvm.C(); + scalarField& fld = cellDistance.internalField(); + + forAll(cc, i) + { + // Signed distance + fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal(); + } + } + + // Patch fields + { + forAll(fvm.C().boundaryField(), patchI) + { + const pointField& cc = fvm.C().boundaryField()[patchI]; + fvPatchScalarField& fld = cellDistance.boundaryField()[patchI]; + + forAll(fld, i) + { + fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal(); + } + } + } + + + // On processor patches the mesh.C() will already be the cell centre + // on the opposite side so no need to swap cellDistance. + + + // Distance to points + pointDistance_.setSize(fvm.nPoints()); + { + const pointField& pts = fvm.points(); + + forAll(pointDistance_, i) + { + pointDistance_[i] = (pts[i] - plane_.refPoint()) & plane_.normal(); + } + } + + + if (debug) + { + Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl; + cellDistance.write(); + pointScalarField pDist + ( + IOobject + ( + "pointDistance", + fvm.time().timeName(), + fvm.time(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + pointMesh::New(fvm), + dimensionedScalar("zero", dimLength, 0) + ); + pDist.internalField() = pointDistance_; + + Pout<< "Writing point distance:" << pDist.objectPath() << endl; + pDist.write(); + } + + + //- Direct from cell field and point field. + isoSurfPtr_.reset + ( + new isoSurface + ( + cellDistance, + pointDistance_, + 0.0, + regularise_ + ) + ); + + if (debug) + { + print(Pout); + Pout<< endl; + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::sampledCuttingPlane::sampledCuttingPlane +( + const word& name, + const polyMesh& mesh, + const dictionary& dict +) +: + sampledSurface(name, mesh, dict), + plane_(dict), + mergeTol_(dict.lookupOrDefault("mergeTol", 1E-6)), + regularise_(dict.lookupOrDefault("regularise", true)), + zoneID_(dict.lookupOrDefault("zone", word::null), mesh.cellZones()), + exposedPatchName_(word::null), + needsUpdate_(true), + subMeshPtr_(NULL), + cellDistancePtr_(NULL), + isoSurfPtr_(NULL), + facesPtr_(NULL) +{ + if (zoneID_.index() != -1) + { + dict.lookup("exposedPatchName") >> exposedPatchName_; + + if (mesh.boundaryMesh().findPatchID(exposedPatchName_) == -1) + { + FatalErrorIn + ( + "sampledCuttingPlane::sampledCuttingPlane" + "(const word&, const polyMesh&, const dictionary&)" + ) << "Cannot find patch " << exposedPatchName_ + << " in which to put exposed faces." << endl + << "Valid patches are " << mesh.boundaryMesh().names() + << exit(FatalError); + } + + if (debug && zoneID_.index() != -1) + { + Info<< "Restricting to cellZone " << zoneID_.name() + << " with exposed internal faces into patch " + << exposedPatchName_ << endl; + } + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::sampledCuttingPlane::~sampledCuttingPlane() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +bool Foam::sampledCuttingPlane::needsUpdate() const +{ + return needsUpdate_; +} + + +bool Foam::sampledCuttingPlane::expire() +{ + if (debug) + { + Pout<< "sampledCuttingPlane::expire :" + << " have-facesPtr_:" << facesPtr_.valid() + << " needsUpdate_:" << needsUpdate_ << endl; + } + + // Clear any stored topologies + facesPtr_.clear(); + + // already marked as expired + if (needsUpdate_) + { + return false; + } + + needsUpdate_ = true; + return true; +} + + +bool Foam::sampledCuttingPlane::update() +{ + if (debug) + { + Pout<< "sampledCuttingPlane::update :" + << " have-facesPtr_:" << facesPtr_.valid() + << " needsUpdate_:" << needsUpdate_ << endl; + } + + if (!needsUpdate_) + { + return false; + } + + createGeometry(); + + needsUpdate_ = false; + return true; +} + + +Foam::tmp +Foam::sampledCuttingPlane::sample +( + const volScalarField& vField +) const +{ + return sampleField(vField); +} + + +Foam::tmp +Foam::sampledCuttingPlane::sample +( + const volVectorField& vField +) const +{ + return sampleField(vField); +} + + +Foam::tmp +Foam::sampledCuttingPlane::sample +( + const volSphericalTensorField& vField +) const +{ + return sampleField(vField); +} + + +Foam::tmp +Foam::sampledCuttingPlane::sample +( + const volSymmTensorField& vField +) const +{ + return sampleField(vField); +} + + +Foam::tmp +Foam::sampledCuttingPlane::sample +( + const volTensorField& vField +) const +{ + return sampleField(vField); +} + + +Foam::tmp +Foam::sampledCuttingPlane::interpolate +( + const interpolation& interpolator +) const +{ + return interpolateField(interpolator); +} + + +Foam::tmp +Foam::sampledCuttingPlane::interpolate +( + const interpolation& interpolator +) const +{ + return interpolateField(interpolator); +} + +Foam::tmp +Foam::sampledCuttingPlane::interpolate +( + const interpolation& interpolator +) const +{ + return interpolateField(interpolator); +} + + +Foam::tmp +Foam::sampledCuttingPlane::interpolate +( + const interpolation& interpolator +) const +{ + return interpolateField(interpolator); +} + + +Foam::tmp +Foam::sampledCuttingPlane::interpolate +( + const interpolation& interpolator +) const +{ + return interpolateField(interpolator); +} + + +void Foam::sampledCuttingPlane::print(Ostream& os) const +{ + os << "sampledCuttingPlane: " << name() << " :" + << " plane:" << plane_ + << " faces:" << faces().size() + << " points:" << points().size(); +} + + +// ************************************************************************* // diff --git a/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.H b/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.H new file mode 100644 index 0000000000..fcf14db8bd --- /dev/null +++ b/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.H @@ -0,0 +1,258 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 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::sampledCuttingPlane + +Description + A sampledSurface defined by a plane + +SourceFiles + sampledCuttingPlane.C + +\*---------------------------------------------------------------------------*/ + +#ifndef sampledCuttingPlane_H +#define sampledCuttingPlane_H + +#include "sampledSurface.H" +#include "isoSurface.H" +#include "plane.H" +#include "ZoneIDs.H" +#include "fvMeshSubset.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class sampledCuttingPlane Declaration +\*---------------------------------------------------------------------------*/ + +class sampledCuttingPlane +: + public sampledSurface +{ + // Private data + + //- Plane + const plane plane_; + + //- Merge tolerance + const scalar mergeTol_; + + //- Whether to coarsen + const Switch regularise_; + + //- zone name/index (if restricted to zones) + mutable cellZoneID zoneID_; + + //- for zones: patch to put exposed faces into + mutable word exposedPatchName_; + + //- Track if the surface needs an update + mutable bool needsUpdate_; + + + //- Optional subsetted mesh + autoPtr subMeshPtr_; + + //- Distance to cell centres + autoPtr cellDistancePtr_; + + //- Distance to points + scalarField pointDistance_; + + //- Constructed iso surface + autoPtr isoSurfPtr_; + + //- triangles converted to faceList + mutable autoPtr facesPtr_; + + + // Private Member Functions + + //- Create iso surface + void createGeometry(); + + //- sample field on faces + template + tmp > sampleField + ( + const GeometricField& vField + ) const; + + + template + tmp > + interpolateField(const interpolation&) const; + + +public: + + //- Runtime type information + TypeName("sampledCuttingPlane"); + + + // Constructors + + //- Construct from dictionary + sampledCuttingPlane + ( + const word& name, + const polyMesh& mesh, + const dictionary& dict + ); + + + // Destructor + + virtual ~sampledCuttingPlane(); + + + // 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 surface().points(); + } + + //- Faces of surface + virtual const faceList& faces() const + { + if (facesPtr_.empty()) + { + const triSurface& s = surface(); + + facesPtr_.reset(new faceList(s.size())); + + forAll(s, i) + { + facesPtr_()[i] = s[i].triFaceFace(); + } + } + return facesPtr_; + } + + + const isoSurface& surface() const + { + return isoSurfPtr_(); + } + + //- sample field on surface + virtual tmp sample + ( + const volScalarField& + ) const; + + //- sample field on surface + virtual tmp sample + ( + const volVectorField& + ) const; + + //- sample field on surface + virtual tmp sample + ( + const volSphericalTensorField& + ) const; + + //- sample field on surface + virtual tmp sample + ( + const volSymmTensorField& + ) const; + + //- sample field on surface + virtual tmp sample + ( + const volTensorField& + ) const; + + + //- interpolate field on surface + virtual tmp interpolate + ( + const interpolation& + ) const; + + //- interpolate field on surface + virtual tmp interpolate + ( + const interpolation& + ) const; + + //- interpolate field on surface + virtual tmp interpolate + ( + const interpolation& + ) const; + + //- interpolate field on surface + virtual tmp interpolate + ( + const interpolation& + ) const; + + //- interpolate field on surface + virtual tmp interpolate + ( + const interpolation& + ) const; + + //- Write + virtual void print(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "sampledCuttingPlaneTemplates.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlaneTemplates.C b/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlaneTemplates.C new file mode 100644 index 0000000000..f83aeca4e5 --- /dev/null +++ b/src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlaneTemplates.C @@ -0,0 +1,97 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 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 "volPointInterpolation.H" +#include "sampledCuttingPlane.H" +#include "isoSurface.H" +#include "volFieldsFwd.H" +#include "pointFields.H" +#include "volPointInterpolation.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +template +Foam::tmp > +Foam::sampledCuttingPlane::sampleField +( + const GeometricField& vField +) const +{ + return tmp >(new Field(vField, surface().meshCells())); +} + + +template +Foam::tmp > +Foam::sampledCuttingPlane::interpolateField +( + const interpolation& interpolator +) const +{ + // Get fields to sample. Assume volPointInterpolation! + const GeometricField& volFld = + interpolator.psi(); + + if (subMeshPtr_.valid()) + { + tmp > tvolSubFld = + subMeshPtr_().interpolate(volFld); + + const GeometricField& volSubFld = + tvolSubFld(); + + tmp > tpointSubFld = + volPointInterpolation::New(volSubFld.mesh()).interpolate(volSubFld); + + // Sample. + return surface().interpolate + ( + cellDistancePtr_(), + pointDistance_, + volSubFld, + tpointSubFld() + ); + } + else + { + tmp > tpointFld + ( + volPointInterpolation::New(volFld.mesh()).interpolate(volFld) + ); + + // Sample. + return surface().interpolate + ( + cellDistancePtr_(), + pointDistance_, + volFld, + tpointFld() + ); + } +} + + +// ************************************************************************* //