Add an alternative of Mpzf using boost cpp_int (#7191)

## TODO:
- [x] branch size
- [x] boost backend should not be the default in 5.6
This commit is contained in:
Sebastien Loriot 2023-05-22 09:23:35 +02:00 committed by GitHub
commit 11b92e94f8
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
17 changed files with 15481 additions and 104 deletions

View File

@ -139,7 +139,7 @@ typedef unspecified_type Is_numerical_sensitive;
This type specifies the return type of the predicates provided
by this traits. The type must be convertible to `bool` and
typically the type indeed maps to `bool`. However, there are also
cases such as interval arithmetic, in which it is `Uncertain<bool>`
cases such as interval arithmetic, in which it is `CGAL::Uncertain<bool>`
or some similar type.
*/
@ -300,4 +300,3 @@ typedef unspecified_type Root_of;
/// @}
}; /* end AlgebraicStructureTraits */

14550
Data/data/points_3/ocean_r.xyz Normal file

File diff suppressed because it is too large Load Diff

View File

@ -88,6 +88,8 @@ public:
result_type
operator()(const Args&... args) const
{
#ifndef CGAL_EPICK_NO_INTERVALS
CGAL_BRANCH_PROFILER(std::string(" failures/calls to : ") + std::string(CGAL_PRETTY_FUNCTION), tmp);
// Protection is outside the try block as VC8 has the CGAL_CFG_FPU_ROUNDING_MODE_UNWINDING_VC_BUG
{
@ -103,6 +105,7 @@ public:
CGAL_BRANCH_PROFILER_BRANCH(tmp);
Protect_FPU_rounding<!Protection> p(CGAL_FE_TONEAREST);
CGAL_expensive_assertion(FPU_get_cw() == CGAL_FE_TONEAREST);
#endif // CGAL_EPICK_NO_INTERVALS
return ep(c2e(args)...);
}
};
@ -154,6 +157,7 @@ public:
result_type
operator()(const Args&... args) const
{
#ifndef CGAL_EPICK_NO_INTERVALS
CGAL_BRANCH_PROFILER(std::string(" failures/calls to : ") + std::string(CGAL_PRETTY_FUNCTION), tmp);
// Protection is outside the try block as VC8 has the CGAL_CFG_FPU_ROUNDING_MODE_UNWINDING_VC_BUG
{
@ -169,7 +173,7 @@ public:
CGAL_BRANCH_PROFILER_BRANCH(tmp);
Protect_FPU_rounding<!Protection> p(CGAL_FE_TONEAREST);
CGAL_expensive_assertion(FPU_get_cw() == CGAL_FE_TONEAREST);
#endif // CGAL_EPICK_NO_INTERVALS
return call(args...);
}
};

View File

@ -27,11 +27,9 @@
#if defined(__has_include) && ( ! defined _MSC_VER || _MSC_VER > 1900)
# if CGAL_USE_GMP && ! __has_include(<gmp.h>)
# pragma CGAL_WARNING("<gmp.h> cannot be found. Less efficient number types will be used instead. Define CGAL_NO_GMP=1 if that is on purpose.")
# undef CGAL_USE_GMP
# undef CGAL_USE_MPFR
# elif CGAL_USE_MPFR && ! __has_include(<mpfr.h>)
# pragma CGAL_WARNING("<mpfr.h> cannot be found and the GMP support in CGAL requires it. Less efficient number types will be used instead. Define CGAL_NO_GMP=1 if that is on purpose.")
# undef CGAL_USE_GMP
# undef CGAL_USE_MPFR
# endif // CGAL_USE_MPFR and no <mpfr.h>

View File

@ -196,3 +196,17 @@ if (NOT TARGET CGAL::CGAL_Basic_viewer)
INTERFACE_LINK_LIBRARIES CGAL::CGAL_Qt5)
endif()
#warning: the order in this list has to match the enum in Exact_type_selector
set(CGAL_CMAKE_EXACT_NT_BACKEND_OPTIONS GMP_BACKEND GMPXX_BACKEND BOOST_GMP_BACKEND BOOST_BACKEND LEDA_BACKEND MP_FLOAT_BACKEND Default)
set(CGAL_CMAKE_EXACT_NT_BACKEND "Default" CACHE STRING "Setting for advanced users that what to change the default number types used in filtered kernels. Some options might not be working depending on how you configured your build.")
set_property(CACHE CGAL_CMAKE_EXACT_NT_BACKEND PROPERTY STRINGS ${CGAL_CMAKE_EXACT_NT_BACKEND_OPTIONS})
if ( NOT "${CGAL_CMAKE_EXACT_NT_BACKEND}" STREQUAL "Default" )
list(FIND CGAL_CMAKE_EXACT_NT_BACKEND_OPTIONS ${CGAL_CMAKE_EXACT_NT_BACKEND} DEB_VAL)
set_property(
TARGET CGAL
APPEND PROPERTY
INTERFACE_COMPILE_DEFINITIONS "CMAKE_OVERRIDDEN_DEFAULT_ENT_BACKEND=${DEB_VAL}"
) # do not use set_target_properties to avoid overwritting
endif()

