#include <CGAL/Simple_cartesian.h>
#include <CGAL/Interval_nt.h>
#include <CGAL/leda_real.h>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
using EK = CGAL::Simple_cartesian<leda::real>;
using FK = CGAL::Simple_cartesian<CGAL::Interval_nt<true> >;
namespace CGAL {
template < typename EFT, typename FFT>
bool is_exact_in_interval(const EFT e, const FFT f)
{
if(!is_finite(f.inf()) && !is_finite(f.sup()))
{
// std::cerr << "infinite interval..." << std::endl;
return true; // not true if both are -inf or +inf, but doesn't matter we're never there with the data
}
// std::cout << "finiteness of the interval: " << is_finite(f.inf()) << " " << is_finite(f.sup()) << std::endl;
if(!is_finite(f.inf()))
{
// std::cout << "infinite inf: " << (e <= EFT(f.sup())) << std::endl;
return (e <= EFT(f.sup()));
}
else if(!is_finite(f.sup()))
{
// std::cout << "infinite sup: " << (EFT(f.inf()) <= e) << std::endl;
return (EFT(f.inf()) <= e);
}
else
{
// std::cout << "general case: " << (EFT(f.inf()) <= e) << " " << (e <= EFT(f.sup())) << std::endl;
return (EFT(f.inf()) <= e) && (e <= EFT(f.sup()));
}
}
} // namespace CGAL
int main()
{
const double px = -1.7976931348622648e+308;
const double py = 5.3779408224704862e-299;
const double pz = 2.5394974196240072e-321;
const double qx = 5.3779407578037265e-299;
const double qy = 8.7171284671510145e-121;
const double qz = 5.3779407512700166e-299;
const double rx = 5.3779414276905704e-299;
const double ry = 5.3779407512682172e-299;
const double rz = 7.1237266591579577e-307;
const double sx = 9.4175561178482311e+298;
const double sy = 5.3779428873390216e-299;
const double sz = 5.3779407512681174e-299;
EK::Point_3 ep = { px, py, pz };
EK::Point_3 eq = { qx, qy, qz };
EK::Point_3 er = { rx, ry, rz };
EK::Point_3 es = { sx, sy, sz };
FK::Point_3 fp = { px, py, pz };
FK::Point_3 fq = { qx, qy, qz };
FK::Point_3 fr = { rx, ry, rz };
FK::Point_3 fs = { sx, sy, sz };
CGAL_assertion(CGAL::is_exact_in_interval(ep.x(), fp.x()));
CGAL_assertion(CGAL::is_exact_in_interval(ep.y(), fp.y()));
CGAL_assertion(CGAL::is_exact_in_interval(ep.z(), fp.z()));
CGAL_assertion(CGAL::is_exact_in_interval(eq.x(), fq.x()));
CGAL_assertion(CGAL::is_exact_in_interval(eq.y(), fq.y()));
CGAL_assertion(CGAL::is_exact_in_interval(eq.z(), fq.z()));
CGAL_assertion(CGAL::is_exact_in_interval(er.x(), fr.x()));
CGAL_assertion(CGAL::is_exact_in_interval(er.y(), fr.y()));
CGAL_assertion(CGAL::is_exact_in_interval(er.z(), fr.z()));
CGAL_assertion(CGAL::is_exact_in_interval(es.x(), fs.x()));
CGAL_assertion(CGAL::is_exact_in_interval(es.y(), fs.y()));
CGAL_assertion(CGAL::is_exact_in_interval(es.z(), fs.z()));
// CGAL::Protect_FPU_rounding<true> protector;
std::cout << "---- a00" << std::endl;
EK::FT e_a00 = eq.x() - ep.x();
FK::FT f_a00 = fq.x() - fp.x();
CGAL_assertion(CGAL::is_exact_in_interval(e_a00, f_a00));
std::cout << "---- a10" << std::endl;
EK::FT e_a10 = eq.y() - ep.y();
FK::FT f_a10 = fq.y() - fp.y();
CGAL_assertion(CGAL::is_exact_in_interval(e_a10, f_a10));
std::cout << "---- a20" << std::endl;
EK::FT e_a20 = eq.z() - ep.z();
FK::FT f_a20 = fq.z() - fp.z();
CGAL_assertion(CGAL::is_exact_in_interval(e_a20, f_a20));
std::cout << "---- a01" << std::endl;
EK::FT e_a01 = er.x() - ep.x();
FK::FT f_a01 = fr.x() - fp.x();
CGAL_assertion(CGAL::is_exact_in_interval(e_a01, f_a01));
std::cout << "---- a11" << std::endl;
EK::FT e_a11 = er.y() - ep.y();
FK::FT f_a11 = fr.y() - fp.y();
CGAL_assertion(CGAL::is_exact_in_interval(e_a11, f_a11));
std::cout << "---- a21" << std::endl;
EK::FT e_a21 = er.z() - ep.z();
FK::FT f_a21 = fr.z() - fp.z();
CGAL_assertion(CGAL::is_exact_in_interval(e_a21, f_a21));
std::cout << "---- a02" << std::endl;
EK::FT e_a02 = es.x() - ep.x();
FK::FT f_a02 = fs.x() - fp.x();
CGAL_assertion(CGAL::is_exact_in_interval(e_a02, f_a02));
std::cout << "---- a12" << std::endl;
EK::FT e_a12 = es.y() - ep.y();
FK::FT f_a12 = fs.y() - fp.y();
CGAL_assertion(CGAL::is_exact_in_interval(e_a12, f_a12));
std::cout << "---- a22" << std::endl;
EK::FT e_a22 = es.z() - ep.z();
FK::FT f_a22 = fs.z() - fp.z();
CGAL_assertion(CGAL::is_exact_in_interval(e_a22, f_a22));
std::cout << "---- m01" << std::endl;
const EK::FT e_m01 = e_a00 * e_a11 - e_a10 * e_a01;
const FK::FT f_m01 = f_a00 * f_a11 - f_a10 * f_a01;
CGAL_assertion(CGAL::is_exact_in_interval(e_m01, f_m01));
std::cout << "---- m02" << std::endl;
const EK::FT e_m02 = e_a00 * e_a21 - e_a20 * e_a01;
const FK::FT f_m02 = f_a00 * f_a21 - f_a20 * f_a01;
CGAL_assertion(CGAL::is_exact_in_interval(e_m02, f_m02));
std::cout << "---- m12" << std::endl;
const EK::FT e_m12 = e_a10 * e_a21 - e_a20 * e_a11;
const FK::FT f_m12 = f_a10 * f_a21 - f_a20 * f_a11;
CGAL_assertion(CGAL::is_exact_in_interval(e_m12, f_m12));
std::cout << "---- part_m012_1" << std::endl;
const EK::FT e_m012_1 = e_m01 * e_a22;
const FK::FT f_m012_1 = f_m01 * f_a22;
CGAL_assertion(CGAL::is_exact_in_interval(e_m012_1, f_m012_1));
std::cout << "---- part_m012_2" << std::endl;
const EK::FT e_m012_2 = e_m02 * e_a12;
const FK::FT f_m012_2 = f_m02 * f_a12;
CGAL_assertion(CGAL::is_exact_in_interval(e_m012_2, f_m012_2));
std::cout << "---- part_m012_3" << std::endl;
const EK::FT e_m012_3 = e_m12 * e_a02;
const FK::FT f_m012_3 = f_m12 * f_a02;
std::cout << "e_m12: " << e_m12 << std::endl;
std::cout << "f_m12: " << f_m12 << std::endl;
std::cout << "e_a02: " << e_a02 << std::endl;
std::cout << "f_a02: " << f_a02 << std::endl;
std::cout << "e_m12 * e_a02: " << e_m12 * e_a02 << std::endl;
std::cout << "f_m12 * f_a02: " << f_m12 * f_a02 << std::endl;
CGAL_assertion(CGAL::is_exact_in_interval(e_m012_3, f_m012_3));
// std::cout << "---- sum_l" << std::endl;
// const EK::FT e_m012_l = e_m01 * e_a22 - e_m02 * e_a12;
// const FK::FT f_m012_l = f_m01 * f_a22 - f_m02 * f_a12;
// CGAL_assertion(CGAL::is_exact_in_interval(e_m012_l, f_m012_l));
// std::cout << "---- sum_r" << std::endl;
// const EK::FT e_m012_r = - e_m02 * e_a12 + e_m12 * e_a02;
// const FK::FT f_m012_r = - f_m02 * f_a12 + f_m12 * f_a02;
// CGAL_assertion(CGAL::is_exact_in_interval(e_m012_r, f_m012_r));
// std::cout << "---- m012" << std::endl;
// const EK::FT e_m012 = e_m01 * e_a22 - e_m02 * e_a12 + e_m12 * e_a02;
// const FK::FT f_m012 = f_m01 * f_a22 - f_m02 * f_a12 + f_m12 * f_a02;
// CGAL_assertion(CGAL::is_exact_in_interval(e_m012, f_m012));
}