Polynomial Class Reference

#include <Poly.h>

Inheritance diagram for Polynomial:

Collaboration diagram for Polynomial:

Classes
struct	PolynomialDataTag

Public Member Functions
	Polynomial ()
	Polynomial (unsigned Dim, unsigned deg=0, double *data=0)
	Polynomial (unsigned Dim, double val)
	Polynomial (MultInd &)
	Polynomial (char *name)
unsigned	dim ()
unsigned	deg ()
unsigned	sz ()
	operator double * () const
	~Polynomial ()
	Polynomial (const Polynomial &A)
Polynomial &	operator= (const Polynomial &A)
Polynomial	clone ()
void	copyFrom (Polynomial a)
Polynomial	operator* (const double)
Polynomial	operator/ (const double)
Polynomial	operator+ (Polynomial)
Polynomial	operator- (Polynomial)
Polynomial	operator- (void)
Polynomial	operator+ (void)
Polynomial	operator+= (Polynomial)
Polynomial	operator-= (Polynomial)
Polynomial	operator*= (const double)
Polynomial	operator/= (const double)
double	shiftedEval (Vector Point, double minusVal)
double	operator() (Vector)
Polynomial	derivate (int i)
void	gradient (Vector P, Vector G)
void	gradientHessian (Vector P, Vector G, Matrix H)
void	translate (Vector translation)
int	operator== (const Polynomial q)
int	equals (Polynomial q)
void	print ()
void	save (char *name)
void	setFlag (unsigned int val)
void	unsetFlag (unsigned int val)
unsigned	queryFlag (unsigned int val)

Static Public Attributes
static const unsigned int	NicePrint = 1
static const unsigned int	Warning = 2
static const unsigned int	Normalized = 4
static unsigned int	flags = Polynomial::Warning||Polynomial::NicePrint
static Polynomial	emptyPolynomial

Protected Types
typedef struct Polynomial::PolynomialDataTag	PolynomialData

Protected Member Functions
void	init (int _dim, int _deg, double *data=NULL)
void	destroyCurrentBuffer ()

Protected Attributes
PolynomialData *	d

Detailed Description

Definition at line 38 of file Poly.h.

Member Typedef Documentation

PolynomialData

typedef struct Polynomial::PolynomialDataTag Polynomial::PolynomialData

protected

Constructor & Destructor Documentation

Polynomial() [1/6]

Polynomial::Polynomial ( )

inline

Definition at line 54 of file Poly.h.

{ init(0,0); };

Polynomial() [2/6]

Polynomial::Polynomial	(	unsigned	Dim,
		unsigned	deg = 0,
		double *	data = 0
	)

Definition at line 66 of file Poly.cpp.

{
    init(Dim,Deg,data);
}

Polynomial() [3/6]

Polynomial::Polynomial	(	unsigned	Dim,
		double	val
	)

Definition at line 71 of file Poly.cpp.

{ 
    init(Dim,0,&val); 
}

Polynomial() [4/6]

Polynomial::Polynomial

(

MultInd &

I

)

Definition at line 76 of file Poly.cpp.

{
    init(I.dim,I.len());
    d->coeff[I.index()] = 1;
}

Polynomial() [5/6]

Polynomial::Polynomial

(

char *

name

)

Definition at line 588 of file Poly.cpp.

{
    unsigned _dim,_deg;
    std::ifstream ifp(name, std::ios::in | std::ios::binary);
    ifp.read(reinterpret_cast<char*>(&_dim), sizeof(int));
    ifp.read(reinterpret_cast<char*>(&_deg), sizeof(int));
    init(_dim,_deg);
    ifp.read(reinterpret_cast<char*>(d->coeff), d->n*sizeof(double));
    ifp.close();
}

~Polynomial()

Polynomial::~Polynomial

(

)

Definition at line 82 of file Poly.cpp.

{
    destroyCurrentBuffer();
}

Polynomial() [6/6]

Polynomial::Polynomial

(

const Polynomial &

A

)

Definition at line 110 of file Poly.cpp.

{
    // shallow copy
    d=A.d;
    (d->ref_count)++ ;
}

Member Function Documentation

clone()

Polynomial Polynomial::clone

(

)

Definition at line 117 of file Poly.cpp.

{
    // a deep copy
    Polynomial m(d->dim,d->deg);
    m.copyFrom(*this);
    return m;
}

