Incorrect result of the function bg::intersection()

[linshu77]

This is my code with the define BOOST_GEOMETRY_NO_ROBUSTNESS，
It can be compiled by DEV-C + + 5.11,Boost version: 1.72.
The answer should be 0, but it's 27.
However, if I delete the "#define BOOST_GEOMETRY_NO_ROBUSTNESS",
the answer will be correct.

#define BOOST_GEOMETRY_NO_ROBUSTNESS
#include<iostream>
#include <boost/assign.hpp>
#include <boost/assign/list_of.hpp>
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/polygon.hpp>
using namespace std;
namespace bg = boost::geometry;
int main()
{
	typedef bg::model::d2::point_xy<double> DPoint;
	typedef bg::model::polygon<DPoint, false> DPolygon;
	typedef bg::model::multi_polygon<DPolygon> MPolygon;
	DPolygon poly1;
	DPolygon poly2;
	MPolygon poly;
	std::list<DPoint> lstOf = boost::assign::list_of
		(DPoint(140.10975222826562, 135.99989776999999))
		(DPoint(140.10975222826562, 126.97244753050222))
		(DPoint(140.16255693232961, 135.99989776999999))
		(DPoint(140.10975222826562, 135.99989776999999));
	poly1.outer().assign(lstOf.begin(), lstOf.end());
	lstOf = boost::assign::list_of
		(DPoint(151.99995129999999, 123.99989777000000))
		(DPoint(139.99994018999999, 126.99990054000000))
		(DPoint(129.99995129999999, 124.99989777000000))
		(DPoint(151.99995129999999, 123.99989777000000));
	poly2.outer().assign(lstOf.begin(), lstOf.end());
	bg::correct(poly1);
	bg::correct(poly2);
	bg::intersection(poly1, poly2, poly);
	double Area = bg::area(poly);
	cout<<Area;
	return 1;
}

[mloskot]

What compiler and what version of the compiler is this DEV-C + + 5.11 thing?

What compilation flags did you use, especially what optimisations?

[linshu77]

I also compiled with MSVC2012 Ver11.0.51106.01.
what optimisations? I don't know.

[mloskot]

It's quite important to know at least if you do non-optimised debug build or you do optimised build (Release in Visual Studio speak).

[linshu77]

DEV-C + + 5.11 : TDM-GCC 4.9.2 32-bit Release.
MSVC2012: Release and Debug
The answers are all the same. They're all Incorrect.
However, if I delete the "#define BOOST_GEOMETRY_NO_ROBUSTNESS",
the answer will be correct.

[mloskot]

Thanks for checking it's not optimisation related.

if I delete the "#define BOOST_GEOMETRY_NO_ROBUSTNESS", the answer will be correct.

For some/many cases, that is what you should do, because without that macro defined, the algorithms will perform rescaling to integer for more numerically robust solution (default behaviour).

geometry/include/boost/geometry/core/config.hpp

Lines 25 to 29 in b930d33

	// Rescaling is turned on, unless NO_ROBUSTNESS is defined
	// In future versions of Boost.Geometry, it will be turned off by default
	#if ! defined(BOOST_GEOMETRY_NO_ROBUSTNESS)
	#define BOOST_GEOMETRY_USE_RESCALING
	#endif

You may still need to define the macro in order to disable the rescaling for some/edge cases as you discovered in the previous issue you reported #662.

I think this case can be closed.

[linshu77]

I have two code.
If I add the "#define BOOST_GEOMETRY_NO_ROBUSTNESS",
the code1 will be incorrect, the code2 will be correct.
but if I delete the "#define BOOST_GEOMETRY_NO_ROBUSTNESS",
the code1 will be correct, the code2 will be incorrect.
what can I do? thanks.

code1

