realroot_doc 0.1.1
Functions
/Users/mourrain/Devel/mmx/realroot/include/realroot/bernstein_eenv_multiply.hpp File Reference

Go to the source code of this file.

Functions

Function Documentation

static void _mvrcvloop_ ( bernstein::eenv *  oenv,
real_t *  dst,
bernstein::eenv *  aenv,
real_t *  sra,
bernstein::eenv *  benv,
real_t *  srb,
unsigned *  asupp = 0,
unsigned  nas = 0,
unsigned *  oasup = 0,
unsigned *  bsupp = 0,
unsigned  nbs = 0,
unsigned *  obsup = 0 
) [inline, static]

operations sous forme monomiale

Definition at line 90 of file bernstein_eenv_multiply.hpp.

References assert, and oaddress().

Referenced by mcrossp().

{
  /* calcul des addresses de sorties            */
  unsigned * loasup, * lobsup; loasup = lobsup = 0;
  /* calcul des addresses de sorties de a sur o */
  if ( !oasup ) { if ( !asupp ) nas = aenv->data_size(); loasup = new unsigned[ nas ]; oasup = loasup; oaddress( oenv, oasup, aenv, asupp, nas ); };
  /* calcul des addresses de sorties de b sur o */
  if ( !obsup ) { if ( !bsupp ) nbs = benv->data_size(); lobsup = new unsigned[ nbs ]; obsup = lobsup; oaddress( oenv, obsup, benv, bsupp, nas ); };

  unsigned ia, ib, out; real_t   ca;
  switch( (asupp!=0)<<1|(bsupp!=0) )
    {
    case 0:
      assert(nas==aenv->data_size() && nbs == benv->data_size()); 
      if ( op == '+' )
        {
          for ( ia = 0; ia < nas; ia ++ )
            for ( ib = 0; ib < nbs;  ib ++ )
              dst[oasup[ia]+obsup[ib]] += sra[ia]*srb[ib];
          break;
        }
      if ( op == '-' )
        {
          for ( ia = 0; ia < nas; ia ++ )
            for ( ib = 0; ib < nbs; out += obsup[ib], ib ++ )
              dst[oasup[ia]+obsup[ib]] -= sra[ia]*srb[ib];
          break;
        };
      break;
    case 1:
      assert(nas == aenv->data_size() && bsupp != 0);
      if ( op == '+' )
        {
          assert(nas==aenv->data_size() && nbs == benv->data_size());
          for ( ia = 0; ia < nas; ia ++ )
            for (  ib = 0; ib < nbs; ib ++ )
              dst[oasup[ia]+obsup[ib]] += sra[ia]*srb[bsupp[ib]];
          break;
        }
      if ( op == '-' )
        {
          for ( ia = 0; ia < nas; ia ++ )
            for ( ib = 0; ib < nbs; ib ++ )
              dst[oasup[ia]+obsup[ib]] -= sra[ia]*srb[bsupp[ib]];
          break;
        };
      break;
    case 2:
      assert(nbs == benv->data_size() && asupp != 0);
      if ( op == '+' )
        {
          for ( ia = 0; ia < nas; ia ++ )
            for ( ib = 0; ib < nbs; ib ++ )
              dst[oasup[ia]+obsup[ib]] += sra[asupp[ia]]*srb[ib];
          break;
        }
      if ( op == '-' )
        {
          for ( ia = 0; ia < nas; ia ++ )
            for ( ib = 0; ib < nbs;  ib ++ )
              dst[oasup[ia]+obsup[ib]] -= sra[asupp[ia]]*srb[ib];
          break;
        };
      break;
    case 3:
      assert(asupp != 0 && bsupp != 0);
      if ( op == '+' )
        {
          for ( ia = 0; ia < nas; ia ++ )
            for ( ib = 0; ib < nbs;  ib ++ )
              dst[oasup[ia]+obsup[ib]] += sra[asupp[ia]]*srb[bsupp[ib]];
          break;
        }
      if ( op == '-' )
        {
          assert(nas==aenv->data_size() && nbs == benv->data_size());
          for ( ia = 0; ia < nas; ia ++ )
            for (  ib = 0; ib < nbs;  ib ++ )
              dst[oasup[ia]+obsup[ib]] -= sra[asupp[ia]]*srb[asupp[ib]];
          break;
        };
      break;
    };
  /* liberation des addresses de sorties   */
  if ( loasup ) delete[] loasup; if ( lobsup ) delete[] lobsup;
};
real_t density ( real_t *  src,
const real_t  prec = 1e-12,
unsigned  nsmp = 0 
)

