songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

rm file Modular_type 

added test for class Modular
added separate test for class Modular_traits
This commit is contained in:
Michael Hemmer 2007-03-02 12:02:31 +00:00
parent 4eb4184763
commit 60ef82986d
8 changed files with 417 additions and 348 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
.gitattributes vendored
View File

@ -1388,6 +1388,7 @@ Modifier/doc_tex/Modifier/idraw/modifier.eps -text svneol=unset#application/post
Modifier/doc_tex/Modifier/idraw/modifier.pdf -text svneol=unset#application/pdf
Modifier/doc_tex/Modifier/modifier.gif -text svneol=unset#image/gif
Modifier/doc_tex/Modifier/modifier_small.gif -text svneol=unset#image/gif
Modular_arithmetic/test/Modular_arithmetic/Modular_traits.C -text
Nef_2/demo/Nef_2/filtered_homogeneous_data/complex.nef -text svneol=native#application/octet-stream
Nef_2/demo/Nef_2/filtered_homogeneous_data/symmdif.nef -text svneol=native#application/octet-stream
Nef_2/demo/Nef_2/help/index.html svneol=native#text/html

285
Modular_arithmetic/include/CGAL/Modular.h
View File

@ -12,14 +12,282 @@
#define CGAL_MODULAR_H 1
#include <CGAL/basic.h>
#include <CGAL/Modular_type.h>
#include <CGAL/Modular_traits.h>
#include <cfloat>
namespace CGAL {
// I/O
class Modular;
Modular operator + (const Modular&);
Modular operator - (const Modular&);
Modular operator + (const Modular&, const Modular&);
Modular operator - (const Modular&, const Modular&);
Modular operator * (const Modular&, const Modular&);
Modular operator / (const Modular&, const Modular&);
std::ostream& operator << (std::ostream& os, const Modular& p);
std::istream& operator >> (std::istream& is, Modular& p);
/*! \ingroup CGAL_Modular_traits
* \brief This class represents the Field Z mod p.
*
* This class uses the type double for representation.
* Therefore the value of p is restricted to primes less than 2^26.
* By default p is set to 67111067.
*
* It provides the standard operators +,-,*,/ as well as in&output.
*
* \see Modular_traits
*/
class Modular{
public:
typedef Modular Self;
typedef Modular NT;
private:
static const double CST_CUT = ((3.0*(1<<30))*(1<<21));
private:
static int prime_int;
static double prime;
static double prime_inv;
/* Quick integer rounding, valid if a<2^51. for double */
static inline
double MOD_round (double a){
#ifdef LiS_HAVE_LEDA
return ( (a + CST_CUT) - CST_CUT);
#else
// TODO:
// In order to get rid of the volatile double
// one should call:
// CGAL/FPU.h : inline void force_ieee_double_precision()
// the problem is where and when it should be called ?
// and whether on should restore the old behaviour
// since it changes the global behaviour of doubles.
// Note that this code works if LEDA is present, since leda automatically
// changes this behaviour in the desired way.
volatile double b = (a + CST_CUT);
return b - CST_CUT;
#endif
}
/* Big modular reduction (e.g. after multiplication) */
static inline
double MOD_reduce (double a){
return a - prime * MOD_round(a * prime_inv);
}
/* Little modular reduction (e.g. after a simple addition). */
static inline
double MOD_soft_reduce (double a){
double b = 2*a;
return (b>prime) ? a-prime :
((b<-prime) ? a+prime : a);
}
/* -a */
static inline
double MOD_negate(double a){
return MOD_soft_reduce(-a);
}
/* a*b */
static inline
double MOD_mul (double a, double b){
double c = a*b;
return MOD_reduce(c);
}
/* a+b */
static inline
double MOD_add (double a, double b){
double c = a+b;
return MOD_soft_reduce(c);
}
/* a^-1, using Bezout (extended Euclidian algorithm). */
static inline
double MOD_inv (double ri1){
double bi = 0.0;
double bi1 = 1.0;
double ri = prime;
double p, tmp, tmp2;
Real_embeddable_traits<double>::Abs double_abs;
while (double_abs(ri1) != 1.0)
{
p = MOD_round(ri/ri1);
tmp = bi - p * bi1;
tmp2 = ri - p * ri1;
bi = bi1;
ri = ri1;
bi1 = tmp;
ri1 = tmp2;
};
return ri1 * MOD_soft_reduce(bi1); /* Quicker !!!! */
}
/* a/b */
static inline
double MOD_div (double a, double b){
return MOD_mul(a, MOD_inv(b));
}
public:
/*! \brief sets the current prime.
*
* Note that you are going to change a static member!
* \pre p is prime, but we abstained from such a test.
* \pre 0 < p < 2^26
*
*/
static int
set_current_prime(int p){
int old_prime = prime_int;
prime_int = p;
prime = (double)p;
prime_inv = (double)1/prime;
return old_prime;
}
/*! \brief return the current prime. */
static int get_current_prime(){
return prime_int;
}
int get_value() const{
return int(x_);
}
private:
double x_;
public:
//! Explicit constructor of Modular, from int
Modular(int n = 0){
x_= MOD_reduce(n);
}
//! Explicit constructor of Modular, from long
Modular(long n){
x_= MOD_reduce(n);
}
//! Access operator for x, \c const
const double& x() const { return x_; }
//! Access operator for x
double& x() { return x_; }
Self& operator += (const Self& p2) {
x() = MOD_add(x(),p2.x());
return (*this);
}
Self& operator -= (const Self& p2){
x() = MOD_add(x(),MOD_negate(p2.x()));
return (*this);
}
Self& operator *= (const Self& p2){
x() = MOD_mul(x(),p2.x());
return (*this);
}
Self& operator /= (const Self& p2) {
x() = MOD_div(x(),p2.x());
return (*this);
}
friend Self operator + (const Self&);
friend Self operator - (const Self&);
friend Self operator + (const Self&, const Self&);
friend Self operator - (const Self&, const Self&);
friend Self operator * (const Self&, const Self&);
friend Self operator / (const Self& p1, const Self& p2);
};
inline Modular operator + (const Modular& p1)
{ return p1; }
inline Modular operator - (const Modular& p1){
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_negate(p1.x());
return r;
}
inline Modular operator + (const Modular& p1,const Modular& p2) {
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_add(p1.x(),p2.x());
return r;
}
inline Modular operator - (const Modular& p1, const Modular& p2) {
return p1+(-p2);
}
inline Modular operator * (const Modular& p1, const Modular& p2) {
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_mul(p1.x(),p2.x());
return r;
}
inline Modular operator / (const Modular& p1, const Modular& p2) {
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_div(p1.x(),p2.x());
return r;
}
inline bool operator == (const Modular& p1, const Modular& p2)
{ return ( p1.x()==p2.x() ); }
inline bool operator != (const Modular& p1, const Modular& p2)
{ return ( p1.x()!=p2.x() ); }
// left hand side
inline bool operator == (int num, const Modular& p)
{ return ( Modular(num) == p );}
inline bool operator != (int num, const Modular& p)
{ return ( Modular(num) != p );}
// right hand side
inline bool operator == (const Modular& p, int num)
{ return ( Modular(num) == p );}
inline bool operator != (const Modular& p, int num)
{ return ( Modular(num) != p );}
// left hand side
inline Modular operator + (int num, const Modular& p2)
{ return (Modular(num) + p2); }
inline Modular operator - (int num, const Modular& p2)
{ return (Modular(num) - p2); }
inline Modular operator * (int num, const Modular& p2)
{ return (Modular(num) * p2); }
inline Modular operator / (int num, const Modular& p2)
{ return (Modular(num)/p2); }
// right hand side
inline Modular operator + (const Modular& p1, int num)
{ return (p1 + Modular(num)); }
inline Modular operator - (const Modular& p1, int num)
{ return (p1 - Modular(num)); }
inline Modular operator * (const Modular& p1, int num)
{ return (p1 * Modular(num)); }
inline Modular operator / (const Modular& p1, int num)
{ return (p1 / Modular(num)); }
// I/O
inline std::ostream& operator << (std::ostream& os, const Modular& p) {
typedef Modular MOD;
os <<"("<< p.x()<<"%"<<MOD::get_current_prime()<<")";
@ -44,12 +312,15 @@ inline std::istream& operator >> (std::istream& is, Modular& p) {
* \ingroup CGAL_NT_traits_spec
*/
template <>
class Algebraic_structure_traits<Modular>
: public Algebraic_structure_traits_base< Modular ,Field_tag >{};
struct Algebraic_structure_traits<Modular>
: public Algebraic_structure_traits_base< Modular ,Field_tag >{
typedef CGAL::Tag_true Is_exact;
};
}///namespace CGAL
//TODO: move to src/
#include <./src_Modular.C>
#endif //#ifnedef CGAL_MODULAR_H 1

15
Modular_arithmetic/include/CGAL/Modular_traits.h
View File

@ -3,7 +3,8 @@
#ifndef CGAL_MODULAR_TRAITS_H
#define CGAL_MODULAR_TRAITS_H 1
#include <CGAL/basic.h>
#include <CGAL/Modular.h>
#include <CGAL/leda_integer.h>
#include <CGAL/Sqrt_extension.h>
#include <vector>
@ -189,6 +190,18 @@ public:
};
};
// TODO: put this into Modular_arithmetic/src/
int primes[64] = {
67089287,67089299,67089329,67089377,67089461,67089469,67089479,67089511,
67089527,67089541,67089577,67089587,67089619,67089683,67089697,67089707,
67089721,67089733,67089739,67089751,67089793,67089809,67089811,67089829,
67089839,67089857,67089877,67089907,67089943,67089949,67089989,67090013,
67090027,67090031,67090033,67090043,67090061,67090073,67090091,67090099,
67090117,67090129,67090151,67090171,67090189,67090207,67090217,67090223,
67090229,67090237,67090259,67090271,67090307,67090321,67090343,67090351,
67090399,67090403,67090411,67090433,67090451,67090459,67090489,67090519
};
}///namespace CGAL
#endif //#ifnedef CGAL_MODULAR_TRAITS_H 1

