-inline
-RS_AK1::Algebraic_1 min
-BOOST_PREVENT_MACRO_SUBSTITUTION(RS_AK1::Algebraic_1
a,
- RS_AK1::Algebraic_1
b){
- return(a
-inline
-RS_AK1::Algebraic_1
max
-BOOST_PREVENT_MACRO_SUBSTITUTION(RS_AK1::Algebraic_1
a,
- RS_AK1::Algebraic_1
b){
- return(a>b?a:b);
-}
-
-template
-inline
-std::ostream& operator<<(std::ostream &o,
- const RS_AK1::Algebraic_1 &a){
- return(o<<'['<
-inline
-std::istream& operator>>(std::istream &i,
- RS_AK1::Algebraic_1 &a){
- std::istream::int_type c;
- P pol;
- B lb,rb;
- c=i.get();
- if(c!='['){
- CGAL_error_msg("error reading istream, \'[\' expected");
- return i;
- }
- i>>pol;
- c=i.get();
- if(c!=','){
- CGAL_error_msg("error reading istream, \',\' expected");
- return i;
- }
- i>>lb;
- c=i.get();
- if(c!=','){
- CGAL_error_msg("error reading istream, \',\' expected");
- return i;
- }
- i>>rb;
- c=i.get();
- if(c!=']'){
- CGAL_error_msg("error reading istream, \']\' expected");
- return i;
- }
- a=RS_AK1::Algebraic_1
(pol,lb,rb);
- return i;
-}
-
-} // namespace CGAL
-
-#endif // CGAL_RS_ALGEBRAIC_1_H
diff --git a/Algebraic_kernel_d/include/CGAL/RS/algebraic_z_1.h b/Algebraic_kernel_d/include/CGAL/RS/algebraic_z_1.h
deleted file mode 100644
index 1b0e3b0f081a..000000000000
--- a/Algebraic_kernel_d/include/CGAL/RS/algebraic_z_1.h
+++ /dev/null
@@ -1,351 +0,0 @@
-// Copyright (c) 2006-2013 INRIA Nancy-Grand Est (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: Luis Peñaranda
-
-#ifndef CGAL_RS_ALGEBRAIC_Z_1_H
-#define CGAL_RS_ALGEBRAIC_Z_1_H
-
-#include "algebraic_1.h"
-
-namespace CGAL{
-namespace RS_AK1{
-
-// This class represents an algebraic number storing two polynomials of
-// which it is root: one of them is the given polynomial; the other one is
-// an integer polynomial, which is used to perform the operations such as
-// refinements. Since RS works only with integer polynomials, this
-// architecture permits to convert the input polynomials only once.
-
-template
-class Algebraic_z_1:
-public Algebraic_1,
-boost::totally_ordered,
- double>{
- protected:
- typedef Polynomial_ Polynomial;
- typedef ZPolynomial_ ZPolynomial;
- typedef Bound_ Bound;
- typedef ZRefiner_ ZRefiner;
- typedef ZComparator_ ZComparator;
- typedef ZPtraits_ ZPtraits;
- typedef typename ZPtraits::Coefficient_type ZCoefficient;
- typedef typename ZPtraits::Scale ZScale;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Coefficient_type Coefficient;
- typedef typename Ptraits::Scale Scale;
- typedef Algebraic_1
- Algebraic_base;
- typedef Algebraic_z_1
- Algebraic;
- ZPolynomial zpol;
-
- public:
- Algebraic_z_1(){};
- Algebraic_z_1(const Polynomial &p,
- const ZPolynomial &zp,
- const Bound &l,
- const Bound &r):Algebraic_base(p,l,r),zpol(zp){
- CGAL_assertion(l<=r);
- }
- Algebraic_z_1(const Algebraic &a):
- Algebraic_base(a.pol,a.left,a.right),zpol(a.zpol){}
- Algebraic_z_1(const Bound &b){
- typedef typename Ptraits::Shift Shift;
- typedef typename ZPtraits::Shift ZShift;
- this->left=b;
- this->right=b;
- Gmpq q(b);
- this->pol=Coefficient(mpq_denref(q.mpq()))*
- Shift()(Polynomial(1),1,0)-
- Coefficient(mpq_numref(q.mpq()));
- zpol=ZCoefficient(mpq_denref(q.mpq()))*
- ZShift()(ZPolynomial(1),1,0)-
- ZCoefficient(mpq_numref(q.mpq()));
- CGAL_assertion(this->left==this->right);
- CGAL_assertion(this->left==b);
- }
- template
- Algebraic_z_1(const T &t){
- typedef typename Ptraits::Shift Shift;
- typedef typename ZPtraits::Shift ZShift;
- CGAL::Gmpq q(t);
- this->pol=Coefficient(mpq_denref(q.mpq()))*
- Shift()(Polynomial(1),1,0)-
- Coefficient(mpq_numref(q.mpq()));
- zpol=ZCoefficient(mpq_denref(q.mpq()))*
- ZShift()(ZPolynomial(1),1,0)-
- ZCoefficient(mpq_numref(q.mpq()));
- this->left=Bound(t,std::round_toward_neg_infinity);
- this->right=Bound(t,std::round_toward_infinity);
- CGAL_assertion(this->left<=t&&this->right>=t);
- }
- Algebraic_z_1(const CGAL::Gmpq &q){
- typedef typename Ptraits::Shift Shift;
- typedef typename ZPtraits::Shift ZShift;
- this->pol=Coefficient(mpq_denref(q.mpq()))*
- Shift()(Polynomial(1),1,0)-
- Coefficient(mpq_numref(q.mpq()));
- zpol=ZCoefficient(mpq_denref(q.mpq()))*
- ZShift()(ZPolynomial(1),1,0)-
- ZCoefficient(mpq_numref(q.mpq()));
- this->left=Bound();
- this->right=Bound();
- mpfr_t b;
- mpfr_init(b);
- mpfr_set_q(b,q.mpq(),GMP_RNDD);
- mpfr_swap(b,this->left.fr());
- mpfr_set_q(b,q.mpq(),GMP_RNDU);
- mpfr_swap(b,this->right.fr());
- mpfr_clear(b);
- CGAL_assertion(this->left<=q&&this->right>=q);
- }
- ~Algebraic_z_1(){}
-
- ZPolynomial get_zpol()const{return zpol;}
- void set_zpol(const ZPolynomial &p){zpol=p;}
-
- Algebraic operator-()const{
- return Algebraic(Scale()(this->get_pol(),Coefficient(-1)),
- ZScale()(get_zpol(),ZCoefficient(-1)),
- -this->right,
- -this->left);
- }
-
-#define CGAL_RS_COMPARE_ALGEBRAIC_Z(_a) \
- (ZComparator()(get_zpol(),this->get_left(),this->get_right(), \
- (_a).get_zpol(),(_a).get_left(),(_a).get_right()))
-
- Comparison_result compare(Algebraic a)const{
- return CGAL_RS_COMPARE_ALGEBRAIC_Z(a);
- };
-
-#define CGAL_RS_COMPARE_ALGEBRAIC_Z_TYPE(_t) \
- bool operator<(_t t)const \
- {Algebraic a(t);return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)==CGAL::SMALLER;} \
- bool operator>(_t t)const \
- {Algebraic a(t);return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)==CGAL::LARGER;} \
- bool operator==(_t t)const \
- {Algebraic a(t);return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)==CGAL::EQUAL;}
-
- bool operator==(Algebraic a)const
- {return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)==CGAL::EQUAL;}
- bool operator!=(Algebraic a)const
- {return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)!=CGAL::EQUAL;}
- bool operator<(Algebraic a)const
- {return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)==CGAL::SMALLER;}
- bool operator<=(Algebraic a)const
- {return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)!=CGAL::LARGER;}
- bool operator>(Algebraic a)const
- {return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)==CGAL::LARGER;}
- bool operator>=(Algebraic a)const
- {return CGAL_RS_COMPARE_ALGEBRAIC_Z(a)!=CGAL::SMALLER;}
-
- CGAL_RS_COMPARE_ALGEBRAIC_Z_TYPE(double)
-
-#undef CGAL_RS_COMPARE_ALGEBRAIC_Z_TYPE
-#undef CGAL_RS_COMPARE_ALGEBRAIC_Z
-
-#ifdef IEEE_DBL_MANT_DIG
-#define CGAL_RS_DBL_PREC IEEE_DBL_MANT_DIG
-#else
-#define CGAL_RS_DBL_PREC 53
-#endif
- double to_double()const{
- typedef Real_embeddable_traits RT;
- typedef typename RT::To_double TD;
- ZRefiner()(get_zpol(),
- this->get_left(),
- this->get_right(),
- CGAL_RS_DBL_PREC);
- CGAL_assertion(TD()(this->get_left())==TD()(this->get_right()));
- return TD()(this->get_left());
- }
- std::pair to_interval()const{
- typedef Real_embeddable_traits RT;
- typedef typename RT::To_interval TI;
- return std::make_pair(TI()(this->get_left()).first,
- TI()(this->get_right()).second);
- }
-#undef CGAL_RS_DBL_PREC
-}; // class Algebraic_z_1
-
-} // namespace RS_AK1
-
-// We define Algebraic_z_1 as real-embeddable
-template
-class Real_embeddable_traits >:
-public INTERN_RET::Real_embeddable_traits_base<
- RS_AK1::Algebraic_z_1,
- CGAL::Tag_true>{
- typedef Polynomial_ P;
- typedef ZPolynomial_ ZP;
- typedef Bound_ B;
- typedef ZRefiner_ R;
- typedef ZComparator_ C;
- typedef Ptraits_ T;
- typedef ZPtraits_ ZT;
-
- public:
-
- typedef RS_AK1::Algebraic_z_1 Type;
- typedef CGAL::Tag_true Is_real_embeddable;
- typedef bool Boolean;
- typedef CGAL::Sign Sign;
- typedef CGAL::Comparison_result Comparison_result;
-
- typedef INTERN_RET::Real_embeddable_traits_base
- Base;
- typedef typename Base::Compare Compare;
-
- class Sgn:public CGAL::cpp98::unary_function{
- public:
- CGAL::Sign operator()(const Type &a)const{
- return Compare()(a,Type(0));
- }
- };
-
- class To_double:public CGAL::cpp98::unary_function{
- public:
- double operator()(const Type &a)const{return a.to_double();}
- };
-
- class To_interval:
- public CGAL::cpp98::unary_function >{
- public:
- std::pair operator()(const Type &a)const{
- return a.to_interval();
- }
- };
-
- class Is_zero:public CGAL::cpp98::unary_function{
- public:
- bool operator()(const Type &a)const{
- return Sgn()(a)==CGAL::ZERO;
- }
- };
-
- class Is_finite:public CGAL::cpp98::unary_function{
- public:
- bool operator()(const Type&)const{return true;}
- };
-
- class Abs:public CGAL::cpp98::unary_function{
- public:
- Type operator()(const Type &a)const{
- return Sgn()(a)==CGAL::NEGATIVE?-a:a;
- }
- };
-};
-
-template
-inline
-RS_AK1::Algebraic_z_1 min
-BOOST_PREVENT_MACRO_SUBSTITUTION(RS_AK1::Algebraic_z_1
a,
- RS_AK1::Algebraic_z_1
b){
- return(a
-inline
-RS_AK1::Algebraic_z_1
max
-BOOST_PREVENT_MACRO_SUBSTITUTION(RS_AK1::Algebraic_z_1
a,
- RS_AK1::Algebraic_z_1
b){
- return(a>b?a:b);
-}
-
-template
-inline
-std::ostream& operator<<(std::ostream &o,
- const RS_AK1::Algebraic_z_1 &a){
- return(o<<'['<
-inline
-std::istream& operator>>(std::istream &i,
- RS_AK1::Algebraic_z_1 &a){
- std::istream::int_type c;
- P pol;
- ZP zpol;
- B lb,rb;
- c=i.get();
- if(c!='['){
- CGAL_error_msg("error reading istream, \'[\' expected");
- return i;
- }
- i>>pol;
- c=i.get();
- if(c!=','){
- CGAL_error_msg("error reading istream, \',\' expected");
- return i;
- }
- i>>zpol;
- c=i.get();
- if(c!=','){
- CGAL_error_msg("error reading istream, \',\' expected");
- return i;
- }
- i>>lb;
- c=i.get();
- if(c!=','){
- CGAL_error_msg("error reading istream, \',\' expected");
- return i;
- }
- i>>rb;
- c=i.get();
- if(c!=']'){
- CGAL_error_msg("error reading istream, \']\' expected");
- return i;
- }
- a=RS_AK1::Algebraic_z_1
(pol,zpol,lb,rb);
- return i;
-}
-
-} // namespace CGAL
-
-#endif // CGAL_RS_ALGEBRAIC_Z_1_H
diff --git a/Algebraic_kernel_d/include/CGAL/RS/bisection_refiner_1.h b/Algebraic_kernel_d/include/CGAL/RS/bisection_refiner_1.h
deleted file mode 100644
index 28c655a8314e..000000000000
--- a/Algebraic_kernel_d/include/CGAL/RS/bisection_refiner_1.h
+++ /dev/null
@@ -1,165 +0,0 @@
-// Copyright (c) 2006-2013 INRIA Nancy-Grand Est (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: Luis Peñaranda
-
-// This file contains the simplest refiner, that bisects the interval a given
-// number of times.
-
-#ifndef CGAL_RS_BISECTION_REFINER_1_H
-#define CGAL_RS_BISECTION_REFINER_1_H
-
-#include
-#include "signat_1.h"
-#include "Gmpfr_make_unique.h"
-
-namespace CGAL{
-
-template
-struct Bisection_refiner_1{
- typedef CGAL::RS_AK1::Signat_1 Signat;
- void operator()(const Polynomial_&,Bound_&,Bound_&,int);
-}; // class Bisection_refiner_1
-
-// TODO: Write in a generic way, if possible (see next function).
-template
-void
-Bisection_refiner_1::
-operator()(const Polynomial_&,Bound_&,Bound_&,int){
- CGAL_error_msg("bisection refiner not implemented for these types");
- return;
-}
-
-// This works with any type of polynomial, but only for Gmpfr bounds.
-// TODO: Beyond writing generically, optimize this function. This would
-// remove the need for the next function, which is essentially the same.
-template<>
-void
-Bisection_refiner_1,Gmpfr>::
-operator()(const Polynomial &pol,Gmpfr &left,Gmpfr &right,int prec){
- typedef Polynomial Polynomial;
- typedef Polynomial_traits_d Ptraits;
- typedef Ptraits::Make_square_free Sfpart;
- typedef CGAL::RS_AK1::Signat_1
- Signat;
- CGAL_precondition(left<=right);
- //std::cout<<"refining ["<pc?pl:pc)+(mp_prec_t)prec;
- mpfr_init2(center,pc);
- CGAL_assertion_code(int round=)
- mpfr_prec_round(left.fr(),pc,GMP_RNDN);
- CGAL_assertion(!round);
- CGAL_assertion_code(round=)
- mpfr_prec_round(right.fr(),pc,GMP_RNDN);
- CGAL_assertion(!round);
- for(int i=0;i
-void
-Bisection_refiner_1,Gmpfr>::
-operator()(const Polynomial &pol,Gmpfr &left,Gmpfr &right,int prec){
- typedef Polynomial Polynomial;
- typedef Polynomial_traits_d Ptraits;
- typedef Ptraits::Make_square_free Sfpart;
- typedef CGAL::RS_AK1::Signat_1
- Signat;
- CGAL_precondition(left<=right);
- //std::cout<<"refining ["<pc?pl:pc)+(mp_prec_t)prec;
- mpfr_init2(center,pc);
- CGAL_assertion_code(int round=)
- mpfr_prec_round(left.fr(),pc,GMP_RNDN);
- CGAL_assertion(!round);
- CGAL_assertion_code(round=)
- mpfr_prec_round(right.fr(),pc,GMP_RNDN);
- CGAL_assertion(!round);
- for(int i=0;i
-
-#ifndef CGAL_RS_COMPARATOR_1_H
-#define CGAL_RS_COMPARATOR_1_H
-
-namespace CGAL{
-namespace RS_AK1{
-
-template
-struct Simple_comparator_1{
- typedef Polynomial_ Polynomial;
- typedef Bound_ Bound;
- typedef Refiner_ Refiner;
- typedef Signat_ Signat;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Degree Degree;
-
- CGAL::Comparison_result
- operator()(const Polynomial &p1,Bound &l1,Bound &r1,
- const Polynomial &p2,Bound &l2,Bound &r2)const{
- CGAL_precondition(l1<=r1&&l2<=r2);
- if(l1<=l2){
- if(r1l2?l1:l2);
- if(sleft==ZERO)
- return EQUAL;
- CGAL::Sign sright=sg(r1=l2:r2>=l1);
- return (r1
-
-#ifndef CGAL_RS_DYADIC_H
-#define CGAL_RS_DYADIC_H
-
-#include
-#include
-#include
-#include
-#include
-
-// for c++, compile with -lgmpxx
-#ifdef __cplusplus
-#include
-#endif
-
-#define CGALRS_dyadic_struct __mpfr_struct
-#define CGALRS_dyadic_t mpfr_t
-#define CGALRS_dyadic_ptr mpfr_ptr
-#define CGALRS_dyadic_srcptr mpfr_srcptr
-
-// some auxiliary defines
-#define CGALRS_dyadic_set_prec(D,P) \
- ( mpfr_set_prec( (D), (P)>MPFR_PREC_MIN?(P):MPFR_PREC_MIN) )
-
-#define CGALRS_dyadic_prec_round(D,P) \
- ( mpfr_prec_round( (D), (P)>MPFR_PREC_MIN?(P):MPFR_PREC_MIN, GMP_RNDN) )
-
-#define CGALRS_dyadic_set_exp(D,E) \
- ( CGAL_assertion( (E) <= mpfr_get_emax() && \
- (E) >= mpfr_get_emin() ) ,\
- mpfr_set_exp(D,E) )
-
-// init functions
-#define CGALRS_dyadic_init(D) mpfr_init2(D,MPFR_PREC_MIN)
-#define CGALRS_dyadic_init2(D,P) mpfr_init2(D,P)
-#define CGALRS_dyadic_clear(D) mpfr_clear(D)
-
-inline void CGALRS_dyadic_set(CGALRS_dyadic_ptr rop,CGALRS_dyadic_srcptr op){
- if(rop!=op){
- CGALRS_dyadic_set_prec(rop,mpfr_get_prec(op));
- mpfr_set(rop,op,GMP_RNDN);
- }
- CGAL_assertion(mpfr_equal_p(rop,op)!=0);
-}
-
-inline void CGALRS_dyadic_set_z(CGALRS_dyadic_ptr rop,mpz_srcptr z){
- size_t prec;
- prec=mpz_sizeinbase(z,2)-(mpz_tstbit(z,0)?0:mpz_scan1(z,0));
- CGALRS_dyadic_set_prec(rop,prec);
- mpfr_set_z(rop,z,GMP_RNDN);
- CGAL_assertion(!mpfr_cmp_z(rop,z));
-}
-
-inline void CGALRS_dyadic_set_si(CGALRS_dyadic_ptr rop,long s){
- CGALRS_dyadic_set_prec(rop,sizeof(long));
- mpfr_set_si(rop,s,GMP_RNDN);
- CGAL_assertion(!mpfr_cmp_si(rop,s));
-}
-
-inline void CGALRS_dyadic_set_ui(CGALRS_dyadic_ptr rop,unsigned long u){
- CGALRS_dyadic_set_prec(rop,sizeof(unsigned long));
- mpfr_set_ui(rop,u,GMP_RNDN);
- CGAL_assertion(!mpfr_cmp_ui(rop,u));
-}
-
-inline void CGALRS_dyadic_set_fr(CGALRS_dyadic_ptr rop,mpfr_srcptr op){
- if(rop!=op){
- CGALRS_dyadic_set_prec(rop,mpfr_get_prec(op));
- mpfr_set(rop,op,GMP_RNDN);
- CGAL_assertion(mpfr_equal_p(rop,op)!=0);
- }
-}
-
-#define CGALRS_dyadic_init_set(R,D) \
- ( CGALRS_dyadic_init(R), CGALRS_dyadic_set((R), (D)) )
-#define CGALRS_dyadic_init_set_z(R,Z) \
- ( CGALRS_dyadic_init(R), CGALRS_dyadic_set_z((R), (Z)) )
-#define CGALRS_dyadic_init_set_si(R,I) \
- ( CGALRS_dyadic_init(R), CGALRS_dyadic_set_si((R), (I)) )
-#define CGALRS_dyadic_init_set_ui(R,I) \
- ( CGALRS_dyadic_init(R), CGALRS_dyadic_set_ui((R), (I)) )
-#define CGALRS_dyadic_init_set_fr(R,F) \
- ( CGALRS_dyadic_init(R), CGALRS_dyadic_set_fr((R), (F)) )
-
-#define CGALRS_dyadic_get_fr(M,D) mpfr_set(M,D,GMP_RNDN)
-#define CGALRS_dyadic_get_d(D,RM) mpfr_get_d(D,RM)
-inline void CGALRS_dyadic_get_exactfr(mpfr_ptr rop,CGALRS_dyadic_srcptr op){
- if(rop!=op){
- CGALRS_dyadic_set_prec(rop,mpfr_get_prec(op));
- mpfr_set(rop,op,GMP_RNDN);
- CGAL_assertion(mpfr_equal_p(rop,op)!=0);
- }
-}
-
-#define CGALRS_dyadic_canonicalize(D) ()
-
-// comparison functions
-#define CGALRS_dyadic_sgn(D) mpfr_sgn(D)
-#define CGALRS_dyadic_zero(D) mpfr_zero_p(D)
-#define CGALRS_dyadic_cmp(D,E) mpfr_cmp(D,E)
-
-// arithmetic functions
-#define CGALRS_dyadic_add(R,D,E) CGALRS_dyadic_ll_add(R,D,E,0)
-#define CGALRS_dyadic_sub(R,D,E) CGALRS_dyadic_ll_sub(R,D,E,0)
-#define CGALRS_dyadic_mul(R,D,E) CGALRS_dyadic_ll_mul(R,D,E,0)
-
-inline void CGALRS_dyadic_neg(CGALRS_dyadic_ptr rop,CGALRS_dyadic_srcptr op){
- if(rop!=op)
- CGALRS_dyadic_set_prec(rop,mpfr_get_prec(op));
- mpfr_neg(rop,op,GMP_RNDN);
- CGAL_assertion(
- rop==op||
- (!mpfr_cmpabs(rop,op)&&
- ((CGALRS_dyadic_zero(op)&&CGALRS_dyadic_zero(rop))||
- (CGALRS_dyadic_sgn(op)!=CGALRS_dyadic_sgn(rop)))));
-}
-
-// low-level addition:
-// add op1 and op2 and reserve b bits for future lowlevel operations
-inline void CGALRS_dyadic_ll_add(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- CGALRS_dyadic_srcptr op2,
- mp_prec_t b){
- mp_exp_t l,r,temp1,temp2;
- mp_prec_t rop_prec;
- if(mpfr_zero_p(op1)){
- if(rop!=op2)
- CGALRS_dyadic_set(rop,op2);
- return;
- }
- if(mpfr_zero_p(op2)){
- if(rop!=op1)
- CGALRS_dyadic_set(rop,op1);
- return;
- }
- l=mpfr_get_exp(op1)>mpfr_get_exp(op2)?
- mpfr_get_exp(op1):
- mpfr_get_exp(op2);
- temp1=mpfr_get_exp(op1)-(mp_exp_t)mpfr_get_prec(op1);
- temp2=mpfr_get_exp(op2)-(mp_exp_t)mpfr_get_prec(op2);
- r=temp1>temp2?temp2:temp1;
- CGAL_assertion(l>r);
-
- rop_prec=b+1+(mp_prec_t)(l-r);
- CGAL_assertion(rop_prec>=mpfr_get_prec(op1)&&
- rop_prec>=mpfr_get_prec(op2));
- if(rop==op1||rop==op2)
- CGALRS_dyadic_prec_round(rop,rop_prec);
- else
- CGALRS_dyadic_set_prec(rop,rop_prec);
- CGAL_assertion_code(int round=)
- mpfr_add(rop,op1,op2,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_add_z(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- mpz_srcptr z){
- mp_exp_t l,r;
- mp_prec_t rop_prec;
- if(mpfr_zero_p(op1)){
- CGALRS_dyadic_set_z(rop,z);
- return;
- }
- if(!mpz_sgn(z)){
- if(rop!=op1)
- CGALRS_dyadic_set(rop,op1);
- return;
- }
- l=mpfr_get_exp(op1)>(mp_exp_t)mpz_sizeinbase(z,2)?
- mpfr_get_exp(op1):
- mpz_sizeinbase(z,2);
- r=mpfr_get_exp(op1)-(mp_exp_t)mpfr_get_prec(op1)<0?
- mpfr_get_exp(op1)-(mp_exp_t)mpfr_get_prec(op1):
- 0;
- CGAL_assertion(l>r);
-
- rop_prec=1+(mp_prec_t)(l-r);
- CGAL_assertion(rop_prec>=mpfr_get_prec(op1)&&
- rop_prec>=(mp_prec_t)mpz_sizeinbase(z,2));
- if(rop==op1)
- CGALRS_dyadic_prec_round(rop,rop_prec);
- else
- CGALRS_dyadic_set_prec(rop,rop_prec);
- CGAL_assertion_code(int round=)
- mpfr_add_z(rop,op1,z,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-// low-level subtraction:
-// subtract op2 to op1 and reserve b bits for future lowlevel operations
-inline void CGALRS_dyadic_ll_sub(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- CGALRS_dyadic_srcptr op2,
- mp_prec_t b){
- mp_exp_t l,r,temp1,temp2;
- mp_prec_t rop_prec;
- if(mpfr_zero_p(op1)){
- CGALRS_dyadic_neg(rop,op2);
- return;
- }
- if(mpfr_zero_p(op2)){
- if(rop!=op1)
- CGALRS_dyadic_set(rop,op1);
- return;
- }
- l=mpfr_get_exp(op1)>mpfr_get_exp(op2)?
- mpfr_get_exp(op1):
- mpfr_get_exp(op2);
- temp1=mpfr_get_exp(op1)-(mp_exp_t)mpfr_get_prec(op1);
- temp2=mpfr_get_exp(op2)-(mp_exp_t)mpfr_get_prec(op2);
- r=temp1>temp2?temp2:temp1;
- CGAL_assertion(l>r);
-
- rop_prec=b+1+(mp_prec_t)(l-r);
- CGAL_assertion(rop_prec>=mpfr_get_prec(op1)&&
- rop_prec>=mpfr_get_prec(op2));
- if(rop==op1||rop==op2)
- CGALRS_dyadic_prec_round(rop,rop_prec);
- else
- CGALRS_dyadic_set_prec(rop,rop_prec);
- CGAL_assertion_code(int round=)
- mpfr_sub(rop,op1,op2,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-// low-level multiplication:
-// multiply op1 and op2 and reserve b bits for future lowlevel operations
-inline void CGALRS_dyadic_ll_mul(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- CGALRS_dyadic_srcptr op2,
- mp_prec_t b){
- if(rop==op1||rop==op2)
- CGALRS_dyadic_prec_round(rop,b+mpfr_get_prec(op1)+mpfr_get_prec(op2));
- else
- CGALRS_dyadic_set_prec(rop,b+mpfr_get_prec(op1)+mpfr_get_prec(op2));
- CGAL_assertion_code(int round=)
- mpfr_mul(rop,op1,op2,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_mul_z(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- mpz_srcptr z){
- if(rop==op1)
- CGALRS_dyadic_prec_round(
- rop,
- mpfr_get_prec(op1)+mpz_sizeinbase(z,2));
- else
- CGALRS_dyadic_set_prec(
- rop,
- mpfr_get_prec(op1)+mpz_sizeinbase(z,2));
- CGAL_assertion_code(int round=)
- mpfr_mul_z(rop,op1,z,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_mul_si(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- long s){
- if(rop==op1)
- CGALRS_dyadic_prec_round(rop,mpfr_get_prec(op1)+sizeof(long));
- else
- CGALRS_dyadic_set_prec(rop,mpfr_get_prec(op1)+sizeof(long));
- CGAL_assertion_code(int round=)
- mpfr_mul_si(rop,op1,s,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_mul_ui(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- unsigned long u){
- if(rop==op1)
- CGALRS_dyadic_prec_round(
- rop,
- mpfr_get_prec(op1)+sizeof(unsigned long));
- else
- CGALRS_dyadic_set_prec(
- rop,
- mpfr_get_prec(op1)+sizeof(unsigned long));
- CGAL_assertion_code(int round=)
- mpfr_mul_ui(rop,op1,u,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_pow_ui(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- unsigned long u){
- if(!u){
- CGAL_assertion_msg(!mpfr_zero_p(op1),"0^0");
- CGALRS_dyadic_set_ui(rop,1);
- return;
- }
- if(u==1){
- if(rop!=op1)
- CGALRS_dyadic_set(rop,op1);
- return;
- }
- if(mpfr_zero_p(op1)){
- CGAL_assertion_msg(u!=0,"0^0");
- CGALRS_dyadic_set_ui(rop,0);
- return;
- }
- if(!mpfr_cmp_ui(op1,1)){
- if(rop!=op1)
- CGALRS_dyadic_set(rop,op1);
- return;
- }
- if(rop==op1)
- CGALRS_dyadic_prec_round(rop,u*mpfr_get_prec(op1));
- else
- CGALRS_dyadic_set_prec(rop,u*mpfr_get_prec(op1));
- CGAL_assertion_code(int round=)
- mpfr_pow_ui(rop,op1,u,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_addmul(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- CGALRS_dyadic_srcptr op2){
- CGALRS_dyadic_t temp;
- CGALRS_dyadic_init(temp);
- CGALRS_dyadic_mul(temp,op1,op2);
- CGALRS_dyadic_add(rop,rop,temp);
- CGALRS_dyadic_clear(temp);
-}
-
-inline void CGALRS_dyadic_addmul_si(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- long op2){
- CGALRS_dyadic_t temp;
- CGALRS_dyadic_init(temp);
- CGALRS_dyadic_mul_si(temp,op1,op2);
- CGALRS_dyadic_add(rop,rop,temp);
- CGALRS_dyadic_clear(temp);
-}
-
-inline void CGALRS_dyadic_addmul_ui(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- unsigned long u){
- //CGALRS_dyadic_t temp;
- //CGALRS_dyadic_init(temp);
- //CGALRS_dyadic_mul_ui(temp,op1,u);
- //CGALRS_dyadic_add(rop,rop,temp);
- //CGALRS_dyadic_clear(temp);
- CGALRS_dyadic_t temp;
- mp_exp_t l,r,temp1,temp2;
- mp_prec_t rop_prec;
- if(u==0||mpfr_zero_p(op1))
- return;
- if(u==1){
- CGALRS_dyadic_add(rop,rop,op1);
- return;
- }
- // TODO: if(op1==1)
- // calculate temp=op1*u
- mpfr_init2(temp,mpfr_get_prec(op1)+sizeof(unsigned int));
- CGAL_assertion_code(int round1=)
- mpfr_mul_ui(temp,op1,u,GMP_RNDN);
- CGAL_assertion(!round1);
- // calculate the precision needed for rop
- l=mpfr_get_exp(op1)>0?mpfr_get_exp(op1):0;
- temp1=mpfr_get_exp(op1)-(mp_exp_t)mpfr_get_prec(op1);
- temp2=sizeof(unsigned long);
- r=temp1>temp2?temp2:temp1;
- CGAL_assertion(l>r);
- rop_prec=sizeof(unsigned long)+1+(mp_prec_t)(l-r);
- CGAL_assertion(rop_prec>=mpfr_get_prec(op1)&&
- rop_prec>=mpfr_get_prec(rop));
- // set precision and add
- CGALRS_dyadic_prec_round(rop,rop_prec);
- CGAL_assertion_code(int round2=)
- mpfr_add(rop,rop,temp,GMP_RNDN);
- CGAL_assertion(!round2);
-}
-
-inline void CGALRS_dyadic_mul_2exp(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- unsigned long ui){
- // mpfr_mul_2ui does change the mantissa!
- if(rop==op1)
- CGALRS_dyadic_prec_round(
- rop,
- sizeof(unsigned long)+mpfr_get_prec(op1));
- else
- CGALRS_dyadic_set_prec(
- rop,
- sizeof(unsigned long)+mpfr_get_prec(op1));
- CGAL_assertion_code(int round=)
- mpfr_mul_2ui(rop,op1,ui,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-inline void CGALRS_dyadic_div_2exp(CGALRS_dyadic_ptr rop,
- CGALRS_dyadic_srcptr op1,
- unsigned long ui){
- // mpfr_div_2ui does not change the mantissa... am I sure?
- CGAL_assertion_code(int round=)
- mpfr_div_2ui(rop,op1,ui,GMP_RNDN);
- CGAL_assertion(!round);
-}
-
-// miscellaneous functions
-#define CGALRS_dyadic_midpoint(R,D,E) \
- ( CGALRS_dyadic_ll_add(R,D,E,1) , mpfr_div_2ui(R,R,1,GMP_RNDN) )
-#define CGALRS_dyadic_swap(D,E) mpfr_swap(D,E)
-
-// I/O functions
-#define CGALRS_dyadic_out_str(F,D) mpfr_out_str(F,10,0,D,GMP_RNDN)
-#ifdef __cplusplus
-inline std::ostream& operator<<(std::ostream &s,CGALRS_dyadic_srcptr op){
- mp_exp_t exponent;
- mpz_t mantissa;
- mpz_init(mantissa);
- exponent=mpfr_get_z_exp(mantissa,op);
- s<<"["<
-
-#ifndef CGAL_RS_EXACT_SIGNAT_1_H
-#define CGAL_RS_EXACT_SIGNAT_1_H
-
-#include
-#include "dyadic.h"
-
-namespace CGAL{
-namespace RS_AK1{
-
-template
-struct ExactSignat_1{
- typedef Polynomial_ Polynomial;
- typedef Bound_ Bound;
- typedef CGAL::Polynomial_traits_d PT;
- typedef typename PT::Degree Degree;
- Polynomial pol;
- ExactSignat_1(const Polynomial &p):pol(p){};
- CGAL::Sign operator()(const Bound&)const;
-}; // struct ExactSignat_1
-
-template <>
-inline CGAL::Sign
-ExactSignat_1,Gmpfr>::operator()(const Gmpfr &x)const{
- int d=Degree()(pol);
- if(d==0)
- return pol[0].sign();
- // Construct a Gmpfr containing exactly the leading coefficient.
- Gmpfr h(pol[d],pol[d].bit_size());
- CGAL_assertion(h==pol[d]);
- // Horner's evaluation.
- for(int i=1;i
-
-#ifndef CGAL_RS_FUNCTORS_1_H
-#define CGAL_RS_FUNCTORS_1_H
-
-#include
-#include
-
-namespace CGAL{
-namespace RS_AK1{
-
-template
-struct Construct_algebraic_real_1{
- typedef Polynomial_ Polynomial;
- typedef Algebraic_ Algebraic;
- typedef Bound_ Bound;
- typedef Coefficient_ Coefficient;
- typedef Isolator_ Isolator;
- typedef Algebraic result_type;
-
- template
- Algebraic operator()(const T &a)const{
- return Algebraic(a);
- }
-
- Algebraic operator()(const Polynomial &p,size_t i)const{
- Isolator isol(p);
- return Algebraic(p,isol.left_bound(i),isol.right_bound(i));
- }
-
- Algebraic operator()(const Polynomial &p,
- const Bound &l,
- const Bound &r)const{
- return Algebraic(p,l,r);
- }
-
-}; // struct Construct_algebraic_1
-
-template
-struct Compute_polynomial_1:
-public CGAL::cpp98::unary_function{
- typedef Polynomial_ Polynomial;
- typedef Algebraic_ Algebraic;
- Polynomial operator()(const Algebraic &x)const{
- return x.get_pol();
- }
-}; // struct Compute_polynomial_1
-
-template
-struct Is_coprime_1:
-public CGAL::cpp98::binary_function{
- typedef Polynomial_ Polynomial;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Degree Degree;
- inline bool operator()(const Polynomial &p1,const Polynomial &p2)const{
- return Degree()(Gcd()(p1,p2))==0;
- }
-}; // struct Is_coprime_1
-
-template
-struct Make_coprime_1{
- typedef Polynomial_ Polynomial;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Degree Degree;
- typedef typename Ptraits::Integral_division_up_to_constant_factor
- IDiv;
- bool operator()(const Polynomial &p1,
- const Polynomial &p2,
- Polynomial &g,
- Polynomial &q1,
- Polynomial &q2)const{
- g=Gcd()(p1,p2);
- q1=IDiv()(p1,g);
- q2=IDiv()(p2,g);
- return Degree()(Gcd()(p1,p2))==0;
- }
-}; // struct Make_coprime_1
-
-template
-struct Solve_1{
- typedef Polynomial_ Polynomial_1;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Isolator_ Isolator;
- typedef Signat_ Signat;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Square_free_factorize_up_to_constant_factor
- Sqfr;
- typedef typename Ptraits::Degree Degree;
- typedef typename Ptraits::Make_square_free Sfpart;
-
- template
- OutputIterator operator()(const Polynomial_1 &p,
- OutputIterator res)const{
- typedef std::pair polmult;
- typedef std::vector sqvec;
-
- Polynomial_1 sfp=Sfpart()(p);
- sqvec sfv;
- Sqfr()(p,std::back_inserter(sfv));
- Isolator isol(sfp);
- int *m=(int*)calloc(isol.number_of_real_roots(),sizeof(int));
- for(typename sqvec::iterator i=sfv.begin();i!=sfv.end();++i){
- int k=Degree()(i->first);
- Signat signof(i->first);
- for(int j=0;k&&jsecond;
- --k;
- }
- }
- }
- }
- for(int l=0;l
- OutputIterator operator()(const Polynomial_1 &p,
- bool,
- OutputIterator res)const{
- Isolator isol(p);
- for(int l=0;l
- OutputIterator operator()(const Polynomial_1 &p,
- const Bound &l,
- const Bound &u,
- OutputIterator res)const{
- typedef std::vector RV;
- typedef std::pair PM;
- typedef std::vector PMV;
- typedef typename PMV::iterator PMVI;
- CGAL_precondition_msg(l<=u,
- "left bound must be <= right bound");
- RV roots; // all roots of the polynomial
- this->operator()(p,false,std::back_inserter(roots));
- size_t nb_roots=roots.size();
- // indices of the first and last roots to be reported:
- size_t index_l=0,index_u;
- while(index_l and match the
- // roots in the interval [index_l,index_u)
- for(PMVI i=sfv.begin();i!=sfv.end();++i){
- int k=Degree()(i->first);
- Signat signof(i->first);
- for(size_t j=index_l;k&&jsecond;
- --k;
- }
- }
- }
- }
- for(size_t l=index_l;l
- OutputIterator operator()(const Polynomial_1 &p,
- bool known_to_be_square_free,
- const Bound &l,
- const Bound &u,
- OutputIterator res)const{
- typedef std::vector RV;
- typedef typename RV::iterator RVI;
- CGAL_precondition_msg(l<=u,
- "left bound must be <= right bound");
- RV roots;
- this->operator()(p,
- known_to_be_square_free,
- std::back_inserter(roots));
- for(RVI it=roots.begin();it!=roots.end();it++)
- if(*it>=l&&*it<=u)
- *res++=*it;
- return res;
- }
-
-}; // struct Solve_1
-
-template
-class Sign_at_1:
-public CGAL::cpp98::binary_function{
- // This implementation will work with any polynomial type whose
- // coefficient type is explicit interoperable with Gmpfi.
- // TODO: Make this function generic.
- public:
- typedef Polynomial_ Polynomial_1;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Refiner_ Refiner;
- typedef Signat_ Signat;
- typedef Ptraits_ Ptraits;
-
- private:
- CGAL::Uncertain eval_interv(const Polynomial_1 &p,
- const Bound &l,
- const Bound &r)const{
- typedef typename Ptraits::Substitute Subst;
- std::vector substitutions;
- substitutions.push_back(CGAL::Gmpfi(l,r));
- CGAL::Gmpfi eval=Subst()(p,
- substitutions.begin(),
- substitutions.end());
- return eval.sign();
- }
-
- // This function assumes that the sign of the evaluation is not zero,
- // it just refines x until having the correct sign.
- CGAL::Sign refine_and_return(const Polynomial_1 &p,Algebraic x)const{
- CGAL::Gmpfr xl(x.get_left());
- CGAL::Gmpfr xr(x.get_right());
- CGAL::Uncertain s;
- for(;;){
- Refiner()(x.get_pol(),
- xl,
- xr,
- 2*CGAL::max(xl.get_precision(),
- xr.get_precision()));
- s=eval_interv(p,xl,xr);
- if(!s.is_same(Uncertain::indeterminate())){
- x.set_left(xl);
- x.set_right(xr);
- return s;
- }
- }
- }
-
- public:
- CGAL::Sign operator()(const Polynomial_1 &p,Algebraic x)const{
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Make_square_free Sfpart;
- typedef typename Ptraits::Degree Degree;
- typedef typename Ptraits::Differentiate Deriv;
- CGAL::Uncertain unknown=
- Uncertain::indeterminate();
- CGAL::Uncertain s=eval_interv(p,
- x.get_left(),
- x.get_right());
- if(!s.is_same(unknown))
- return s;
- // We are not sure about the sign. We calculate the gcd in
- // order to know if both polynomials have common roots.
- Polynomial_1 sfpp=Sfpart()(p);
- Polynomial_1 gcd=Gcd()(sfpp,Sfpart()(x.get_pol()));
- if(Degree()(gcd)==0)
- return refine_and_return(p,x);
-
- // At this point, gcd is not 1; we proceed as follows:
- // -we refine x until having p monotonic in x's interval (to be
- // sure that p has at most one root on that interval),
- // -if the gcd has a root on this interval, both roots are
- // equal (we return 0), otherwise, we refine until having a
- // result.
-
- // How to assure that p is monotonic in an interval: when its
- // derivative is never zero in that interval.
- Polynomial_1 dsfpp=Deriv()(sfpp);
- CGAL::Gmpfr xl(x.get_left());
- CGAL::Gmpfr xr(x.get_right());
- while(eval_interv(dsfpp,xl,xr).is_same(unknown)){
- Refiner()(x.get_pol(),
- xl,
- xr,
- 2*CGAL::max(xl.get_precision(),
- xr.get_precision()));
- }
- x.set_left(xl);
- x.set_right(xr);
-
- // How to know that the gcd has a root: evaluate endpoints.
- CGAL::Sign sleft,sright;
- Signat sign_at_gcd(gcd);
- if((sleft=sign_at_gcd(x.get_left()))==CGAL::ZERO||
- (sright=sign_at_gcd(x.get_right()))==CGAL::ZERO||
- (sleft!=sright))
- return CGAL::ZERO;
- return refine_and_return(p,x);
- }
-}; // struct Sign_at_1
-
-template
-class Is_zero_at_1:
-public CGAL::cpp98::binary_function{
- // This implementation will work with any polynomial type whose
- // coefficient type is explicit interoperable with Gmpfi.
- // TODO: Make this function generic.
- public:
- typedef Polynomial_ Polynomial_1;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Refiner_ Refiner;
- typedef Signat_ Signat;
- typedef Ptraits_ Ptraits;
-
- private:
- CGAL::Uncertain eval_interv(const Polynomial_1 &p,
- const Bound &l,
- const Bound &r)const{
- typedef typename Ptraits::Substitute Subst;
- std::vector substitutions;
- substitutions.push_back(CGAL::Gmpfi(l,r));
- CGAL::Gmpfi eval=Subst()(p,
- substitutions.begin(),
- substitutions.end());
- return eval.sign();
- }
-
- public:
- bool operator()(const Polynomial_1 &p,Algebraic x)const{
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Make_square_free Sfpart;
- typedef typename Ptraits::Degree Degree;
- typedef typename Ptraits::Differentiate Deriv;
- CGAL::Uncertain unknown=
- Uncertain::indeterminate();
- CGAL::Uncertain s=eval_interv(p,
- x.get_left(),
- x.get_right());
- if(!s.is_same(unknown))
- return (s==CGAL::ZERO);
- // We are not sure about the sign. We calculate the gcd in
- // order to know if both polynomials have common roots.
- Polynomial_1 sfpp=Sfpart()(p);
- Polynomial_1 gcd=Gcd()(sfpp,Sfpart()(x.get_pol()));
- if(Degree()(gcd)==0)
- return false;
-
- // At this point, gcd is not 1; we proceed as follows:
- // -we refine x until having p monotonic in x's interval (to be
- // sure that p has at most one root on that interval),
- // -if the gcd has a root on this interval, both roots are
- // equal (we return 0), otherwise, we refine until having a
- // result.
-
- // How to assure that p is monotonic in an interval: when its
- // derivative is never zero in that interval.
- Polynomial_1 dsfpp=Deriv()(sfpp);
- CGAL::Gmpfr xl(x.get_left());
- CGAL::Gmpfr xr(x.get_right());
- while(eval_interv(dsfpp,xl,xr).is_same(unknown)){
- Refiner()(x.get_pol(),
- xl,
- xr,
- 2*CGAL::max(xl.get_precision(),
- xr.get_precision()));
- }
- x.set_left(xl);
- x.set_right(xr);
-
- // How to know that the gcd has a root: evaluate endpoints.
- CGAL::Sign sleft,sright;
- Signat sign_at_gcd(gcd);
- return((sleft=sign_at_gcd(x.get_left()))==CGAL::ZERO||
- (sright=sign_at_gcd(x.get_right()))==CGAL::ZERO||
- (sleft!=sright));
- }
-}; // class Is_zero_at_1
-
-// TODO: it says in the manual that this should return a size_type, but test
-// programs assume that this is equal to int
-template
-struct Number_of_solutions_1:
-public CGAL::cpp98::unary_function{
- typedef Polynomial_ Polynomial_1;
- typedef Isolator_ Isolator;
- size_t operator()(const Polynomial_1 &p)const{
- // TODO: make sure that p is square free (precondition)
- Isolator isol(p);
- return isol.number_of_real_roots();
- }
-}; // struct Number_of_solutions_1
-
-// This functor not only compares two algebraic numbers. In case they are
-// different, it refines them until they do not overlap.
-template
-struct Compare_1:
-public CGAL::cpp98::binary_function{
- typedef Algebraic_ Algebraic;
- typedef Bound_ Bound;
- typedef Comparator_ Comparator;
-
- CGAL::Comparison_result operator()(Algebraic a,Algebraic b)const{
- Bound al=a.get_left();
- Bound ar=a.get_right();
- Bound bl=b.get_left();
- Bound br=b.get_right();
- CGAL::Comparison_result c=Comparator()(a.get_pol(),al,ar,
- b.get_pol(),bl,br);
- a.set_left(al);
- a.set_right(ar);
- b.set_left(bl);
- b.set_right(br);
- return c;
- }
-
- CGAL::Comparison_result operator()(Algebraic a,const Bound &b)const{
- Bound al=a.get_left();
- Bound ar=a.get_right();
- Algebraic balg(b);
- CGAL::Comparison_result c=Comparator()(a.get_pol(),al,ar,
- balg.get_pol(),b,b);
- a.set_left(al);
- a.set_right(ar);
- return c;
- }
-
- template
- CGAL::Comparison_result operator()(Algebraic a,const T &b)const{
- Bound al=a.get_left();
- Bound ar=a.get_right();
- Algebraic balg(b);
- CGAL::Comparison_result c=Comparator()(a.get_pol(),
- al,
- ar,
- balg.get_pol(),
- balg.get_left(),
- balg.get_right());
- a.set_left(al);
- a.set_right(ar);
- return c;
- }
-
-}; // class Compare_1
-
-template
-struct Bound_between_1:
-public CGAL::cpp98::binary_function{
- typedef Algebraic_ Algebraic;
- typedef Bound_ Bound;
- typedef Comparator_ Comparator;
-
- Bound operator()(Algebraic a,Algebraic b)const{
- typedef Compare_1 Compare;
- typename Bound::Precision_type prec;
- switch(Compare()(a,b)){
- case CGAL::LARGER:
- CGAL_assertion(b.get_right()
-struct Isolate_1:
-public CGAL::cpp98::binary_function >{
- typedef Polynomial_ Polynomial_1;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Isolator_ Isolator;
- typedef Comparator_ Comparator;
- typedef Signat_ Signat;
- typedef Ptraits_ Ptraits;
-
- std::pair
- operator()(const Algebraic &a,const Polynomial_1 &p)const{
- std::vector roots;
- std::back_insert_iterator > rit(roots);
- typedef Solve_1 Solve;
- typedef Compare_1 Compare;
- Solve()(p,false,rit);
- for(typename std::vector::size_type i=0;
- i
-struct Approximate_absolute_1:
-public CGAL::cpp98::binary_function >{
- typedef Polynomial_ Polynomial_1;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Refiner_ Refiner;
-
- // This implementation assumes that Bound is Gmpfr.
- // TODO: make generic.
- std::pair operator()(const Algebraic &x,const int &a)const{
- Bound xl(x.get_left()),xr(x.get_right());
- // refsteps=log2(xl-xr)
- mpfr_t temp;
- mpfr_init(temp);
- mpfr_sub(temp,xr.fr(),xl.fr(),GMP_RNDU);
- mpfr_log2(temp,temp,GMP_RNDU);
- long refsteps=mpfr_get_si(temp,GMP_RNDU);
- mpfr_clear(temp);
- Refiner()(x.get_pol(),xl,xr,CGAL::abs(refsteps+a));
- x.set_left(xl);
- x.set_right(xr);
- CGAL_assertion(a>0?
- (xr-xl)*CGAL::ipower(Bound(2),a)<=Bound(1):
- (xr-xl)<=CGAL::ipower(Bound(2),-a));
- return std::make_pair(xl,xr);
- }
-}; // struct Approximate_absolute_1
-
-template
-struct Approximate_relative_1:
-public CGAL::cpp98::binary_function >{
- typedef Polynomial_ Polynomial_1;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Refiner_ Refiner;
-
- std::pair operator()(const Algebraic &x,const int &a)const{
- if(CGAL::is_zero(x))
- return std::make_pair(Bound(0),Bound(0));
- Bound error=CGAL::ipower(Bound(2),CGAL::abs(a));
- Bound xl(x.get_left()),xr(x.get_right());
- Bound max_b=(CGAL::max)(CGAL::abs(xr),CGAL::abs(xl));
- while(a>0?(xr-xl)*error>max_b:(xr-xl)>error*max_b){
- Refiner()(x.get_pol(),
- xl,
- xr,
- std::max(
- CGAL::abs(a),
- CGAL::max(xl.get_precision(),
- xr.get_precision())));
- max_b=(CGAL::max)(CGAL::abs(xr),CGAL::abs(xl));
- }
- x.set_left(xl);
- x.set_right(xr);
- CGAL_assertion(
- a>0?
- (xr-xl)*CGAL::ipower(Bound(2),a)<=max_b:
- (xr-xl)<=CGAL::ipower(Bound(2),-a)*max_b);
- return std::make_pair(xl,xr);
- }
-}; // struct Approximate_relative_1
-
-} // namespace RS_AK1
-} // namespace CGAL
-
-#endif // CGAL_RS_FUNCTORS_1_H
diff --git a/Algebraic_kernel_d/include/CGAL/RS/functors_z_1.h b/Algebraic_kernel_d/include/CGAL/RS/functors_z_1.h
deleted file mode 100644
index 56961c92e4f0..000000000000
--- a/Algebraic_kernel_d/include/CGAL/RS/functors_z_1.h
+++ /dev/null
@@ -1,699 +0,0 @@
-// Copyright (c) 2006-2013 INRIA Nancy-Grand Est (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: Luis Peñaranda
-
-#ifndef CGAL_RS_FUNCTORS_Z_1_H
-#define CGAL_RS_FUNCTORS_Z_1_H
-
-#include
-#include
-
-namespace CGAL{
-namespace RS_AK1{
-
-template
-struct Construct_algebraic_real_z_1{
- typedef Polynomial_ Polynomial;
- typedef ZPolynomial_ ZPolynomial;
- typedef PolConverter_ PolConverter;
- typedef Algebraic_ Algebraic;
- typedef Bound_ Bound;
- typedef Coefficient_ Coefficient;
- typedef Isolator_ Isolator;
- typedef Algebraic result_type;
-
- template
- Algebraic operator()(const T &a)const{
- return Algebraic(a);
- }
-
- Algebraic operator()(const Polynomial &p,size_t i)const{
- ZPolynomial zp=PolConverter()(p);
- Isolator isol(zp);
- return Algebraic(p,
- zp,
- isol.left_bound(i),
- isol.right_bound(i));
- }
-
- Algebraic operator()(const Polynomial &p,
- const Bound &l,
- const Bound &r)const{
- return Algebraic(p,PolConverter()(p),l,r);
- }
-
-}; // struct Construct_algebraic_real_z_1
-
-template
-struct Compute_polynomial_z_1:
-public CGAL::cpp98::unary_function{
- typedef Polynomial_ Polynomial;
- typedef Algebraic_ Algebraic;
- Polynomial operator()(const Algebraic &x)const{
- return x.get_pol();
- }
-}; // struct Compute_polynomial_z_1
-
-template
-struct Is_coprime_z_1:
-public CGAL::cpp98::binary_function{
- typedef Polynomial_ Polynomial;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Degree Degree;
- inline bool operator()(const Polynomial &p1,const Polynomial &p2)const{
- return Degree()(Gcd()(p1,p2))==0;
- }
-}; // struct Is_coprime_z_1
-
-template
-struct Make_coprime_z_1{
- typedef Polynomial_ Polynomial;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Degree Degree;
- typedef typename Ptraits::Integral_division_up_to_constant_factor
- IDiv;
- bool operator()(const Polynomial &p1,
- const Polynomial &p2,
- Polynomial &g,
- Polynomial &q1,
- Polynomial &q2)const{
- g=Gcd()(p1,p2);
- q1=IDiv()(p1,g);
- q2=IDiv()(p2,g);
- return Degree()(Gcd()(p1,p2))==0;
- }
-}; // struct Make_coprime_z_1
-
-template
-struct Solve_z_1{
- typedef Polynomial_ Polynomial_1;
- typedef ZPolynomial_ ZPolynomial_1;
- typedef PolConverter_ PolConverter;
- typedef Bound_ Bound;
- typedef Algebraic_ Algebraic;
- typedef Isolator_ Isolator;
- typedef Signat_ ZSignat;
- typedef Ptraits_ Ptraits;
- typedef typename Ptraits::Gcd_up_to_constant_factor Gcd;
- typedef typename Ptraits::Square_free_factorize_up_to_constant_factor
- Sqfr;
- typedef typename Ptraits::Degree Degree;
- typedef typename Ptraits::Make_square_free Sfpart;
- typedef ZPtraits_ ZPtraits;
- typedef typename ZPtraits::Gcd_up_to_constant_factor ZGcd;
- typedef typename ZPtraits::Square_free_factorize_up_to_constant_factor
- ZSqfr;
- typedef typename ZPtraits::Degree ZDegree;
- typedef typename ZPtraits::Make_square_free ZSfpart;
-
- template
- OutputIterator operator()(const Polynomial_1 &p,
- OutputIterator res)const{
- typedef std::pair zpolmult;
- typedef std::vector zsqvec;
-
- ZPolynomial_1 zp=PolConverter()(p);
- Polynomial_1 sfp=Sfpart()(p);
- ZPolynomial_1 zsfp=PolConverter()(sfp);
- zsqvec zsfv;
- ZSqfr()(zp,std::back_inserter(zsfv));
- Isolator isol(zsfp);
- int *m=(int*)calloc(isol.number_of_real_roots(),sizeof(int));
- for(typename zsqvec::iterator i=zsfv.begin();i!=zsfv.end();++i){
- int k=ZDegree()(i->first);
- ZSignat signof(i->first);
- for(int j=0;k&&jsecond;
- --k;
- }
- }
- }
- }
- for(int l=0;l
- OutputIterator operator()(const Polynomial_1 &p,
- bool,
- OutputIterator res)const{
- ZPolynomial_1 zp=PolConverter()(p);
- Isolator isol(zp);
- for(int l=0;l