Definition at line 8 of file bernstein_eenv_multiply.hpp.

References mmx::abs().

{
  unsigned c = 0;
  if ( nsmp == 0 )
    {
      nsmp = data_size();
      for ( unsigned i = 0; i < nsmp; i ++ )
        if ( std::abs(src[i]) > prec ) c++;
    }
  else
    {
      for ( unsigned i = 0; i < nsmp; i ++ )
        if ( std::abs(src[rand()%data_size()] > prec) ) c++;
    };
  return ((real_t)c)/nsmp;
};
static void mcrossp ( eenv *  oenv,
real_t **  dst,
eenv *  aenv,
real_t **  src_a,
eenv *  benv,
real_t **  src_b,
bool  clear = true,
unsigned *  asupp = 0,
unsigned  nas = 0,
unsigned *  oasup = 0,
unsigned *  bsupp = 0,
unsigned  nbs = 0,
unsigned *  obsup = 0 
) [static]

Definition at line 265 of file bernstein_eenv_multiply.hpp.

References _mvrcvloop_(), mmx::tensor::clear(), mmx::vct::fill(), and oaddress().

Referenced by mcrossp(), and monoms_crossprod().

{
  if ( clear ) for (int d=0;d<3;std::fill(dst[d],dst[d]+oenv->data_size(),(real_t)0.0),d++);
  unsigned * loasup, * lobsup; loasup = lobsup = 0;
  /* calcul des addresses de sorties de a sur o */
  if ( !oasup ) { if ( !asupp ) nas = aenv->data_size(); loasup = new unsigned[ nas ]; oasup = loasup; oaddress( oenv, oasup, aenv, asupp, nas ); };
  /* calcul des addresses de sorties de b sur o */
  if ( !obsup ) { if ( !bsupp ) nbs = benv->data_size(); lobsup = new unsigned[ nbs ]; obsup = lobsup; oaddress( oenv, obsup, benv, bsupp, nbs ); };
  _mvrcvloop_<real_t,'+'>(oenv,dst[2],aenv,src_a[0],benv,src_b[1],asupp,nas,oasup,bsupp,nbs,obsup);
  _mvrcvloop_<real_t,'-'>(oenv,dst[2],aenv,src_a[1],benv,src_b[0],asupp,nas,oasup,bsupp,nbs,obsup);
  _mvrcvloop_<real_t,'+'>(oenv,dst[0],aenv,src_a[1],benv,src_b[2],asupp,nas,oasup,bsupp,nbs,obsup);
  _mvrcvloop_<real_t,'-'>(oenv,dst[0],aenv,src_a[2],benv,src_b[1],asupp,nas,oasup,bsupp,nbs,obsup);
  _mvrcvloop_<real_t,'+'>(oenv,dst[1],aenv,src_a[2],benv,src_b[0],asupp,nas,oasup,bsupp,nbs,obsup);
  _mvrcvloop_<real_t,'-'>(oenv,dst[1],aenv,src_a[0],benv,src_b[2],asupp,nas,oasup,bsupp,nbs,obsup);
  if ( loasup ) delete[] loasup; if ( lobsup ) delete[] lobsup;
};
static void mcrossp ( eenv *  res,
eenv *  a,
eenv *  b,
real_t **  dst,
real_t **  src_a,
real_t **  src_b 
) [static]

