mirror of https://github.com/CGAL/cgal
Removed files with outdated algorithms/data structures, as pre-decided with Michael Hemmer.
This commit is contained in:
Sebastian Limbach
parent
a907ead17b
commit
1793da5074
|
|
@ -1,125 +0,0 @@
|
|||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
#ifndef CGAL_CHINESE_REMAINDER_TRAITS_H
|
||||
#define CGAL_CHINESE_REMAINDER_TRAITS_H 1
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/chinese_remainder.h>
|
||||
#include <CGAL/Algebraic_structure_traits.h>
|
||||
#include <CGAL/Sqrt_extension.h>
|
||||
#include <CGAL/Polynomial.h>
|
||||
#include <vector>
|
||||
|
||||
namespace CGAL{
|
||||
|
||||
//TODO: 'm' is recomputed again and again in the current scheme.
|
||||
|
||||
template <class T> class Chinese_remainder_traits;
|
||||
template <class T, class TAG> class Chinese_remainder_traits_base;
|
||||
|
||||
|
||||
template <class T> class Chinese_remainder_traits
|
||||
:public Chinese_remainder_traits_base<T,
|
||||
typename Algebraic_structure_traits<T>::Algebraic_category>{};
|
||||
|
||||
template <class T_>
|
||||
struct Chinese_remainder_traits_base<T_,Euclidean_ring_tag>{
|
||||
typedef T_ T;
|
||||
typedef T_ Scalar_type;
|
||||
|
||||
struct Chinese_remainder{
|
||||
void operator() (
|
||||
Scalar_type m1, T u1,
|
||||
Scalar_type m2, T u2,
|
||||
Scalar_type& m, T& u){
|
||||
CGAL::chinese_remainder(m1,u1,m2,u2,m,u);
|
||||
}
|
||||
};
|
||||
};
|
||||
|
||||
template <class T_, class TAG>
|
||||
class Chinese_remainder_traits_base{
|
||||
typedef T_ T;
|
||||
typedef void Scalar_type;
|
||||
typedef Null_functor Chinese_remainder;
|
||||
};
|
||||
|
||||
|
||||
// Spec for Sqrt_extension
|
||||
// TODO mv to Sqrt_extension.h
|
||||
|
||||
template <class NT, class ROOT> class Sqrt_extension;
|
||||
template <class NT, class ROOT>
|
||||
struct Chinese_remainder_traits<Sqrt_extension<NT,ROOT> >{
|
||||
typedef Sqrt_extension<NT,ROOT> T;
|
||||
typedef Chinese_remainder_traits<NT> CRT_NT;
|
||||
typedef Chinese_remainder_traits<ROOT> CRT_ROOT;
|
||||
|
||||
// SAME AS CRT_ROOT::Scalar_type
|
||||
typedef typename CRT_NT::Scalar_type Scalar_type;
|
||||
|
||||
struct Chinese_remainder{
|
||||
void operator() (
|
||||
Scalar_type m1, T u1,
|
||||
Scalar_type m2, T u2,
|
||||
Scalar_type& m, T& u){
|
||||
|
||||
NT a0,a1;
|
||||
ROOT root;
|
||||
|
||||
typename CRT_NT::Chinese_remainder chinese_remainder_nt;
|
||||
chinese_remainder_nt(m1,u1.a0(),m2,u2.a0(),m,a0);
|
||||
if(u1.is_extended() || u2.is_extended()){
|
||||
chinese_remainder_nt(m1,u1.a1(),m2,u2.a1(),m,a1);
|
||||
typename CRT_ROOT::Chinese_remainder chinese_remainder_root;
|
||||
chinese_remainder_root(m1,u1.root(),m2,u2.root(),m,root);
|
||||
u=T(a0,a1,root);
|
||||
}else{
|
||||
u=T(a0);
|
||||
}
|
||||
}
|
||||
};
|
||||
|
||||
};
|
||||
|
||||
|
||||
|
||||
// Spec for Polynomial
|
||||
// TODO mv to Polynomial.h
|
||||
|
||||
template <class NT> class Polynomial;
|
||||
template <class NT>
|
||||
struct Chinese_remainder_traits<Polynomial<NT> >{
|
||||
typedef Polynomial<NT> T;
|
||||
typedef Chinese_remainder_traits<NT> CRT_NT;
|
||||
|
||||
typedef typename CRT_NT::Scalar_type Scalar_type;
|
||||
|
||||
struct Chinese_remainder{
|
||||
void operator() (
|
||||
Scalar_type m1, T u1,
|
||||
Scalar_type m2, T u2,
|
||||
Scalar_type& m, T& u){
|
||||
|
||||
typename CRT_NT::Chinese_remainder chinese_remainder_nt;
|
||||
|
||||
CGAL_precondition(u1.degree() == u2.degree());
|
||||
|
||||
std::vector<NT> coeffs;
|
||||
coeffs.reserve(u1.degree()+1);
|
||||
for(int i = 0; i <= u1.degree(); i++){
|
||||
NT c;
|
||||
chinese_remainder_nt(m1,u1[i],m2,u2[i],m,c);
|
||||
coeffs.push_back(c);
|
||||
}
|
||||
u = Polynomial<NT>(coeffs.begin(),coeffs.end());
|
||||
}
|
||||
};
|
||||
};
|
||||
|
||||
|
||||
|
||||
} // namespace CGAL
|
||||
|
||||
#endif // CGAL_CHINESE_REMAINDER_TRAITS_H //
|
||||
|
||||
|
|
@ -1,56 +0,0 @@
|
|||
//Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
#ifndef CGAL_CHINESE_REMAINDER_H
|
||||
#define CGAL_CHINESE_REMAINDER_H 1
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/extended_euclidean_algorithm.h>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
// this is just the version for 'integers'
|
||||
// NT must be model of RealEmbeddable
|
||||
// NT must be model of EuclideanRing
|
||||
template <class NT>
|
||||
void chinese_remainder(
|
||||
NT m1, NT u1,
|
||||
NT m2, NT u2,
|
||||
NT& m , NT& u ){
|
||||
typedef Algebraic_structure_traits<NT> AST;
|
||||
typename AST::Mod mod;
|
||||
//typename AST::Unit_part unit_part;
|
||||
typename AST::Integral_division idiv;
|
||||
|
||||
if(u1 < NT(0) ) u1 += m1;
|
||||
if(u2 < NT(0) ) u2 += m2;
|
||||
|
||||
CGAL_precondition(0 < m1);
|
||||
CGAL_precondition(u1 < m1);
|
||||
CGAL_precondition(u1 >= NT(0));
|
||||
|
||||
CGAL_precondition(0 < m2);
|
||||
CGAL_precondition(u2 < m2);
|
||||
CGAL_precondition(u2 >= NT(0));
|
||||
|
||||
|
||||
NT tmp,c,dummy;
|
||||
tmp = CGAL::extended_euclidean_algorithm(m1,m2,c,dummy);
|
||||
CGAL_postcondition(tmp == NT(1));
|
||||
CGAL_postcondition(m1*c + m2*dummy == NT(1));
|
||||
|
||||
|
||||
m = m1*m2;
|
||||
NT v = mod(c*(u2-u1),m2);
|
||||
u = m1*v + u1;
|
||||
|
||||
// u is not unique yet!
|
||||
NT m_half = idiv(m-mod(m,NT(2)),NT(2));
|
||||
if (u > m_half) u -= m ;
|
||||
if (u <= -m_half) u += m ;
|
||||
|
||||
}
|
||||
|
||||
}///namespace CGAL
|
||||
|
||||
#endif //#ifnedef CGAL_CHINESE_REMAINDER_H 1
|
||||
|
||||
|
|
@ -1,104 +0,0 @@
|
|||
|
||||
// Author(s) : Lutz Kettner <kettner@mpi-inf.mpg.de>
|
||||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
|
||||
/*! \file CGAL/euclidean_algorithm.h
|
||||
\brief Defines funciton related to euclids algorithm.
|
||||
*/
|
||||
|
||||
#ifndef CGAL_EUCLIDEAN_ALGORITHM_H
|
||||
#define CGAL_EUCLIDEAN_ALGORITHM_H 1
|
||||
|
||||
// This forward declaration is required to resolve the circular dependency
|
||||
// between euclidean_algorithm and the partial specializations of NT_Traits
|
||||
// for the built-in number types.
|
||||
|
||||
namespace CGAL {
|
||||
template <class NT>
|
||||
NT euclidean_algorithm( const NT& a, const NT& b);
|
||||
}
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Algebraic_structure_traits.h>
|
||||
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
// We have a circular header file inclusion dependency with
|
||||
// CGAL/Algebraic_structure_traits.h.
|
||||
// As a consequence, we might not get the declaration for
|
||||
// CGAL::Algebraic_structure_traits
|
||||
// although we include the header file above. We repeat the declaration
|
||||
// here. We still include the header file to hide this dependency from users
|
||||
// such that they get the full Algebraic_structure_traits declaration after
|
||||
// including euclid_algorithm.h.
|
||||
|
||||
template <class NT> class Algebraic_structure_traits;
|
||||
|
||||
|
||||
/*! \brief generic Euclids algorithm, returns the unit
|
||||
normal greatest common devisor (gcd) of \a a and \a b.
|
||||
|
||||
Requires the number type \c NT to be a model of the concepts
|
||||
\c EuclideanRing, however, it uses only
|
||||
the functors \c Mod and \c Unit_part from the \c
|
||||
Algebraic_structure_traits<NT>, and the equality comparison operator.
|
||||
The implementation uses loop unrolling to avoid swapping the local
|
||||
variables all the time.
|
||||
*/
|
||||
template <class NT>
|
||||
NT euclidean_algorithm( const NT& a, const NT& b) {
|
||||
typedef Algebraic_structure_traits<NT> AST;
|
||||
typename AST::Mod mod;
|
||||
typename AST::Unit_part unit_part;
|
||||
typename AST::Integral_division idiv;
|
||||
// First: the extreme cases and negative sign corrections.
|
||||
if (a == NT(0)) {
|
||||
if (b == NT(0))
|
||||
return NT(0);
|
||||
return idiv(b,unit_part(b));
|
||||
}
|
||||
if (b == NT(0))
|
||||
return idiv(a,unit_part(a));
|
||||
NT u = idiv(a,unit_part(a));
|
||||
NT v = idiv(b,unit_part(b));
|
||||
// Second: assuming mod is the most expensive op here, we don't compute it
|
||||
// unnecessarily if u < v
|
||||
if (u < v) {
|
||||
v = mod(v,u);
|
||||
// maintain invariant of v > 0 for the loop below
|
||||
if ( v == 0)
|
||||
return idiv(u,unit_part(u));
|
||||
}
|
||||
// Third: generic case of two positive integer values and u >= v.
|
||||
// The standard loop would be:
|
||||
// while ( v != 0) {
|
||||
// int tmp = mod(u,v);
|
||||
// u = v;
|
||||
// v = tmp;
|
||||
// }
|
||||
// return u;
|
||||
//
|
||||
// But we want to save us all the variable assignments and unroll
|
||||
// the loop. Before that, we transform it into a do {...} while()
|
||||
// loop to reduce branching statements.
|
||||
NT w;
|
||||
do {
|
||||
w = mod(u,v);
|
||||
if ( w == 0)
|
||||
return idiv(v,unit_part(v));
|
||||
u = mod(v,w);
|
||||
if ( u == 0)
|
||||
return idiv(w,unit_part(w));
|
||||
v = mod(w,u);
|
||||
} while (v != 0);
|
||||
return idiv(u,unit_part(u));;
|
||||
}
|
||||
|
||||
// TODO: do we need a variant for unit normal inputs?
|
||||
|
||||
} // namespace CGAL
|
||||
|
||||
#endif // CGAL_EUCLIDEAN_ALGORITHM_H //
|
||||
// EOF
|
||||
|
|
@ -1,59 +0,0 @@
|
|||
//Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
#ifndef CGAL_EXTENDED_EUCLIDEAN_ALGORITHM_H
|
||||
#define CGAL_EXTENDED_EUCLIDEAN_ALGORITHM_H 1
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
template< class NT >
|
||||
NT extended_euclidean_algorithm(const NT& a_, const NT& b_, NT& u, NT& v){
|
||||
|
||||
typedef Algebraic_structure_traits<NT> AST;
|
||||
typename AST::Div_mod div_mod;
|
||||
typename AST::Unit_part unit_part;
|
||||
typename AST::Integral_division idiv;
|
||||
|
||||
NT unit_part_a(unit_part(a_));
|
||||
NT unit_part_b(unit_part(b_));
|
||||
|
||||
NT a(idiv(a_,unit_part_a));
|
||||
NT b(idiv(b_,unit_part_b));
|
||||
|
||||
|
||||
NT x(0),y(1),last_x(1),last_y(0);
|
||||
NT temp, quotient;
|
||||
// typename AST::Div div;
|
||||
// typename AST::Mod mod;
|
||||
|
||||
//TODO: unroll to avoid swapping
|
||||
while (b != 0){
|
||||
temp = b;
|
||||
div_mod(a,b,quotient,b);
|
||||
a = temp;
|
||||
temp = x;
|
||||
x = last_x-quotient*x;
|
||||
last_x = temp;
|
||||
|
||||
temp = y;
|
||||
y = last_y-quotient*y;
|
||||
last_y = temp;
|
||||
}
|
||||
u = last_x * unit_part_a;
|
||||
v = last_y * unit_part_b;
|
||||
|
||||
// std::cout <<"a_: "<<a_ <<std::endl;
|
||||
// std::cout <<"b_: "<<b_ <<std::endl;
|
||||
// std::cout <<"gcd: "<<a <<std::endl;
|
||||
// std::cout <<"u: "<<u <<std::endl;
|
||||
// std::cout <<"v: "<<v <<std::endl;
|
||||
// std::cout <<std::endl;
|
||||
CGAL_precondition(unit_part(a) == NT(1));
|
||||
CGAL_precondition(a == a_*u + b_*v);
|
||||
return a;
|
||||
}
|
||||
|
||||
} // namespace CGAL
|
||||
|
||||
#endif // CGAL_EXTENDED_EUCLIDEAN_ALGORITHM_H //
|
||||
|
|
@ -1,427 +0,0 @@
|
|||
//Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
/*! \file CGAL/modular_gcd.h
|
||||
provides gcd for Polynomials, based on Modular arithmetic.
|
||||
*/
|
||||
|
||||
|
||||
#ifndef CGAL_MODULAR_GCD_H
|
||||
#define CGAL_MODULAR_GCD_H 1
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Modular_traits.h>
|
||||
#include <CGAL/Polynomial.h>
|
||||
#include <CGAL/Polynomial_traits_d.h>
|
||||
#include <CGAL/Scalar_factor_traits.h>
|
||||
#include <CGAL/Chinese_remainder_traits.h>
|
||||
|
||||
//#include <CGAL/Polynomial_traits_d_d.h>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
template <class NT>
|
||||
typename Scalar_factor_traits<NT>::Scalar
|
||||
scalar_factor(const NT& x){
|
||||
typename Scalar_factor_traits<NT>::Scalar_factor scalar_factor;
|
||||
return scalar_factor(x);
|
||||
}
|
||||
template <class NT>
|
||||
typename Scalar_factor_traits<NT>::Scalar
|
||||
scalar_factor(const NT& x,const typename Scalar_factor_traits<NT>::Scalar& d){
|
||||
typename Scalar_factor_traits<NT>::Scalar_factor scalar_factor;
|
||||
return scalar_factor(x,d);
|
||||
}
|
||||
|
||||
template <class NT>
|
||||
typename Modular_traits<NT>::Modular_NT
|
||||
modular_image(const NT& x){
|
||||
typename Modular_traits<NT>::Modular_image modular_image;
|
||||
return modular_image(x);
|
||||
}
|
||||
|
||||
template <int> class MY_INT_TAG{};
|
||||
|
||||
template <class T>
|
||||
bool operator < (const std::vector<T>& a, const std::vector<T>& b){
|
||||
for(unsigned int i = 0; i < a.size(); i++){
|
||||
if (a[i] < b[i]) return true;
|
||||
}
|
||||
return false;
|
||||
}
|
||||
|
||||
template <class T>
|
||||
std::vector<T> min(const std::vector<T>& a, const std::vector<T>& b){
|
||||
return (a < b)?a:b;
|
||||
}
|
||||
|
||||
//ALGORITHM P (TODO)
|
||||
template <class Coeff, class TAG >
|
||||
Polynomial<Coeff> algorithm_x(
|
||||
const Polynomial <Coeff>& p1, const Polynomial <Coeff>& p2, TAG){
|
||||
CGAL_precondition(Polynomial_traits_d< Polynomial<Coeff> >::d > 1);
|
||||
|
||||
typedef Polynomial<Coeff> Poly;
|
||||
typedef Polynomial_traits_d<Poly> PT;
|
||||
typedef typename PT::Innermost_coefficient IC;
|
||||
|
||||
const int num_of_vars = PT::d;
|
||||
typename PT::Innermost_leading_coefficient ilcoeff;
|
||||
typename PT::Degree_vector degree_vector;
|
||||
|
||||
// will play the role of content
|
||||
typedef typename Scalar_factor_traits<Poly>::Scalar Scalar;
|
||||
|
||||
typedef typename Modular_traits<Poly>::Modular_NT MPoly;
|
||||
typename Polynomial_traits_d<MPoly>::Degree_vector mdegree_vector;
|
||||
typedef typename Modular_traits<Scalar>::Modular_NT MScalar;
|
||||
|
||||
typedef Chinese_remainder_traits<Poly> CRT;
|
||||
typename CRT::Chinese_remainder chinese_remainder;
|
||||
|
||||
typename Polynomial_traits_d<Poly>::Canonicalize canonicalize;
|
||||
Poly F1 = canonicalize(p1);
|
||||
Poly F2 = canonicalize(p2);
|
||||
|
||||
//std::cout <<" F1 : " << F1 <<std::endl;
|
||||
//std::cout <<" F2 : " << F2 <<std::endl;
|
||||
{
|
||||
// this part is needed for algebraic extensions e.g. Sqrt_extesnion
|
||||
// We have to ensure that G,H1,H2 can be expressed in terms of algebraic integers
|
||||
// Therefore we multiply F1 and F2 by denominatior for algebraic integer.
|
||||
//typename PT::Innermost_coefficient_to_polynomial ictp;
|
||||
typename PT::Innermost_coefficient_begin begin;
|
||||
typename PT::Innermost_coefficient_end end;
|
||||
typename Algebraic_extension_traits<IC>::Denominator_for_algebraic_integers dfai;
|
||||
typename Algebraic_extension_traits<IC>::Normalization_factor nfac;
|
||||
|
||||
// in case IC is an algebriac extension it may happen, that
|
||||
// Fx=G*Hx is not possible if the coefficients are algebraic integers
|
||||
Poly tmp = F1+F2;
|
||||
|
||||
IC denom = dfai(begin(tmp),end(tmp)); // TODO use this
|
||||
//IC denom = dfai(tmp.begin(),tmp.end());
|
||||
denom *= nfac(denom);
|
||||
tmp = Poly(denom);
|
||||
F1 *=tmp;
|
||||
F2 *=tmp;
|
||||
}
|
||||
|
||||
//std::cout <<" F1*denom*nafc: " << F1 <<std::endl;
|
||||
//std::cout <<" F2*denom*nfac: " << F2 <<std::endl;
|
||||
|
||||
Scalar f1 = CGAL::scalar_factor(ilcoeff(F1)); // ilcoeff(F1)
|
||||
Scalar f2 = CGAL::scalar_factor(ilcoeff(F2)); // ilcoeff(F2)
|
||||
Scalar g_ = CGAL::scalar_factor(f1,f2);
|
||||
|
||||
Poly F1_ = F1*Poly(g_);
|
||||
Poly F2_ = F2*Poly(g_);
|
||||
|
||||
//std::cout <<" g_ : "<< g_ << std::endl;
|
||||
//std::cout <<" F1*denom*nafc*g_: " << F1_ <<std::endl;
|
||||
//std::cout <<" F2*denom*nfac*g_: " << F2_ <<std::endl;
|
||||
|
||||
|
||||
bool solved = false;
|
||||
int prime_index = -1;
|
||||
int n = 0; // number of lucky primes
|
||||
std::vector<int> dv_F1 = degree_vector(F1);
|
||||
std::vector<int> dv_F2 = degree_vector(F2);
|
||||
std::vector<int> dv_e = min(dv_F1,dv_F2);;
|
||||
|
||||
MScalar mg_;
|
||||
MPoly mF1,mF2,mG_,mH1,mH2;
|
||||
|
||||
typename CRT::Scalar_type p,q,pq;
|
||||
Poly Gs,H1s,H2s; // s =^ star
|
||||
while(!solved){
|
||||
do{
|
||||
//---------------------------------------
|
||||
//choose prime not deviding f1 or f2
|
||||
do{
|
||||
prime_index++;
|
||||
CGAL_precondition(0<= prime_index && prime_index < 64);
|
||||
int current_prime = primes[prime_index];
|
||||
Modular::set_current_prime(current_prime);
|
||||
}
|
||||
while(!(( modular_image(f1) != 0 ) && ( modular_image(f2) != 0 )));
|
||||
// --------------------------------------
|
||||
// invoke gcd for current prime
|
||||
mg_ = CGAL::modular_image(g_);
|
||||
mF1 = CGAL::modular_image(F1_);
|
||||
mF2 = CGAL::modular_image(F2_);
|
||||
// replace mG_ = gcd (mF1,mF2)*MPoly(mg_); for multivariat
|
||||
mG_ = algorithm_x(mF1,mF2,MY_INT_TAG<num_of_vars>())*MPoly(mg_);
|
||||
mH1 = CGAL::integral_division(mF1,mG_);
|
||||
mH2 = CGAL::integral_division(mF2,mG_);
|
||||
//---------------------------------------
|
||||
// return if G is constant
|
||||
if (mG_ == MPoly(1)) return Poly(1);
|
||||
// --------------------------------------
|
||||
}// repeat until mG_ degree is less equal the known bound
|
||||
// check prime
|
||||
while( mdegree_vector(mG_) > dv_e);
|
||||
|
||||
if(mdegree_vector(mG_) < dv_e ){
|
||||
// restart chinese remainder
|
||||
// ignore previous unlucky primes
|
||||
n=1;
|
||||
dv_e= mdegree_vector(mG_);
|
||||
}else{
|
||||
CGAL_postcondition( mdegree_vector(mG_)== dv_e);
|
||||
n++; // increase number of lucky primes
|
||||
}
|
||||
|
||||
// --------------------------------------
|
||||
// try chinese remainder
|
||||
|
||||
//std::cout <<" chinese remainder round :" << n << std::endl;
|
||||
typename Modular_traits<Poly>::Modular_image_inv inv_map;
|
||||
if(n == 1){
|
||||
// init chinese remainder
|
||||
q = Modular::get_current_prime(); // implicit !
|
||||
Gs = inv_map(mG_);
|
||||
H1s = inv_map(mH1);
|
||||
H2s = inv_map(mH2);
|
||||
}else{
|
||||
// continue chinese remainder
|
||||
|
||||
int p = Modular::get_current_prime(); // implicit!
|
||||
//std::cout <<" p: "<< p<<std::endl;
|
||||
//std::cout <<" q: "<< q<<std::endl;
|
||||
//std::cout <<" gcd(p,q): "<< gcd(p,q)<<std::endl;
|
||||
chinese_remainder(q,Gs ,p,inv_map(mG_),pq,Gs);
|
||||
chinese_remainder(q,H1s,p,inv_map(mH1),pq,H1s);
|
||||
chinese_remainder(q,H2s,p,inv_map(mH2),pq,H2s);
|
||||
q=pq;
|
||||
}
|
||||
|
||||
// std::cout << "Gs: "<< Gs << std::endl;
|
||||
// std::cout << "H1s: "<< H1s << std::endl;
|
||||
// std::cout << "H2s: "<< H2s << std::endl;
|
||||
// std::cout <<std::endl;
|
||||
// std::cout << "F1s: "<<Gs*H1s<< std::endl;
|
||||
// std::cout << "F1 : "<<F1_<< std::endl;
|
||||
// std::cout << "diff : "<<F1_-Gs*H1s<< std::endl;
|
||||
// std::cout <<std::endl;
|
||||
// std::cout << "F2s: "<<Gs*H2s<< std::endl;
|
||||
// std::cout << "F2 : "<<F2_<< std::endl;
|
||||
// std::cout << "diff : "<<F2_-Gs*H2s<< std::endl;
|
||||
|
||||
try{// This is a HACK!!!!
|
||||
// TODO: in case of Sqrt_extension it may happen that the disr (root)
|
||||
// is not correct, in this case the behavior of the code is unclear
|
||||
// if CGAL is in debug mode it throws an error
|
||||
if( Gs*H1s == F1_ && Gs*H2s == F2_ ){
|
||||
solved = true;
|
||||
}
|
||||
}catch(...){}
|
||||
|
||||
//std::cout << "Gs: " << CGAL::canonicalize_polynomial(Gs)<<std::endl;
|
||||
// std::cout << "canonical(Gs): " << CGAL::canonicalize_polynomial(Gs)<<std::endl;
|
||||
|
||||
//std::cout << std::endl;
|
||||
|
||||
|
||||
}
|
||||
|
||||
//std::cout << "G: " << CGAL::canonicalize_polynomial(gcd_utcf(F1,F2)) << std::endl;
|
||||
|
||||
return canonicalize(Gs);
|
||||
|
||||
}
|
||||
|
||||
// ALGORITHM U (done)
|
||||
template <class Field>
|
||||
Polynomial<Field> algorithm_x(
|
||||
const Polynomial <Field>& p1, const Polynomial <Field>& p2, MY_INT_TAG<1> ){
|
||||
typedef Polynomial<Field> Poly;
|
||||
BOOST_STATIC_ASSERT(Polynomial_traits_d<Poly>::d == 1);
|
||||
typedef Algebraic_structure_traits<Field> AST;
|
||||
typedef typename AST::Algebraic_category TAG;
|
||||
BOOST_STATIC_ASSERT((boost::is_same<TAG, Field_tag>::value));
|
||||
return gcd(p1,p2);
|
||||
}
|
||||
|
||||
// TODO: ALGORITHM M
|
||||
template <class NT>
|
||||
Polynomial<NT> modular_gcd_utcf(
|
||||
const Polynomial<NT>& FF1 ,
|
||||
const Polynomial<NT>& FF2 ){
|
||||
CGAL_precondition(Polynomial_traits_d<Polynomial<NT> >::d == 1);
|
||||
|
||||
typedef Polynomial<NT> Poly;
|
||||
typedef Polynomial_traits_d<Poly> PT;
|
||||
|
||||
const int num_of_vars = PT::d;
|
||||
typedef typename PT::Innermost_coefficient IC;
|
||||
typename PT::Innermost_leading_coefficient ilcoeff;
|
||||
typename PT::Degree_vector degree_vector;
|
||||
|
||||
// will paly the role of content
|
||||
typedef typename Scalar_factor_traits<Poly>::Scalar Scalar;
|
||||
|
||||
typedef typename Modular_traits<Poly>::Modular_NT MPoly;
|
||||
typename Polynomial_traits_d<MPoly>::Degree_vector mdegree_vector;
|
||||
typedef typename Modular_traits<Scalar>::Modular_NT MScalar;
|
||||
|
||||
typedef Chinese_remainder_traits<Poly> CRT;
|
||||
typename CRT::Chinese_remainder chinese_remainder;
|
||||
|
||||
typename Polynomial_traits_d<Poly>::Canonicalize canonicalize;
|
||||
Poly F1 = canonicalize(FF1);
|
||||
Poly F2 = canonicalize(FF2);
|
||||
|
||||
//std::cout <<" F1 : " << F1 <<std::endl;
|
||||
//std::cout <<" F2 : " << F2 <<std::endl;
|
||||
{
|
||||
// this part is needed for algebraic extensions e.g. Sqrt_extesnion
|
||||
// We have to ensure that G,H1,H2 can be expressed in terms of algebraic integers
|
||||
// Therefore we multiply F1 and F2 by denominatior for algebraic integer.
|
||||
//typedef Polynomial<NT> POLY;
|
||||
//typename Polynomial_traits_d<POLY>::Innermost_coefficient_to_polynomial ictp;
|
||||
//typename Polynomial_traits_d<POLY>::Innermost_coefficient_begin begin;
|
||||
//typename Polynomial_traits_d<POLY>::Innermost_coefficient_end end;
|
||||
typename Algebraic_extension_traits<IC>::Denominator_for_algebraic_integers dfai;
|
||||
typename Algebraic_extension_traits<IC>::Normalization_factor nfac;
|
||||
|
||||
// in case IC is an algebriac extension it may happen, that
|
||||
// Fx=G*Hx is not possible if the coefficients are algebraic integers
|
||||
Poly tmp = F1+F2;
|
||||
|
||||
//IC denom = dfai(begin(tmp),end(tmp)); // TODO use this
|
||||
IC denom = dfai(tmp.begin(),tmp.end());
|
||||
denom *= nfac(denom);
|
||||
tmp = Poly(denom);
|
||||
F1 *=tmp;
|
||||
F2 *=tmp;
|
||||
}
|
||||
|
||||
//std::cout <<" F1*denom*nafc: " << F1 <<std::endl;
|
||||
//std::cout <<" F2*denom*nfac: " << F2 <<std::endl;
|
||||
|
||||
Scalar f1 = CGAL::scalar_factor(ilcoeff(F1)); // ilcoeff(F1)
|
||||
Scalar f2 = CGAL::scalar_factor(ilcoeff(F2)); // ilcoeff(F2)
|
||||
Scalar g_ = CGAL::scalar_factor(f1,f2);
|
||||
|
||||
Poly F1_ = F1*Poly(g_);
|
||||
Poly F2_ = F2*Poly(g_);
|
||||
|
||||
//std::cout <<" g_ : "<< g_ << std::endl;
|
||||
//std::cout <<" F1*denom*nafc*g_: " << F1_ <<std::endl;
|
||||
//std::cout <<" F2*denom*nfac*g_: " << F2_ <<std::endl;
|
||||
|
||||
|
||||
bool solved = false;
|
||||
int prime_index = -1;
|
||||
int n = 0; // number of lucky primes
|
||||
std::vector<int> dv_F1 = degree_vector(F1);
|
||||
std::vector<int> dv_F2 = degree_vector(F1);
|
||||
std::vector<int> dv_e = min(dv_F1,dv_F2);;
|
||||
|
||||
MScalar mg_;
|
||||
MPoly mF1,mF2,mG_,mH1,mH2;
|
||||
|
||||
typename CRT::Scalar_type p,q,pq;
|
||||
Poly Gs,H1s,H2s; // s =^ star
|
||||
while(!solved){
|
||||
do{
|
||||
//---------------------------------------
|
||||
//choose prime not deviding f1 or f2
|
||||
do{
|
||||
prime_index++;
|
||||
CGAL_precondition(0<= prime_index && prime_index < 64);
|
||||
int current_prime = primes[prime_index];
|
||||
Modular::set_current_prime(current_prime);
|
||||
}
|
||||
while(!(( modular_image(f1) != 0 ) && ( modular_image(f2) != 0 )));
|
||||
// --------------------------------------
|
||||
// invoke gcd for current prime
|
||||
mg_ = CGAL::modular_image(g_);
|
||||
mF1 = CGAL::modular_image(F1_);
|
||||
mF2 = CGAL::modular_image(F2_);
|
||||
// replace mG_ = gcd (mF1,mF2)*MPoly(mg_); for multivariat
|
||||
mG_ = algorithm_x(mF1,mF2,MY_INT_TAG<num_of_vars>())*MPoly(mg_);
|
||||
mH1 = CGAL::integral_division(mF1,mG_);
|
||||
mH2 = CGAL::integral_division(mF2,mG_);
|
||||
//---------------------------------------
|
||||
// return if G is constant
|
||||
if (mG_ == MPoly(1)) return Poly(1);
|
||||
// --------------------------------------
|
||||
}// repeat until mG_ degree is less equal the known bound
|
||||
// check prime
|
||||
while( mdegree_vector(mG_) > dv_e);
|
||||
|
||||
if(mdegree_vector(mG_) < dv_e ){
|
||||
// restart chinese remainder
|
||||
// ignore previous unlucky primes
|
||||
n=1;
|
||||
dv_e= mdegree_vector(mG_);
|
||||
}else{
|
||||
CGAL_postcondition( mdegree_vector(mG_)== dv_e);
|
||||
n++; // increase number of lucky primes
|
||||
}
|
||||
|
||||
// --------------------------------------
|
||||
// try chinese remainder
|
||||
|
||||
//std::cout <<" chinese remainder round :" << n << std::endl;
|
||||
typename Modular_traits<Poly>::Modular_image_inv inv_map;
|
||||
if(n == 1){
|
||||
// init chinese remainder
|
||||
q = Modular::get_current_prime(); // implicit !
|
||||
Gs = inv_map(mG_);
|
||||
H1s = inv_map(mH1);
|
||||
H2s = inv_map(mH2);
|
||||
}else{
|
||||
// continue chinese remainder
|
||||
|
||||
int p = Modular::get_current_prime(); // implicit!
|
||||
//std::cout <<" p: "<< p<<std::endl;
|
||||
//std::cout <<" q: "<< q<<std::endl;
|
||||
//std::cout <<" gcd(p,q): "<< gcd(p,q)<<std::endl;
|
||||
chinese_remainder(q,Gs ,p,inv_map(mG_),pq,Gs);
|
||||
chinese_remainder(q,H1s,p,inv_map(mH1),pq,H1s);
|
||||
chinese_remainder(q,H2s,p,inv_map(mH2),pq,H2s);
|
||||
q=pq;
|
||||
}
|
||||
|
||||
// std::cout << "Gs: "<< Gs << std::endl;
|
||||
// std::cout << "H1s: "<< H1s << std::endl;
|
||||
// std::cout << "H2s: "<< H2s << std::endl;
|
||||
// std::cout <<std::endl;
|
||||
// std::cout << "F1s: "<<Gs*H1s<< std::endl;
|
||||
// std::cout << "F1 : "<<F1_<< std::endl;
|
||||
// std::cout << "diff : "<<F1_-Gs*H1s<< std::endl;
|
||||
// std::cout <<std::endl;
|
||||
// std::cout << "F2s: "<<Gs*H2s<< std::endl;
|
||||
// std::cout << "F2 : "<<F2_<< std::endl;
|
||||
// std::cout << "diff : "<<F2_-Gs*H2s<< std::endl;
|
||||
|
||||
try{// This is a HACK!!!!
|
||||
// TODO: in case of Sqrt_extension it may happen that the disr (root)
|
||||
// is not correct, in this case the behavior of the code is unclear
|
||||
// if CGAL is in debug mode it throws an error
|
||||
if( Gs*H1s == F1_ && Gs*H2s == F2_ ){
|
||||
solved = true;
|
||||
}
|
||||
}catch(...){}
|
||||
|
||||
//std::cout << "Gs: " << CGAL::canonicalize_polynomial(Gs)<<std::endl;
|
||||
// std::cout << "canonical(Gs): " << CGAL::canonicalize_polynomial(Gs)<<std::endl;
|
||||
|
||||
//std::cout << std::endl;
|
||||
|
||||
|
||||
}
|
||||
|
||||
//std::cout << "G: " << CGAL::canonicalize_polynomial(gcd_utcf(F1,F2)) << std::endl;
|
||||
|
||||
return canonicalize(Gs);
|
||||
|
||||
}
|
||||
|
||||
|
||||
}///namespace CGAL
|
||||
|
||||
#endif //#ifnedef CGAL_MODULAR_GCD_H 1
|
||||
|
||||
|
|
@ -1,69 +0,0 @@
|
|||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Testsuite/assert.h>
|
||||
#include <CGAL/Chinese_remainder_traits.h>
|
||||
|
||||
#ifdef CGAL_USE_LEDA
|
||||
#include <CGAL/leda_integer.h>
|
||||
#endif
|
||||
|
||||
#ifdef CGAL_USE_CORE
|
||||
#include <CGAL/CORE_BigInt.h>
|
||||
#endif
|
||||
|
||||
template <class CRT>
|
||||
void test_chinese_remainder_traits(){
|
||||
typedef typename CRT::T T;
|
||||
typedef typename CRT::Scalar_type Scalar_type;
|
||||
typedef typename CRT::Chinese_remainder Chinese_remainder;
|
||||
Chinese_remainder chinese_remainder;
|
||||
|
||||
T x(-574);
|
||||
Scalar_type m1 = Scalar_type(23);
|
||||
Scalar_type m2 = Scalar_type(17);
|
||||
Scalar_type m3 = Scalar_type(29);
|
||||
T u1 = -T(22);
|
||||
T u2 = -T(13);
|
||||
T u3 = -T(23);
|
||||
Scalar_type m;
|
||||
T u;
|
||||
|
||||
chinese_remainder(m1,u1,m2,u2,m,u);
|
||||
chinese_remainder(m ,u ,m3,u3,m,u);
|
||||
|
||||
CGAL_test_assert( m == m1*m2*m3 );
|
||||
CGAL_test_assert( x == u );
|
||||
|
||||
}
|
||||
|
||||
|
||||
int main(){
|
||||
|
||||
test_chinese_remainder_traits<CGAL::Chinese_remainder_traits<int> >();
|
||||
typedef CGAL::Sqrt_extension<int,int> Extn_1;
|
||||
typedef CGAL::Sqrt_extension<Extn_1,int> Extn_2;
|
||||
typedef CGAL::Sqrt_extension<Extn_1,Extn_1> Extn_n2;
|
||||
|
||||
test_chinese_remainder_traits<CGAL::Chinese_remainder_traits<Extn_1 > >();
|
||||
test_chinese_remainder_traits<CGAL::Chinese_remainder_traits<Extn_2 > >();
|
||||
test_chinese_remainder_traits<CGAL::Chinese_remainder_traits<Extn_n2 > >();
|
||||
|
||||
typedef CGAL::Polynomial<int> Poly_1;
|
||||
typedef CGAL::Polynomial<Poly_1> Poly_2;
|
||||
test_chinese_remainder_traits<CGAL::Chinese_remainder_traits<Poly_1 > >();
|
||||
test_chinese_remainder_traits<CGAL::Chinese_remainder_traits<Poly_2 > >();
|
||||
|
||||
|
||||
#ifdef CGAL_USE_CORE
|
||||
test_chinese_remainder_traits<
|
||||
CGAL::Chinese_remainder_traits<CORE::BigInt> >();
|
||||
#endif
|
||||
|
||||
#ifdef CGAL_USE_LEDA
|
||||
test_chinese_remainder_traits<
|
||||
CGAL::Chinese_remainder_traits<leda::integer> >();
|
||||
#endif
|
||||
|
||||
|
||||
}
|
||||
|
|
@ -1,62 +0,0 @@
|
|||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
|
||||
#include <CGAL/Testsuite/assert.h>
|
||||
#include <CGAL/chinese_remainder.h>
|
||||
#include <cstdlib>
|
||||
|
||||
#ifdef CGAL_USE_LEDA
|
||||
#include <CGAL/leda_integer.h>
|
||||
#endif
|
||||
|
||||
#ifdef CGAL_USE_CORE
|
||||
#include <CGAL/CORE_BigInt.h>
|
||||
#endif
|
||||
|
||||
template <class NT>
|
||||
void test_chinese_remainder(NT x){
|
||||
typedef CGAL::Algebraic_structure_traits<NT> AST;
|
||||
typename AST::Mod mod;
|
||||
|
||||
NT m1 = 23;
|
||||
NT m2 = 17;
|
||||
NT m3 = 29;
|
||||
NT u1 = mod(x,m1);
|
||||
NT u2 = mod(x,m2);
|
||||
NT u3 = mod(x,m3);
|
||||
NT m,u;
|
||||
|
||||
CGAL::chinese_remainder(m1,u1,m2,u2,m,u);
|
||||
CGAL::chinese_remainder(m ,u ,m3,u3,m,u);
|
||||
|
||||
CGAL_test_assert( m == m1*m2*m3 );
|
||||
CGAL_test_assert( x == u );
|
||||
}
|
||||
|
||||
template <class NT>
|
||||
void test_chinese_remainder(){
|
||||
test_chinese_remainder(NT(0));
|
||||
test_chinese_remainder(NT(1));
|
||||
test_chinese_remainder(NT(-1));
|
||||
test_chinese_remainder(NT(23));
|
||||
test_chinese_remainder(NT(17));
|
||||
test_chinese_remainder(NT(-29));
|
||||
test_chinese_remainder(NT(2456));
|
||||
test_chinese_remainder(NT(-2456));
|
||||
}
|
||||
|
||||
int main(){
|
||||
test_chinese_remainder<int>();
|
||||
|
||||
#ifdef CGAL_USE_LEDA
|
||||
test_chinese_remainder<leda_integer>();
|
||||
#endif
|
||||
|
||||
#ifdef CGAL_USE_CORE
|
||||
test_chinese_remainder<CORE::BigInt>();
|
||||
#endif
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -1,77 +0,0 @@
|
|||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
|
||||
#include <CGAL/Testsuite/assert.h>
|
||||
#include <CGAL/euclidean_algorithm.h>
|
||||
|
||||
#ifdef CGAL_USE_LEDA
|
||||
#include <CGAL/leda_integer.h>
|
||||
#endif
|
||||
|
||||
#ifdef CGAL_USE_CORE
|
||||
#include <CGAL/CORE_BigInt.h>
|
||||
#endif
|
||||
|
||||
#include <cstdlib>
|
||||
|
||||
|
||||
template <class NT>
|
||||
void test_euclidean_algorithm(){
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 0),NT(0)) == NT( 0));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 7),NT(0)) == NT( 7));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT(-7),NT(0)) == NT( 7));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 0),NT(7)) == NT( 7));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT(0),NT(-7)) == NT( 7));
|
||||
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 1),NT(1)) == NT( 1));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 7),NT(1)) == NT( 1));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT(-7),NT(1)) == NT( 1));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 1),NT(7)) == NT( 1));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT(1),NT(-7)) == NT( 1));
|
||||
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 7),NT(7)) == NT( 7));
|
||||
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 3),NT(1)) == NT( 1));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 3),NT(6)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 6),NT(9)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 9),NT(15)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 15),NT(24)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 24),NT(39)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 39),NT(63)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 6),NT(3)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 9),NT(6)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 15),NT(9)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 24),NT(15)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 39),NT(24)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 63),NT(39)) == NT( 3));
|
||||
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -3),NT(6)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -6),NT(9)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -9),NT(15)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -15),NT(24)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -6),NT(3)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -9),NT(6)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -15),NT(9)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -24),NT(15)) == NT( 3));
|
||||
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 3),NT(-1)) == NT( 1));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 3),NT(-6)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 6),NT(-9)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 6),NT(-3)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( 9),NT(-6)) == NT( 3));
|
||||
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -3),NT(-6)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -6),NT(-9)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -6),NT(-3)) == NT( 3));
|
||||
CGAL_test_assert( CGAL::euclidean_algorithm(NT( -9),NT(-6)) == NT( 3));
|
||||
}
|
||||
|
||||
|
||||
int main(){
|
||||
test_euclidean_algorithm<int>();
|
||||
test_euclidean_algorithm<leda_integer>();
|
||||
test_euclidean_algorithm<CORE::BigInt>();
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -1,77 +0,0 @@
|
|||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
|
||||
#include <CGAL/Testsuite/assert.h>
|
||||
#include <CGAL/extended_euclidean_algorithm.h>
|
||||
#include <cstdlib>
|
||||
|
||||
|
||||
template <class NT>
|
||||
void test_extended_euclidean_algorithm
|
||||
(const NT& a, const NT& b, const NT& g_){
|
||||
NT u,v;
|
||||
NT g = CGAL::extended_euclidean_algorithm(a,b,u,v);
|
||||
CGAL_test_assert(g_ == g) ;
|
||||
CGAL_test_assert(g_ == u*a+v*b);
|
||||
}
|
||||
|
||||
|
||||
template <class NT>
|
||||
void test_extended_euclidean_algorithm(){
|
||||
|
||||
test_extended_euclidean_algorithm(NT( 0),NT(0) , NT( 0));
|
||||
test_extended_euclidean_algorithm(NT( 7),NT(0) , NT( 7));
|
||||
test_extended_euclidean_algorithm(NT(-7),NT(0) , NT( 7));
|
||||
test_extended_euclidean_algorithm(NT( 0),NT(7) , NT( 7));
|
||||
test_extended_euclidean_algorithm(NT(0),NT(-7) , NT( 7));
|
||||
|
||||
test_extended_euclidean_algorithm(NT( 1),NT(1) , NT( 1));
|
||||
test_extended_euclidean_algorithm(NT( 7),NT(1) , NT( 1));
|
||||
test_extended_euclidean_algorithm(NT(-7),NT(1) , NT( 1));
|
||||
test_extended_euclidean_algorithm(NT( 1),NT(7) , NT( 1));
|
||||
test_extended_euclidean_algorithm(NT(1),NT(-7) , NT( 1));
|
||||
|
||||
test_extended_euclidean_algorithm(NT( 7),NT(7) , NT( 7));
|
||||
|
||||
test_extended_euclidean_algorithm(NT( 3),NT(1) , NT( 1));
|
||||
test_extended_euclidean_algorithm(NT( 3),NT(6) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 6),NT(9) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 9),NT(15) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 15),NT(24) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 24),NT(39) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 39),NT(63) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 6),NT(3) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 9),NT(6) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 15),NT(9) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 24),NT(15) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 39),NT(24) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 63),NT(39) , NT( 3));
|
||||
|
||||
test_extended_euclidean_algorithm(NT( -3),NT(6) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -6),NT(9) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -9),NT(15) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -15),NT(24) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -6),NT(3) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -9),NT(6) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -15),NT(9) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -24),NT(15) , NT( 3));
|
||||
|
||||
test_extended_euclidean_algorithm(NT( 3),NT(-1) , NT( 1));
|
||||
test_extended_euclidean_algorithm(NT( 3),NT(-6) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 6),NT(-9) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 6),NT(-3) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( 9),NT(-6) , NT( 3));
|
||||
|
||||
test_extended_euclidean_algorithm(NT( -3),NT(-6) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -6),NT(-9) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -6),NT(-3) , NT( 3));
|
||||
test_extended_euclidean_algorithm(NT( -9),NT(-6), NT( 3));
|
||||
}
|
||||
|
||||
|
||||
int main(){
|
||||
test_extended_euclidean_algorithm<int>();
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -15,7 +15,6 @@ include $(CGAL_MAKEFILE)
|
|||
|
||||
CXXFLAGS = \
|
||||
-I../../include \
|
||||
-I../../../Polynomial/include \
|
||||
-I../../../Number_types/test/Number_types/include \
|
||||
$(CGAL_CXXFLAGS) \
|
||||
$(LONG_NAME_PROBLEM_CXXFLAGS)
|
||||
|
|
@ -36,44 +35,18 @@ LDFLAGS = \
|
|||
#---------------------------------------------------------------------#
|
||||
|
||||
all: \
|
||||
chinese_remainder$(EXE_EXT) \
|
||||
Chinese_remainder_traits$(EXE_EXT) \
|
||||
euclidean_algorithm$(EXE_EXT) \
|
||||
extended_euclidean_algorithm$(EXE_EXT) \
|
||||
Modular$(EXE_EXT) \
|
||||
modular_gcd$(EXE_EXT) \
|
||||
Modular_traits$(EXE_EXT)
|
||||
|
||||
chinese_remainder$(EXE_EXT): chinese_remainder$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)chinese_remainder chinese_remainder$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
Chinese_remainder_traits$(EXE_EXT): Chinese_remainder_traits$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)Chinese_remainder_traits Chinese_remainder_traits$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
euclidean_algorithm$(EXE_EXT): euclidean_algorithm$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)euclidean_algorithm euclidean_algorithm$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
extended_euclidean_algorithm$(EXE_EXT): extended_euclidean_algorithm$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)extended_euclidean_algorithm extended_euclidean_algorithm$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
Modular$(EXE_EXT): Modular$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)Modular Modular$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
modular_gcd$(EXE_EXT): modular_gcd$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)modular_gcd modular_gcd$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
Modular_traits$(EXE_EXT): Modular_traits$(OBJ_EXT)
|
||||
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)Modular_traits Modular_traits$(OBJ_EXT) $(LDFLAGS)
|
||||
|
||||
clean: \
|
||||
chinese_remainder.clean \
|
||||
Chinese_remainder_traits.clean \
|
||||
euclidean_algorithm.clean \
|
||||
extended_euclidean_algorithm.clean \
|
||||
Modular.clean \
|
||||
modular_gcd.clean \
|
||||
Modular_traits.clean \
|
||||
src_Modular.clean
|
||||
Modular_traits.clean
|
||||
|
||||
#---------------------------------------------------------------------#
|
||||
# suffix rules
|
||||
|
|
|
|||
|
|
@ -1,128 +0,0 @@
|
|||
// Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
|
||||
|
||||
/*! \file CGAL/Modular.C
|
||||
test for number type modul
|
||||
*/
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Testsuite/assert.h>
|
||||
#include <CGAL/Modular.h>
|
||||
#include <CGAL/modular_gcd.h>
|
||||
#include <CGAL/Polynomial.h>
|
||||
|
||||
#ifdef CGAL_USE_LEDA
|
||||
#include <CGAL/leda_integer.h>
|
||||
#endif // CGAL_USE_LEDA
|
||||
|
||||
#include <cstdlib>
|
||||
|
||||
#include <boost/type_traits.hpp>
|
||||
|
||||
int main()
|
||||
{
|
||||
{
|
||||
typedef leda::integer Integer;
|
||||
typedef CGAL::Polynomial<Integer> Polynomial;
|
||||
Polynomial p1(123,431,2134);
|
||||
Polynomial p2(123,421,234);
|
||||
Polynomial g(1);
|
||||
Polynomial result = modular_gcd_utcf(p1,p2);
|
||||
//std::cout <<" result : " << result <<std::endl;
|
||||
//std::cout <<" true gcd : " << g <<std::endl;
|
||||
CGAL_test_assert(result == g);
|
||||
|
||||
}
|
||||
{ // unlucky prime test
|
||||
typedef leda::integer Integer;
|
||||
typedef CGAL::Polynomial<Integer> Polynomial;
|
||||
|
||||
Polynomial f1(5,234,445);
|
||||
Polynomial f2(12,-234,345);
|
||||
f1 *= Polynomial(Integer(CGAL::primes[0]+3),Integer(1));
|
||||
f2 *= Polynomial(Integer(3),Integer(1));
|
||||
Polynomial g(13,96,2345);
|
||||
Polynomial p1 = f1*g;
|
||||
Polynomial p2 = f2*g;
|
||||
|
||||
Polynomial result = modular_gcd_utcf(p1,p2);
|
||||
//std::cout <<" result : " << result <<std::endl;
|
||||
//std::cout <<" true gcd : " << g <<std::endl;
|
||||
CGAL_test_assert(result == g);
|
||||
|
||||
}
|
||||
|
||||
{
|
||||
typedef leda::integer Integer;
|
||||
typedef CGAL::Polynomial<Integer> Polynomial;
|
||||
|
||||
Polynomial f1(5,234,-26,243,745);
|
||||
Polynomial f2(12,-234,26,243,-731);
|
||||
Polynomial g(13,-5676,234,96);
|
||||
Polynomial p1 = Polynomial(8)*f1*f1*g;
|
||||
Polynomial p2 = Polynomial(5)*f2*f2*g;
|
||||
|
||||
Polynomial result = modular_gcd_utcf(p1,p2);
|
||||
//std::cout <<" result : " << result <<std::endl;
|
||||
//std::cout <<" true gcd : " << g <<std::endl;
|
||||
CGAL_test_assert(result == g);
|
||||
|
||||
}
|
||||
{ // test for sqrt
|
||||
typedef leda::integer Integer;
|
||||
typedef CGAL::Sqrt_extension<Integer,Integer> EXT;
|
||||
typedef CGAL::Polynomial<EXT> Polynomial;
|
||||
|
||||
Integer root(Integer(789234));
|
||||
Polynomial f1(EXT(235143,-2234,root),EXT(232543,-2334,root),EXT(235403,-2394,root),EXT(235483,-2364,root),EXT(223443,-2234,root));
|
||||
Polynomial f2(EXT(25143,-2134,root),EXT(212543,-2315,root),EXT(255453,-5394,root),EXT(535483,-2354,root),EXT(22333,-2214,root));
|
||||
Polynomial g(EXT(215143,-2134,root),EXT(2122422543,-2115,root),EXT(255453,-1394,root),EXT(135483,-2354,root),EXT(7));
|
||||
g=g*g;g=g*g;g=g*g;
|
||||
|
||||
Polynomial p1 = Polynomial(8)*f1*f1*g;
|
||||
Polynomial p2 = Polynomial(5)*f2*f2*g;
|
||||
|
||||
Polynomial result = modular_gcd_utcf(p1,p2);
|
||||
//std::cout <<" result : " << result <<std::endl;
|
||||
//std::cout <<" true gcd : " << g <<std::endl;
|
||||
CGAL_test_assert(result == g);
|
||||
}
|
||||
|
||||
{ // test for sqrt / denom for algebraic integer
|
||||
CGAL::set_pretty_mode(std::cout);
|
||||
typedef leda::integer Integer;
|
||||
typedef CGAL::Sqrt_extension<Integer,Integer> EXT;
|
||||
typedef CGAL::Polynomial<EXT> Polynomial;
|
||||
|
||||
Integer root(4*3);
|
||||
Polynomial f1(EXT(0,1,root),EXT(2));
|
||||
Polynomial f2(EXT(0,1,root),EXT(4));
|
||||
Polynomial g(EXT(0,-1,root),EXT(2));
|
||||
|
||||
|
||||
//std::cout <<" f1 : " << f1 <<std::endl;
|
||||
//std::cout <<" f2 : " << f2 <<std::endl;
|
||||
//std::cout <<" g : " << g <<std::endl;
|
||||
|
||||
Polynomial p1 = f1*g;
|
||||
Polynomial p2 = f2*g;
|
||||
|
||||
|
||||
//std::cout <<" p1 : " << p1 <<std::endl;
|
||||
//std::cout <<" p2 : " << p2 <<std::endl;
|
||||
//std::cout << std::endl;
|
||||
|
||||
Polynomial result = modular_gcd_utcf(p1,p2);
|
||||
//std::cout <<" result : " << result <<std::endl;
|
||||
//std::cout <<" true gcd : " << g <<std::endl;
|
||||
|
||||
CGAL_test_assert(result == g);
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
return 0 ;
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -1,5 +0,0 @@
|
|||
namespace CGAL{
|
||||
int Modular::prime_int = 67111067;
|
||||
double Modular::prime =67111067.0;
|
||||
double Modular::prime_inv =1/67111067.0;
|
||||
} // namespace CGAL
|
||||
Loading…
Reference in New Issue