diff --git a/src/sampling/Make/files b/src/sampling/Make/files index a2cca21562..c6191009cb 100644 --- a/src/sampling/Make/files +++ b/src/sampling/Make/files @@ -22,6 +22,8 @@ sampledSet/shortestPath/shortestPathSet.C surface/cutting/cuttingPlane.C surface/cutting/cuttingPlaneCuts.C surface/cutting/cuttingPlaneSelection.C +surface/cutting/cuttingSurface.C +surface/cutting/cuttingSurfaceCuts.C surface/cutting/cuttingSurfaceBase.C surface/cutting/cuttingSurfaceBaseSelection.C surface/distanceSurface/distanceSurface.C @@ -43,6 +45,7 @@ sampledSurface/isoSurface/sampledIsoSurface.C sampledSurface/isoSurface/sampledIsoSurfaceCell.C sampledSurface/distanceSurface/sampledDistanceSurface.C sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C +sampledSurface/sampledCuttingSurface/sampledCuttingSurface.C sampledSurface/sampledSurface/sampledSurface.C sampledSurface/sampledSurfaces/sampledSurfaces.C sampledSurface/sampledSurfaces/sampledSurfacesGrouping.C diff --git a/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurface.C b/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurface.C new file mode 100644 index 0000000000..f97e86473b --- /dev/null +++ b/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurface.C @@ -0,0 +1,241 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2018 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 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "sampledCuttingSurface.H" +#include "dictionary.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(sampledCuttingSurface, 0); + addToRunTimeSelectionTable + ( + sampledSurface, + sampledCuttingSurface, + word + ); +} + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +Foam::bitSet Foam::sampledCuttingSurface::cellSelection(const bool warn) const +{ + boundBox meshBounds; + + bitSet cellsToSelect = + cuttingSurface::cellSelection + ( + mesh(), bounds_, zoneNames_, meshBounds + ); + + if (warn) + { + cuttingSurface::checkOverlap(name(), meshBounds, bounds_); + } + + return cellsToSelect; +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::sampledCuttingSurface::sampledCuttingSurface +( + const polyMesh& mesh, + const word& surfaceType, + const word& surfaceName, + const bool triangulate, + const boundBox& bounds +) +: + sampledSurface(surfaceName, mesh), + cuttingSurface(mesh, surfaceType, surfaceName), + zoneNames_(), + bounds_(bounds), + triangulate_(triangulate), + needsUpdate_(true) +{} + + +Foam::sampledCuttingSurface::sampledCuttingSurface +( + const word& defaultSurfaceName, + const polyMesh& mesh, + const dictionary& dict +) +: + sampledSurface(defaultSurfaceName, mesh, dict), + cuttingSurface(defaultSurfaceName, mesh, dict), + zoneNames_(), + bounds_(dict.lookupOrDefault("bounds", boundBox::invertedBox)), + triangulate_(dict.lookupOrDefault("triangulate", true)), + needsUpdate_(true) +{ + if (!dict.readIfPresent("zones", zoneNames_) && dict.found("zone")) + { + zoneNames_.resize(1); + dict.readEntry("zone", zoneNames_.first()); + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +bool Foam::sampledCuttingSurface::needsUpdate() const +{ + return needsUpdate_; +} + + +bool Foam::sampledCuttingSurface::expire() +{ + // Already marked as expired + if (needsUpdate_) + { + return false; + } + + sampledSurface::clearGeom(); + + needsUpdate_ = true; + return true; +} + + +bool Foam::sampledCuttingSurface::update() +{ + if (!needsUpdate_) + { + return false; + } + + sampledSurface::clearGeom(); + + performCut(mesh(), triangulate_, cellSelection(true)); + + if (debug) + { + cuttingSurface::print(Pout); + Pout<< endl; + } + + needsUpdate_ = false; + return true; +} + + +Foam::tmp Foam::sampledCuttingSurface::sample +( + const interpolation& sampler +) const +{ + return sampleOnFaces(sampler); +} + + +Foam::tmp Foam::sampledCuttingSurface::sample +( + const interpolation& sampler +) const +{ + return sampleOnFaces(sampler); +} + + +Foam::tmp Foam::sampledCuttingSurface::sample +( + const interpolation& sampler +) const +{ + return sampleOnFaces(sampler); +} + + +Foam::tmp Foam::sampledCuttingSurface::sample +( + const interpolation& sampler +) const +{ + return sampleOnFaces(sampler); +} + + +Foam::tmp Foam::sampledCuttingSurface::sample +( + const interpolation& sampler +) const +{ + return sampleOnFaces(sampler); +} + + +Foam::tmp Foam::sampledCuttingSurface::interpolate +( + const interpolation& interpolator +) const +{ + return sampleOnPoints(interpolator); +} + + +Foam::tmp Foam::sampledCuttingSurface::interpolate +( + const interpolation& interpolator +) const +{ + return sampleOnPoints(interpolator); +} + +Foam::tmp Foam::sampledCuttingSurface::interpolate +( + const interpolation& interpolator +) const +{ + return sampleOnPoints(interpolator); +} + + +Foam::tmp Foam::sampledCuttingSurface::interpolate +( + const interpolation& interpolator +) const +{ + return sampleOnPoints(interpolator); +} + + +Foam::tmp Foam::sampledCuttingSurface::interpolate +( + const interpolation& interpolator +) const +{ + return sampleOnPoints(interpolator); +} + + +// ************************************************************************* // diff --git a/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurface.H b/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurface.H new file mode 100644 index 0000000000..575ee6936c --- /dev/null +++ b/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurface.H @@ -0,0 +1,286 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2018 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 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::sampledCuttingSurface + +Description + A surface define by using an input surface to cut the mesh cells + +Usage + Dictionary controls: + \table + Property | Description | Required | Default + type | surfaceCut | yes | + surfaceType | type of surface | yes | + surfaceName | name of surface in \c triSurface/ | no | dict name + triangulate | triangulate faces | no | true + bounds | limit with bounding box | no | + zone | limit to cell zone (name or regex) | no | + zones | limit to cell zones (names, regexs) | no | + \endtable + +Note + No attempt at resolving degenerate cases. + Since the cut faces can be quite ugly, they will often be triangulated. + +SourceFiles + sampledCuttingSurface.C + sampledCuttingSurfaceTemplates.C + +\*---------------------------------------------------------------------------*/ + +#ifndef sampledCuttingSurface_H +#define sampledCuttingSurface_H + +#include "sampledSurface.H" +#include "cuttingSurface.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class sampledCuttingSurface Declaration +\*---------------------------------------------------------------------------*/ + +class sampledCuttingSurface +: + public sampledSurface, + public cuttingSurface +{ + // Private Data + + //- The zone or zones in which cutting is to occur + wordRes zoneNames_; + + //- Optional bounding box to trim against + const boundBox bounds_; + + //- Triangulate faces or not + const bool triangulate_; + + //- Track if the surface needs an update + mutable bool needsUpdate_; + + + // Private Member Functions + + //- Define cell selection from zones and bounding box. + // Optionally check and warn if the bounding box + // does not overlap with the mesh (or submesh) + bitSet cellSelection(const bool warn=false) const; + + + //- Sample volume field onto surface faces + template + tmp> sampleOnFaces + ( + const interpolation& sampler + ) const; + + //- Interpolate volume field onto surface points + template + tmp> sampleOnPoints + ( + const interpolation& interpolator + ) const; + + +public: + + //- Runtime type information + TypeName("surfaceCut"); + + + // Constructors + + //- Construct from components + sampledCuttingSurface + ( + const polyMesh& mesh, + const word& surfaceType, + const word& surfaceName, + const bool triangulate = true, + const boundBox& bounds = boundBox::invertedBox + ); + + //- Construct from dictionary + sampledCuttingSurface + ( + const word& defaultSurfaceName, + const polyMesh& mesh, + const dictionary& dict + ); + + + //- Destructor + virtual ~sampledCuttingSurface() = default; + + + // 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 cuttingSurface::points(); + } + + //- Faces of surface + virtual const faceList& faces() const + { + return cuttingSurface::surfFaces(); + } + + //- Const access to per-face zone/region information + // Could instead return meshCells or cellZoneId of the meshCells. + virtual const labelList& zoneIds() const + { + return Foam::emptyLabelList; + } + + //- Face area magnitudes + virtual const vectorField& Sf() const + { + return cuttingSurface::Sf(); + } + + //- Face area magnitudes + virtual const scalarField& magSf() const + { + return cuttingSurface::magSf(); + } + + //- Face centres + virtual const vectorField& Cf() const + { + return cuttingSurface::Cf(); + } + + //- For each face, the original cell in mesh + using cuttingSurface::meshCells; + + + // Sample + + //- Sample volume field onto surface faces + virtual tmp sample + ( + const interpolation& sampler + ) const; + + //- Sample volume field onto surface faces + virtual tmp sample + ( + const interpolation& sampler + ) const; + + //- Sample volume field onto surface faces + virtual tmp sample + ( + const interpolation& sampler + ) const; + + //- Sample volume field onto surface faces + virtual tmp sample + ( + const interpolation& sampler + ) const; + + //- Sample volume field onto surface faces + virtual tmp sample + ( + const interpolation& sampler + ) const; + + + // Interpolate + + //- Interpolate volume field onto surface points + virtual tmp interpolate + ( + const interpolation& interpolator + ) const; + + //- Interpolate volume field onto surface points + virtual tmp interpolate + ( + const interpolation& interpolator + ) const; + + //- Interpolate volume field onto surface points + virtual tmp interpolate + ( + const interpolation& interpolator + ) const; + + //- Interpolate volume field onto surface points + virtual tmp interpolate + ( + const interpolation& interpolator + ) const; + + //- Interpolate volume field onto surface points + virtual tmp interpolate + ( + const interpolation& interpolator + ) const; + + + // Output + +// //- Print information +// void print(Ostream& os) const; + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository + #include "sampledCuttingSurfaceTemplates.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurfaceTemplates.C b/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurfaceTemplates.C new file mode 100644 index 0000000000..a3fae55bd1 --- /dev/null +++ b/src/sampling/sampledSurface/sampledCuttingSurface/sampledCuttingSurfaceTemplates.C @@ -0,0 +1,89 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2018 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 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "sampledCuttingSurface.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +template +Foam::tmp> +Foam::sampledCuttingSurface::sampleOnFaces +( + const interpolation& sampler +) const +{ + return sampledSurface::sampleOnFaces + ( + sampler, + meshCells(), + faces(), + points() + ); +} + + +template +Foam::tmp> +Foam::sampledCuttingSurface::sampleOnPoints +( + const interpolation& interpolator +) const +{ + // elements to sample + const labelList& elements = meshCells(); + + // One value per point. + // Initialize with Zero to handle missed/degenerate faces + + auto tvalues = tmp>::New(points().size(), Zero); + auto& values = tvalues.ref(); + + bitSet pointDone(points().size()); + + const faceList& fcs = faces(); + + forAll(fcs, facei) + { + const face& f = fcs[facei]; + const label celli = elements[facei]; + + for (const label pointi : f) + { + if (pointDone.set(pointi)) + { + values[pointi] = interpolator.interpolate + ( + points()[pointi], + celli + ); + } + } + } + + return tvalues; +} + + +// ************************************************************************* // diff --git a/src/sampling/surface/cutting/cuttingSurface.C b/src/sampling/surface/cutting/cuttingSurface.C new file mode 100644 index 0000000000..862fb53779 --- /dev/null +++ b/src/sampling/surface/cutting/cuttingSurface.C @@ -0,0 +1,176 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2018 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 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "cuttingSurface.H" +#include "dictionary.H" +#include "fvMesh.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::cuttingSurface::cuttingSurface +( + const polyMesh& mesh, + const word& surfaceType, + const word& surfaceName +) +: + cuttingSurfaceBase(), + surfPtr_ + ( + searchableSurface::New + ( + surfaceType, + IOobject + ( + surfaceName, // name + mesh.time().constant(), // directory + "triSurface", // instance + mesh.time(), // registry + IOobject::MUST_READ, + IOobject::NO_WRITE + ), + dictionary() + ) + ) +{} + + +Foam::cuttingSurface::cuttingSurface +( + const word& defaultSurfaceName, + const polyMesh& mesh, + const dictionary& dict +) +: + cuttingSurfaceBase(), + surfPtr_ + ( + searchableSurface::New + ( + dict.get("surfaceType"), + IOobject + ( + dict.lookupOrDefault("surfaceName", defaultSurfaceName), + mesh.time().constant(), // directory + "triSurface", // instance + mesh.time(), // registry + IOobject::MUST_READ, + IOobject::NO_WRITE + ), + dict + ) + ) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::cuttingSurface::performCut +( + const primitiveMesh& mesh, + const bool triangulate, + bitSet&& cellIdLabels +) +{ + const fvMesh& fvm = static_cast(mesh); + + MeshStorage::clear(); + meshCells_.clear(); + + // Pre-populate with restriction + bitSet cellCuts(std::move(cellIdLabels)); + + if (cellCuts.size()) + { + cellCuts.resize(mesh.nCells()); + } + + scalarField pointDist; + calcCellCuts(fvm, pointDist, cellCuts); + + + // Walk cell cuts to create faces + + // Action #1: + // - Orient edge so it points in the positive gradient direction. + // - Edge intersection when it spans across point-distance == 0. + const auto edgeOrientIntersect = + [=](edge& e) -> bool + { + if (pointDist[e.last()] < pointDist[e.first()]) + { + e.flip(); + } + + const scalar s0 = pointDist[e.first()]; + const scalar s1 = pointDist[e.last()]; + + if + ( + s0 > ROOTVSMALL // Edge is all positive + || s1 < ROOTVSMALL // Edge is all negative + || Foam::mag(s1-s0) < ROOTVSMALL + ) + { + return false; + } + + return true; + }; + + // Action #2: + // - The edge intersection alpha for point-distance == 0 + // Return -1 for error. + // This is like the iso-fraction for an iso-surface. + const auto edgeAlphaIntersect = + [=](const edge& e) -> scalar + { + const scalar s0 = pointDist[e.first()]; + const scalar s1 = pointDist[e.last()]; + const scalar d = s1-s0; + + return Foam::mag(d) < ROOTVSMALL ? -1 : (-s0/d); + }; + + walkCellCuts + ( + mesh, + cellCuts, + edgeOrientIntersect, + edgeAlphaIntersect, + triangulate + ); +} + + +void Foam::cuttingSurface::print(Ostream& os) const +{ + os << " surface:" << surfaceName() + << " faces:" << MeshStorage::surfFaces().size() + << " points:" << MeshStorage::points().size(); +} + + +// ************************************************************************* // diff --git a/src/sampling/surface/cutting/cuttingSurface.H b/src/sampling/surface/cutting/cuttingSurface.H new file mode 100644 index 0000000000..d122734606 --- /dev/null +++ b/src/sampling/surface/cutting/cuttingSurface.H @@ -0,0 +1,147 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2018 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 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::cuttingSurface + +Description + Constructs a cutting surface through a mesh. + + No attempt at resolving degenerate cases. + Since the cut faces can be quite ugly, they will often be triangulated. + +SourceFiles + cuttingSurface.C + +\*---------------------------------------------------------------------------*/ + +#ifndef cuttingSurface_H +#define cuttingSurface_H + +#include "searchableSurface.H" +#include "cuttingSurfaceBase.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class cuttingSurface Declaration +\*---------------------------------------------------------------------------*/ + +class cuttingSurface +: + public cuttingSurfaceBase +{ + // Private Data + + //- The surface + const autoPtr surfPtr_; + + + // Private Member Functions + + //- Determine the cut cells, possibly restricted to a list of cells + // + // \param pointDist [out] For each mesh point, the (signed) distance + // to the surface. + // \param cellCuts [in,out] On input an empty set (ie, no restriction) + // or subsetted cells. On output, the cells cut according to the + // planeSides detection. + void calcCellCuts + ( + const fvMesh& fvm, + scalarField& pointDist, + bitSet& cellCuts + ); + + +protected: + + // Protected Member Functions + + //- Cut mesh, restricted to a list of cells + using cuttingSurfaceBase::performCut; + + //- Cut mesh, restricted to a list of cells + // Reclaim memory for cellSelectionMask + virtual void performCut + ( + const primitiveMesh& mesh, + const bool triangulate, + bitSet&& cellIdLabels + ); + + +public: + + // Constructors + + //- Construct from components + cuttingSurface + ( + const polyMesh& mesh, + const word& surfaceType, + const word& surfaceName + ); + + //- Construct from dictionary + cuttingSurface + ( + const word& defaultSurfaceName, + const polyMesh& mesh, + const dictionary& dict + ); + + + //- Destructor + virtual ~cuttingSurface() = default; + + + // Member Functions + + //- The name of the underlying searchableSurface + const word& surfaceName() const + { + return surfPtr_->name(); + } + + + // Output + + //- Print information + void print(Ostream& os) const; + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/sampling/surface/cutting/cuttingSurfaceCuts.C b/src/sampling/surface/cutting/cuttingSurfaceCuts.C new file mode 100644 index 0000000000..54d65bb880 --- /dev/null +++ b/src/sampling/surface/cutting/cuttingSurfaceCuts.C @@ -0,0 +1,173 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2018 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 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "cuttingSurface.H" +#include "fvMesh.H" +#include "volFields.H" +#include "pointFields.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::cuttingSurface::calcCellCuts +( + const fvMesh& fvm, + scalarField& pointDist, + bitSet& cellCuts +) +{ + const pointField& cc = fvm.C(); + + const faceList& faces = fvm.faces(); + const cellList& cells = fvm.cells(); + const pointField& pts = fvm.points(); + + const label nCells = fvm.nCells(); + const label nPoints = fvm.nPoints(); + + // The point distances + + List nearest; + surfPtr_().findNearest + ( + pts, + scalarField(nPoints, GREAT), + nearest + ); + + vectorField norms; + surfPtr_().getNormal(nearest, norms); + + pointDist.resize(nPoints); + + for (label i=0; i < nPoints; ++i) + { + const point diff(pts[i] - nearest[i].hitPoint()); + + // Normal distance + pointDist[i] = (diff & norms[i]); + } + + + if (cellCuts.size()) + { + cellCuts.resize(nCells); // safety + } + else + { + cellCuts.resize(nCells, true); + } + + + // Check based on cell distance + + surfPtr_().findNearest + ( + cc, + scalarField(cc.size(), GREAT), + nearest + ); + + + boundBox cellBb; + labelHashSet cpts; + + for (const label celli : cellCuts) + { + cellBb.clear(); + + if (fvm.hasCellPoints()) + { + cellBb.add(pts, fvm.cellPoints(celli)); + } + else + { + // DIY version of cellPoints + cpts.clear(); + + const cell& cFaces = cells[celli]; + + for (const label facei : cFaces) + { + cpts.insert(faces[facei]); + } + + cellBb.add(pts, cpts); + } + + if (!cellBb.contains(nearest[celli].hitPoint())) + { + cellCuts.unset(celli); + } + } + + + if (debug) + { + volScalarField cCuts + ( + IOobject + ( + "cuttingSurface.cellCuts", + fvm.time().timeName(), + fvm.time(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + fvm, + dimensionedScalar(dimless, Zero) + ); + + for (const label celli : cellCuts) + { + cCuts[celli] = 1; + } + + Pout<< "Writing cellCuts:" << cCuts.objectPath() << endl; + cCuts.write(); + + pointScalarField pDist + ( + IOobject + ( + "cuttingSurface.pointDistance", + fvm.time().timeName(), + fvm.time(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + pointMesh::New(fvm), + dimensionedScalar(dimLength, Zero) + ); + pDist.primitiveFieldRef() = pointDist; + + Pout<< "Writing point distance:" << pDist.objectPath() << endl; + pDist.write(); + } +} + + +// ************************************************************************* // diff --git a/tutorials/compressible/rhoSimpleFoam/squareBend/constant/triSurface/angledPlane1.obj b/tutorials/compressible/rhoSimpleFoam/squareBend/constant/triSurface/angledPlane1.obj new file mode 100644 index 0000000000..9b1b81d49b --- /dev/null +++ b/tutorials/compressible/rhoSimpleFoam/squareBend/constant/triSurface/angledPlane1.obj @@ -0,0 +1,25 @@ +# An angled plane with (1,1,1) normal +o angledPlane1 + +# Full extent +# v -0.15 -0.025 0.025 +# v -0.10 -0.025 -0.025 +# v -0.05 -0.075 -0.025 +# v -0.10 -0.075 0.025 + +# Embedded, ending at cell +# v -0.135 -0.0325 0.0175 +# v -0.100 -0.0325 -0.0175 +# v -0.065 -0.0675 -0.0175 +# v -0.100 -0.0675 0.0175 + +# Embedded, ending mid-cell + +v -0.1425 -0.02875 0.02125 +v -0.100 -0.02875 -0.02125 +v -0.0575 -0.07125 -0.02125 +v -0.100 -0.07125 0.02125 + +f 1 2 4 +f 2 3 4 +# EOF diff --git a/tutorials/compressible/rhoSimpleFoam/squareBend/system/samplingDebug b/tutorials/compressible/rhoSimpleFoam/squareBend/system/samplingDebug index bc86fdcfb4..67ba6478ec 100644 --- a/tutorials/compressible/rhoSimpleFoam/squareBend/system/samplingDebug +++ b/tutorials/compressible/rhoSimpleFoam/squareBend/system/samplingDebug @@ -62,6 +62,32 @@ debug surfaceName angledPlane.obj; } + angledPlaneCut + { + type surfaceCut; + triangulate false; + surfaceType triSurfaceMesh; + surfaceName angledPlane.obj; + } + + angledPlane1 + { + type distanceSurface; + distance 0; + signed true; + regularise true; + surfaceType triSurfaceMesh; + surfaceName angledPlane1.obj; + } + + angledPlane1Cut + { + type surfaceCut; + triangulate false; + surfaceType triSurfaceMesh; + surfaceName angledPlane1.obj; + } + disk1 { ${_disk1}