Definition at line 313 of file bernstein_eenv_multiply.hpp.

References mcrossp().

{ mcrossp( res, dst, a, src_a, b, src_b ); };
static void mdotp ( eenv *  oenv,
real_t *  dst,
eenv *  aenv,
real_t **  src_a,
eenv *  benv,
real_t **  src_b,
int  n,
bool  clear = true,
unsigned *  asupp = 0,
unsigned  nas = 0,
unsigned *  oasup = 0,
unsigned *  bsupp = 0,
unsigned  nbs = 0,
unsigned *  obsup = 0 
) [static]

Definition at line 328 of file bernstein_eenv_multiply.hpp.

References mmx::tensor::clear(), mmx::vct::fill(), and oaddress().

Referenced by monoms_dotprod().

{
  if ( clear ) std::fill(dst,dst+oenv->data_size(),(real_t)0.0);
  unsigned * loasup, * lobsup; loasup = lobsup = 0;
  /* calcul des addresses de sorties de a sur o */
  if ( !oasup ) { if ( !asupp ) nas = aenv->data_size(); loasup = new unsigned[ nas ]; oasup = loasup; oaddress( oenv, oasup, aenv, asupp, nas ); };
  /* calcul des addresses de sorties de b sur o */
  if ( !obsup ) { if ( !bsupp ) nbs = benv->data_size(); lobsup = new unsigned[ nbs ]; obsup = lobsup; oaddress( oenv, obsup, benv, bsupp, nas ); };
  for ( int i = 0; i < n; i ++ )
    _mvrcvloop_<real_t,'+'>(oenv,dst,aenv,src_a[i],benv,src_b[i],asupp,nas,oasup,bsupp,nbs,obsup);
  if ( loasup ) delete[] loasup; if ( lobsup ) delete[] lobsup;
};
static void mmul ( eenv *  oenv,
real_t *  dst,
eenv *  aenv,
real_t *  src_a,
eenv *  benv,
real_t *  src_b,
bool  clear = true,
unsigned *  asupp = 0,
unsigned  nas = 0,
unsigned *  oasupp = 0,
unsigned *  bsupp = 0,
unsigned  nbs = 0,
unsigned *  obsupp = 0 
) [inline, static]

Definition at line 231 of file bernstein_eenv_multiply.hpp.

References mmx::tensor::clear(), and pmmul().

Referenced by spmmul().

{ pmmul( oenv, &dst, aenv, &src_a, benv, &src_b, 1, clear, asupp, nas, oasupp, bsupp, nbs, obsupp ); };
static void monoms_crossprod ( eenv *  res,
eenv *  a,
eenv *  b,
real_t **  dst,
real_t **  src_a,
real_t **  src_b 
) [static]

Definition at line 317 of file bernstein_eenv_multiply.hpp.

References mcrossp().

{
  res->op_mul(a,b);
  for ( int d = 0; d < 3; d ++ ) dst[d] = new real_t[ res->data_size() ];
  mcrossp( res, dst, a, src_a, b, src_b );
};
static void monoms_dotprod ( eenv *  res,
eenv *  a,
eenv *  b,
real_t *&  dst,
real_t **  src_a,
real_t **  src_b,
int  n 
) [static]

Definition at line 349 of file bernstein_eenv_multiply.hpp.

References mdotp().

{
  res->op_mul(a,b);
  dst = new real_t[res->data_size()];
  mdotp( res, dst, a, src_a, b, src_b, n );
};
static void monoms_multiply ( bernstein::eenv *  res,
bernstein::eenv *  a,
bernstein::eenv *  b,
real_t *  dst,
real_t *  src_a,
real_t *  src_b 
) [static]

Definition at line 206 of file bernstein_eenv_multiply.hpp.

References monoms_multiply_loop().

{ res->op_mul(a,b); monoms_multiply_loop( res, a, b, dst, src_a, src_b  ); };
static void monoms_multiply_loop ( eenv *  o,
eenv *  a,
eenv *  b,
real_t *  dst,
real_t *  src_a,
real_t *  src_b 
) [static]