#define BOOST_GEOMETRY_NO_ROBUSTNESS
#include<iostream>
#include <boost/assign.hpp>
#include <boost/assign/list_of.hpp>
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/polygon.hpp>
using namespace std;
namespace bg = boost::geometry;
int main()
{
	typedef bg::model::d2::point_xy<double> DPoint;
	typedef bg::model::polygon<DPoint, false> DPolygon;
	typedef bg::model::multi_polygon<DPolygon> MPolygon;
	DPolygon poly1;
	DPolygon poly2;
	MPolygon poly;
	std::list<DPoint> lstOf = boost::assign::list_of
		(DPoint(140.10975222826562, 135.99989776999999))
		(DPoint(140.10975222826562, 126.97244753050222))
		(DPoint(140.16255693232961, 135.99989776999999))
		(DPoint(140.10975222826562, 135.99989776999999));
	poly1.outer().assign(lstOf.begin(), lstOf.end());
	lstOf = boost::assign::list_of
		(DPoint(151.99995129999999, 123.99989777000000))
		(DPoint(139.99994018999999, 126.99990054000000))
		(DPoint(129.99995129999999, 124.99989777000000))
		(DPoint(151.99995129999999, 123.99989777000000));
	poly2.outer().assign(lstOf.begin(), lstOf.end());
	bg::correct(poly1);
	bg::correct(poly2);
	bg::intersection(poly1, poly2, poly);
	double Area = bg::area(poly);
	cout<<Area;
	return 1;
}

code2:

#define BOOST_GEOMETRY_NO_ROBUSTNESS
#include<iostream>
#include <boost/assign.hpp>
#include <boost/assign/list_of.hpp>
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/polygon.hpp>
using namespace std;
namespace bg = boost::geometry;
int main()
{
	typedef bg::model::d2::point_xy<double> DPoint;
	typedef bg::model::polygon<DPoint, false> DPolygon;
	typedef bg::model::multi_polygon<DPolygon> MPolygon;
	DPolygon poly1;
	DPolygon poly2;
	MPolygon poly;
	std::list<DPoint> lstOf = boost::assign::list_of
		(DPoint(171.60480351137778, 162.32096650635697))
		(DPoint(172.21795926332567, 162.99736073096460))
		(DPoint(172.21795926332567, 170.00000000000000))
		(DPoint(171.60480351137778, 170.00000000000000))
		(DPoint(171.60480351137778, 162.32096650635697));
	poly1.outer().assign(lstOf.begin(), lstOf.end());
	lstOf = boost::assign::list_of
		(DPoint(200.00000000000000, 166.00002000000001))
		(DPoint(184.00004999999999, 166.00002000000001))
		(DPoint(180.00004999999999, 164.00002000000001))
		(DPoint(169.99998889000000, 162.00000277999999))
		(DPoint(182.00000000000000, 159.00000000000000))
		(DPoint(200.00000000000000, 159.00000000000000))
		(DPoint(200.00000000000000, 166.00002000000001));
	poly2.outer().assign(lstOf.begin(), lstOf.end());
	bg::correct(poly1);
	bg::correct(poly2);
	bg::intersection(poly1, poly2, poly);
	double Area = bg::area(poly);
	cout<<Area;
	return 1;
}

[awulkiew]

The first case looks like that:

red - poly1
green - poly2
blue - poly (result)

With BOOST_GEOMETRY_NO_ROBUSTNESS the result is equal to the poly2 (with additional point). Without the result is empty. Even if the result was non-empty the intersection should not be equal to poly2.

[0] {140.10975222826562, 126.97244753050222}
[1] {139.99994018999999, 126.99990054000000}
[2] {129.99995129999999, 124.99989777000000}
[3] {151.99995129999999, 123.99989777000000}
[4] {140.10975222826562, 126.97244753050222}

The CW/CCW point order does not affect the result. Tested with VS2017.

[awulkiew]

The second case:

As @linshu77 said it's the other way around WRT BOOST_GEOMETRY_NO_ROBUSTNESS. However this time the poly1 is the result.

{171.60480351137778, 162.32096650635697}
{172.21795926332567, 162.99736073096460}
{172.21795926332567, 170.00000000000000}
{171.60480351137778, 170.00000000000000}
{171.60480351137778, 162.32096650635697}