diff --git a/dune/tectonic/polyhedrondistance.hh b/dune/tectonic/polyhedrondistance.hh index a133f296a2841cb8bd8ee4c13683d701bab5e722..e45c7f036d93022a525a2be82000dd46c3e03d6b 100644 --- a/dune/tectonic/polyhedrondistance.hh +++ b/dune/tectonic/polyhedrondistance.hh @@ -188,4 +188,32 @@ double distance(const ConvexPolyhedron<Coordinate> &s1, assert(false); } +// Point-to-Geometry convenience method +template <class Geometry> +double distance(typename Geometry::GlobalCoordinate const &point, + Geometry const &geo, double valueCorrectionTolerance) { + using Coordinate = typename Geometry::GlobalCoordinate; + + ConvexPolyhedron<Coordinate> polyhedron; + polyhedron.vertices.resize(geo.corners()); + for (size_t i = 0; i < polyhedron.vertices.size(); ++i) + polyhedron.vertices[i] = geo.corner(i); + + return distance(point, polyhedron, valueCorrectionTolerance); +} + +// Polyhedron-to-Geometry convenience method +template <class Geometry> +double distance(ConvexPolyhedron<typename Geometry::GlobalCoordinate> const &p1, + Geometry const &geo, double valueCorrectionTolerance) { + using Coordinate = typename Geometry::GlobalCoordinate; + + ConvexPolyhedron<Coordinate> p2; + p2.vertices.resize(geo.corners()); + for (size_t i = 0; i < p2.vertices.size(); ++i) + p2.vertices[i] = geo.corner(i); + + return distance(p2, p1, valueCorrectionTolerance); +} + #endif diff --git a/src/diameter.hh b/src/diameter.hh new file mode 100644 index 0000000000000000000000000000000000000000..45b446f0e56f81bbdfaf793d48c5d55c23a7bc38 --- /dev/null +++ b/src/diameter.hh @@ -0,0 +1,16 @@ +#ifndef SRC_DIAMETER_HH +#define SRC_DIAMETER_HH + +template <class Geometry> double diameter(Geometry const &geometry) { + auto const numCorners = geometry.corners(); + std::vector<typename Geometry::GlobalCoordinate> corners(numCorners); + + double diameter = 0.0; + for (int i = 0; i < numCorners; ++i) { + corners[i] = geometry.corner(i); + for (int j = 0; j < i; ++j) + diameter = std::max(diameter, (corners[i] - corners[j]).two_norm()); + } + return diameter; +} +#endif diff --git a/src/distances.hh b/src/distances.hh deleted file mode 100644 index ffe7dd514a37de6c6c370fb9978284683b6986eb..0000000000000000000000000000000000000000 --- a/src/distances.hh +++ /dev/null @@ -1,33 +0,0 @@ -#ifndef SRC_DISTANCES_HH -#define SRC_DISTANCES_HH - -#include <dune/tectonic/polyhedrondistance.hh> -#include "sand-wedge-data/mygeometry.hh" - -template <class Geometry> double diameter(Geometry const &geometry) { - auto const numCorners = geometry.corners(); - std::vector<typename Geometry::GlobalCoordinate> corners(numCorners); - - double diameter = 0.0; - for (int i = 0; i < numCorners; ++i) { - corners[i] = geometry.corner(i); - for (int j = 0; j < i; ++j) - diameter = std::max(diameter, (corners[i] - corners[j]).two_norm()); - } - return diameter; -} - -template <class Geometry> -double distance( - Geometry const &g, - ConvexPolyhedron<typename Geometry::GlobalCoordinate> const &weakPatch) { - using Coordinate = typename Geometry::GlobalCoordinate; - - ConvexPolyhedron<Coordinate> bsg; - bsg.vertices.resize(g.corners()); - for (size_t i = 0; i < bsg.vertices.size(); ++i) - bsg.vertices[i] = g.corner(i); - - return distance(bsg, weakPatch, 1e-6 * MyGeometry::lengthScale); -} -#endif diff --git a/src/sand-wedge-data/mygrid.cc b/src/sand-wedge-data/mygrid.cc index 93e3a883936cbed8f0eeb8d69ebbb703042924ab..0e7d6e391c1f8461b2532ffc2db57469d153ad65 100644 --- a/src/sand-wedge-data/mygrid.cc +++ b/src/sand-wedge-data/mygrid.cc @@ -2,9 +2,11 @@ #include "config.h" #endif +#include <dune/tectonic/polyhedrondistance.hh> + #include "mygrid.hh" #include "midpoint.hh" -#include "../distances.hh" +#include "../diameter.hh" #if MY_DIM == 3 SimplexManager::SimplexManager(unsigned int shift) : shift_(shift) {} @@ -203,7 +205,8 @@ void refine(Grid &grid, ConvexPolyhedron<LocalVector> const &weakPatch, it != grid.template leafend<0>(); ++it) { auto const geometry = it->geometry(); - auto const weakeningRegionDistance = distance(geometry, weakPatch); + auto const weakeningRegionDistance = + distance(weakPatch, geometry, 1e-6 * MyGeometry::lengthScale); auto const admissibleDiameter = computeAdmissibleDiameter(weakeningRegionDistance, smallestDiameter); diff --git a/src/sand-wedge-data/patchfunction.hh b/src/sand-wedge-data/patchfunction.hh index 4fcc1eb1f423db84c05a30480ec39083831ae04a..c3b7b9c51401fbf14415adb709439de54b9ccd67 100644 --- a/src/sand-wedge-data/patchfunction.hh +++ b/src/sand-wedge-data/patchfunction.hh @@ -5,7 +5,7 @@ #include <dune/common/fvector.hh> #include <dune/common/parametertree.hh> -#include <dune/tectonic/convexpolyhedron.hh> +#include <dune/tectonic/polyhedrondistance.hh> class PatchFunction : public Dune::VirtualFunction<Dune::FieldVector<double, MY_DIM>, diff --git a/src/sand-wedge.cc b/src/sand-wedge.cc index f6ed06934ca0bc0b882e1ee603e600caae7d756c..3b7aa548b7a879dbeae11c300b41b67f4646e4dd 100644 --- a/src/sand-wedge.cc +++ b/src/sand-wedge.cc @@ -55,7 +55,7 @@ #include "adaptivetimestepper.hh" #include "assemblers.hh" -#include "distances.hh" +#include "diameter.hh" #include "enumparser.hh" #include "enums.hh" #include "friction_writer.hh"