numerix_doc 0.4
modular_test.cpp
/******************************************************************************
* MODULE     : modular_test.cpp
* DESCRIPTION: Test modulars
* COPYRIGHT  : (C) 2008  Gregoire Lecerf
*******************************************************************************
* This software falls under the GNU general public license and comes WITHOUT
* ANY WARRANTY WHATSOEVER. See the file $TEXMACS_PATH/LICENSE for more details.
* If you don't have this file, write to the Free Software Foundation, Inc.,
* 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
******************************************************************************/

#include <stdlib.h>
#include <typeinfo>
#include <numerix/modular_int.hpp>
#include <numerix/modular_integer.hpp>

using namespace mmx;

long long unsigned int long_random (nat n) {
  long long unsigned int a = random ();
  a <<= 4*sizeof(long long int);
  a += rand ();
  return (a >> (8*sizeof(long long int) - n)) | (1ull << (n-1));
}

// Test operations modulo p
template<typename C, typename V>
static void test_to_integer (const C& a, const C& b, const C& p) {
  typedef modulus<C,V> M0;
  typedef modular<M0> R0;
  typedef modulus<integer> M1;
  typedef modular<M1> R1;
  
  M0 p0 (p);
  static nat m = V::template maximum_size_helper<C>::value;
  M1 p1 (p == 0 ? integer(1) << m : p);
  
  R0::set_modulus (p0);
  R1::set_modulus (p1);

  R0 x0 (a), y0 (b);
  R1 x1 (a), y1 (b);

  ASSERT (*(-x0) == *(-x1),         "error in neg");
  ASSERT (*(x0 + y0) == *(x1 + y1), "error in add");
  ASSERT (*(x0 - y0) == *(x1 - y1), "error in sub");
  ASSERT (*(x0 * y0) == *(x1 * y1), "error in mul");

  bool invertible = true;
  try { x0 / y0; x1 / y1; } catch (...) { invertible = false; }
  if (invertible && p != 0)
    ASSERT (*(x0 / y0) == *(x1 / y1), "error in div");
}

template<typename C, typename V>
static void test_random_int (C p) {
  nat bit_size= is_signed_helper<C>::value ? 8*sizeof (C) - 1 : 8*sizeof (C);
  C a (long_random (bit_size)), b (long_random (bit_size));
  try {
    test_to_integer<C,V> (a, b, p);
  }
  catch (exception e) {
    V v;
    mmout << e << "\n";
    mmout << "with arguments " << a << " and "
          << b << " modulo " << p << "\n";
    mmout << "base type C = " << cpp_demangle (typeid(a).name ()) << "\n";
    mmout << "variant V   = " << cpp_demangle (typeid(v).name ()) << "\n";
    exit (-1);
  }
}

template<typename C, nat nail>
static void test_modulus_int_naive () {
  static const nat m = is_signed_helper<C>::value ?
    8*sizeof(C) - 1 - nail : 8*sizeof(C) - nail;
  C p;
  // Random prime
  for (nat i = 1; i <= m; i++) {
    p = long_random (i);
    if (p == 0) p = 1;
    test_random_int<C, modulus_int_naive<m> > (p);
  }

  // Powers of 2
  for (nat i = 0; i<= m; i++) {
    p = ((C) 1) << i;
    test_random_int<C, modulus_int_naive<m> > (p);
  }
}

template<nat nail>
static void test_all_modulus_int_naive () {

  test_modulus_int_naive<unsigned char, nail> ();
  test_modulus_int_naive<short unsigned int, nail> ();
  test_modulus_int_naive<unsigned int, nail> ();
  test_modulus_int_naive<long unsigned int, nail> ();
  test_modulus_int_naive<long long unsigned int, nail> ();

  test_modulus_int_naive<signed char, nail> ();
  test_modulus_int_naive<short int, nail> ();
  test_modulus_int_naive<int, nail> ();
  test_modulus_int_naive<long int, nail> ();
  test_modulus_int_naive<long long int, nail> ();
}

template<typename C, nat m, C n> 
struct test_modulus_int_preinverse_fixed {