295
Modular_arithmetic/include/CGAL/Modular_type.h
View File

@ -1,295 +0,0 @@
//Author(s) : Michael Hemmer <mhemmer@uni-mainz.de>
// Sylvain Pion
// defines the pure type Modular, i.e. no CGAL support
#ifndef CGAL_MODULAR_TYPE_H
#define CGAL_MODULAR_TYPE_H 1
#include <CGAL/basic.h>
namespace CGAL {
int primes[] = { // 64 primes
67089287,67089299,67089329,67089377,67089461,67089469,67089479,67089511,
67089527,67089541,67089577,67089587,67089619,67089683,67089697,67089707,
67089721,67089733,67089739,67089751,67089793,67089809,67089811,67089829,
67089839,67089857,67089877,67089907,67089943,67089949,67089989,67090013,
67090027,67090031,67090033,67090043,67090061,67090073,67090091,67090099,
67090117,67090129,67090151,67090171,67090189,67090207,67090217,67090223,
67090229,67090237,67090259,67090271,67090307,67090321,67090343,67090351,
67090399,67090403,67090411,67090433,67090451,67090459,67090489,67090519
};
class Modular;
Modular operator + (const Modular&);
Modular operator - (const Modular&);
Modular operator + (const Modular&, const Modular&);
Modular operator - (const Modular&, const Modular&);
Modular operator * (const Modular&, const Modular&);
Modular operator / (const Modular&, const Modular&);
std::ostream& operator << (std::ostream& os, const Modular& p);
std::istream& operator >> (std::istream& is, Modular& p);
/*! \ingroup CGAL_Modular_traits
* \brief This class represents the Field Z mod p.
*
* This class uses the type double for representation.
* Therefore the value of p is restricted to primes less than 2^26.
* By default p is set to 67111067.
*
* It provides the standard operators +,-,*,/ as well as in&output.
*
* \see Modular_traits
*/
class Modular{
public:
typedef Modular Self;
typedef Modular NT;
private:
static const double CST_CUT = ((3.0*(1<<30))*(1<<21));
private:
static int prime_int;
static double prime;
static double prime_inv;
/* Quick integer rounding, valid if a<2^51. for double */
static inline
double MOD_round (double a){
#ifdef LiS_HAVE_LEDA
return ( (a + CST_CUT) - CST_CUT);
#else
// TODO:
// In order to get rid of the volatile double
// one should call:
// CGAL/FPU.h : inline void force_ieee_double_precision()
// the problem is where and when it should be called ?
// and whether on should restore the old behaviour
// since it changes the global behaviour of doubles.
// Note that this code works if LEDA is present, since leda automatically
// changes this behaviour in the desired way.
volatile double b = (a + CST_CUT);
return b - CST_CUT;
#endif
}
/* Big modular reduction (e.g. after multiplication) */
static inline
double MOD_reduce (double a){
return a - prime * MOD_round(a * prime_inv);
}
/* Little modular reduction (e.g. after a simple addition). */
static inline
double MOD_soft_reduce (double a){
double b = 2*a;
return (b>prime) ? a-prime :
((b<-prime) ? a+prime : a);
}
/* -a */
static inline
double MOD_negate(double a){
return MOD_soft_reduce(-a);
}
/* a*b */
static inline
double MOD_mul (double a, double b){
double c = a*b;
return MOD_reduce(c);
}
/* a+b */
static inline
double MOD_add (double a, double b){
double c = a+b;
return MOD_soft_reduce(c);
}
/* a^-1, using Bezout (extended Euclidian algorithm). */
static inline
double MOD_inv (double ri1){
double bi = 0.0;
double bi1 = 1.0;
double ri = prime;
double p, tmp, tmp2;
Real_embeddable_traits<double>::Abs double_abs;
while (double_abs(ri1) != 1.0)
{
p = MOD_round(ri/ri1);
tmp = bi - p * bi1;
tmp2 = ri - p * ri1;
bi = bi1;
ri = ri1;
bi1 = tmp;
ri1 = tmp2;
};
return ri1 * MOD_soft_reduce(bi1); /* Quicker !!!! */
}
/* a/b */
static inline
double MOD_div (double a, double b){
return MOD_mul(a, MOD_inv(b));
}
public:
/*! \brief sets the current prime.
*
* Note that you are going to change a static member!
* \pre p is prime, but we abstained from such a test.
* \pre 0 < p < 2^26
*
*/
static void
set_current_prime(int p){
prime_int=p;
prime = (double)p;
prime_inv = (double)1/prime;
}
/*! \brief return the current prime. */
static int get_current_prime(){
return prime_int;
}
int get_value() const{
return int(x_);
}
private:
double x_;
public:
//! Explicit constructor of Modular, from int
Modular(int n = 0){
x_= MOD_reduce(n);
}
//! Explicit constructor of Modular, from long
Modular(long n){
x_= MOD_reduce(n);
}
//! Access operator for x, \c const
const double& x() const { return x_; }
//! Access operator for x
double& x() { return x_; }
Self& operator += (const Self& p2) {
x() = MOD_add(x(),p2.x());
return (*this);
}
Self& operator -= (const Self& p2){
x() = MOD_add(x(),MOD_negate(p2.x()));
return (*this);
}
Self& operator *= (const Self& p2){
x() = MOD_mul(x(),p2.x());
return (*this);
}
Self& operator /= (const Self& p2) {
x() = MOD_div(x(),p2.x());
return (*this);
}
friend Self operator + (const Self&);
friend Self operator - (const Self&);
friend Self operator + (const Self&, const Self&);
friend Self operator - (const Self&, const Self&);
friend Self operator * (const Self&, const Self&);
friend Self operator / (const Self& p1, const Self& p2);
};
inline Modular operator + (const Modular& p1)
{ return p1; }
inline Modular operator - (const Modular& p1){
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_negate(p1.x());
return r;
}
inline Modular operator + (const Modular& p1,const Modular& p2) {
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_add(p1.x(),p2.x());
return r;
}
inline Modular operator - (const Modular& p1, const Modular& p2) {
return p1+(-p2);
}
inline Modular operator * (const Modular& p1, const Modular& p2) {
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_mul(p1.x(),p2.x());
return r;
}
inline Modular operator / (const Modular& p1, const Modular& p2) {
typedef Modular MOD;
Modular r;
r.x() = MOD::MOD_div(p1.x(),p2.x());
return r;
}
inline bool operator == (const Modular& p1, const Modular& p2)
{ return ( p1.x()==p2.x() ); }
inline bool operator != (const Modular& p1, const Modular& p2)
{ return ( p1.x()!=p2.x() ); }
// left hand side
inline bool operator == (int num, const Modular& p)
{ return ( Modular(num) == p );}
inline bool operator != (int num, const Modular& p)
{ return ( Modular(num) != p );}
// right hand side
inline bool operator == (const Modular& p, int num)
{ return ( Modular(num) == p );}
inline bool operator != (const Modular& p, int num)
{ return ( Modular(num) != p );}
// left hand side
inline Modular operator + (int num, const Modular& p2)
{ return (Modular(num) + p2); }
inline Modular operator - (int num, const Modular& p2)
{ return (Modular(num) - p2); }
inline Modular operator * (int num, const Modular& p2)
{ return (Modular(num) * p2); }
inline Modular operator / (int num, const Modular& p2)
{ return (Modular(num)/p2); }
// right hand side
inline Modular operator + (const Modular& p1, int num)
{ return (p1 + Modular(num)); }
inline Modular operator - (const Modular& p1, int num)
{ return (p1 - Modular(num)); }
inline Modular operator * (const Modular& p1, int num)
{ return (p1 * Modular(num)); }
inline Modular operator / (const Modular& p1, int num)
{ return (p1 / Modular(num)); }
}// namespace CGAL
#endif // CGAL_MODULAR_H 1