View File

@ -30,6 +30,7 @@ template <class R_>
class Iso_cuboid_3 : public R_::Kernel_base::Iso_cuboid_3
{
typedef typename R_::RT RT;
typedef typename R_::FT FT;
typedef typename R_::Point_3 Point_3;
typedef typename R_::Aff_transformation_3 Aff_transformation_3;
@ -84,8 +85,9 @@ public:
max_hx, max_hy, max_hz)) {}
Iso_cuboid_3(const Bbox_3& bbox)
: Rep(typename R::Construct_iso_cuboid_3()(Return_base_tag(), bbox.xmin(), bbox.ymin(), bbox.zmin(),
bbox.xmax(), bbox.ymax(), bbox.zmax())) {}
: Rep(typename R::Construct_iso_cuboid_3()(Return_base_tag(),
typename R::Construct_point_3()(FT(bbox.xmin()), FT(bbox.ymin()), FT(bbox.zmin())),
typename R::Construct_point_3()(FT(bbox.xmax()), FT(bbox.ymax()), FT(bbox.zmax())))) {}
decltype(auto)
min BOOST_PREVENT_MACRO_SUBSTITUTION () const

View File

@ -78,7 +78,9 @@ public:
: Rep(typename R::Construct_iso_rectangle_2()(Return_base_tag(), min_hx, min_hy, max_hx, max_hy, hw)) {}
Iso_rectangle_2(const Bbox_2& bbox)
: Rep(typename R::Construct_iso_rectangle_2()(Return_base_tag(), bbox.xmin(), bbox.ymin(), bbox.xmax(), bbox.ymax())) {}
: Rep(typename R::Construct_iso_rectangle_2()(Return_base_tag(),
typename R::Construct_point_2()(FT(bbox.xmin()), FT(bbox.ymin())),
typename R::Construct_point_2()(FT(bbox.xmax()), FT(bbox.ymax())))) {}
decltype(auto)
min BOOST_PREVENT_MACRO_SUBSTITUTION () const

View File

@ -18,25 +18,10 @@
#define CGAL_PRECISE_NUMBERS_H
#include <CGAL/config.h>
#if defined CGAL_USE_GMPXX
# include <CGAL/gmpxx.h>
typedef mpz_class Precise_integer;
typedef mpq_class Precise_rational;
#elif defined CGAL_USE_LEDA
# include <CGAL/leda_integer.h>
# include <CGAL/leda_rational.h>
typedef leda_integer Precise_integer;
typedef leda_rational Precise_rational;
#elif defined CGAL_USE_GMP
# include <CGAL/Gmpz.h>
# include <CGAL/Gmpq.h>
typedef CGAL::Gmpz Precise_integer;
typedef CGAL::Gmpq Precise_rational;
#else
# include <CGAL/MP_Float.h>
# include <CGAL/Quotient.h>
typedef CGAL::MP_Float Precise_integer;
typedef CGAL::Quotient<CGAL::MP_Float> Precise_rational;
#endif
#include <CGAL/Exact_integer.h>
#include <CGAL/Exact_rational.h>
using Precise_integer = CGAL::Exact_integer;
using Precise_rational = CGAL::Exact_rational;
#endif // CGAL_PRECISE_NUMBERS_H

View File

@ -14,20 +14,10 @@
//
// Author(s) : Laurent Rineau
#include <CGAL/config.h>
#include <CGAL/boost_mp.h>
#if CGAL_USE_GMPXX
# include <CGAL/gmpxx.h>
#elif CGAL_USE_GMP
# include <CGAL/Gmpz.h>
#elif CGAL_USE_LEDA
# include <CGAL/leda_integer.h>
#elif CGAL_USE_CORE
# include <CGAL/CORE_BigInt.h>
#elif defined CGAL_USE_BOOST_MP
#else
# error CGAL is configured with none of GMP, LEDA, Boost.Multiprecision and CORE. <CGAL/Exact_integer.h> cannot be used.
#endif
#ifndef CGAL_EXACT_INTEGER_H
#define CGAL_EXACT_INTEGER_H
#include <CGAL/Number_types/internal/Exact_type_selector.h>
namespace CGAL {
@ -50,29 +40,10 @@ typedef unspecified_type Exact_integer;
#else // not DOXYGEN_RUNNING
#if BOOST_VERSION > 107900 && defined(CGAL_USE_BOOST_MP)
// use boost-mp by default
typedef BOOST_cpp_arithmetic_kernel::Integer Exact_integer;
#else // BOOST_VERSION > 107900
#ifdef CGAL_USE_GMPXX
typedef mpz_class Exact_integer;
#elif defined(CGAL_USE_GMP)
#if defined(CGAL_USE_BOOST_MP)
typedef BOOST_gmp_arithmetic_kernel::Integer Exact_integer;
#else
typedef Gmpz Exact_integer;
#endif
#elif defined(CGAL_USE_LEDA)
typedef leda_integer Exact_integer;
#elif defined(CGAL_USE_BOOST_MP)
typedef BOOST_cpp_arithmetic_kernel::Integer Exact_integer;
#elif defined(CGAL_USE_CORE)
typedef CORE::BigInt Exact_integer;
#else
#error "ERROR: Cannot determine a BigInt type!"
#endif // CGAL_USE_CORE
#endif // BOOST_VERSION > 107800
typedef internal::Exact_ring_selector<int>::Type Exact_integer;
#endif // not DOXYGEN_RUNNING
} /* end namespace CGAL */
#endif // CGAL_EXACT_INTEGER_H