copyFrom()

void Polynomial::copyFrom

(

Polynomial

a

)

Definition at line 125 of file Poly.cpp.

{
    if (m.d->dim!=d->dim)
    {
        printf("poly: copyFrom: dim do not agree");
        getchar(); exit(254);
    }

    d->deg = mmax(d->deg,m.d->deg); // New degree
    unsigned N1=sz(), N2=m.sz();
    if (N1!=N2) 
    {
        d->coeff=(double*)realloc(d->coeff,N2*sizeof(double));
        d->n=m.d->n;
    }
    memcpy((*this),m,N2*sizeof(double));
}

deg()

unsigned Polynomial::deg ( )

inline

Definition at line 63 of file Poly.h.

{ return d->deg; };

derivate()

Polynomial Polynomial::derivate

(

int

i

)

Definition at line 469 of file Poly.cpp.

{
    unsigned dim=d->dim, deg=d->deg;
    if (deg<1) return Polynomial(dim,0.0);

    Polynomial r(dim, deg-1);
    MultInd I( dim );
    MultInd J( dim );
    double *tS=(*this), *tD=r;
    unsigned j=sz(), k, *cc, sum,
             *allExpo=(unsigned*)I, *expo=allExpo+i, *firstOfJ=(unsigned*)J;

    while (j--) 
    { 
        if (*expo)
        {
            (*expo)--;

            sum=0; cc=allExpo; k=dim;
            while (k--) sum+=*(cc++);
            if (sum) k=choose( sum-1+dim, dim ); else k=0;
            J.resetCounter(); *firstOfJ=sum;
            while (!(J==I)) { k++; J++; }

            (*expo)++;
            tD[k]=(*tS) * (double)*expo;
        }
        tS++;
        I++;
    }
    return r;
}

destroyCurrentBuffer()

void Polynomial::destroyCurrentBuffer ( )

protected

Definition at line 87 of file Poly.cpp.

{
    if (!d) return;
    (d->ref_count) --;
    if (d->ref_count==0)
    {
        if (d->coeff) free(d->coeff);
        free(d);
    }
}

dim()

unsigned Polynomial::dim ( )

inline

Definition at line 62 of file Poly.h.

{ return d->dim; };

equals()

int Polynomial::equals

(

Polynomial

q

)

Definition at line 305 of file Poly.cpp.

{
    if (d==q.d) return 1;
    if ( (d->deg != q.d->deg) || (d->dim != q.d->dim) ) return 0;

    unsigned N = sz();
    double  *tp = (*this),*tq = q;

    while (N--)
        if ( *(tp++) != *(tq++) ) return 0;

    return 1;
}

gradient()

void Polynomial::gradient	(	Vector	P,
		Vector	G
	)

Definition at line 502 of file Poly.cpp.

{
    unsigned i=d->dim;
    G.setSize(i);
    double *r=G;
    if (P.equals(Vector::emptyVector))
    {
        memcpy(r,d->coeff+1,i*sizeof(double));
        return;
    }
    while (i--) r[i]=(derivate(i))(P);
}

gradientHessian()

void Polynomial::gradientHessian	(	Vector	P,
		Vector	G,
		Matrix	H
	)

Definition at line 515 of file Poly.cpp.

{
    unsigned dim=d->dim;
    G.setSize(dim);
    H.setSize(dim,dim);
    double *r=G, **h=H;
    unsigned i,j;

    if (d->deg==2)
    {        
        double *c=d->coeff+1;
        memcpy(r,c,dim*sizeof(double));
        c+=dim;
        for (i=0; i<dim; i++)
        {
            h[i][i]=2* *(c++);
            for (j=i+1; j<dim; j++)
                h[i][j]=h[j][i]=*(c++);
        }
        if (P.equals(Vector::emptyVector)) return;
        G+=H.multiply(P);
        return;
    }

    Polynomial *tmp=new Polynomial[dim], a;
    i=dim;
    while (i--)
    {
        tmp[i]=derivate(i);
        r[i]=(tmp[i])(P);
    }

    i=dim;
    while (i--)
    {
        j=i+1;
        while (j--)
        {
            a=tmp[i].derivate(j);
            h[i][j]=h[j][i]=a(P);
        }
    }

//    _CrtCheckMemory(); 

    delete []tmp;
}