1
Modular_arithmetic/include/CGAL/modular_gcd.h
View File

@ -9,7 +9,6 @@
#define CGAL_MODULAR_GCD_H 1
#include <CGAL/basic.h>
#include <CGAL/Modular_type.h>
#include <CGAL/Modular_traits.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>

109
Modular_arithmetic/test/Modular_arithmetic/Modular.C
View File

@ -5,55 +5,76 @@
*/
#include <CGAL/basic.h>
#include <CGAL/Testsuite/assert.h>
#include <CGAL/Modular.h>
#include <./src_Modular.C>
#include <CGAL/_test_algebraic_structure.h>
#ifdef CGAL_USE_LEDA
#include <CGAL/leda_integer.h>
#include <CGAL/leda_rational.h>
#endif // CGAL_USE_LEDA
#ifdef CGAL_USE_CORE
#include <CGAL/CORE_BigInt.h>
#include <CGAL/CORE_BigRat.h>
#endif // CGAL_USE_CORE
#include <cstdlib>
#include <boost/type_traits.hpp>
template <class TESTT>
void test_modular_traits(){
typedef CGAL::Modular Modular;
typedef CGAL::Modular_traits<TESTT> MT;
typedef typename MT::Modular_NT Modular_NT;
typedef typename MT::Modular_image Modular_image;
typedef typename MT::Is_convertible Is_convertible;
typedef typename MT::NT NT;
CGAL_test_assert(
!(::boost::is_same<CGAL::Null_functor,Modular_image>::value));
CGAL_test_assert(
(::boost::is_same<CGAL::Tag_true,Is_convertible>::value));
CGAL_test_assert(
(::boost::is_same<TESTT,NT>::value));
Modular::set_current_prime(7);
Modular_image modular_image;
CGAL_test_assert(modular_image(TESTT(21)) == Modular_NT(0));
CGAL_test_assert(modular_image(TESTT(22)) == Modular_NT(1));
CGAL_test_assert(modular_image(TESTT(777777722)) == Modular_NT(1));
}
int main()
{
test_modular_traits<int>();
// test_modular_traits<long>();
for (int i = 0 ; i < 64 ; i++){
std::cout <<i<< ": "<<CGAL::primes[i] << std::endl;
}
typedef CGAL::Modular NT;
typedef CGAL::Field_tag Tag;
typedef CGAL::Tag_true Is_exact;
CGAL::test_algebraic_structure<NT,Tag, Is_exact>();
int old_prime = NT::get_current_prime();
CGAL_test_assert(old_prime == NT::set_current_prime(7));
CGAL_test_assert(7 == NT::get_current_prime());
NT x(4),y(5),z(12),t;
// operator ==
CGAL_test_assert(!(x==y));
CGAL_test_assert(y==z);
// operator !=
CGAL_test_assert(x!=y);
CGAL_test_assert(!(z!=y));
// constructor
CGAL_test_assert(NT(2)==NT(2-5*NT::get_current_prime()));
CGAL_test_assert(NT(2)==NT(2+5*NT::get_current_prime()));
// operator unary +
CGAL_test_assert((+x)==x);
// operator unary -
CGAL_test_assert(-x==x*NT(-1));
// operator binary +
CGAL_test_assert((x+y)==NT(2));
// operator binary -
CGAL_test_assert((x-y)==NT(6));
// operator *
CGAL_test_assert((x*y)==NT(6));
// operator /
CGAL_test_assert((x/y)==NT(5));
// operator +=
t=x; CGAL_test_assert((x+y)==(t+=y));
// operator -=
t=x; CGAL_test_assert((x-y)==(t-=y));
// operator *=
t=x; CGAL_test_assert((x*y)==(t*=y));
// operator /=
t=x; CGAL_test_assert((x/y)==(t/=y));
// left/right Hand
// operator ==
CGAL_test_assert(x==4);
CGAL_test_assert(5==y);
// operator !=
CGAL_test_assert(x!=5);
CGAL_test_assert(4!=y);
// operator +
t=x; CGAL_test_assert((x+5)==(5+x));
// operator -
t=x; CGAL_test_assert((x-5)==(4-y));
// operator *
t=x; CGAL_test_assert((x*5)==(5*x));
// operator =
t=x; CGAL_test_assert((x/5)==(4/y));
//cout << x << endl;
}

