diff --git a/test/common/geometry/test_makegeometry.cc b/test/common/geometry/test_makegeometry.cc index b1ec9d08e386a2c032f93e6a51ea680c730b2e95..718db663d5e55b65479866a245d60e9f7a5cb06e 100644 --- a/test/common/geometry/test_makegeometry.cc +++ b/test/common/geometry/test_makegeometry.cc @@ -30,6 +30,7 @@ #include <dune/common/float_cmp.hh> #include <dumux/common/geometry/intersectspointgeometry.hh> #include <dumux/common/geometry/makegeometry.hh> +#include <dumux/common/geometry/grahamconvexhull.hh> //! test if the ordering of the input points has an effect on the resulting geometry template<class GlobalPosition> @@ -40,60 +41,73 @@ void permutatePointsAndTest(const std::vector<GlobalPosition>& cornerPoints, { std::array<int, 4> s = {0,1,2,3}; + if (Dumux::grahamConvexHull<2>(cornerPoints).size() != 4) + { + if (verbose) + std::cout << "Input points do not span a strictly convex polygon. Skipping." << std::endl; + return; + } + const auto area = Dumux::makeDuneQuadrilaterial(cornerPoints).volume(); do { - std::vector<GlobalPosition> tmp = {cornerPoints[s[0]], cornerPoints[s[1]], cornerPoints[s[2]], cornerPoints[s[3]]}; + const std::vector<GlobalPosition> corners = {cornerPoints[s[0]], cornerPoints[s[1]], cornerPoints[s[2]], cornerPoints[s[3]]}; - if(verbose) + if (Dumux::grahamConvexHull<2>(corners).size() != 4) { - std::cout << "input corner points: " << std::endl; - for(auto&& i : tmp) - std::cout << "(" << i << ") "; - std::cout << std::endl; + if (verbose) + std::cout << "Input points do not span a strictly convex polygon. Skipping." << std::endl; + continue; } - const auto quad = Dumux::makeDuneQuadrilaterial(tmp); + const auto quad = Dumux::makeDuneQuadrilaterial(corners); - if(verbose) + auto printCorners = [&quad]() + { + std::ostringstream tmp; + for (int i = 0; i < quad.corners(); ++i) + tmp << "(" << quad.corner(i) << ") "; + return tmp.str(); + }; + + if (verbose) { std::cout << "resulting corners \n"; - for(int i = 0; i < quad.corners(); ++i) - std::cout << "(" << quad.corner(i) << ") "; + std::cout << printCorners(); std::cout << std::endl; std::cout << "actual area: " << area << ", area after permuation: " << quad.volume() << std::endl; } - if(!Dune::FloatCmp::eq(quad.volume(), area)) + if (!Dune::FloatCmp::eq(quad.volume(), area)) DUNE_THROW(Dune::InvalidStateException, "Area of quadrilateral after permuation of input points is wrong"); - if(verbose) + if (verbose) std::cout << "checking point intersection: \n"; - for(const auto& p: pointsWithinGeometry) + for (const auto& p: pointsWithinGeometry) { - if(Dumux::intersectsPointGeometry(p, quad)) + if (Dumux::intersectsPointGeometry(p, quad)) { if(verbose) std::cout << "point " << p << " lies within the quadrilateral" << std::endl; } else - DUNE_THROW(Dune::InvalidStateException, "False negative: Check for point " << p << " which is inside the geometry failed"); + DUNE_THROW(Dune::InvalidStateException, "False negative: Check for point " << p << " which is inside the geometry " << printCorners() << " failed"); } - for(const auto& p : pointsOutsideGeometry) + for (const auto& p : pointsOutsideGeometry) { - if(!Dumux::intersectsPointGeometry(p, quad)) + if (!Dumux::intersectsPointGeometry(p, quad)) { - if(verbose) + if (verbose) std::cout << "point " << p << " lies outside of the quadrilateral" << std::endl; } else - DUNE_THROW(Dune::InvalidStateException, "False positive: Check for point " << p << " which is outside the geometry failed"); + DUNE_THROW(Dune::InvalidStateException, "False positive: Check for point " << p << " which is outside the geometry " << printCorners() << " failed"); } - } while(std::next_permutation(s.begin(), s.end())); + } while (std::next_permutation(s.begin(), s.end())); } template<class GlobalPosition>