bernstein_bzenv.hpp

Go to the documentation of this file.

/********************************************************************
 *   This file is part of the source code of the realroot library.
 *   Author(s): J.P. Pavone, GALAAD, INRIA
 *   $Id: bzenv.hpp,v 1.1 2005/07/11 10:03:55 jppavone Exp $
 ********************************************************************/
#ifndef realroot_SOLVE_SBDSLV_BERNSTEIN_BZENV_H
#define realroot_SOLVE_SBDSLV_BERNSTEIN_BZENV_H
//--------------------------------------------------------------------
#include <realroot/loops_vctops.hpp>
#include <realroot/bernstein_binomials.hpp>
//====================================================================
namespace mmx {
//--------------------------------------------------------------------
namespace bernstein
{
  template<class X>
  struct bzenv
  {
    typedef int size_t;
    binomials<X> m_bins;
    static bzenv * const  _default_ ;//= (bzenv*)(&(binomials<X>::_default_));
    
    bzenv( int n = 100 ) : m_bins(n) {};
    
    template< typename T >
    void scale( T * bzc, size_t sz, int st = 1 )
    {
      X * bins = m_bins.get(sz-1);
      int p = 0;
      for ( size_t i = 0; i < sz; i++, p += st ) bzc[p] *= bins[i]; 
    };
    
    template< typename T >
    void uscale( T * bzc, size_t sz, int st = 1 ) 
    {
      X * bins = m_bins.get(sz-1);
      int p = 0;
      for ( int i = 0; i < sz; i++, p+= st ) bzc[p] /= bins[i];
    };
    
    
    /* degree elevation in scaled bernstein basis */
    template < typename T >
    void scaled_elevate( T * dst, T * src, size_t sz, int sts, int stdst, size_t n )
    {
      X * bins = m_bins.get(n);
      vctops::convolution( dst, src, bins, sz, n+1, sts, 1, stdst );
    };
    
    /* degree elevation in bernstein basis        */
    template < typename T >
    void elevate( T * dst, T * src, size_t sz, int sts, int stdst, size_t n )
    {
      scale(src,sz,sts);
      scaled_elevate(dst,src,sz,sts,stdst,n);
      uscale(src,sz,sts);
      uscale(dst,sz+n,stdst);
    };
    
    /* bernstein basis -> monomial basis conversion */
    template< typename T >
    void toMonoms ( T * bzc, size_t sz, int st = 1 )
    {
      int i,k;
      for ( i = st; i != sz*st; i += st )
    for ( k = (sz-1)*st; k != i-st; k -= st )
      bzc[k] -= bzc[k-st];
      scale(bzc,sz,st);
    };
    
    /* monomials basis -> bernstein basis conversion */
    template< typename T > inline
    void fromMonomsToScaled( T * bzc, size_t sz, int st )
    {
      T tmp[sz];
      int i,k,p,l;
      X * bin;
      for ( p = 0, i = 0; i < sz; i++, p += st ) 
    { tmp[i] = bzc[p]; bzc[p] = 0; };
      for ( p = 0, i = 0; i < sz; i++, p += st )
    {
      bin = m_bins.get( sz-i-1 );
      for ( l = p, k = 0; k < sz-i; k++, l += st )
        bzc[l] += tmp[i]*bin[k];
    };
    };
    
    template< typename T >  inline
    void fromMonoms( T * bzc, size_t sz, int st = 1 )
    {
      fromMonomsToScaled( bzc, sz, st );
      uscale(bzc,sz,st);
    };
    X * get( size_t i ) { return m_bins.get(i); };
    const X& binomial( size_t n, size_t i ) { return m_bins.binomial(n,i); };
  };
  template<class X>
  bzenv<X> * const  bzenv<X>::_default_ = (bzenv<X>*)(&(binomials<X>::_default_));
};
//--------------------------------------------------------------------
} //namespace mmx
/********************************************************************/
#endif //