init()

void Polynomial::init ( int  _dim, int  _deg, double *  data = NULL  )

protected

Definition at line 43 of file Poly.cpp.

{
    int n;
    d=(PolynomialData*)malloc(sizeof(PolynomialData));
    if (_dim) n=d->n=choose( _dim+_deg, _dim );
    else n=d->n=0;

    d->dim=_dim;
    d->deg=_deg;
    d->ref_count=1;

    if (n==0) { d->coeff=NULL; return; };

    d->coeff=(double*)malloc(n*sizeof(double)); 
    if (d->coeff==NULL) { printf("memory allocation error\n"); getchar(); exit(253); }

    if (data) memcpy(d->coeff, data, d->n*sizeof(double));
    else memset(d->coeff, 0, d->n*sizeof(double));
}

operator double *()

Polynomial::operator double * ( ) const

inline

Definition at line 65 of file Poly.h.

{ return d->coeff; };

operator()()

double Polynomial::operator()

(

Vector

Point

)

Definition at line 393 of file Poly.cpp.

{
// I didn't notice any difference in precision:
//  return simpleEval(P);

  unsigned dim=d->dim, deg=d->deg;
  double r,r0;                             // no static here because of the 2 threads !
  double rbuf[100]; // That should suffice // no static here because of the 2 threads !
  double *rbufp = rbuf;
  unsigned lsize = 100;
  double *rptr;
  int i,j;

  if (Point==Vector::emptyVector) return *d->coeff;
  if (deg==0) return *d->coeff;

  if ( dim != (unsigned)Point.sz() )
  {
    printf( "Polynomial::operator()( Vector& ) : Improper size\n");
    getchar(); exit(-1);
  }

  if ( !sz() )
  {
    if ( queryFlag( Warning ) )
    {
      printf( "Polynomial::operator()( Vector& ) : evaluating void polynomial\n");
    }
    return 0;
  }

  if ( dim > lsize )    // Someone must be crazy !!!
  {
    if ( queryFlag( Warning ) )
    {
    printf( "Polynomial::operator()( Vector& ) : Warning -> 100 variables\n");
    }

    lsize=dim;
    rbufp = (double*)malloc(lsize*sizeof(double));  // So be it ...

    if ( !rbufp )
    {
      printf( "Polynomial::operator()( Vector& ) : Cannot allocate <rbufp>\n");
      getchar(); exit( -1 );
    }
  }

  // Initialize
  MultInd *mic=cacheMultInd.get( dim, deg );
  unsigned *nextI=mic->indexesOfCoefInLexOrder(),
           *lcI=mic->lastChanges();
  double *cc = (*this), *P=Point;
  unsigned nxt, lc;

  // Empty buffer (all registers = 0)
  memset(rbufp,0,dim*sizeof(double));

  r0=cc[*(nextI++)];
  i=sz()-1;
  while (i--)
  {
    nxt= *(nextI++);
    lc = *(lcI++);

    r=r0; rptr=rbufp+lc; j=dim-lc;
    while (j--) { r+=*rptr; *(rptr++)=0; }
    rbufp[lc]=P[lc]*r; 
    r0=cc[nxt];
  }
  r=r0; rptr=rbufp; i=(int)dim;
  while (i--) r+=*(rptr++);

  return r;
}

operator*()

Polynomial Polynomial::operator*

(

const double

t

)

Definition at line 143 of file Poly.cpp.

{
    int i=sz();
    Polynomial q( d->dim, d->deg );
    double *tq = q.d->coeff, *tp = d->coeff;

    while (i--) *(tq++) = *(tp++) * t;
    return q;
}

operator*=()

Polynomial Polynomial::operator*=

(

const double

t

)

Definition at line 280 of file Poly.cpp.

{
    int i=sz();
    double *tp = (*this);

    while (i--) *(tp++) *=t;
    return *this;
}

operator+() [1/2]

Polynomial Polynomial::operator+

(

Polynomial

q

)

Definition at line 169 of file Poly.cpp.