58
Modular_arithmetic/test/Modular_arithmetic/Modular_traits.C Normal file
View File

@ -0,0 +1,58 @@
// 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_traits.h>
#include <./src_Modular.C>
#ifdef CGAL_USE_LEDA
#include <CGAL/leda_integer.h>
#include <CGAL/leda_rational.h>
#endif // CGAL_USE_LEDA
#ifdef CGAL_USE_CORE
#include <CGAL/CORE_BigInt.h>
#include <CGAL/CORE_BigRat.h>
#endif // CGAL_USE_CORE
#include <cstdlib>
#include <boost/type_traits.hpp>
template <class TESTT>
void test_modular_traits(){
typedef CGAL::Modular Modular;
typedef CGAL::Modular_traits<TESTT> MT;
typedef typename MT::Modular_NT Modular_NT;
typedef typename MT::Modular_image Modular_image;
typedef typename MT::Is_convertible Is_convertible;
typedef typename MT::NT NT;
CGAL_test_assert(
!(::boost::is_same<CGAL::Null_functor,Modular_image>::value));
CGAL_test_assert(
(::boost::is_same<CGAL::Tag_true,Is_convertible>::value));
CGAL_test_assert(
(::boost::is_same<TESTT,NT>::value));
Modular::set_current_prime(7);
Modular_image modular_image;
CGAL_test_assert(modular_image(TESTT(21)) == Modular_NT(0));
CGAL_test_assert(modular_image(TESTT(22)) == Modular_NT(1));
CGAL_test_assert(modular_image(TESTT(777777722)) == Modular_NT(1));
}
int main()
{
test_modular_traits<int>();
// test_modular_traits<long>();
for (int i = 0 ; i < 64 ; i++){
std::cout <<i<< ": "<<CGAL::primes[i] << std::endl;
}
}

1
Modular_arithmetic/test/Modular_arithmetic/makefile
View File

@ -17,6 +17,7 @@ CXXFLAGS = \
-I. \
-I../../include \
-I../../../Polynomial/include \
-I../../../Number_types/test/Number_types/include/ \
$(CGAL_CXXFLAGS) \
$(LONG_NAME_PROBLEM_CXXFLAGS)