[BACK]Return to poly2.c CVS log [TXT][DIR] Up to [local] / OpenXM / src / kan96xx / Kan

Diff for /OpenXM/src/kan96xx/Kan/poly2.c between version 1.2 and 1.7

version 1.2, 2000/01/16 07:55:40 version 1.7, 2005/07/03 11:08:54
Line 1 
Line 1 
 /* $OpenXM$ */  /* $OpenXM: OpenXM/src/kan96xx/Kan/poly2.c,v 1.6 2005/06/16 05:07:23 takayama Exp $ */
 #include <stdio.h>  #include <stdio.h>
   #include <stdlib.h>
 #include "datatype.h"  #include "datatype.h"
 #include "stackm.h"  #include "stackm.h"
 #include "extern.h"  #include "extern.h"
Line 8 
Line 9 
 static POLY mapZmonom(POLY f,struct ring *ringp);  static POLY mapZmonom(POLY f,struct ring *ringp);
   
 POLY ppAdd(f,g)  POLY ppAdd(f,g)
 POLY f; POLY g;  /* The result is read only. */           POLY f; POLY g;  /* The result is read only. */
 {  {
   POLY node;    POLY node;
   struct listPoly nod;    struct listPoly nod;
Line 33  POLY f; POLY g;  /* The result is read only. */
Line 34  POLY f; POLY g;  /* The result is read only. */
       h = h->next;        h = h->next;
       f = f->next;        f = f->next;
       if (f == POLYNULL) {        if (f == POLYNULL) {
         h->next = g;                  h->next = g;
         return(node->next);                  return(node->next);
       }        }
       break;        break;
     case 0: /* f < g */      case 0: /* f < g */
Line 42  POLY f; POLY g;  /* The result is read only. */
Line 43  POLY f; POLY g;  /* The result is read only. */
       h = h->next;        h = h->next;
       g = g->next;        g = g->next;
       if (g == POLYNULL) {        if (g == POLYNULL) {
         h->next = f;                  h->next = f;
         return(node->next);                  return(node->next);
       }        }
       break;        break;
     case 2:/* f == g */      case 2:/* f == g */
       c = coeffCopy(f->coeffp);        c = coeffCopy(f->coeffp);
       Cadd(c,c,g->coeffp);        Cadd(c,c,g->coeffp);
       if (!isZero(c)) {        if (!isZero(c)) {
         h->next = newCell(c,f->m);                  h->next = newCell(c,f->m);
         h = h->next;                  h = h->next;
         f = f->next;                  f = f->next;
         g = g->next;                  g = g->next;
         if (f == POLYNULL) {                  if (f == POLYNULL) {
           h->next = g;                    h->next = g;
           return(node->next);                    return(node->next);
         }                  }
         if (g == POLYNULL) {                  if (g == POLYNULL) {
           h->next = f;                    h->next = f;
           return(node->next);                    return(node->next);
         }                  }
       }else{        }else{
         f = f->next;                  f = f->next;
         g = g->next;                  g = g->next;
         if (f == POLYNULL) {                  if (f == POLYNULL) {
           h->next = g;                    h->next = g;
           return(node->next);                    return(node->next);
         }                  }
         if (g == POLYNULL) {                  if (g == POLYNULL) {
           h->next = f;                    h->next = f;
           return(node->next);                    return(node->next);
         }                  }
       }        }
       break;        break;
     default:      default:
Line 84  POLY f; POLY g;  /* The result is read only. */
Line 85  POLY f; POLY g;  /* The result is read only. */
 }  }
   
 POLY ppSub(f,g)  POLY ppSub(f,g)
 POLY f; POLY g;  /* The result is read only. */           POLY f; POLY g;  /* The result is read only. */
 {  {
   POLY h;    POLY h;
   struct coeff *c;    struct coeff *c;
Line 99  POLY f; POLY g;  /* The result is read only. */
Line 100  POLY f; POLY g;  /* The result is read only. */
   
   
 POLY cpMult(c,f)  POLY cpMult(c,f)
 struct coeff *c;           struct coeff *c;
 POLY f;           POLY f;
 {  {
   POLY node;    POLY node;
   struct listPoly nod;    struct listPoly nod;
Line 125  POLY f;
Line 126  POLY f;
 }  }
   
 MONOMIAL monomialAdd_poly(m,m2)  MONOMIAL monomialAdd_poly(m,m2)
 MONOMIAL m,m2;           MONOMIAL m,m2;
 {  {
   extern int Msize;    extern int Msize;
   MONOMIAL f;    MONOMIAL f;
Line 145  MONOMIAL m,m2;
Line 146  MONOMIAL m,m2;
   
 /* Note that mpMult_poly is called from mmLarger_tower! */  /* Note that mpMult_poly is called from mmLarger_tower! */
 POLY mpMult_poly(f,g)  POLY mpMult_poly(f,g)
 POLY f;           POLY f;
 POLY g;           POLY g;
 {  {
   POLY node;    POLY node;
   struct listPoly nod;    struct listPoly nod;
Line 176  POLY g;
Line 177  POLY g;
 }  }
   
 POLY ppMult_old(f,g)  POLY ppMult_old(f,g)
 POLY f,g;           POLY f,g;
 {  {
   POLY r;    POLY r;
   POLY tmp;    POLY tmp;
Line 190  POLY f,g;
Line 191  POLY f,g;
 }  }
   
 POLY ppAddv(f,g)  POLY ppAddv(f,g)
 POLY f; POLY g;  /* It breaks f and g. Use it just after calling mpMult() */           POLY f; POLY g;  /* It breaks f and g. Use it just after calling mpMult() */
 {  {
   POLY node;    POLY node;
   struct listPoly nod;    struct listPoly nod;
Line 215  POLY f; POLY g;  /* It breaks f and g. Use it just aft
Line 216  POLY f; POLY g;  /* It breaks f and g. Use it just aft
       h->next = f;        h->next = f;
       h = h->next; f = f->next;;        h = h->next; f = f->next;;
       if (f == POLYNULL) {        if (f == POLYNULL) {
         h->next = g;                  h->next = g;
         return(node->next);                  return(node->next);
       }        }
       break;        break;
     case 0: /* f < g */      case 0: /* f < g */
       h->next =  g;        h->next =  g;
       h = h->next; g = g->next;        h = h->next; g = g->next;
       if (g == POLYNULL) {        if (g == POLYNULL) {
         h->next = f;                  h->next = f;
         return(node->next);                  return(node->next);
       }        }
       break;        break;
     case 2:/* f == g */      case 2:/* f == g */
       c = f->coeffp;        c = f->coeffp;
       Cadd(c,c,g->coeffp);        Cadd(c,c,g->coeffp);
       if (!isZero(c)) {        if (!isZero(c)) {
         h->next = f;                  h->next = f;
         h = h->next; f = f->next;;                  h = h->next; f = f->next;;
         g = g->next;                  g = g->next;
         if (f == POLYNULL) {                  if (f == POLYNULL) {
           h->next = g;                    h->next = g;
           return(node->next);                    return(node->next);
         }                  }
         if (g == POLYNULL) {                  if (g == POLYNULL) {
           h->next = f;                    h->next = f;
           return(node->next);                    return(node->next);
         }                  }
       }else{        }else{
         f = f->next;                  f = f->next;
         g = g->next;                  g = g->next;
         if (f == POLYNULL) {                  if (f == POLYNULL) {
           h->next = g;                    h->next = g;
           return(node->next);                    return(node->next);
         }                  }
         if (g == POLYNULL) {                  if (g == POLYNULL) {
           h->next = f;                    h->next = f;
           return(node->next);                    return(node->next);
         }                  }
       }        }
       break;        break;
     default:      default:
Line 264  POLY f; POLY g;  /* It breaks f and g. Use it just aft
Line 265  POLY f; POLY g;  /* It breaks f and g. Use it just aft
 }  }
   
 POLY pPower(f,k)  POLY pPower(f,k)
 POLY f;           POLY f;
 int k;           int k;
 {  {
   POLY r;    POLY r;
   int i,n;    int i,n;
Line 296  int k;
Line 297  int k;
 }  }
   
 POLY pPower_poly(f,k)  POLY pPower_poly(f,k)
 POLY f;           POLY f;
 int k;           int k;
 {  {
   POLY r;    POLY r;
   int i,n;    int i,n;
Line 328  int k;
Line 329  int k;
 }  }
   
 POLY modulop_trash(f,ringp)  POLY modulop_trash(f,ringp)
 POLY f;           POLY f;
 struct ring *ringp;           struct ring *ringp;
 {  {
   int p;    int p;
   POLY h;    POLY h;
Line 382  struct ring *ringp;
Line 383  struct ring *ringp;
       mpz_mod_ui(c,f->coeffp->val.bigp,(unsigned long int)p);        mpz_mod_ui(c,f->coeffp->val.bigp,(unsigned long int)p);
       cc = (int) mpz_get_si(c);        cc = (int) mpz_get_si(c);
       if (cc != 0) {        if (cc != 0) {
         h->next = newCell(newCoeff(),monomialCopy(f->m));                  h->next = newCell(newCoeff(),monomialCopy(f->m));
         h = h->next;                  h = h->next;
         h->m->ringp = ringp;                  h->m->ringp = ringp;
         h->coeffp->tag = INTEGER;                  h->coeffp->tag = INTEGER;
         h->coeffp->p = p;                  h->coeffp->p = p;
         h->coeffp->val.i = cc;                  h->coeffp->val.i = cc;
       }        }
       f = f->next;        f = f->next;
     }      }
Line 404  struct ring *ringp;
Line 405  struct ring *ringp;
 }  }
   
 POLY modulop(f,ringp)  POLY modulop(f,ringp)
 POLY f;           POLY f;
 struct ring *ringp;           struct ring *ringp;
 /* Z[x] ---> R[x] where R=Z, Z/Zp, ringp->next. */           /* Z[x] ---> R[x] where R=Z, Z/Zp, ringp->next. */
 {  {
   int p;    int p;
   POLY h;    POLY h;
Line 447  struct ring *ringp;
Line 448  struct ring *ringp;
       mpz_mod_ui(c,f->coeffp->val.bigp,(unsigned long int)p);        mpz_mod_ui(c,f->coeffp->val.bigp,(unsigned long int)p);
       cc = (int) mpz_get_si(c);        cc = (int) mpz_get_si(c);
       if (cc != 0) {        if (cc != 0) {
         h->next = newCell(newCoeff(),monomialCopy(f->m));                  h->next = newCell(newCoeff(),monomialCopy(f->m));
         h = h->next;                  h = h->next;
         h->m->ringp = ringp;                  h->m->ringp = ringp;
         h->coeffp->tag = INTEGER;                  h->coeffp->tag = INTEGER;
         h->coeffp->p = p;                  h->coeffp->p = p;
         h->coeffp->val.i = cc;                  h->coeffp->val.i = cc;
       }        }
       f = f->next;        f = f->next;
     }      }
Line 471  struct ring *ringp;
Line 472  struct ring *ringp;
 }  }
   
 POLY modulopZ(f,pcoeff)  POLY modulopZ(f,pcoeff)
 POLY f;           POLY f;
 struct coeff *pcoeff;           struct coeff *pcoeff;
 /* Z[x] ---> Z[x] , f ---> f mod pcoeff*/           /* Z[x] ---> Z[x] , f ---> f mod pcoeff*/
 {  {
   int p;    int p;
   POLY h;    POLY h;
Line 535  struct coeff *pcoeff;
Line 536  struct coeff *pcoeff;
 }  }
   
 struct pairOfPOLY quotientByNumber(f,pcoeff)  struct pairOfPOLY quotientByNumber(f,pcoeff)
 POLY f;           POLY f;
 struct coeff *pcoeff;           struct coeff *pcoeff;
 /* Z[x] ---> Z[x],Z[x] ,  f = first*pcoeff + second */           /* Z[x] ---> Z[x],Z[x] ,  f = first*pcoeff + second */
 {  {
   int p;    int p;
   POLY h;    POLY h;
Line 631  struct coeff *pcoeff;
Line 632  struct coeff *pcoeff;
   
   
 POLY modulo0(f,ringp)  POLY modulo0(f,ringp)
 POLY f;           POLY f;
 struct ring *ringp;           struct ring *ringp;
 {  {
   int p;    int p;
   POLY h;    POLY h;
Line 650  struct ring *ringp;
Line 651  struct ring *ringp;
       node = pcmCopy(f);        node = pcmCopy(f);
       f = node;        f = node;
       while (f != POLYNULL) {        while (f != POLYNULL) {
         f->m->ringp = ringp; /* Touch the monomial "ringp" field. */                  f->m->ringp = ringp; /* Touch the monomial "ringp" field. */
         f = f->next;                  f = f->next;
       }        }
       return(node);        return(node);
     }      }
Line 692  struct ring *ringp;
Line 693  struct ring *ringp;
   
   
 struct object test(ob)  /* test3 */  struct object test(ob)  /* test3 */
 struct object ob;           struct object ob;
 {  {
   struct object rob;    struct object rob = OINIT;
   int k;    int k;
   static POLY f0;    static POLY f0;
   static POLY f1;    static POLY f1;
Line 714  struct object ob;
Line 715  struct object ob;
     for (i=k; i>=0; i--) {      for (i=k; i>=0; i--) {
       f0->next = bxx(BiiPower(-k,i),0,i,CurrentRingp);        f0->next = bxx(BiiPower(-k,i),0,i,CurrentRingp);
       if (f0->next != POLYNULL) {        if (f0->next != POLYNULL) {
         f0 = f0->next;                  f0 = f0->next;
       }        }
     }      }
     f0 = addNode->next;      f0 = addNode->next;
Line 726  struct object ob;
Line 727  struct object ob;
     for (i=k; i>=0; i--) {      for (i=k; i>=0; i--) {
       f1->next = bxx(BiiPower(k,i),0,i,CurrentRingp);        f1->next = bxx(BiiPower(k,i),0,i,CurrentRingp);
       if (f1->next != POLYNULL) {        if (f1->next != POLYNULL) {
         f1 = f1->next;                  f1 = f1->next;
       }        }
     }      }
     f1 = addNode->next;      f1 = addNode->next;
Line 746  struct object ob;
Line 747  struct object ob;
   
   
 int pLength(f)  int pLength(f)
 POLY f;           POLY f;
 {  {
   int c=0;    int c=0;
   if (f ISZERO) return(0);    if (f ISZERO) return(0);
Line 759  POLY f;
Line 760  POLY f;
   
   
 POLY ppAddv2(f,g,top,nexttop)  POLY ppAddv2(f,g,top,nexttop)
 POLY f; POLY g;  /* It breaks f and g. Use it just after calling mpMult() */           POLY f; POLY g;  /* It breaks f and g. Use it just after calling mpMult() */
 POLY top;           POLY top;
 POLY *nexttop;           POLY *nexttop;
 /* top is the starting address in the list f.           /* top is the starting address in the list f.
    if top == POLYNULL, start from f.                  if top == POLYNULL, start from f.
   
    *nexttop == 0                  *nexttop == 0
             == g                  == g
             == h or 0                  == h or 0
   
   It must be called as r = ppAddv2(r,g,...);                  It must be called as r = ppAddv2(r,g,...);
 */  */
 {  {
   POLY node;    POLY node;
Line 793  POLY *nexttop;
Line 794  POLY *nexttop;
   if (top != POLYNULL) {    if (top != POLYNULL) {
     while (f != top) {      while (f != top) {
       if (f == POLYNULL) {        if (f == POLYNULL) {
         fprintf(stderr,"\nppAddv2(): Internal error.\n");fflush(stderr);                  fprintf(stderr,"\nppAddv2(): Internal error.\n");fflush(stderr);
         fprintf(stderr,"f = %s\n",POLYToString(f0,'*',0));                  fprintf(stderr,"f = %s\n",POLYToString(f0,'*',0));
         fprintf(stderr,"g = %s\n",POLYToString(g0,'*',0));                  fprintf(stderr,"g = %s\n",POLYToString(g0,'*',0));
         fprintf(stderr,"top=%s\n",POLYToString(top,'*',0));                  fprintf(stderr,"top=%s\n",POLYToString(top,'*',0));
         errorPoly("ppAddv2(). Internal error=1.");                  errorPoly("ppAddv2(). Internal error=1.");
       }        }
       h->next = f;        h->next = f;
       h = h->next;        h = h->next;
Line 816  POLY *nexttop;
Line 817  POLY *nexttop;
       h->next = f;        h->next = f;
       h = h->next; f = f->next;;        h = h->next; f = f->next;;
       if (f == POLYNULL) {        if (f == POLYNULL) {
         h->next = g;                  h->next = g;
         return(node->next);                  return(node->next);
       }        }
       break;        break;
     case 0: /* f < g */      case 0: /* f < g */
       h->next =  g;        h->next =  g;
       h = h->next; g = g->next;        h = h->next; g = g->next;
       if (g == POLYNULL) {        if (g == POLYNULL) {
         h->next = f;                  h->next = f;
         return(node->next);                  return(node->next);
       }        }
       break;        break;
     case 2:/* f == g */      case 2:/* f == g */
       c = g->coeffp;        c = g->coeffp;
       Cadd(c,f->coeffp,c);        Cadd(c,f->coeffp,c);
       if (!isZero(c)) {        if (!isZero(c)) {
         h->next = g;                  h->next = g;
         h = h->next;                  h = h->next;
         f = f->next;;                  f = f->next;;
         g = g->next;                  g = g->next;
         if (f == POLYNULL) {                  if (f == POLYNULL) {
           h->next = g;                    h->next = g;
           return(node->next);                    return(node->next);
         }                  }
         if (g == POLYNULL) {                  if (g == POLYNULL) {
           h->next = f;                    h->next = f;
           return(node->next);                    return(node->next);
         }                  }
       }else{        }else{
         if (g == g0) {                  if (g == g0) {
           if (h != node) {                    if (h != node) {
             *nexttop = h;                          *nexttop = h;
           }else{                    }else{
             *nexttop = POLYNULL;                          *nexttop = POLYNULL;
           }                    }
         }                  }
   
         f = f->next;                  f = f->next;
         g = g->next;                  g = g->next;
   
         if (f == POLYNULL) {                  if (f == POLYNULL) {
           h->next = g;                    h->next = g;
           return(node->next);                    return(node->next);
         }                  }
         if (g == POLYNULL) {                  if (g == POLYNULL) {
           h->next = f;                    h->next = f;
           return(node->next);                    return(node->next);
         }                  }
       }        }
       break;        break;
     default:      default:
Line 875  POLY *nexttop;
Line 876  POLY *nexttop;
 }  }
   
 POLY ppMult(f,g)  POLY ppMult(f,g)
 POLY f,g;           POLY f,g;
 {  {
   POLY r;    POLY r;
   POLY tmp;    POLY tmp;
Line 894  POLY f,g;
Line 895  POLY f,g;
 }  }
   
 POLY ppMult_poly(f,g)  POLY ppMult_poly(f,g)
 POLY f,g;           POLY f,g;
 {  {
   POLY r;    POLY r;
   POLY tmp;    POLY tmp;
Line 911  POLY f,g;
Line 912  POLY f,g;
 }  }
   
 POLY mapZmonom(f,ringp)  POLY mapZmonom(f,ringp)
 POLY f; /* assumes monomial. f \in Z[x] */           POLY f; /* assumes monomial. f \in Z[x] */
 struct ring *ringp;  /* R[x] */           struct ring *ringp;  /* R[x] */
 {  {
   struct ring *nextRing;    struct ring *nextRing;
   struct ring nextRing0;    struct ring nextRing0;
Line 947  struct ring *ringp;  /* R[x] */
Line 948  struct ring *ringp;  /* R[x] */
   node->coeffp->val.f = ff;    node->coeffp->val.f = ff;
   return(node);    return(node);
 }  }
   
   POLY reduceContentOfPoly(POLY f,struct coeff **contp) {
     struct coeff *cont;
     struct coeff *cOne = NULL;
     extern struct ring *SmallRingp;
     if (cOne == NULL) cOne = intToCoeff(1,SmallRingp);
     *contp = cOne;
   
     if (f == POLYNULL) return f;
     if (f->m->ringp->p != 0) return f;
     if (f->coeffp->tag != MP_INTEGER) return f;
     cont = gcdOfCoeff(f);
     *contp = cont;
     if (coeffGreater(cont,cOne)) {
           f = quotientByNumber(f,cont).first;
     }
     return f;
   }
   
   int coeffSizeMin(POLY f) {
     int size;
     int t;
     if (f == POLYNULL) return 0;
     if (f->m->ringp->p != 0) return 0;
     if (f->coeffp->tag != MP_INTEGER) return 0;
     size = mpz_size(f->coeffp->val.bigp);
     while (f != POLYNULL) {
           t = mpz_size(f->coeffp->val.bigp);
           if (t < size)  size = t;
           if (size == 1) return size;
           f = f->next;
     }
   }
   
   struct coeff *gcdOfCoeff(POLY f) {
     extern struct ring *SmallRingp;
     struct coeff *t;
     MP_INT *tmp;
     MP_INT *tmp2;
     static MP_INT *cOne = NULL;
     if (cOne == NULL) {
           cOne = newMP_INT();
           mpz_set_si(cOne,(long) 1);
     }
     if (f == POLYNULL) return intToCoeff(0,SmallRingp);
     if (f->m->ringp->p != 0) return intToCoeff(0,SmallRingp);
     if (f->coeffp->tag != MP_INTEGER) return intToCoeff(0,SmallRingp);
     tmp = f->coeffp->val.bigp;
     tmp2 = newMP_INT();
     while (f != POLYNULL) {
       mpz_gcd(tmp2,tmp,f->coeffp->val.bigp); /* tmp = tmp2 OK? */
           tmp = tmp2;
           if (mpz_cmp(tmp,cOne)==0) return intToCoeff(1,SmallRingp);
           f = f->next;
     }
     return mpintToCoeff(tmp,SmallRingp);
   
   }
   
   int shouldReduceContent(POLY f,int ss) {
     extern DoCancel;
     static int prevSize = 1;
     int size;
     if (f == POLYNULL) return 0;
     if (f->m->ringp->p != 0) return 0;
     if (f->coeffp->tag != MP_INTEGER) return 0;
     if (DoCancel & 2) return 1;
     /* Apply the Noro strategy to reduce content. */
     size = mpz_size(f->coeffp->val.bigp);
     if (ss > 0) {
           prevSize = size;
           return 0;
     }
     if (size > 2*prevSize) {
           return 1;
     }else{
           return 0;
     }
   }
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.7
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>