View File

@ -25,6 +25,7 @@
#include <CGAL/Lazy_exact_nt.h>
#include <CGAL/boost_mp.h>
#ifdef CGAL_USE_GMP
# include <CGAL/Gmpz.h>
# include <CGAL/Gmpq.h>
@ -49,28 +50,59 @@ class Expr;
namespace CGAL { namespace internal {
// Two classes which tell the preferred "exact number types" corresponding to a type.
// Exact_ring_selector<double> and Exact_field_selector<double> are used by EPICK as exact number type
// to answer predicates at the end of the filtering chain of predicates and EPECK uses
// Exact_field_selector<double> for as its exact number type.
// The default template chooses mpq_class, Gmpq, leda_rational, or Quotient<MP_Float>.
// It should support the built-in types.
template < typename >
struct Exact_field_selector
// Warning, the order in this list must match the one in Installation/lib/cmake/CGALConfig.cmake
enum ENT_backend_choice
{
GMP_BACKEND,
GMPXX_BACKEND,
BOOST_GMP_BACKEND,
BOOST_BACKEND,
LEDA_BACKEND,
MP_FLOAT_BACKEND
};
#if BOOST_VERSION > 107900 && defined(CGAL_USE_BOOST_MP)
// use boost-mp by default
// Boost
{ typedef BOOST_cpp_arithmetic_kernel::Rational Type; };
#else // BOOST_VERSION > 107900
#ifdef CGAL_USE_GMPXX
{ typedef mpq_class Type; };
#elif defined(CGAL_USE_GMP)
#if defined(CGAL_USE_BOOST_MP)
{ typedef BOOST_gmp_arithmetic_kernel::Rational Type; };
template <ENT_backend_choice>
struct Exact_NT_backend;
#ifdef CGAL_USE_GMP
template <>
struct Exact_NT_backend<GMP_BACKEND>
: public GMP_arithmetic_kernel
{
#ifdef CGAL_HAS_MPZF
typedef Mpzf Ring_for_float;
#else
{ typedef Gmpq Type; };
typedef Gmpzf Ring_for_float;
#endif
#elif defined(CGAL_USE_LEDA)
{ typedef leda_rational Type; };
#elif defined(CGAL_USE_BOOST_MP)
};
#endif
#ifdef CGAL_USE_GMPXX
template <>
struct Exact_NT_backend<GMPXX_BACKEND>
: public GMPXX_arithmetic_kernel
{
typedef Exact_NT_backend<GMP_BACKEND>::Ring_for_float Ring_for_float;
};
#endif
#if defined (CGAL_USE_BOOST_MP) && defined(CGAL_USE_GMP)
template <>
struct Exact_NT_backend<BOOST_GMP_BACKEND>
: public BOOST_gmp_arithmetic_kernel
{
typedef Exact_NT_backend<GMP_BACKEND>::Ring_for_float Ring_for_float;
};
#endif
#ifdef CGAL_USE_BOOST_MP
template <>
struct Exact_NT_backend<BOOST_BACKEND>
{
// See the discussion in https://github.com/CGAL/cgal/pull/3614
// This is disabled for now because cpp_rational is even slower than Quotient<MP_Float>. Quotient<cpp_int> will be a good candidate after some polishing.
// In fact, the new version of cpp_rational from here: https://github.com/boostorg/multiprecision/pull/366
@ -78,28 +110,69 @@ struct Exact_field_selector
// while Quotient does not. Though, we can still use it if needed.
#if BOOST_VERSION <= 107800
// See this comment: https://github.com/CGAL/cgal/pull/5937#discussion_r721533675
{ typedef Quotient<boost::multiprecision::cpp_int> Type; };
typedef Quotient<boost::multiprecision::cpp_int> Rational;
#else
{ typedef BOOST_cpp_arithmetic_kernel::Rational Type; };
typedef BOOST_cpp_arithmetic_kernel::Rational Rational;
#endif
#else
{ typedef Quotient<MP_Float> Type; };
typedef boost::multiprecision::cpp_int Integer;
typedef cpp_float Ring_for_float;
};
#endif
#endif // BOOST_VERSION > 107900
// By default, a field is a safe choice of ring.
template < typename T >
struct Exact_ring_selector : Exact_field_selector < T > { };
#ifdef CGAL_USE_LEDA
template <>
struct Exact_NT_backend<LEDA_BACKEND>
: public LEDA_arithmetic_kernel
{
typedef leda_rational Ring_for_float;
};
#endif
template <>
struct Exact_NT_backend<MP_FLOAT_BACKEND>
: public MP_Float_arithmetic_kernel
{
typedef MP_Float Ring_for_float;
};
#ifndef CMAKE_OVERRIDDEN_DEFAULT_ENT_BACKEND
constexpr ENT_backend_choice Default_exact_nt_backend =
#ifdef CGAL_USE_GMPXX
GMPXX_BACKEND;
#elif defined(CGAL_USE_GMP)
#if defined(CGAL_USE_BOOST_MP)
BOOST_GMP_BACKEND;
#else
GMP_BACKEND;
#endif
#elif BOOST_VERSION > 107900 && defined(CGAL_USE_BOOST_MP)
BOOST_BACKEND;
#elif defined(CGAL_USE_LEDA)
LEDA_BACKEND;
#else
MP_FLOAT_BACKEND;
#endif
#else
constexpr ENT_backend_choice Default_exact_nt_backend = static_cast<ENT_backend_choice>(CMAKE_OVERRIDDEN_DEFAULT_ENT_BACKEND);
#endif
template < typename >
struct Exact_field_selector
{
using Type = typename Exact_NT_backend<Default_exact_nt_backend>::Rational;
};
template < typename >
struct Exact_ring_selector
{
using Type = typename Exact_NT_backend<Default_exact_nt_backend>::Integer;
};
template <>
struct Exact_ring_selector<double>
#ifdef CGAL_HAS_MPZF
{ typedef Mpzf Type; };
#elif defined(CGAL_HAS_THREADS) || !defined(CGAL_USE_GMP)
{ typedef MP_Float Type; };
#else
{ typedef Gmpzf Type; };
#endif
{
using Type = typename Exact_NT_backend<Default_exact_nt_backend>::Ring_for_float;
};
template <>
struct Exact_ring_selector<float> : Exact_ring_selector<double> { };
@ -133,6 +206,10 @@ struct Exact_ring_selector<Gmpzf>
template <>
struct Exact_field_selector<Gmpq>
{ typedef Gmpq Type; };
template <>
struct Exact_ring_selector<Gmpq>
{ typedef Gmpq Type; };
#endif
#ifdef CGAL_USE_GMPXX
@ -147,6 +224,10 @@ struct Exact_ring_selector< ::mpz_class>
template <>
struct Exact_field_selector< ::mpq_class>
{ typedef ::mpq_class Type; };
template <>
struct Exact_ring_selector< ::mpq_class>
{ typedef ::mpq_class Type; };
#endif
#ifdef CGAL_USE_LEDA
@ -162,6 +243,10 @@ template <>
struct Exact_field_selector<leda_rational>
{ typedef leda_rational Type; };
template <>
struct Exact_ring_selector<leda_rational>
{ typedef leda_rational Type; };
template <>
struct Exact_field_selector<leda_real>
{ typedef leda_real Type; };
@ -171,6 +256,28 @@ struct Exact_field_selector<leda_real>
template <>
struct Exact_field_selector<CORE::Expr>
{ typedef CORE::Expr Type; };
template <>
struct Exact_ring_selector<CORE::Expr>
{ typedef CORE::Expr Type; };
#endif
#ifdef CGAL_USE_BOOST_MP
template <>
struct Exact_field_selector<Exact_NT_backend<BOOST_BACKEND>::Integer>
{ typedef Exact_NT_backend<BOOST_BACKEND>::Rational Type; };
template <>
struct Exact_ring_selector<Exact_NT_backend<BOOST_BACKEND>::Integer>
{ typedef Exact_NT_backend<BOOST_BACKEND>::Integer Type; };
template <>
struct Exact_field_selector<Exact_NT_backend<BOOST_BACKEND>::Rational>
{ typedef Exact_NT_backend<BOOST_BACKEND>::Rational Type; };
template <>
struct Exact_ring_selector<Exact_NT_backend<BOOST_BACKEND>::Rational>
{ typedef Exact_NT_backend<BOOST_BACKEND>::Rational Type; };
#endif
template < typename ET >

View File

@ -897,6 +897,7 @@ template< > class Real_embeddable_traits< Quotient<boost::multiprecision::cpp_in
} //namespace CGAL
#include <CGAL/BOOST_MP_arithmetic_kernel.h>
#include <CGAL/cpp_float.h>
#endif // BOOST_VERSION
#endif

View File

@ -0,0 +1,581 @@
// Copyright (c) 2023 GeometryFactory (France), INRIA Saclay - Ile de France (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org)
//
// $URL$
// $Id$
// SPDX-License-Identifier: LGPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Andreas Fabri, Marc Glisse
#ifndef CGAL_CPP_FLOAT_H
#define CGAL_CPP_FLOAT_H
//#define CGAL_DEBUG_CPPF
#include <CGAL/boost_mp.h>
#include <CGAL/assertions.h>
#include <boost/multiprecision/cpp_int.hpp>
#include <iostream>
namespace CGAL {
namespace internal {
// Only used with an argument known not to be 0.
inline int low_bit (boost::uint64_t x) {
#if defined(_MSC_VER)
unsigned long ret;
_BitScanForward64(&ret, x);
return (int)ret;
#elif defined(__xlC__)
return __cnttz8 (x);
#else
// Assume long long is 64 bits
return __builtin_ctzll (x);
#endif
}
inline int high_bit (boost::uint64_t x) {
#if defined(_MSC_VER)
unsigned long ret;
_BitScanReverse64(&ret, x);
return (int)ret;
#elif defined(__xlC__)
// Macro supposedly not defined on z/OS.
return 63 - __cntlz8 (x);
#else
return 63 - __builtin_clzll (x);
#endif
}
} // namespace internal
#if 0 // needs C++20
template <class T>
void fmt(const T& t)
{
std::cout << std::format("{:b}", t) << std::endl;
}
#endif
// It would have made sense to make this
// boost::multiprecision::number<some_new_backend>, but we keep that
// for later when we contribute to boost::mp
class cpp_float {
#ifdef CGAL_DEBUG_CPPF
boost::multiprecision::cpp_rational rat;
#endif
typedef boost::multiprecision::number<boost::multiprecision::cpp_int_backend<512> > Mantissa;
Mantissa man;
int exp; /* The number man (an integer) * 2 ^ exp */
public:
cpp_float()
:
#ifdef CGAL_DEBUG_CPPF
rat(),
#endif
man(), exp()
{}
cpp_float(short i)
:
#ifdef CGAL_DEBUG_CPPF
rat(i),
#endif
man(i),exp(0)
{}
cpp_float(int i)
:
#ifdef CGAL_DEBUG_CPPF
rat(i),
#endif
man(i),exp(0)
{}
cpp_float(long i)
:
#ifdef CGAL_DEBUG_CPPF
rat(i),
#endif
man(i),exp(0)
{}
#ifdef CGAL_DEBUG_CPPF
cpp_float(const Mantissa& man, int exp, const boost::multiprecision::cpp_rational& rat)
: rat(rat), man(man),exp(exp)
{}
#else
cpp_float(const Mantissa& man, int exp)
: man(man), exp(exp)
{
CGAL_HISTOGRAM_PROFILER("size (man/exp)", static_cast<int>(man.backend().size()));
}
#ifndef CGAL_DEBUG_CPPF
template <typename Expression>
cpp_float(const Expression& man, int exp)
:man(man), exp(exp)
{
CGAL_HISTOGRAM_PROFILER("size (expression/exp)", static_cast<int>(this->man.backend().size()));
}
#endif
#endif
cpp_float(double d)
#ifdef CGAL_DEBUG_CPPF
: rat(d)
#endif
{
// std::cout << "\ndouble = " << d << std::endl;
using boost::uint64_t;
union {
#ifdef CGAL_LITTLE_ENDIAN
struct { uint64_t man:52; uint64_t exp:11; uint64_t sig:1; } s;
#else /* CGAL_BIG_ENDIAN */
//WARNING: untested!
struct { uint64_t sig:1; uint64_t exp:11; uint64_t man:52; } s;
#endif
double d;
} u;
u.d = d;
uint64_t m;
uint64_t dexp = u.s.exp;
CGAL_assertion_msg(dexp != 2047, "Creating an cpp_float from infinity or NaN.");
if (dexp == 0) {
if (d == 0) { exp=0; return; }
else { // denormal number
m = u.s.man;
++dexp;
}
} else {
m = (1LL << 52) | u.s.man;
}
int idexp = (int)dexp;
idexp -= 1023;
// std::cout << "m = " << m << std::endl;
// std::cout << "idexp = " << idexp << std::endl;
int shifted = internal::low_bit(m);
m >>= shifted;
int nbits = internal::high_bit(m);
// std::cout << "nbits = " << nbits << std::endl;
exp = idexp - nbits;
man = m;
if(u.s.sig){
man.backend().negate();
}
#ifdef CGAL_DEBUG_CPPF
assert(rat.sign() == man.sign());
#endif
// std::cout << "m = " << m << " * 2^" << exp << std::endl;
// fmt(m);
CGAL_HISTOGRAM_PROFILER("size when constructed from double", static_cast<int>(man.backend().size()));
}
friend std::ostream& operator<<(std::ostream& os, const cpp_float& m)
{
return os << m.to_double();
#if 0 // dehug output
return os << m.man << " * 2 ^ " << m.exp << " ( " << m.to_double() << ") "
#ifdef CGAL_DEBUG_CPPF
<< " " << m.rat
#endif
;
#endif
}
friend std::istream& operator>>(std::istream& is, cpp_float& m)
{
double d;
is >> d;
m = cpp_float(d);
return is;
}
friend cpp_float operator-(cpp_float const&x)
{
#ifdef CGAL_DEBUG_CPPF
return cpp_float(-x.man,x.exp, -x.rat);
#else
return cpp_float(-x.man,x.exp);
#endif
}
cpp_float& operator*=(const cpp_float& other)
{
#ifdef CGAL_DEBUG_CPPF
rat *= other.rat;
#endif
man *= other.man;
exp += other.exp;
return *this;
}
friend
cpp_float operator*(const cpp_float& a, const cpp_float&b){
#ifdef CGAL_DEBUG_CPPF
return cpp_float(a.man*b.man, a.exp+b.exp, a.rat * b.rat);
#else
return cpp_float(a.man*b.man, a.exp+b.exp);
#endif
}
// Marc Glisse commented on github:
// We can sometimes end up with a mantissa that has quite a few zeros at the end,
// but the cases where the mantissa is too long by more than 1 limb should be negligible,
// and normalizing so the mantissa is always odd would have a cost.
cpp_float operator+=(const cpp_float& other)
{
#ifdef CGAL_DEBUG_CPPF
rat += other.rat;
#endif
int shift = exp - other.exp;
if(shift > 0){
man = (man << shift) + other.man;
exp = other.exp;
}else if(shift < 0){
man = man + (other.man << -shift);
}else{
man += other.man;
}
return *this;
}
#ifdef CGAL_DEBUG_CPPF
friend
cpp_float operator+(const cpp_float& a, const cpp_float&b){
int shift = a.exp - b.exp;
if(shift > 0){
return cpp_float((a.man << shift) + b.man, b.exp, a.rat+b.rat);
}else if(shift < 0){
return cpp_float(a.man + (b.man << -shift), a.exp, a.rat+b.rat);
}
return cpp_float(a.man + b.man, a.exp, a.rat+b.rat);
}
#else
friend
cpp_float operator+(const cpp_float& a, const cpp_float&b){
int shift = a.exp - b.exp;
CGAL_HISTOGRAM_PROFILER("shift+", CGAL::abs(shift));
if(shift > 0){
return cpp_float((a.man << shift) + b.man, b.exp);
}else if(shift < 0){
return cpp_float(a.man + (b.man << -shift), a.exp);
}
return cpp_float(a.man + b.man, a.exp);
}
#endif
cpp_float operator-=(const cpp_float& other)
{
#ifdef CGAL_DEBUG_CPPF
rat -= other.rat;
#endif
int shift = exp - other.exp;
if(shift > 0){
man <<= shift;
man -= other.man;
exp = other.exp;
}else if(shift < 0){
man -= (other.man << -shift);
}else{
man -= other.man;
}
return *this;
}
#ifdef CGAL_DEBUG_CPPF
friend
cpp_float operator-(const cpp_float& a, const cpp_float&b){
int shift = a.exp - b.exp;
if(shift > 0){
return cpp_float((a.man << shift) - b.man, b.exp, a.rat-b.rat);
}else if(shift < 0){
return cpp_float(a.man - (b.man << -shift), a.exp, a.rat-b.rat);
}
return cpp_float(a.man - b.man, a.exp, a.rat-b.rat);
}
#else
friend
cpp_float operator-(const cpp_float& a, const cpp_float&b){
int shift = a.exp - b.exp;
CGAL_HISTOGRAM_PROFILER("shift-", CGAL::abs(shift));
if(shift > 0){
return cpp_float((a.man << shift) - b.man, b.exp);
}else if(shift < 0){
return cpp_float(a.man - (b.man << -shift), a.exp);
}
return cpp_float(a.man - b.man, a.exp);
}
#endif
bool is_positive() const
{
return CGAL::is_positive(man);
}
bool is_negative() const
{
return CGAL::is_negative(man);
}
// Would it make sense to compare the sign of the exponent?
// to distinguish the interval between ]-1,1[ from the rest.
friend bool operator<(const cpp_float& a, const cpp_float& b)
{
if(((! a.is_positive()) && b.is_positive())
|| (a.is_negative() && b.is_zero()))return true;
if(((! b.is_positive()) && a.is_positive())
|| (b.is_negative() && a.is_zero()))return false;
#ifdef CGAL_DEBUG_CPPF
bool qres = a.rat < b.rat;
#endif
cpp_float d = b-a;
#ifdef CGAL_DEBUG_CPPF
assert(qres == d.is_positive());
#endif
return d.is_positive();
}
friend bool operator>(cpp_float const&a, cpp_float const&b){
return b<a;
}
friend bool operator>=(cpp_float const&a, cpp_float const&b){
return !(a<b);
}
friend bool operator<=(cpp_float const&a, cpp_float const&b){
return !(a>b);
}
friend bool operator==(cpp_float const&a, cpp_float const&b){
#ifdef CGAL_DEBUG_CPPF
bool qres = a.rat == b.rat;
#endif
int shift = a.exp - b.exp;
CGAL_HISTOGRAM_PROFILER("shift==", CGAL::abs(shift));
if(shift > 0){
Mantissa ac = a.man << shift;
#ifdef CGAL_DEBUG_CPPF
assert( qres == (ac == b.man));
#endif
return ac == b.man;
}else if(shift < 0){
Mantissa bc = b.man << -shift;
#ifdef CGAL_DEBUG_CPPF
assert(qres == (a.man == bc));
#endif
return a.man == bc;
}
#ifdef CGAL_DEBUG_CPPF
assert(qres == (a.man == b.man));
#endif
return a.man==b.man;
}
Comparison_result compare(const cpp_float& other) const
{
if(*this < other) return SMALLER;
if(*this > other) return LARGER;
return EQUAL;
}
friend bool operator!=(cpp_float const&a, cpp_float const&b){
return !(a==b);
}
double to_double() const
{
if(exp == 0){
return CGAL::to_double(man);
}
if(exp > 0){
Mantissa as(man);
as <<= exp;
return CGAL::to_double(as);
}
Mantissa pow(1);
pow <<= -exp;
boost::multiprecision::cpp_rational r(man, pow);
return CGAL::to_double(r);
}
std::pair<double,double> to_interval() const
{
if(exp == 0){
return CGAL::to_interval(man);
}
if(exp > 0){
Mantissa as = man << exp;
return CGAL::to_interval(as);
}
Mantissa pow(1);
pow <<= -exp;
boost::multiprecision::cpp_rational r(man, pow);
return CGAL::to_interval(r);
}
bool is_zero () const {
return CGAL::is_zero(man);
}
bool is_one () const {
if(! is_positive()) return false;
int msb = static_cast<int>(boost::multiprecision::msb(man));
if (msb != -exp) return false;
int lsb = static_cast<int>(boost::multiprecision::lsb(man));
return (msb == lsb);
}
CGAL::Sign sign () const
{
return CGAL::sign(man);
}
std::size_t size() const
{
return man.backend().size();
}
};
template <> struct Algebraic_structure_traits< cpp_float >
: public Algebraic_structure_traits_base< cpp_float, Integral_domain_without_division_tag > {
typedef Tag_true Is_exact;
typedef Tag_false Is_numerical_sensitive;
struct Is_zero
: public CGAL::cpp98::unary_function< Type, bool > {
bool operator()( const Type& x ) const {
return x.is_zero();
}
};
struct Is_one
: public CGAL::cpp98::unary_function< Type, bool > {
bool operator()( const Type& x ) const {
return x.is_one();
}
};
struct Gcd
: public CGAL::cpp98::binary_function< Type, Type, Type > {
Type operator()(
const Type& /* x */,
const Type& /* y */ ) const {
assert(false);
return Type(); // cpp_float_gcd(x, y);
}
};
struct Square
: public CGAL::cpp98::unary_function< Type, Type > {
Type operator()( const Type& x ) const {
return x*x ; // cpp_float_square(x);
}
};
struct Integral_division
: public CGAL::cpp98::binary_function< Type, Type, Type > {
Type operator()(
const Type& /* x */,
const Type& /* y */ ) const {
assert(false);
return Type(); // x / y;
}
};
struct Sqrt
: public CGAL::cpp98::unary_function< Type, Type > {
Type operator()( const Type& /* x */) const {
assert(false);
return Type(); // cpp_float_sqrt(x);
}
};
struct Is_square
: public CGAL::cpp98::binary_function< Type, Type&, bool > {
bool operator()( const Type& /* x */, Type& /* y */ ) const {
// TODO: avoid doing 2 calls.
assert(false);
return true;
}
bool operator()( const Type& /* x */) const {
assert(false);
return true;
}
};
};
template <> struct Real_embeddable_traits< cpp_float >
: public INTERN_RET::Real_embeddable_traits_base< cpp_float , CGAL::Tag_true > {
struct Sgn
: public CGAL::cpp98::unary_function< Type, ::CGAL::Sign > {
::CGAL::Sign operator()( const Type& x ) const {
return x.sign();
}
};
struct To_double
: public CGAL::cpp98::unary_function< Type, double > {
double operator()( const Type& x ) const {
return x.to_double();
}
};
struct Compare
: public CGAL::cpp98::binary_function< Type, Type, Comparison_result > {
Comparison_result operator()(
const Type& x,
const Type& y ) const {
return x.compare(y);
}
};
struct To_interval
: public CGAL::cpp98::unary_function< Type, std::pair< double, double > > {
std::pair<double, double> operator()( const Type& x ) const {
return x.to_interval();
}
};
};
CGAL_DEFINE_COERCION_TRAITS_FOR_SELF(cpp_float)
CGAL_DEFINE_COERCION_TRAITS_FROM_TO(short ,cpp_float)
CGAL_DEFINE_COERCION_TRAITS_FROM_TO(int ,cpp_float)
CGAL_DEFINE_COERCION_TRAITS_FROM_TO(long ,cpp_float)
CGAL_DEFINE_COERCION_TRAITS_FROM_TO(float ,cpp_float)
CGAL_DEFINE_COERCION_TRAITS_FROM_TO(double ,cpp_float)
} // namespace CGAL
#endif // CGAL_CPP_FLOAT_H

View File

@ -11,6 +11,7 @@ find_package(CGAL REQUIRED COMPONENTS Core)
include_directories(BEFORE include)
create_single_source_cgal_program("bench_interval.cpp")
create_single_source_cgal_program("cpp_float.cpp")
create_single_source_cgal_program("constant.cpp")
create_single_source_cgal_program("CORE_BigFloat.cpp")
create_single_source_cgal_program("CORE_BigInt.cpp")

View File

@ -0,0 +1,116 @@
#ifdef CGAL_DO_NOT_USE_BOOST_MP
#include <iostream>
int main()
{
std::cout << "The class CGAL::cpp_float is not tested on this platform" << std::endl;
return 0;
}
#else
#include <CGAL/cpp_float.h>
template<class NT>
void test1(){
NT z;
NT a=3;
NT b=4.5;
NT c=2*(a+b)+-a*5;
assert(CGAL::sign(c)==0);
NT e=.0003;
NT f=1e-90;
assert(CGAL::to_double(b) == 4.5);
std::pair<double,double> p=CGAL::to_interval(b);
assert(p.first<=4.5 && p.second >= 4.5);
assert(a<b && CGAL::compare(b,a)>0);
assert(z<f && CGAL::compare(z,f)<0 && CGAL::compare(f,z)>0);
assert(z==z && CGAL::compare(z,z)==0);
assert(CGAL::square(b)*4==81);
assert(CGAL::is_zero(c));
assert(!CGAL::is_zero(a));
assert(!CGAL::is_one(a));
assert(CGAL::is_one(a-2));
assert(e-e==0);
}
int main()
{
test1<CGAL::cpp_float>();
double d = -0;
CGAL::cpp_float zero(d);
assert(! CGAL::is_positive(zero));
assert(! CGAL::is_negative(zero));
assert(CGAL::is_zero(zero));
CGAL::cpp_float m(0);
std::cout << m << std::endl;
CGAL::cpp_float m0(23.0);
CGAL::cpp_float m1(5.0);
CGAL::cpp_float m2(5.125);
CGAL::cpp_float m3(2.5);
CGAL::cpp_float m4(0.625);
CGAL::cpp_float m5(0.5);
CGAL::cpp_float m6(-0.625);
CGAL::is_positive(m5);
CGAL::is_negative(m6);
assert(m < m5);
assert(m6 < m);
assert(! (m5 < m));
assert(! (m < m6));
assert(! (m6 < m6));
assert(m4 > m5);
assert(-m4 < -m5);
assert(-m4 == -m4);
assert(-m4 != -m5);
assert(-m5 != -m4);
CGAL::cpp_float m15 = m1 + m5;
m15 = m15 - m5;
assert(m15 == m1);
assert(! (m15 < m1));
assert(! (m1 < m15));
m15 = m15 - m15;
std::cout << m15 << std::endl;
assert(m15.is_zero());
m1 *= m5;
m0 += m4;
std::cout << m0 << std::endl;
CGAL::cpp_float one(1);
assert(one.is_one());
one += m4;
one -= m4;
std::cout << one << std::endl;
assert(one.is_one());
return 0;
}
#endif

View File

@ -8,6 +8,26 @@ project(Triangulation_3)
# CGAL and its components
find_package(CGAL REQUIRED)
# Boost and its components
find_package(Boost REQUIRED)
if(NOT Boost_FOUND)
message(
STATUS "This project requires the Boost library, and will not be compiled.")
return()
endif()
# include for local directory
# include for local package
# Creating entries for all C++ files with "main" routine
# ##########################################################
create_single_source_cgal_program("ocean.cpp")
create_single_source_cgal_program("incident_edges.cpp")
create_single_source_cgal_program("simple_2.cpp")
create_single_source_cgal_program("simple.cpp")

View File

@ -0,0 +1,20 @@
#include <fstream>
#include <CGAL/Random.h>
int main()
{
int N=100;
std::cout.precision(17);
CGAL::Random rng;
for(int i = 0; i < N; ++i){
double x = rng.get_double(-1.0, 1.0);
double y = rng.get_double(-1.0, 1.0);
for(int j = 0; j < N; j++){
std::cout << x << " " << y << " " << rng.get_double(-1.0, 1.0) << std::endl;
}
}
return 0;
}

View File

@ -1,6 +1,6 @@
//#define CGAL_PROFILE
#define CGAL_USE_SSE2_FABS
#define CGAL_USE_SSE2_MAX
//#define CGAL_USE_SSE2_FABS
//#define CGAL_USE_SSE2_MAX
//#define CGAL_MSVC_USE_STD_FABS // use this one with precise
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
@ -11,25 +11,31 @@
#include <CGAL/Timer.h>
#include <iostream>
#include <string>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Delaunay_triangulation_3<K> DT;
typedef DT::Point Point_3;
typedef CGAL::Timer Timer;
int main()
int main(int argc, char* argv[])
{
const std::string filename = (argc > 1) ? argv[1] : CGAL::data_file_path("points_3/ocean_r.xyz");
std::ifstream in(filename.c_str());
std::vector<Point_3> points;
Point_3 p;
Point_3 p, q;
while(std::cin >> p){
while(in >> p ){
points.push_back(p);
}
std::cout << points.size() << " points read\n";
Timer timer;
timer.start();
size_t N = 0;
for(int i = 0; i < 5; i++){
for(int i = 0; i < 1; i++){
DT dt;
dt.insert(points.begin(), points.end());
N += dt.number_of_cells();