00001
00002
00003
00004
00005
00006 #ifndef realroot_SOLVE_SBDSLV_BERNSTEIN_BZENV_H
00007 #define realroot_SOLVE_SBDSLV_BERNSTEIN_BZENV_H
00008
00009 #include <realroot/loops_vctops.hpp>
00010 #include <realroot/bernstein_binomials.hpp>
00011
00012 namespace mmx {
00013
00014 namespace bernstein
00015 {
00016 template<class X>
00017 struct bzenv
00018 {
00019 typedef int size_t;
00020 binomials<X> m_bins;
00021 static bzenv * const _default_ ;
00022
00023 bzenv( int n = 100 ) : m_bins(n) {};
00024
00026 template< typename T >
00027 void scale( T * bzc, size_t sz, int st = 1 )
00028 {
00029 X * bins = m_bins.get(sz-1);
00030 int p = 0;
00031 for ( size_t i = 0; i < sz; i++, p += st ) bzc[p] *= bins[i];
00032 };
00033
00035 template< typename T >
00036 void uscale( T * bzc, size_t sz, int st = 1 )
00037 {
00038 X * bins = m_bins.get(sz-1);
00039 int p = 0;
00040 for ( int i = 0; i < sz; i++, p+= st ) bzc[p] /= bins[i];
00041 };
00042
00043
00044
00045 template < typename T >
00046 void scaled_elevate( T * dst, T * src, size_t sz, int sts, int stdst, size_t n )
00047 {
00048 X * bins = m_bins.get(n);
00049 vctops::convolution( dst, src, bins, sz, n+1, sts, 1, stdst );
00050 };
00051
00052
00053 template < typename T >
00054 void elevate( T * dst, T * src, size_t sz, int sts, int stdst, size_t n )
00055 {
00056 scale(src,sz,sts);
00057 scaled_elevate(dst,src,sz,sts,stdst,n);
00058 uscale(src,sz,sts);
00059 uscale(dst,sz+n,stdst);
00060 };
00061
00062
00063 template< typename T >
00064 void toMonoms ( T * bzc, size_t sz, int st = 1 )
00065 {
00066 int i,k;
00067 for ( i = st; i != sz*st; i += st )
00068 for ( k = (sz-1)*st; k != i-st; k -= st )
00069 bzc[k] -= bzc[k-st];
00070 scale(bzc,sz,st);
00071 };
00072
00073
00074 template< typename T > inline
00075 void fromMonomsToScaled( T * bzc, size_t sz, int st )
00076 {
00077 T tmp[sz];
00078 int i,k,p,l;
00079 X * bin;
00080 for ( p = 0, i = 0; i < sz; i++, p += st )
00081 { tmp[i] = bzc[p]; bzc[p] = 0; };
00082 for ( p = 0, i = 0; i < sz; i++, p += st )
00083 {
00084 bin = m_bins.get( sz-i-1 );
00085 for ( l = p, k = 0; k < sz-i; k++, l += st )
00086 bzc[l] += tmp[i]*bin[k];
00087 };
00088 };
00089
00090 template< typename T > inline
00091 void fromMonoms( T * bzc, size_t sz, int st = 1 )
00092 {
00093 fromMonomsToScaled( bzc, sz, st );
00094 uscale(bzc,sz,st);
00095 };
00096 X * get( size_t i ) { return m_bins.get(i); };
00097 const X& binomial( size_t n, size_t i ) { return m_bins.binomial(n,i); };
00098 };
00099 template<class X>
00100 bzenv<X> * const bzenv<X>::_default_ = (bzenv<X>*)(&(binomials<X>::_default_));
00101 };
00102
00103 }
00104
00105 #endif //