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

Diff for /OpenXM/src/kan96xx/Kan/poly4.c between version 1.5 and 1.12

version 1.5, 2003/07/10 08:20:04 version 1.12, 2003/08/24 05:19:42
Line 1 
Line 1 
 /* $OpenXM: OpenXM/src/kan96xx/Kan/poly4.c,v 1.4 2002/09/08 10:49:50 takayama Exp $ */  /* $OpenXM: OpenXM/src/kan96xx/Kan/poly4.c,v 1.11 2003/08/23 02:28:38 takayama Exp $ */
 #include <stdio.h>  #include <stdio.h>
 #include "datatype.h"  #include "datatype.h"
 #include "stackm.h"  #include "stackm.h"
Line 8 
Line 8 
 static void shell(int v[],int n);  static void shell(int v[],int n);
 static int degreeOfPrincipalPart(POLY f);  static int degreeOfPrincipalPart(POLY f);
 static int degreeOfInitW(POLY f,int w[]);  static int degreeOfInitW(POLY f,int w[]);
   static int degreeOfInitWS(POLY f,int w[],int s[]);
   
   
 static void shell(v,n)  static void shell(v,n)
Line 437  POLY POLYToInitW(f,w)
Line 438  POLY POLYToInitW(f,w)
      POLY f;       POLY f;
      int w[]; /* weight vector */       int w[]; /* weight vector */
 {  {
   POLY node;  
   struct listPoly nod;  
   POLY h;    POLY h;
   POLY g;    POLY g;
   int maxd;    int maxd;
   int dd;    int dd;
   node = &nod; node->next = POLYNULL; h = node;    h = POLYNULL;
   
     /*printf("1:%s\n",POLYToString(f,'*',1));*/
   if (f ISZERO) return(f);    if (f ISZERO) return(f);
   maxd = degreeOfInitW(f,w);    maxd = degreeOfInitW(f,w);
   g = pCopy(f); /* shallow copy */    g = f;
   while (!(f ISZERO)) {    while (!(f ISZERO)) {
     dd = degreeOfInitW(f,w);      dd = degreeOfInitW(f,w);
     if (dd > maxd) maxd = dd;      if (dd > maxd) maxd = dd;
Line 456  POLY POLYToInitW(f,w)
Line 456  POLY POLYToInitW(f,w)
   while (!(g ISZERO)) {    while (!(g ISZERO)) {
     dd = degreeOfInitW(g,w);      dd = degreeOfInitW(g,w);
     if (dd == maxd) {      if (dd == maxd) {
       h->next = g;        h = ppAdd(h,newCell(g->coeffp,g->m)); /* it might be slow. */
       h = h->next;  
     }      }
     g = g->next;      g = g->next;
   }    }
   h->next = POLYNULL;    /*printf("2:%s\n",POLYToString(h,'*',1));*/
   return(node->next);    return(h);
 }  }
   
   static int degreeOfInitWS(f,w,s)
        POLY f;
        int w[];
            int s[];
   {
     int n,i,dd;
     if (f ISZERO) {
       errorPoly("degreeOfInitWS(0,w) ");
     }
     if (s == (int *) NULL) return degreeOfInitW(f,w);
     n = f->m->ringp->n; dd = 0;
     for (i=0; i<n-1; i++) {
       dd += (f->m->e[i].D)*w[n+i];
       dd += (f->m->e[i].x)*w[i];
     }
     dd += s[(f->m->e[n-1].x)];
     return(dd);
   }
   
   POLY POLYToInitWS(f,w,s)
        POLY f;
        int w[]; /* weight vector */
            int s[]; /* shift vector */
   {
     POLY h;
     POLY g;
     int maxd;
     int dd;
     h = POLYNULL;
   
     /*printf("1s:%s\n",POLYToString(f,'*',1));*/
     if (f ISZERO) return(f);
     maxd = degreeOfInitWS(f,w,s);
     g = f;
     while (!(f ISZERO)) {
       dd = degreeOfInitWS(f,w,s);
       if (dd > maxd) maxd = dd;
       f = f->next;
     }
     while (!(g ISZERO)) {
       dd = degreeOfInitWS(g,w,s);
       if (dd == maxd) {
         h = ppAdd(h,newCell(g->coeffp,g->m)); /* it might be slow. */
       }
       g = g->next;
     }
     /*printf("2s:%s\n",POLYToString(h,'*',1));*/
     return(h);
   }
   
   int ordWsAll(f,w,s)
        POLY f;
        int w[]; /* weight vector */
            int s[]; /* shift vector */
   {
     int maxd;
     int dd;
   
     if (f ISZERO)  errorPoly("ordWsAll(0,w,s) ");
     maxd = degreeOfInitWS(f,w,s);
     while (!(f ISZERO)) {
       dd = degreeOfInitWS(f,w,s);
       if (dd > maxd) maxd = dd;
       f = f->next;
     }
     return maxd;
   }
   
   
 /*  /*
 1.The substitution  "ringp->multiplication = ...." is allowed only in  1.The substitution  "ringp->multiplication = ...." is allowed only in
   KsetUpRing(), so the check in KswitchFunction is not necessary.    KsetUpRing(), so the check in KswitchFunction is not necessary.
Line 508  int isTheSameRing(struct ring *rstack[],int rp, struct
Line 575  int isTheSameRing(struct ring *rstack[],int rp, struct
     /* if (rrr->schreyer != newRingp->schreyer) { a=16; goto bbb ; }*/      /* if (rrr->schreyer != newRingp->schreyer) { a=16; goto bbb ; }*/
     if (newRingp->schreyer == 1) { a=16; goto bbb; }      if (newRingp->schreyer == 1) { a=16; goto bbb; }
     if (rrr->weightedHomogenization != newRingp->weightedHomogenization) { a=16; goto bbb; }      if (rrr->weightedHomogenization != newRingp->weightedHomogenization) { a=16; goto bbb; }
       if (rrr->degreeShiftSize != newRingp->degreeShiftSize) {
         a = 17; goto bbb;
       }
       if (rrr->degreeShiftN != newRingp->degreeShiftN) {
         a = 17; goto bbb;
       }
       for (i=0; i < rrr->degreeShiftSize; i++) {
         for (j=0; j< rrr->degreeShiftN; j++) {
           if (rrr->degreeShift[i*(rrr->degreeShiftN)+j] !=
               newRingp->degreeShift[i*(rrr->degreeShiftN)+j]) {
             a = 17; goto bbb;
           }
         }
       }
   
     /* The following fields are ignored.      /* The following fields are ignored.
        void *gbListTower;         void *gbListTower;
        int *outputOrder;         int *outputOrder;
Line 521  int isTheSameRing(struct ring *rstack[],int rp, struct
Line 603  int isTheSameRing(struct ring *rstack[],int rp, struct
   }    }
   return(-1);    return(-1);
 }  }
   
   /* s->1 */
   POLY goDeHomogenizeS(POLY f) {
     POLY lRule[1];
     POLY rRule[1];
     struct ring *rp;
     POLY ans;
     /* printf("1:[%s]\n",POLYToString(f,'*',1)); */
     if (f == POLYNULL) return f;
     rp = f->m->ringp;
     if (rp->next == NULL) {
           lRule[0] = cxx(1,0,1,rp);
           rRule[0] = cxx(1,0,0,rp);
           ans=replace(f,lRule,rRule,1);
     }else{
           struct coeff *cp;
           POLY t;
           POLY nc;
       ans = POLYNULL;
           while (f != POLYNULL) {
             cp = f->coeffp;
             if (cp->tag == POLY_COEFF) {
                   t = goDeHomogenizeS((cp->val).f);
                   nc = newCell(polyToCoeff(t,f->m->ringp),monomialCopy(f->m));
                   ans = ppAddv(ans,nc);
                   f = f->next;
             }else{
                   ans = f; break;
             }
           }
     }
     /* printf("2:[%s]\n",POLYToString(ans,'*',1)); */
     return ans;
   }
   
   POLY goDeHomogenizeS_buggy(POLY f) {
     POLY node;
     POLY lastf;
     struct listPoly nod;
     POLY h;
     POLY tf;
     int gt,first;
   
     printf("1:[%s]\n",POLYToString(f,'*',1));
     if (f == POLYNULL) return(POLYNULL);
     node = &nod;
     node->next = POLYNULL;
     lastf = POLYNULL;
     first = 1;
     while (f != POLYNULL) {
       tf = newCell(f->coeffp,monomialCopy(f->m));
       tf->m->e[0].x = 0;  /* H, s variable in the G-O paper. */
       if (first) {
         node->next = tf;
         lastf = tf;
         first = 0;
       }else{
         gt = (*mmLarger)(lastf,tf);
         if (gt == 1) {
           lastf->next = tf;
           lastf = tf;
         }else{
           h = node->next;
           h = ppAddv(h,tf);
           node->next = h;
           lastf = h;
           while (lastf->next != POLYNULL) {
             lastf = lastf->next;
           }
         }
       }
       f = f->next;
     }
     printf("2:[%s]\n",POLYToString(node->next,'*',1));
     return (node->next);
   }
   
 /* Granger-Oaku's homogenization for the ecart tangent cone.  /* Granger-Oaku's homogenization for the ecart tangent cone.
    Note: 2003.07.10.     Note: 2003.07.10.
    ds[] is the degree shift.     ds[] is the degree shift.
    ei ( element index ). If it is < 0, then e[n-1]->x will be used,     ei ( element index ). If it is < 0, then e[n-1]->x will be used,
                          else ei is used.                           else ei is used.
      if onlyS is set to 1, then input is assumed to be (u,v)-h-homogeneous.
 */  */
 POLY goHomogenize(POLY f,int u[],int v[],int ds[],int dssize,int ei)  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int dssize,int ei,int onlyS)
 {  {
   POLY node;    POLY node;
   POLY lastf;    POLY lastf;
Line 536  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
Line 695  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
   POLY h;    POLY h;
   POLY tf;    POLY tf;
   int gt,first,m,mp,t,tp,dsIdx,message;    int gt,first,m,mp,t,tp,dsIdx,message;
     struct ring *rp;
   
   message = 1;    message = 1;
   if (f == POLYNULL) return(POLYNULL);    if (f == POLYNULL) return(POLYNULL);
     rp = f->m->ringp;
     /*
     if ((rp->degreeShiftSize == 0) && (dssize > 0)) {
           warningPoly("You are trying to homogenize a polynomial with degree shift. However, the polynomial belongs to the ring without degreeShift option. It may cause a trouble in comparison in free module.\n");
     }
     */
   node = &nod;    node = &nod;
   node->next = POLYNULL;    node->next = POLYNULL;
   lastf = POLYNULL;    lastf = POLYNULL;
Line 560  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
Line 726  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
       /*go-debug fprintf(stderr,"Automatic dehomogenize and homogenize.\n"); */        /*go-debug fprintf(stderr,"Automatic dehomogenize and homogenize.\n"); */
           message = 0;            message = 0;
         }          }
     tf->m->e[0].D = -t;  /* h */          if (!onlyS) {
             tf->m->e[0].D = -t;  /* h */
           }
     tf->m->e[0].x = tp;  /* H, s variable in the G-O paper. */      tf->m->e[0].x = tp;  /* H, s variable in the G-O paper. */
         /*go-debug printf("t(h)=%d, tp(uv+ds)=%d\n",t,tp); */          /*go-debug printf("t(h)=%d, tp(uv+ds)=%d\n",t,tp); */
     if (first) {      if (first) {
Line 588  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
Line 756  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
   h = node->next;    h = node->next;
   /*go-debug printf("m=%d, mp=%d\n",m,mp); */    /*go-debug printf("m=%d, mp=%d\n",m,mp); */
   while (h != POLYNULL) {    while (h != POLYNULL) {
     /*go-debug printf("Old: h=%d, s=%d\n",h->m->e[0].D,h->m->e[0].x); */      /*go-debug printf("Old: h=%d, s=%d\n",h->m->e[0].D,h->m->e[0].x);  */
     h->m->e[0].D += m;   /* h */      if (!onlyS) h->m->e[0].D += m;   /* h */
     h->m->e[0].x += -mp; /* H, s*/      h->m->e[0].x += -mp; /* H, s*/
     /*go-debug printf("New: h=%d, s=%d\n",h->m->e[0].D,h->m->e[0].x); */      /*go-debug printf("New: h=%d, s=%d\n",h->m->e[0].D,h->m->e[0].x); */
     h = h->next;      h = h->next;
Line 598  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
Line 766  POLY goHomogenize(POLY f,int u[],int v[],int ds[],int 
 }  }
   
 /* u[] = -1, v[] = 1 */  /* u[] = -1, v[] = 1 */
 POLY goHomogenize11(POLY f,int ds[],int dssize,int ei)  POLY goHomogenize11(POLY f,int ds[],int dssize,int ei,int onlyS)
 {  {
   int r;    int r;
   int i,t,n,m,nn;    int i,t,n,m,nn;
Line 622  POLY goHomogenize11(POLY f,int ds[],int dssize,int ei)
Line 790  POLY goHomogenize11(POLY f,int ds[],int dssize,int ei)
           u[i] = -1; v[i] = 1;            u[i] = -1; v[i] = 1;
         }          }
   }    }
   return(goHomogenize(f,u,v,ds,dssize,ei));    return(goHomogenize(f,u,v,ds,dssize,ei,onlyS));
 }  }
   
 POLY goHomogenize_dsIdx(POLY f,int u[],int v[],int dsIdx,int ei)  POLY goHomogenize_dsIdx(POLY f,int u[],int v[],int dsIdx,int ei,int onlyS)
 {  {
   if (f == POLYNULL) return POLYNULL;    if (f == POLYNULL) return POLYNULL;
 }  }
 POLY goHomogenize11_dsIdx(POLY f,int ds[],int dsIdx,int ei)  POLY goHomogenize11_dsIdx(POLY f,int ds[],int dsIdx,int ei,int onlyS)
 {  {
   if (f == POLYNULL) return POLYNULL;    if (f == POLYNULL) return POLYNULL;
   }
   
   /* cf. KsetUpRing() in kanExport0.c */
   struct ring *newRingOverFp(struct ring *rp,int p) {
     struct ring *newRingp;
     char *ringName = NULL;
     char pstr[64];
     sprintf(pstr,"_%d",p);
     ringName = (char *)sGC_malloc(128);
     newRingp = (struct ring *)sGC_malloc(sizeof(struct ring));
     if (newRingp == NULL) errorPoly("No more memory.\n");
     strcpy(ringName,rp->name);
     strcat(ringName,pstr);
     *newRingp = *rp;
     newRingp->p = p;
     newRingp->name = ringName;
     return newRingp;
   }
   
   /*
      P = 3001;
      L = [ ];
      while (P<10000) {
        L=cons(P,L);
        P = pari(nextprime,P+1);
           }
      print(L);
   */
   #define N799  799
   static int nextPrime(void) {
     static int pt = 0;
     static int tb[N799] =
   {3001,3011,3019,3023,3037,3041,3049,3061,3067,3079,3083,3089,3109,3119,3121,3137,3163,3167,3169,3181,3187,3191,3203,3209,3217,3221,3229,3251,3253,3257,3259,3271,3299,3301,3307,3313,3319,3323,3329,3331,3343,3347,3359,3361,3371,3373,3389,3391,3407,3413,3433,3449,3457,3461,3463,3467,3469,3491,3499,3511,3517,3527,3529,3533,3539,3541,3547,3557,3559,3571,3581,3583,3593,3607,3613,3617,3623,3631,3637,3643,3659,3671,3673,3677,3691,3697,3701,3709,3719,3727,3733,3739,3761,3767,3769,3779,3793,3797,3803,3821,3823,3833,3847,3851,3853,3863,3877,3881,3889,3907,3911,3917,3919,3923,3929,3931,3943,3947,3967,3989,
     4001,4003,4007,4013,4019,4021,4027,4049,4051,4057,4073,4079,4091,4093,4099,4111,4127,4129,4133,4139,4153,4157,4159,4177,4201,4211,4217,4219,4229,4231,4241,4243,4253,4259,4261,4271,4273,4283,4289,4297,4327,4337,4339,4349,4357,4363,4373,4391,4397,4409,4421,4423,4441,4447,4451,4457,4463,4481,4483,4493,4507,4513,4517,4519,4523,4547,4549,4561,4567,4583,4591,4597,4603,4621,4637,4639,4643,4649,4651,4657,4663,4673,4679,4691,4703,4721,4723,4729,4733,4751,4759,4783,4787,4789,4793,4799,4801,4813,4817,4831,4861,4871,4877,4889,4903,4909,4919,4931,4933,4937,4943,4951,4957,4967,4969,4973,4987,4993,4999,
     5003,5009,5011,5021,5023,5039,5051,5059,5077,5081,5087,5099,5101,5107,5113,5119,5147,5153,5167,5171,5179,5189,5197,5209,5227,5231,5233,5237,5261,5273,5279,5281,5297,5303,5309,5323,5333,5347,5351,5381,5387,5393,5399,5407,5413,5417,5419,5431,5437,5441,5443,5449,5471,5477,5479,5483,5501,5503,5507,5519,5521,5527,5531,5557,5563,5569,5573,5581,5591,5623,5639,5641,5647,5651,5653,5657,5659,5669,5683,5689,5693,5701,5711,5717,5737,5741,5743,5749,5779,5783,5791,5801,5807,5813,5821,5827,5839,5843,5849,5851,5857,5861,5867,5869,5879,5881,5897,5903,5923,5927,5939,5953,5981,5987,
     6007,6011,6029,6037,6043,6047,6053,6067,6073,6079,6089,6091,6101,6113,6121,6131,6133,6143,6151,6163,6173,6197,6199,6203,6211,6217,6221,6229,6247,6257,6263,6269,6271,6277,6287,6299,6301,6311,6317,6323,6329,6337,6343,6353,6359,6361,6367,6373,6379,6389,6397,6421,6427,6449,6451,6469,6473,6481,6491,6521,6529,6547,6551,6553,6563,6569,6571,6577,6581,6599,6607,6619,6637,6653,6659,6661,6673,6679,6689,6691,6701,6703,6709,6719,6733,6737,6761,6763,6779,6781,6791,6793,6803,6823,6827,6829,6833,6841,6857,6863,6869,6871,6883,6899,6907,6911,6917,6947,6949,6959,6961,6967,6971,6977,6983,6991,6997,
     7001,7013,7019,7027,7039,7043,7057,7069,7079,7103,7109,7121,7127,7129,7151,7159,7177,7187,7193,7207,7211,7213,7219,7229,7237,7243,7247,7253,7283,7297,7307,7309,7321,7331,7333,7349,7351,7369,7393,7411,7417,7433,7451,7457,7459,7477,7481,7487,7489,7499,7507,7517,7523,7529,7537,7541,7547,7549,7559,7561,7573,7577,7583,7589,7591,7603,7607,7621,7639,7643,7649,7669,7673,7681,7687,7691,7699,7703,7717,7723,7727,7741,7753,7757,7759,7789,7793,7817,7823,7829,7841,7853,7867,7873,7877,7879,7883,7901,7907,7919,7927,7933,7937,7949,7951,7963,7993,
   8009,8011,8017,8039,8053,8059,8069,8081,8087,8089,8093,8101,8111,8117,8123,8147,8161,8167,8171,8179,8191,8209,8219,8221,8231,8233,8237,8243,8263,8269,8273,8287,8291,8293,8297,8311,8317,8329,8353,8363,8369,8377,8387,8389,8419,8423,8429,8431,8443,8447,8461,8467,8501,8513,8521,8527,8537,8539,8543,8563,8573,8581,8597,8599,8609,8623,8627,8629,8641,8647,8663,8669,8677,8681,8689,8693,8699,8707,8713,8719,8731,8737,8741,8747,8753,8761,8779,8783,8803,8807,8819,8821,8831,8837,8839,8849,8861,8863,8867,8887,8893,8923,8929,8933,8941,8951,8963,8969,8971,8999,
     9001,9007,9011,9013,9029,9041,9043,9049,9059,9067,9091,9103,9109,9127,9133,9137,9151,9157,9161,9173,9181,9187,9199,9203,9209,9221,9227,9239,9241,9257,9277,9281,9283,9293,9311,9319,9323,9337,9341,9343,9349,9371,9377,9391,9397,9403,9413,9419,9421,9431,9433,9437,9439,9461,9463,9467,9473,9479,9491,9497,9511,9521,9533,9539,9547,9551,9587,9601,9613,9619,9623,9629,9631,9643,9649,9661,9677,9679,9689,9697,9719,9721,9733,9739,9743,9749,9767,9769,9781,9787,9791,9803,9811,9817,9829,9833,9839,9851,9857,9859,9871,9883,9887,9901,9907,9923,9929,9931,9941,9949,9967,9973};
   
     if (pt <N799) {
           return tb[pt++];
     }else{
           pt = 0;
           return tb[pt++];
     }
   }
   
   int getPrime(int p) {
     int i;
     if (p <= 2) return nextPrime();
     for (i=2; i<p; i++) {
           if (p % i == 0) {
             return nextPrime();
           }
     }
     return p;
 }  }
Legend:
Removed from v.1.5  
changed lines
  Added in v.1.12
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>