Definition at line 202 of file bernstein_eenv_multiply.hpp.

Referenced by monoms_multiply().

{ _mvrcvloop_<real_t,'+'>( o, dst, a, src_a, b, src_b ); };
static void oaddress ( bernstein::eenv *  oenv,
unsigned *  osupp,
bernstein::eenv *  ienv,
unsigned *  isupp = 0,
unsigned  nsp = 0 
) [static]

Definition at line 66 of file bernstein_eenv_multiply.hpp.

References vmap().

Referenced by _mvrcvloop_(), mcrossp(), mdotp(), msimplify(), pmmul(), and vrestrict().

{
  int v;
  unsigned c;
  sz_t vmap_[ ienv->m_nvr ];
  unsigned addi,addo;
  vmap( vmap_, oenv, ienv );
  if ( isupp )
    for ( c = 0; c < nsp; osupp[c] = addo, c ++ )
      for ( addi = isupp[c], addo = 0, v = ienv->m_nvr-1; addi; addi /= ienv->m_szs[v], v-- )
        addo += (addi%ienv->m_szs[v])*oenv->m_str[vmap_[v]];
  else
    {
      nsp = ienv->data_size();
    for ( c = 0; c < nsp; osupp[c] = addo, c ++ )
      for ( addi = c, addo = 0, v = ienv->m_nvr-1; addi; addi /= ienv->m_szs[v], v -- )
        addo += (addi%ienv->m_szs[v])*oenv->m_str[vmap_[v]];
    };
};
void op_mul ( eenv *  a,
eenv *  b 
)

Definition at line 37 of file bernstein_eenv_multiply.hpp.

References mmx::sparse::copy(), mmx::vctops::max(), and mmx::vctops::set_conversion().

{
  if ( nvars() ) this->~eenv();
  _alloc_( a->nvars()+b->nvars() );
  std::copy( a->m_vrs, a->m_vrs+a->nvars(), m_vrs );
  std::copy( b->m_vrs, b->m_vrs+b->nvars(), m_vrs + a->nvars() );
  m_nvr = vctops::set_conversion( m_vrs, m_nvr );
  for ( sz_t v = 0; v < m_nvr; v ++ )
    {
      int av = a->indexofvar(var(v));
      int bv = b->indexofvar(var(v));
      m_szs[ v ] = ((av!=-1)?(a->sz(av)-1):0) + ((bv!=-1)?(b->sz(bv)-1):0) + 1;
    };
  m_mxvr = std::max(a->m_mxvr,b->m_mxvr);
  compute_strides();
};
static void pmmul ( eenv *  oenv,
real_t **  dst,
eenv *  aenv,
real_t **  src_a,
eenv *  benv,
real_t **  src_b,
int  n,
bool  clear = true,
unsigned *  asupp = 0,
unsigned  nas = 0,
unsigned *  oasup = 0,
unsigned *  bsupp = 0,
unsigned  nbs = 0,
unsigned *  obsup = 0 
) [static]

Definition at line 212 of file bernstein_eenv_multiply.hpp.

References mmx::tensor::clear(), mmx::vct::fill(), and oaddress().

Referenced by mmul().

{
  if ( clear ) for (int d=0;d<n;std::fill(dst[d],dst[d]+oenv->data_size(),(real_t)0.0),d++);
  unsigned * loasup, * lobsup; loasup = lobsup = 0;
  /* calcul des addresses de sorties de a sur o */
  if ( !oasup ) { if ( !asupp ) nas = aenv->data_size(); loasup = new unsigned[ nas ]; oasup = loasup; oaddress( oenv, oasup, aenv, asupp, nas ); };
  /* calcul des addresses de sorties de b sur o */
  if ( !obsup ) { if ( !bsupp ) nbs = benv->data_size(); lobsup = new unsigned[ nbs ]; obsup = lobsup; oaddress( oenv, obsup, benv, bsupp, nbs ); };
  for ( unsigned i = 0; i < n; i ++ )
    _mvrcvloop_<real_t,'+'>(oenv,dst[i],aenv,src_a[i],benv,src_b[i],asupp,nas,oasup,bsupp,nbs,obsup);
  if ( loasup ) delete[] loasup; if ( lobsup ) delete[] lobsup;
};
static void rvbinoms ( bernstein::eenv *  ienv,
real_t *  bcff,
unsigned *  isupp = 0,
unsigned  nsp = 0,
bernstein::bzenv< real_t > *  bznv = bzenv<real_t>::_default_ 
) [static]