  static void op () {
    fixed_value<C, n> _c;
    modulus<C, modulus_int_preinverse<m> > mf (_c);
    modulus<C, modulus_int_preinverse<m> > md (n);
    try {
      ASSERT (mf.p == md.p, "error in auxiliaries");
      ASSERT (mf.q == md.q, "error in auxiliaries");
      ASSERT (mf.r == md.r, "error in auxiliaries");
      ASSERT (mf.s == md.s, "error in auxiliaries");
      ASSERT (mf.t == md.t, "error in auxiliaries");
    }
    catch (exception e) {
      modulus_int_preinverse<m> v;
      mmout << e << "\n";
      mmout << "base type C = " << cpp_demangle (typeid(mf.p).name ()) << "\n";
      mmout << "variant V   = " << cpp_demangle (typeid(v).name ()) << "\n";
      mmout << "static auxiliaries  (p, q, r, s, t) = ("
            << mf.p << ", " << mf.q << ", " << mf.r << ", "
            << mf.s << ", " << mf.t << ")\n";
      mmout << "dynamic auxiliaries (p, q, r, s, t) = ("
            << md.p << ", " << md.q << ", " << md.r << ", "
            << md.s << ", " << md.t << ")\n";
      exit (-1);
    }
  }
};

template<typename C, nat nail>
static void test_modulus_int_preinverse () {
  static const nat m = is_signed_helper<C>::value ?
    (8*sizeof(C) - 2 - nail) : (8*sizeof(C) - nail);
  C p;
  // Random prime
  for (nat i = 1; i <= m; i++) {
    p = long_random (i);
    if (p == 0) p = 1;
    test_random_int<C, modulus_int_preinverse<m> > (p);
  }

  // Powers of 2
  for (nat i = 0; i<= m; i++) {
    p = ((C) 1) << i;
    test_random_int<C, modulus_int_preinverse<m> > (p);
  }

  // Test auxiliaries
  test_modulus_int_preinverse_fixed <C, m, (C) 0>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 1>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 2>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 3>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 4>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 5>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 6>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 7>::op ();
  test_modulus_int_preinverse_fixed <C, m, (C) 8>::op ();
  test_modulus_int_preinverse_fixed <C, m, ((C) 1) << (m/2)>::op ();

  // Bug with old version of g++ when -O3 is activated
  static const C x = MMX_SAFE_LEFT_SHIFT_INT(C, 1, m);
  test_modulus_int_preinverse_fixed <C, m, x>::op ();
  test_modulus_int_preinverse_fixed <C, m, x + C(-1)>::op ();
  test_modulus_int_preinverse_fixed <C, m, x + C(-2)>::op ();
  test_modulus_int_preinverse_fixed <C, m, x + C(-3)>::op ();
  test_modulus_int_preinverse_fixed <C, m, x + C(-4)>::op ();
  test_modulus_int_preinverse_fixed <C, m, x + C(-5)>::op ();

  modulus<C,modulus_int_preinverse<m> > dummy;
  mmout << "Modular with " << cpp_demangle (typeid (dummy).name ())
        << " works\n";
}

template<nat nail>
static void test_all_modulus_int_preinverse () {

  test_modulus_int_preinverse<unsigned char, nail> ();
  test_modulus_int_preinverse<short unsigned int, nail> ();
  test_modulus_int_preinverse<unsigned int, nail> ();
  test_modulus_int_preinverse<long unsigned int, nail> ();
  test_modulus_int_preinverse<long long unsigned int, nail> ();

  test_modulus_int_preinverse<signed char, nail> ();
  test_modulus_int_preinverse<short int, nail> ();
  test_modulus_int_preinverse<int, nail> ();
  test_modulus_int_preinverse<long int, nail> ();
  test_modulus_int_preinverse<long long int, nail> ();
}

int
main () {
  test_all_modulus_int_naive<0> ();
  test_all_modulus_int_naive<1> ();
  test_all_modulus_int_naive<2> ();
  test_all_modulus_int_naive<3> ();
  test_all_modulus_int_preinverse<0> ();
  test_all_modulus_int_preinverse<1> ();
  test_all_modulus_int_preinverse<2> ();
  test_all_modulus_int_preinverse<3> ();

  return 0;
}
 All Classes Namespaces Files Functions Variables Typedefs Friends Defines