{
    if (d->dim != q.d->dim)
    {
        printf( "Poly::op+ : Different dimension\n");
        getchar(); exit(-1);
    }

    Polynomial r(d->dim,mmax(d->deg,q.d->deg));
    unsigned N1=sz(), N2=q.sz(), Ni=mmin(N1,N2);
    double  *tr = r, *tp = (*this), *tq = q;
    while (Ni--) *(tr++) = *(tp++) + *(tq++);
    if (N1<N2)
    {
        memcpy(tr,tq,(N2-N1)*sizeof(double));
//        N2-=N1; while (N2--) *(tr++)=*(tq++);
    }
    return r;
}

operator+() [2/2]

Polynomial Polynomial::operator+ ( void  )

inline

Definition at line 84 of file Poly.h.

        { return *this; }

operator+=()

Polynomial Polynomial::operator+=

(

Polynomial

p

)

Definition at line 223 of file Poly.cpp.

{
    if (d->dim != p.d->dim)
    {
        printf("Poly::op+= : Different dimension\n");
        getchar(); exit(-1);
    }

    d->deg = mmax(d->deg,p.d->deg); // New degree
    unsigned N1=sz(), N2=p.sz(), Ni=mmin(N1,N2);
    if (N1<N2) 
    {
        d->coeff=(double*)realloc(d->coeff,N2*sizeof(double));
        d->n=p.d->n;
    }
    double *tt = (*this),*tp = p;

    while (Ni--) *(tt++) += *(tp++);
    
    if (N1<N2) 
    {
        memcpy(tt,tp,(N2-N1)*sizeof(double));
//        N2-=N1; while (N2--) *(tt++)=*(tp++);
    }

    return *this;
}

operator-() [1/2]

Polynomial Polynomial::operator-

(

Polynomial

q

)

Definition at line 190 of file Poly.cpp.

{
    if (d->dim != q.d->dim)
    {
        printf("Poly::op- : Different dimension\n");
        getchar(); exit(-1);
    }

    Polynomial r(d->dim,mmax(d->deg,q.d->deg));
    unsigned N1=sz(), N2=q.sz(), Ni=mmin(N1,N2);
    double  *tr = r, *tp = (*this), *tq = q;
    while (Ni--) *(tr++) = *(tp++) - *(tq++);
    if (N1<N2)
    {
        N2-=N1; while (N2--) *(tr++)=-(*(tq++));
    }
    return r;
}

operator-() [2/2]

Polynomial Polynomial::operator-

(

void

)

Definition at line 209 of file Poly.cpp.

{
    unsigned Ni = sz();
    double *tp = (*this);

    if (!Ni || !tp) return *this;   // Take it like it is ...
    
    Polynomial r(d->dim,d->deg);
    double *tq = (r);
    while( Ni-- ) *(tq++) = -(*(tp++));
    return r;
}

operator-=()

Polynomial Polynomial::operator-=

(

Polynomial

p

)

Definition at line 252 of file Poly.cpp.

{
    if (d->dim != p.d->dim)
    {
        printf( "Poly::op-= : Different dimension\n");
        getchar(); exit(-1);
    }

    d->deg = mmax(d->deg,p.d->deg); // New degree
    unsigned N1=sz(), N2=p.sz(), Ni=mmin(N1,N2);
    if (N1<N2) 
    {
        d->coeff=(double*)realloc(d->coeff,N2*sizeof(double));
        d->n=p.d->n;
    }
    double *tt = (*this),*tp = p;

    while (Ni--) *(tt++) -= *(tp++);
    
    if (N1<N2) 
    {
        N2-=N1; while (N2--) *(tt++)=-(*(tp++));
    }

    return *this;
}

operator/()

Polynomial Polynomial::operator/

(

const double

t

)

Definition at line 154 of file Poly.cpp.

{
    if (t == 0)
    {
        printf( "op/(Poly,double): Division by zero\n");
        getchar(); exit(-1);
    }
    int i=sz();
    Polynomial q( d->dim, d->deg );
    double *tq = q.d->coeff, *tp = d->coeff;

    while (i--) *(tq++) = *(tp++) / t;
    return q;
}

operator/=()

Polynomial Polynomial::operator/=

(

const double

t

)

Definition at line 290 of file Poly.cpp.

{
    if (t == 0)
    {
        printf( "Poly::op/= : Division by zero\n");
        getchar(); exit(-1);
    }

    int i=sz();
    double *tp = (*this);

    while (i--) *(tp++) /=t;
    return *this;
}

operator=()

Polynomial & Polynomial::operator=

(

const Polynomial &

A

)

