include/algebramix/base_blocks.hpp

Go to the documentation of this file.

/******************************************************************************
* MODULE     : base_blocks.hpp
* DESCRIPTION: Implementation change of bases by blocks
* COPYRIGHT  : (C) 2009  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.
******************************************************************************/

#ifndef __MMX__BASE_BLOCKS__HPP
#define __MMX__BASE_BLOCKS__HPP
#include <algebramix/base_naive.hpp>
#include <algebramix/base_dicho.hpp>

namespace mmx {

/******************************************************************************
* Blocks variant
******************************************************************************/

#define Base_blocks_variant(C) base_blocks_variant_helper<C>::BV

template<typename V>
struct base_blocks : public V {};

template<typename F, typename V, typename W>
struct implementation<F,V,base_blocks<W> >:
  public implementation<F,V,W> {};

template<typename C>
struct base_blocks_variant_helper {
  typedef base_blocks<typename base_naive_variant_helper<C>::BV> BV;
};

/******************************************************************************
* Direct and inverse transforms
******************************************************************************/

template<typename V, typename W>
struct implementation<base_transform,V,base_blocks<W> > {

  template<typename C, typename I, typename WL, typename WH> static inline nat
  direct (I* c, nat n, const C& a, WL* low, WH* high, nat s) {
    typedef typename WH::modulus_base J;
    if (a == 0) return 0;
    nat m= (n + s - 1) / s;
    J* aux= mmx_new<J> (m);
    nat oh= direct_base (aux, m, a, * high), ol= 0;
    for (nat i= 0; i < oh; i++, n -= s) {
      ol= direct_base (c + s * i, min (s, n), aux[i], * low);
      for (nat j= ol; j < min (s, n); j++) c[s * i + j]= 0;
    }
    mmx_delete<J> (aux, m);
    return oh == 0 ? 0 : (oh - 1) * s + ol; }

  template<typename C, typename I, typename WL, typename WH> static inline void
  inverse (C& a, const I* c, nat n, WL* low, WH* high, nat s) {
    typedef typename WH::modulus_base J;
    nat m= (n + s - 1) / s;
    J* aux= mmx_new<J> (m);
    for (nat i= 0; i < m; i++, n -= s) 
      inverse_base (aux[i], c + s * i, min (s, n), * low);
    inverse_base (a, aux, m, * high);
    mmx_delete<J> (aux, m); }

};

/******************************************************************************
* Default threshold
******************************************************************************/

struct base_blocks_threshold {};

template<typename C>
struct threshold_helper<C,base_blocks_threshold> {
  typedef fixed_value<nat,32> impl;
};

/******************************************************************************
* The blocks transformer class
******************************************************************************/

template<typename WL, typename WH,
         typename V=typename Base_blocks_variant(typename WH::base) >
struct base_blocks_transformer {
  typedef typename WH::base base;
  typedef typename WL::modulus_base modulus_base;
  typedef typename WL::modulus_base_variant modulus_base_variant;
  typedef typename WH::modulus_base modulus_middle;
  typedef typename WH::modulus_base_variant modulus_middle_variant;

  typedef base C;
  typedef modulus_base I;
  typedef modulus_middle J;
  typedef modulus<I, modulus_base_variant> M;
  typedef modulus<J, modulus_middle_variant> N;
  typedef implementation<base_transform,V> Base;

  M p;
  nat s;
  WL* low;
  WH* high;

public:

  static const nat th= Threshold(typename WH::base,base_blocks_threshold);
  template<typename K> 
  inline base_blocks_transformer (const K& _p, nat _s= th)
    : p(_p), s(_s) {
    ASSERT (_p != 0 && _p != 1, "invalid base");
    low= mmx_new<WL> (1, p);
    J p_s= binpow (J(* p), s);
    high= mmx_new<WH> (1, p_s); }

  inline base_blocks_transformer (WL* _low, nat _s= th)
    : p(_low->p), s(_s), low(_low) {
    J p_s= binpow (J(* p), s);
    high= mmx_new<WH> (1, p_s); }

  inline ~base_blocks_transformer () {
    mmx_delete<WL> (low, 1);
    mmx_delete<WH> (high, 1); }

  inline nat direct_transform (I* c, nat n, const C& a) {
    return Base::direct (c, n, a, low, high, s); }

  inline void inverse_transform (C& a, const I* c, nat n) {
    Base::inverse (a, c, n, low, high, s); }
};

} // namespace mmx
#endif //__MMX__BASE_BLOCKS__HPP

---