operations sous forme de bernstein

Definition at line 360 of file bernstein_eenv_multiply.hpp.

{
  real_t cff;
  unsigned c,addi;
  int v;
  if ( isupp )
    for ( c = 0; c < nsp; bcff[c] = cff, c ++ )
      for ( cff = (real_t)1.0, addi = isupp[c],  v = ienv->m_nvr-1; addi; addi /= ienv->m_szs[v], v-- )
        cff /= bznv->binomial(ienv->m_szs[v]-1,addi % ienv->m_szs[v]);
  else
    for ( c = 0; c < ienv->data_size(); bcff[c] = cff,  c ++ )
      for ( addi = c, v = ienv->m_nvr-1; addi; addi /= ienv->m_szs[v], v -- )
        cff /= bznv->binomial(ienv->m_szs[v]-1,addi % ienv->m_szs[v]);
};
static void spmmul ( eenv *  oenv,
real_t *  dst,
eenv *  aenv,
real_t *  src_a,
eenv *  benv,
real_t *  src_b,
bool  clear = true,
bool  stats = false 
) [static]

Definition at line 258 of file bernstein_eenv_multiply.hpp.

References mmx::tensor::clear(), and spmmul().

{
  spmmul( oenv, &dst, aenv, &src_a, benv, &src_b, 1, clear, stats );
};
static void spmmul ( eenv *  oenv,
real_t **  dst,
eenv *  aenv,
real_t **  src_a,
eenv *  benv,
real_t **  src_b,
int  n,
bool  clear = true,
bool  stats = false 
) [static]

Definition at line 240 of file bernstein_eenv_multiply.hpp.

References assert, and mmul().

Referenced by spmmul().

{
  unsigned nas,nbs;
  unsigned * assup = new unsigned[ aenv->data_size() ];
  unsigned * bssup = new unsigned[ benv->data_size() ];
  assert( assup && bssup );
  for ( int i = 0; i < n; i ++ )
    {
      nas = aenv->support(assup,src_a[i]);
      nbs = benv->support(bssup,src_b[i]);
      mmul( oenv, dst[i], aenv, src_a[i], benv, src_b[i], true, assup, nas, 0, bssup, nbs, 0 );
    };
  delete[] assup;
  delete[] bssup;
};
unsigned support ( unsigned *  supp,
real_t *  src,
const real_t &  prec = 1e-6 
)

Definition at line 27 of file bernstein_eenv_multiply.hpp.

References mmx::abs().

{
//  supp = new unsigned[ data_size() ];
  unsigned      c = 0;
  for ( unsigned i = 0; i < data_size(); i ++ )
    if ( std::abs(src[i]) > prec ) supp[c++] = i;
  return c;
};
static void vmap ( sz_t *  vmap,
bernstein::eenv *  o,
bernstein::eenv *  i 
) [static]

Definition at line 55 of file bernstein_eenv_multiply.hpp.

References assert.

Referenced by oaddress().

{
  for ( sz_t v = 0; v < i->m_nvr; v ++ )
    {
      vmap[v] = -1;
      for ( sz_t k = 0; k < o->m_nvr; k ++ )
        if ( o->m_vrs[k] == i->m_vrs[v] ) vmap[v] = k;
      assert(vmap[v]!=-1);
    };
};