Definition at line 98 of file Poly.cpp.

{
    // shallow copy
    if (this != &A)
    {
        destroyCurrentBuffer();
        d=A.d;
        (d->ref_count) ++ ;
    }
    return *this;
}

operator==()

int Polynomial::operator== ( const Polynomial  q )

inline

Definition at line 106 of file Poly.h.

{ return d==q.d; };

print()

void Polynomial::print

(

)

Definition at line 320 of file Poly.cpp.

{
    MultInd I( d->dim );
    double *tt = (*this);
    unsigned N = sz();
    bool IsFirst=true;
    
    if ( !N || !tt ) { printf("[Void polynomial]\n"); return; }

    if (*tt) { IsFirst=false; printf("%f", *tt); }
    tt++; ++I; 

    for (unsigned i = 1; i < N; i++,tt++,++I)
    {
        if (*tt != 0)
        {
            if (IsFirst)
            {
                if (queryFlag( NicePrint ))
                {
                    if (*tt<0) printf("-");
                    printf("%f x^",condorAbs(*tt)); I.print();
                }
                else
                {
                    printf("+%f x^",*tt); I.print();
                }
                IsFirst = false;
                continue;
            }
            if (queryFlag( NicePrint ))
            {
                if (*tt<0) printf("-"); else printf("+");
                printf("%f x^",condorAbs(*tt)); I.print();
            }
            else
            {
                printf("+%f x^",*tt); I.print();
            }

        }
    }
}

queryFlag()

unsigned Polynomial::queryFlag ( unsigned int  val )

inline

Definition at line 124 of file Poly.h.

{ return flags & val; }

save()

void Polynomial::save

(

char *

name

)

Definition at line 579 of file Poly.cpp.

{
    std::ofstream ofp(name, std::ios::out | std::ios::binary);
    ofp.write(reinterpret_cast<const char*>(&d->dim), sizeof(int));
    ofp.write(reinterpret_cast<const char*>(&d->deg), sizeof(int));
    ofp.write(reinterpret_cast<const char*>(d->coeff), d->n*sizeof(double));
    ofp.close();
}

setFlag()

void Polynomial::setFlag ( unsigned int  val )

inline

Definition at line 122 of file Poly.h.

{ flags |= val; }

shiftedEval()

double Polynomial::shiftedEval	(	Vector	Point,
		double	minusVal
	)

Definition at line 380 of file Poly.cpp.

{
    double tmp1=d->coeff[0], tmp2;
    d->coeff[0]-=minusVal;
    tmp2=(*this)(Point);
    d->coeff[0]=tmp1;
    return tmp2;
}

sz()

unsigned Polynomial::sz ( )

inline

Definition at line 64 of file Poly.h.

{ return d->n; };

translate()

void Polynomial::translate

(

Vector

translation

)

Definition at line 563 of file Poly.cpp.

{
    if (d->deg>2)
    {
        printf("Translation only for polynomial of degree lower than 3.\n");
        getchar(); exit(255);
    }
    d->coeff[0]=(*this)(translation);
    if (d->deg==1) return;
    int dim=d->dim;
    Vector G(dim);
    Matrix H(dim,dim);
    gradientHessian(translation, G, H);
    memcpy(((double*)d->coeff)+1, (double*)G, dim*sizeof(double));
}

unsetFlag()

void Polynomial::unsetFlag ( unsigned int  val )

inline

Definition at line 123 of file Poly.h.

{ flags &= ~val; }

Member Data Documentation

d

PolynomialData* Polynomial::d

protected

Definition at line 49 of file Poly.h.

emptyPolynomial

Polynomial Polynomial::emptyPolynomial

static

Definition at line 126 of file Poly.h.

flags

unsigned int Polynomial::flags = Polynomial::Warning||Polynomial::NicePrint

static

Definition at line 121 of file Poly.h.

NicePrint

const unsigned int Polynomial::NicePrint = 1

static

Definition at line 117 of file Poly.h.

Normalized

const unsigned int Polynomial::Normalized = 4

static

Definition at line 119 of file Poly.h.

Warning

const unsigned int Polynomial::Warning = 2

static

Definition at line 118 of file Poly.h.

---

The documentation for this class was generated from the following files:

• xmipp/external/condor/Poly.h
• xmipp/external/condor/Poly.cpp