#include #include #include #include "symtab.h" #include "poly.h" #define DEBUG 0 #if DEBUG #define TRACE(f,s,d) fprintf(stderr,f,s,d); #define TRACEP(s,p) fprintf(stderr,": %s",s);pprint(stderr,p);fprintf(stderr," :\n"); #else #define TRACE(f,s,d) ; #define TRACEP(s,p) ; #endif static node *variables[MAXEXP]={0}; static int max_var=0; static int pcreatevar(node *n) /* create new variable to be used with polynoms */ { int i; assert(max_vari=i; TRACE("created %s (%d)\n",n->symbol,i) return i; } int pgetvar(char *name) { node *n; n=lookup(name); if (n!=NULL) return n->i; n= insert(name); return pcreatevar(n); } static poly *pcreate1(void) { poly *p; int i; p=malloc(sizeof(*p)); assert(p!=NULL); p->c=1.0; for(i=0;ie[i]=0; p->next=NULL; TRACE("new\n",0,0); return p; } static poly *pfree(poly *p) { if (p!=NULL) free(p); return NULL; } static poly *pfreeall(poly *p) { while (p!=NULL) { poly *next=p->next; free(p); p=next; } return NULL; } static int ecmp(poly *p, poly *q) /* compare two exponents. return 1 if p>q; 0 if p=q; -1 if pe[i]>q->e[i]) return 1; if (p->e[i]e[i]) return -1; } return 0; } #define EPSILON (1e-10) static int tcmp(poly *p, poly *q) /* compare two terms. return 1 if p>q; 0 if p=q; -1 if pc-q->c)< fabs(EPSILON*p->c) ||fabs(p->c-q->c)< fabs(EPSILON*q->c)) return 0; else if (p->c < q->c) return -1; else if (p->c > q->c ) return 1; else return 0; } static int pcmp(poly *p, poly *q) /* compare two terms. return 1 if p>q; 0 if p=q; -1 if pnext; q=q->next; } return 0; } poly * pconst(double c) { poly *p; if (c==0) return NULL; p=pcreate1(); p->c=c; TRACEP("const ",p) return p; } poly * pvar(int i, int e) { poly *p; p=pcreate1(); p->e[i]=e; TRACEP("var ",p) return p; } poly *padd(poly *p, poly *q) /* returns p+q; destroys p and q. */ { if (q==NULL) return p; TRACE("adding\n",NULL,0); while (p!=NULL) { poly *head,*tail, *next; int x; head = NULL; tail = q; while ((x=ecmp(tail,p))>0) { head=tail; tail=tail->next; } next=p->next; if (x==0) /*merge*/ { tail->c+=p->c; if (fabs(tail->c)c)) { poly *n=tail; tail=tail->next; if (head!=NULL) head->next=tail; else q=tail; pfree(n); } pfree(p); } else /*insert*/ { if (head==NULL) /*in front*/ { p->next=q; q=p; } else /* in middle */ { head->next=p; p->next=tail; } } p=next; } TRACEP("add ",q) return q; } poly * pcopy(poly *p) /* return a copy of p; leave p unchanged */ { poly *x=NULL; poly **tail=&x; while (p!=NULL) { poly *n; int i; n=pcreate1(); n->c=p->c; for(i=0;ie[i]=p->e[i]; *tail=n; tail=&(n->next); p=p->next; } return x; } poly *pmulc(double c, poly *p) /* multiply p with c */ { poly *x; if (c==0) { pfreeall(p); return NULL;} x=p; while (p!=NULL) { p->c = p->c*c; p=p->next; } TRACEP("mul c ",x) return x; } static poly *pmulhead(poly *p, poly *q) /* assume multiply head of p with q leave p and q unchanged */ { poly *x=NULL; poly **tail=&x; if (p==NULL) return NULL; while (q!=NULL) { poly *n; int i; n=pcreate1(); n->c=p->c*q->c; for(i=0;ie[i]=p->e[i]+q->e[i]; *tail=n; tail=&(n->next); q=q->next; } TRACEP("mulhead ",x) return x; } static poly *pmulcopy(poly *p, poly *q) /* assume p!=NULL, q!=NULL, multiply head of p with q leave p and q unchanged */ { if (p==NULL) return NULL; if (q==NULL) return NULL; return padd(pmulhead(p,q), pmulcopy(p->next,q)); } poly * pmul(poly *p, poly *q) { poly *x; if (p==NULL) x=NULL; else if (p->next==NULL) x = pmulhead(p,q); else x = padd(pmulhead(p,q), pmulcopy(p->next,q)); pfreeall(p); pfreeall(q); TRACEP("mul all ",x); return x; } void pprint(FILE *f,poly *p) { int c=0; if (p==NULL) { fprintf(f,"0"); TRACE("printing NULL polynom\n",0,0); return; } TRACE("start printing nodes\n",0,0) while (p!=NULL) { int i; TRACE("printing node %d, c=%f\n",++c,p->c); fprintf(f,"%+g",p->c); for(i=0;ie[i]!=0) { TRACE("printing variable %d exponent %d",i,p->e[i]) fprintf(f,"%s",variables[i]->symbol); if (p->e[i]!=1) fprintf(f,"^%d",p->e[i]); } TRACE("end node %d\n",c,0); p=p->next; } } static int phasvar(int i, poly *p) /* tests wheter the polynom contains this variable */ { while (p!=NULL) { if (p->e[i]!=0) return 1; p=p->next; } return 0; } int pverify(poly *p, poly *oop, poly *mem) { int u = pgetvar("\\upsilon"); int m = pgetvar("\\mu"); poly *rt=NULL; if (phasvar(u,p)) { poly *v=pvar(u,1); rt=pmulhead(v,oop); } if (phasvar(m,p)) { poly *v=pvar(m,1); rt = padd(rt,pmulhead(v,mem)); } return pcmp(p,rt)==0; }