/Users/mourrain/Devel/mmx/algebramix/include/algebramix/polynomial_tangent.hpp

Go to the documentation of this file.

/******************************************************************************
* MODULE     : polynomial_tangent.hpp
* DESCRIPTION: Multiplication of polynomials with order one jet coefficients
* COPYRIGHT  : (C) 2008  Joris van der Hoeven
*******************************************************************************
* 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__POLYNOMIAL_TANGENT__HPP
#define __MMX__POLYNOMIAL_TANGENT__HPP
#include <algebramix/polynomial.hpp>
#include <algebramix/polynomial_dicho.hpp>
#include <numerix/tangent.hpp>

namespace mmx {
#define TMPL template<typename C>

/******************************************************************************
* Default variant for tangent coefficients
******************************************************************************/

template<typename V>
struct polynomial_tangent {
  typedef typename V::Vec Vec;
  typedef typename V::Naive Naive;
  typedef polynomial_tangent<typename V::Positive> Positive;
  typedef polynomial_tangent<typename V::No_simd> No_simd;
  typedef polynomial_tangent<typename V::No_thread> No_thread;
  typedef polynomial_tangent<typename V::No_scaled> No_scaled;
};

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

template<typename C>
struct polynomial_variant_helper<tangent<C> > {
  typedef typename Polynomial_variant(C) CPV;
  typedef polynomial_tangent<CPV> PV;
};

/******************************************************************************
* Tangent multiplication
******************************************************************************/

template<typename V, typename CV>
struct implementation<polynomial_multiply,V,polynomial_tangent<CV> >:
  public implementation<polynomial_linear,V>
{
  typedef implementation<vector_linear,CV> Vec_C;
  typedef implementation<polynomial_multiply,CV> Pol_C;

TMPL static inline void
mul (tangent<C>* dest, const tangent<C>* s1, const tangent<C>* s2,
     nat n1, nat n2)
{
  nat nd     = n1 + n2 - 1;
  nat spc1   = aligned_size<C,CV> (n1);
  nat spc2   = aligned_size<C,CV> (n2);
  nat spcd   = aligned_size<C,CV> (nd);
  nat spc    = 5 * (spc1 + spc2);
  C* base_x1 = mmx_new<C> (spc);
  C* slope_x1= base_x1  + spc1;
  C* base_x2 = slope_x1 + spc1;
  C* slope_x2= base_x2  + spc2;
  C* base_d  = slope_x2 + spc2;
  C* slope_d = base_d   + spcd;
  C* extra_d = slope_d  + spcd;

  Vec_C::half_copy (base_x1, (C*) ((void*) s1), n1);
  Vec_C::half_copy (slope_x1, ((C*) ((void*) s1)) + 1, n1);
  Vec_C::half_copy (base_x2, (C*) ((void*) s2), n2);
  Vec_C::half_copy (slope_x2, ((C*) ((void*) s2)) + 1, n2);
  Pol_C::mul (base_d, base_x1, base_x2, n1, n2);
  Pol_C::mul (slope_d, base_x1, slope_x2, n1, n2);
  Pol_C::mul (extra_d, slope_x1, base_x2, n1, n2);
  Pol_C::add (slope_d, extra_d, nd);
  Vec_C::double_copy ((C*) ((void*) dest), base_d, nd);
  Vec_C::double_copy (((C*) ((void*) dest)) + 1, slope_d, nd);

  mmx_delete<C> (base_x1, spc);
}

TMPL static inline void
square (tangent<C>* dest, const tangent<C>* s1, nat n1) {
  nat nd     = n1 + n1 - 1;
  nat spc1   = aligned_size<C,CV> (n1);
  nat spcd   = aligned_size<C,CV> (nd);
  nat spc    = 6 * spc1;
  C* base_x1 = mmx_new<C> (spc);
  C* slope_x1= base_x1  + spc1;
  C* base_d  = slope_x1 + spc1;
  C* slope_d = base_d   + spcd;
  
  Vec_C::half_copy (base_x1, (C*) ((void*) s1), n1);
  Vec_C::half_copy (slope_x1, ((C*) ((void*) s1)) + 1, n1);
  Pol_C::square (base_d, base_x1, n1);
  Pol_C::mul (slope_d, base_x1, slope_x1, n1, n1);
  Pol_C::add (slope_d, slope_d, nd);
  Vec_C::double_copy ((C*) ((void*) dest), base_d, nd);
  Vec_C::double_copy (((C*) ((void*) dest)) + 1, slope_d, nd);

  mmx_delete<C> (base_x1, spc);
}

TMPL static inline void
tmul (C* dest, const C* s1, const C* s2, nat n1, nat n2) {
  (void) dest; (void) s1; (void) s2; (void) n1; (void) n2;
  ERROR ("transposed tangent multiplication not implemented");
}

}; // implementation<polynomial_multiply,V,polynomial_tangent<CV> >

#undef TMPL
} // namespace mmx
#endif //__MMX__POLYNOMIAL_TANGENT__HPP