Changeset a2d993 in git for libpolys/polys/simpleideals.cc


Ignore:
Timestamp:
Apr 15, 2011, 11:32:40 AM (13 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
3d65bed222aced007a5c042b50c9f815f724b2e5
Parents:
2f5547f23ef54968be591558663090e0f192f644
git-author:
Hans Schoenemann <hannes@mathematik.uni-kl.de>2011-04-15 11:32:40+02:00
git-committer:
Mohamed Barakat <mohamed.barakat@rwth-aachen.de>2011-11-09 12:30:33+01:00
Message:
fix simpleideals
File:
1 edited

Legend:

Unmodified
Added
Removed

  • libpolys/polys/simpleideals.cc

    r2f5547 ra2d993  
    1414#include <omalloc/omalloc.h>
    1515#include <polys/monomials/p_polys.h>
     16#include <polys/weight.h>
     17#include <polys/matpol.h>
    1618#include <misc/intvec.h>
    1719#include <polys/simpleideals.h>
     20#include "sbuckets.h"
    1821
    1922static poly * idpower;
    … …  
    790793
    791794/*2
    792 *the minimal index of used variables - 1
    793 */
    794 int p_LowVar (poly p, const ring r)
    795 {
    796   int k,l,lex;
    797 
    798   if (p == NULL) return -1;
    799 
    800   k = 32000;/*a very large dummy value*/
    801   while (p != NULL)
    802   {
    803     l = 1;
    804     lex = p_GetExp(p,l,r);
    805     while ((l < rVar(r)) && (lex == 0))
    806     {
    807       l++;
    808       lex = p_GetExp(p,l,r);
    809     }
    810     l--;
    811     if (l < k) k = l;
    812     pIter(p);
    813   }
    814   return k;
    815 }
    816 
    817 /*2
    818795*initialized a field with r numbers between beg and end for the
    819796*procedure idNextChoise
    … …  
    10851062//  return id;
    10861063//}
    1087 static void idNextPotence(ideal given, ideal result,
    1088   int begin, int end, int deg, int restdeg, poly ap)
     1064static void id_NextPotence(ideal given, ideal result,
     1065  int begin, int end, int deg, int restdeg, poly ap, const ring r)
    10891066{
    10901067  poly p;
    10911068  int i;
    10921069
    1093   p = pPower(pCopy(given->m[begin]),restdeg);
     1070  p = p_Power(p_Copy(given->m[begin],r),restdeg,r);
    10941071  i = result->nrows;
    1095   result->m[i] = pMult(pCopy(ap),p);
     1072  result->m[i] = p_Mult_q(p_Copy(ap,r),p,r);
    10961073//PrintS(".");
    10971074  (result->nrows)++;
    … …  
    11041081  for (i=restdeg-1;i>0;i--)
    11051082  {
    1106     p = pPower(pCopy(given->m[begin]),i);
    1107     p = pMult(pCopy(ap),p);
    1108     idNextPotence(given, result, begin+1, end, deg, restdeg-i, p);
    1109     pDelete(&p);
    1110   }
    1111   idNextPotence(given, result, begin+1, end, deg, restdeg, ap);
     1083    p = p_Power(p_Copy(given->m[begin],r),i,r);
     1084    p = p_Mult_q(p_Copy(ap,r),p,r);
     1085    id_NextPotence(given, result, begin+1, end, deg, restdeg-i, p,r);
     1086    p_Delete(&p,r);
     1087  }
     1088  id_NextPotence(given, result, begin+1, end, deg, restdeg, ap,r);
    11121089}
    11131090
    … …  
    11261103//Print("ideal contains %d elements\n",i);
    11271104  p1=p_One(r);
    1128   idNextPotence(temp,result,0,IDELEMS(temp)-1,exp,exp,p1);
     1105  id_NextPotence(temp,result,0,IDELEMS(temp)-1,exp,exp,p1,r);
    11291106  p_Delete(&p1,r);
    11301107  id_Delete(&temp,r);
    … …  
    11361113
    11371114/*2
    1138 * compute the which-th ar-minor of the matrix a
    1139 */
    1140 poly idMinor(matrix a, int ar, unsigned long which, ideal R)
    1141 {
    1142   int     i,j,k,size;
    1143   unsigned long curr;
    1144   int *rowchoise,*colchoise;
    1145   BOOLEAN rowch,colch;
    1146   ideal result;
    1147   matrix tmp;
    1148   poly p,q;
    1149 
    1150   i = binom(a->rows(),ar);
    1151   j = binom(a->cols(),ar);
    1152 
    1153   rowchoise=(int *)omAlloc(ar*sizeof(int));
    1154   colchoise=(int *)omAlloc(ar*sizeof(int));
    1155   if ((i>512) || (j>512) || (i*j >512)) size=512;
    1156   else size=i*j;
    1157   result=idInit(size,1);
    1158   tmp=mpNew(ar,ar);
    1159   k = 0; /* the index in result*/
    1160   curr = 0; /* index of current minor */
    1161   idInitChoise(ar,1,a->rows(),&rowch,rowchoise);
    1162   while (!rowch)
    1163   {
    1164     idInitChoise(ar,1,a->cols(),&colch,colchoise);
    1165     while (!colch)
    1166     {
    1167       if (curr == which)
    1168       {
    1169         for (i=1; i<=ar; i++)
    1170         {
    1171           for (j=1; j<=ar; j++)
    1172           {
    1173             MATELEM(tmp,i,j) = MATELEM(a,rowchoise[i-1],colchoise[j-1]);
    1174           }
    1175         }
    1176         p = mpDetBareiss(tmp);
    1177         if (p!=NULL)
    1178         {
    1179           if (R!=NULL)
    1180           {
    1181             q = p;
    1182             p = kNF(R,currQuotient,q);
    1183             pDelete(&q);
    1184           }
    1185           /*delete the matrix tmp*/
    1186           for (i=1; i<=ar; i++)
    1187           {
    1188             for (j=1; j<=ar; j++) MATELEM(tmp,i,j) = NULL;
    1189           }
    1190           idDelete((ideal*)&tmp);
    1191           omFreeSize((ADDRESS)rowchoise,ar*sizeof(int));
    1192           omFreeSize((ADDRESS)colchoise,ar*sizeof(int));
    1193           return (p);
    1194         }
    1195       }
    1196       curr++;
    1197       idGetNextChoise(ar,a->cols(),&colch,colchoise);
    1198     }
    1199     idGetNextChoise(ar,a->rows(),&rowch,rowchoise);
    1200   }
    1201   return (poly) 1;
    1202 }
    1203 
    1204 #ifdef WITH_OLD_MINOR
    1205 /*2
    1206 * compute all ar-minors of the matrix a
    1207 */
    1208 ideal idMinors(matrix a, int ar, ideal R)
    1209 {
    1210   int     i,j,k,size;
    1211   int *rowchoise,*colchoise;
    1212   BOOLEAN rowch,colch;
    1213   ideal result;
    1214   matrix tmp;
    1215   poly p,q;
    1216 
    1217   i = binom(a->rows(),ar);
    1218   j = binom(a->cols(),ar);
    1219 
    1220   rowchoise=(int *)omAlloc(ar*sizeof(int));
    1221   colchoise=(int *)omAlloc(ar*sizeof(int));
    1222   if ((i>512) || (j>512) || (i*j >512)) size=512;
    1223   else size=i*j;
    1224   result=idInit(size,1);
    1225   tmp=mpNew(ar,ar);
    1226   k = 0; /* the index in result*/
    1227   idInitChoise(ar,1,a->rows(),&rowch,rowchoise);
    1228   while (!rowch)
    1229   {
    1230     idInitChoise(ar,1,a->cols(),&colch,colchoise);
    1231     while (!colch)
    1232     {
    1233       for (i=1; i<=ar; i++)
    1234       {
    1235         for (j=1; j<=ar; j++)
    1236         {
    1237           MATELEM(tmp,i,j) = MATELEM(a,rowchoise[i-1],colchoise[j-1]);
    1238         }
    1239       }
    1240       p = mpDetBareiss(tmp);
    1241       if (p!=NULL)
    1242       {
    1243         if (R!=NULL)
    1244         {
    1245           q = p;
    1246           p = kNF(R,currQuotient,q);
    1247           pDelete(&q);
    1248         }
    1249         if (p!=NULL)
    1250         {
    1251           if (k>=size)
    1252           {
    1253             pEnlargeSet(&result->m,size,32);
    1254             size += 32;
    1255           }
    1256           result->m[k] = p;
    1257           k++;
    1258         }
    1259       }
    1260       idGetNextChoise(ar,a->cols(),&colch,colchoise);
    1261     }
    1262     idGetNextChoise(ar,a->rows(),&rowch,rowchoise);
    1263   }
    1264   /*delete the matrix tmp*/
    1265   for (i=1; i<=ar; i++)
    1266   {
    1267     for (j=1; j<=ar; j++) MATELEM(tmp,i,j) = NULL;
    1268   }
    1269   idDelete((ideal*)&tmp);
    1270   if (k==0)
    1271   {
    1272     k=1;
    1273     result->m[0]=NULL;
    1274   }
    1275   omFreeSize((ADDRESS)rowchoise,ar*sizeof(int));
    1276   omFreeSize((ADDRESS)colchoise,ar*sizeof(int));
    1277   pEnlargeSet(&result->m,size,k-size);
    1278   IDELEMS(result) = k;
    1279   return (result);
    1280 }
    1281 #else
    1282 /*2
    1283 * compute all ar-minors of the matrix a
    1284 * the caller of mpRecMin
    1285 * the elements of the result are not in R (if R!=NULL)
    1286 */
    1287 ideal idMinors(matrix a, int ar, ideal R)
    1288 {
    1289   int elems=0;
    1290   int r=a->nrows,c=a->ncols;
    1291   int i;
    1292   matrix b;
    1293   ideal result,h;
    1294   ring origR;
    1295   ring tmpR;
    1296   long bound;
    1297 
    1298   if((ar<=0) || (ar>r) || (ar>c))
    1299   {
    1300     Werror("%d-th minor, matrix is %dx%d",ar,r,c);
    1301     return NULL;
    1302   }
    1303   h = idMatrix2Module(mpCopy(a));
    1304   bound = smExpBound(h,c,r,ar);
    1305   idDelete(&h);
    1306   tmpR=smRingChange(&origR,bound);
    1307   b = mpNew(r,c);
    1308   for (i=r*c-1;i>=0;i--)
    1309   {
    1310     if (a->m[i])
    1311       b->m[i] = prCopyR(a->m[i],origR);
    1312   }
    1313   if (R!=NULL)
    1314   {
    1315     R = idrCopyR(R,origR);
    1316     //if (ar>1) // otherwise done in mpMinorToResult
    1317     //{
    1318     //  matrix bb=(matrix)kNF(R,currQuotient,(ideal)b);
    1319     //  bb->rank=b->rank; bb->nrows=b->nrows; bb->ncols=b->ncols;
    1320     //  idDelete((ideal*)&b); b=bb;
    1321     //}
    1322   }
    1323   result=idInit(32,1);
    1324   if(ar>1) mpRecMin(ar-1,result,elems,b,r,c,NULL,R);
    1325   else mpMinorToResult(result,elems,b,r,c,R);
    1326   idDelete((ideal *)&b);
    1327   if (R!=NULL) idDelete(&R);
    1328   idSkipZeroes(result);
    1329   rChangeCurrRing(origR);
    1330   result = idrMoveR(result,tmpR);
    1331   smKillModifiedRing(tmpR);
    1332   idTest(result);
    1333   return result;
    1334 }
    1335 #endif
    1336 
    1337 /*2
    13381115*skips all zeroes and double elements, searches also for units
    13391116*/
    … …  
    13621139
    13631140/*2
    1364 *returns TRUE if id1 is a submodule of id2
    1365 */
    1366 BOOLEAN idIsSubModule(ideal id1,ideal id2)
    1367 {
    1368   int i;
    1369   poly p;
    1370 
    1371   if (idIs0(id1)) return TRUE;
    1372   for (i=0;i<IDELEMS(id1);i++)
    1373   {
    1374     if (id1->m[i] != NULL)
    1375     {
    1376       p = kNF(id2,currQuotient,id1->m[i]);
    1377       if (p != NULL)
    1378       {
    1379         pDelete(&p);
    1380         return FALSE;
    1381       }
    1382     }
    1383   }
    1384   return TRUE;
    1385 }
    1386 
    1387 /*2
    13881141* returns the ideals of initial terms
    13891142*/
    1390 ideal idHead(ideal h)
     1143ideal id_Head(ideal h,const ring r)
    13911144{
    13921145  ideal m = idInit(IDELEMS(h),h->rank);
    … …  
    13951148  for (i=IDELEMS(h)-1;i>=0; i--)
    13961149  {
    1397     if (h->m[i]!=NULL) m->m[i]=pHead(h->m[i]);
     1150    if (h->m[i]!=NULL) m->m[i]=p_Head(h->m[i],r);
    13981151  }
    13991152  return m;
    14001153}
    14011154
    1402 ideal idHomogen(ideal h, int varnum)
     1155ideal id_Homogen(ideal h, int varnum,const ring r)
    14031156{
    14041157  ideal m = idInit(IDELEMS(h),h->rank);
    … …  
    14071160  for (i=IDELEMS(h)-1;i>=0; i--)
    14081161  {
    1409     m->m[i]=pHomogen(h->m[i],varnum);
     1162    m->m[i]=p_Homogen(h->m[i],varnum,r);
    14101163  }
    14111164  return m;
    … …  
    14131166
    14141167/*------------------type conversions----------------*/
    1415 ideal idVec2Ideal(poly vec)
     1168ideal id_Vec2Ideal(poly vec, const ring R)
    14161169{
    14171170   ideal result=idInit(1,1);
    14181171   omFree((ADDRESS)result->m);
    14191172   result->m=NULL; // remove later
    1420    pVec2Polys(vec, &(result->m), &(IDELEMS(result)));
     1173   p_Vec2Polys(vec, &(result->m), &(IDELEMS(result)),R);
    14211174   return result;
    14221175}
    14231176
    1424 #define NEW_STUFF
    1425 #ifndef NEW_STUFF
     1177
    14261178// converts mat to module, destroys mat
    1427 ideal idMatrix2Module(matrix mat)
    1428 {
    1429   int mc=MATCOLS(mat);
    1430   int mr=MATROWS(mat);
    1431   ideal result = idInit(si_max(mc,1),si_max(mr,1));
    1432   int i,j;
    1433   poly h;
    1434 
    1435   for(j=0;j<mc /*MATCOLS(mat)*/;j++) /* j is also index in result->m */
    1436   {
    1437     for (i=1;i<=mr /*MATROWS(mat)*/;i++)
    1438     {
    1439       h = MATELEM(mat,i,j+1);
    1440       if (h!=NULL)
    1441       {
    1442         MATELEM(mat,i,j+1)=NULL;
    1443         pSetCompP(h,i);
    1444         result->m[j] = pAdd(result->m[j],h);
    1445       }
    1446     }
    1447   }
    1448   // obachman: need to clean this up
    1449   idDelete((ideal*) &mat);
    1450   return result;
    1451 }
    1452 #else
    1453 
    1454 #include "sbuckets.h"
    1455 
    1456 // converts mat to module, destroys mat
    1457 ideal idMatrix2Module(matrix mat)
     1179ideal id_Matrix2Module(matrix mat, const ring R)
    14581180{
    14591181  int mc=MATCOLS(mat);
    … …  
    14631185  poly h;
    14641186  poly p;
    1465   sBucket_pt bucket = sBucketCreate(currRing);
     1187  sBucket_pt bucket = sBucketCreate(R);
    14661188
    14671189  for(j=0;j<mc /*MATCOLS(mat)*/;j++) /* j is also index in result->m */
    … …  
    14741196        l=pLength(h);
    14751197        MATELEM(mat,i,j+1)=NULL;
    1476         p_SetCompP(h,i, currRing);
     1198        p_SetCompP(h,i, R);
    14771199        sBucket_Merge_p(bucket, h, l);
    14781200      }
    … …  
    14831205
    14841206  // obachman: need to clean this up
    1485   idDelete((ideal*) &mat);
     1207  id_Delete((ideal*) &mat,R);
    14861208  return result;
    14871209}
    1488 #endif
    14891210
    14901211/*2
    14911212* converts a module into a matrix, destroyes the input
    14921213*/
    1493 matrix idModule2Matrix(ideal mod)
     1214matrix id_Module2Matrix(ideal mod, const ring R)
    14941215{
    14951216  matrix result = mpNew(mod->rank,IDELEMS(mod));
    … …  
    15071228      pNext(h)=NULL;
    15081229//      cp = si_max(1,pGetComp(h));     // if used for ideals too
    1509       cp = pGetComp(h);
    1510       pSetComp(h,0);
    1511       pSetmComp(h);
     1230      cp = p_GetComp(h,R);
     1231      p_SetComp(h,0,R);
     1232      p_SetmComp(h,R);
    15121233#ifdef TEST
    15131234      if (cp>mod->rank)
    … …  
    15251246          }
    15261247        }
    1527         idDelete((ideal *)&result);
     1248        id_Delete((ideal *)&result,R);
    15281249        result=d;
    15291250      }
    15301251#endif
    1531       MATELEM(result,cp,i+1) = pAdd(MATELEM(result,cp,i+1),h);
     1252      MATELEM(result,cp,i+1) = p_Add_q(MATELEM(result,cp,i+1),h,R);
    15321253    }
    15331254  }
    15341255  // obachman 10/99: added the following line, otherwise memory leack!
    1535   idDelete(&mod);
     1256  id_Delete(&mod,R);
    15361257  return result;
    15371258}
    15381259
    1539 matrix idModule2formatedMatrix(ideal mod,int rows, int cols)
     1260matrix id_Module2formatedMatrix(ideal mod,int rows, int cols, const ring R)
    15401261{
    15411262  matrix result = mpNew(rows,cols);
    1542   int i,cp,r=idRankFreeModule(mod),c=IDELEMS(mod);
     1263  int i,cp,r=id_RankFreeModule(mod,R),c=IDELEMS(mod);
    15431264  poly p,h;
    15441265
    … …  
    15541275      pIter(p);
    15551276      pNext(h)=NULL;
    1556       cp = pGetComp(h);
     1277      cp = p_GetComp(h,R);
    15571278      if (cp<=r)
    15581279      {
    1559         pSetComp(h,0);
    1560         pSetmComp(h);
    1561         MATELEM(result,cp,i+1) = pAdd(MATELEM(result,cp,i+1),h);
     1280        p_SetComp(h,0,R);
     1281        p_SetmComp(h,R);
     1282        MATELEM(result,cp,i+1) = p_Add_q(MATELEM(result,cp,i+1),h,R);
    15621283      }
    15631284      else
    1564         pDelete(&h);
    1565     }
    1566   }
    1567   idDelete(&mod);
     1285        p_Delete(&h,R);
     1286    }
     1287  }
     1288  id_Delete(&mod,R);
    15681289  return result;
    15691290}
    … …  
    15731294* destroy id
    15741295*/
    1575 ideal  idSubst(ideal id, int n, poly e)
     1296ideal  id_Subst(ideal id, int n, poly e, const ring r)
    15761297{
    15771298  int k=MATROWS((matrix)id)*MATCOLS((matrix)id);
    … …  
    15811302  for(k--;k>=0;k--)
    15821303  {
    1583     res->m[k]=pSubst(id->m[k],n,e);
     1304    res->m[k]=p_Subst(id->m[k],n,e,r);
    15841305    id->m[k]=NULL;
    15851306  }
    1586   idDelete(&id);
     1307  id_Delete(&id,r);
    15871308  return res;
    15881309}
    15891310
    1590 BOOLEAN idHomModule(ideal m, ideal Q, intvec **w)
     1311BOOLEAN id_HomModule(ideal m, ideal Q, intvec **w, const ring R)
    15911312{
    15921313  if (w!=NULL) *w=NULL;
    1593   if ((Q!=NULL) && (!idHomIdeal(Q,NULL))) return FALSE;
     1314  if ((Q!=NULL) && (!id_HomIdeal(Q,NULL,R))) return FALSE;
    15941315  if (idIs0(m))
    15951316  {
    … …  
    16031324  poly p=NULL;
    16041325  pFDegProc d;
    1605   if (pLexOrder && (currRing->order[0]==ringorder_lp))
     1326  if (R->pLexOrder && (R->order[0]==ringorder_lp))
    16061327     d=p_Totaldegree;
    16071328  else
    16081329     d=pFDeg;
    16091330  int length=IDELEMS(m);
    1610   polyset P=m->m;
    1611   polyset F=(polyset)omAlloc(length*sizeof(poly));
     1331  poly* P=m->m;
     1332  poly* F=(poly*)omAlloc(length*sizeof(poly));
    16121333  for (i=length-1;i>=0;i--)
    16131334  {
    16141335    p=F[i]=P[i];
    1615     cmax=si_max(cmax,(long)pMaxComp(p));
     1336    cmax=si_max(cmax,(long)p_MaxComp(p,R));
    16161337  }
    16171338  cmax++;
    … …  
    16251346    {
    16261347      p=F[i];
    1627       while ((p!=NULL) && (iscom[pGetComp(p)]==0)) pIter(p);
     1348      while ((p!=NULL) && (iscom[p_GetComp(p,R)]==0)) pIter(p);
    16281349    }
    16291350    if ((p==NULL) && (i<length))
    … …  
    16451366      //  order = p->order;
    16461367      //  order = pFDeg(p,currRing);
    1647       order = d(p,currRing) +diff[pGetComp(p)];
     1368      order = d(p,R) +diff[p_GetComp(p,R)];
    16481369      //order += diff[pGetComp(p)];
    16491370      p = F[i];
    … …  
    16541375    while (p!=NULL)
    16551376    {
    1656       if (pLexOrder && (currRing->order[0]==ringorder_lp))
    1657         ord=pTotaldegree(p);
     1377      if (R->pLexOrder && (R->order[0]==ringorder_lp))
     1378        ord=p_Totaldegree(p,R);
    16581379      else
    16591380      //  ord = p->order;
    1660         ord = pFDeg(p,currRing);
    1661       if (iscom[pGetComp(p)]==0)
    1662       {
    1663         diff[pGetComp(p)] = order-ord;
    1664         iscom[pGetComp(p)] = 1;
     1381        ord = pFDeg(p,R);
     1382      if (iscom[p_GetComp(p,R)]==0)
     1383      {
     1384        diff[p_GetComp(p,R)] = order-ord;
     1385        iscom[p_GetComp(p,R)] = 1;
    16651386/*
    16661387*PrintS("new diff: ");
    … …  
    16811402*Print("order %d, ord %d, diff %d\n",order,ord,diff[pGetComp(p)]);
    16821403*/
    1683         if (order != (ord+diff[pGetComp(p)]))
     1404        if (order != (ord+diff[p_GetComp(p,R)]))
    16841405        {
    16851406          omFreeSize((ADDRESS) iscom,cmax*sizeof(int));
    … …  
    17111432}
    17121433
    1713 BOOLEAN idTestHomModule(ideal m, ideal Q, intvec *w)
    1714 {
    1715   if ((Q!=NULL) && (!idHomIdeal(Q,NULL)))  { PrintS(" Q not hom\n"); return FALSE;}
    1716   if (idIs0(m)) return TRUE;
    1717 
    1718   int cmax=-1;
    1719   int i;
    1720   poly p=NULL;
    1721   int length=IDELEMS(m);
    1722   polyset P=m->m;
    1723   for (i=length-1;i>=0;i--)
    1724   {
    1725     p=P[i];
    1726     if (p!=NULL) cmax=si_max(cmax,(int)pMaxComp(p)+1);
    1727   }
    1728   if (w != NULL)
    1729   if (w->length()+1 < cmax)
    1730   {
    1731     // Print("length: %d - %d \n", w->length(),cmax);
    1732     return FALSE;
    1733   }
    1734 
    1735   if(w!=NULL)
    1736     pSetModDeg(w);
    1737 
    1738   for (i=length-1;i>=0;i--)
    1739   {
    1740     p=P[i];
    1741     poly q=p;
    1742     if (p!=NULL)
    1743     {
    1744       int d=pFDeg(p,currRing);
    1745       loop
    1746       {
    1747         pIter(p);
    1748         if (p==NULL) break;
    1749         if (d!=pFDeg(p,currRing))
    1750         {
    1751           //pWrite(q); wrp(p); Print(" -> %d - %d\n",d,pFDeg(p,currRing));
    1752           if(w!=NULL)
    1753             pSetModDeg(NULL);
    1754           return FALSE;
    1755         }
    1756       }
    1757     }
    1758   }
    1759 
    1760   if(w!=NULL)
    1761     pSetModDeg(NULL);
    1762 
    1763   return TRUE;
    1764 }
    1765 
    1766 ideal idJet(ideal i,int d)
     1434// uses glabl vars via pSetModDeg
     1435//BOOLEAN idTestHomModule(ideal m, ideal Q, intvec *w)
     1436//{
     1437//  if ((Q!=NULL) && (!idHomIdeal(Q,NULL)))  { PrintS(" Q not hom\n"); return FALSE;}
     1438//  if (idIs0(m)) return TRUE;
     1439//
     1440//  int cmax=-1;
     1441//  int i;
     1442//  poly p=NULL;
     1443//  int length=IDELEMS(m);
     1444//  poly* P=m->m;
     1445//  for (i=length-1;i>=0;i--)
     1446//  {
     1447//    p=P[i];
     1448//    if (p!=NULL) cmax=si_max(cmax,(int)pMaxComp(p)+1);
     1449//  }
     1450//  if (w != NULL)
     1451//  if (w->length()+1 < cmax)
     1452//  {
     1453//    // Print("length: %d - %d \n", w->length(),cmax);
     1454//    return FALSE;
     1455//  }
     1456//
     1457//  if(w!=NULL)
     1458//    pSetModDeg(w);
     1459//
     1460//  for (i=length-1;i>=0;i--)
     1461//  {
     1462//    p=P[i];
     1463//    poly q=p;
     1464//    if (p!=NULL)
     1465//    {
     1466//      int d=pFDeg(p,currRing);
     1467//      loop
     1468//      {
     1469//        pIter(p);
     1470//        if (p==NULL) break;
     1471//        if (d!=pFDeg(p,currRing))
     1472//        {
     1473//          //pWrite(q); wrp(p); Print(" -> %d - %d\n",d,pFDeg(p,currRing));
     1474//          if(w!=NULL)
     1475//            pSetModDeg(NULL);
     1476//          return FALSE;
     1477//        }
     1478//      }
     1479//    }
     1480//  }
     1481//
     1482//  if(w!=NULL)
     1483//    pSetModDeg(NULL);
     1484//
     1485//  return TRUE;
     1486//}
     1487
     1488ideal id_Jet(ideal i,int d, const ring R)
    17671489{
    17681490  ideal r=idInit((i->nrows)*(i->ncols),i->rank);
    … …  
    17731495  for(k=(i->nrows)*(i->ncols)-1;k>=0; k--)
    17741496  {
    1775     r->m[k]=ppJet(i->m[k],d);
     1497    r->m[k]=pp_Jet(i->m[k],d,R);
    17761498  }
    17771499  return r;
    17781500}
    17791501
    1780 ideal idJetW(ideal i,int d, intvec * iv)
     1502ideal id_JetW(ideal i,int d, intvec * iv, const ring R)
    17811503{
    17821504  ideal r=idInit(IDELEMS(i),i->rank);
    … …  
    17871509  else
    17881510  {
    1789     short *w=iv2array(iv);
     1511    short *w=iv2array(iv,R);
    17901512    int k;
    17911513    for(k=0; k<IDELEMS(i); k++)
    17921514    {
    1793       r->m[k]=ppJetW(i->m[k],d,w);
    1794     }
    1795     omFreeSize((ADDRESS)w,(rVar(r)+1)*sizeof(short));
     1515      r->m[k]=pp_JetW(i->m[k],d,w,R);
     1516    }
     1517    omFreeSize((ADDRESS)w,(rVar(R)+1)*sizeof(short));
    17961518  }
    17971519  return r;
    1798 }
    1799 
    1800 int idMinDegW(ideal M,intvec *w)
    1801 {
    1802   int d=-1;
    1803   for(int i=0;i<IDELEMS(M);i++)
    1804   {
    1805     int d0=pMinDeg(M->m[i],w);
    1806     if(-1<d0&&(d0<d||d==-1))
    1807       d=d0;
    1808   }
    1809   return d;
    1810 }
    1811 
    1812 ideal idSeries(int n,ideal M,matrix U,intvec *w)
    1813 {
    1814   for(int i=IDELEMS(M)-1;i>=0;i--)
    1815   {
    1816     if(U==NULL)
    1817       M->m[i]=pSeries(n,M->m[i],NULL,w);
    1818     else
    1819     {
    1820       M->m[i]=pSeries(n,M->m[i],MATELEM(U,i+1,i+1),w);
    1821       MATELEM(U,i+1,i+1)=NULL;
    1822     }
    1823   }
    1824   if(U!=NULL)
    1825     idDelete((ideal*)&U);
    1826   return M;
    1827 }
    1828 
    1829 matrix idDiff(matrix i, int k)
    1830 {
    1831   int e=MATCOLS(i)*MATROWS(i);
    1832   matrix r=mpNew(MATROWS(i),MATCOLS(i));
    1833   r->rank=i->rank;
    1834   int j;
    1835   for(j=0; j<e; j++)
    1836   {
    1837     r->m[j]=pDiff(i->m[j],k);
    1838   }
    1839   return r;
    1840 }
    1841 
    1842 matrix idDiffOp(ideal I, ideal J,BOOLEAN multiply)
    1843 {
    1844   matrix r=mpNew(IDELEMS(I),IDELEMS(J));
    1845   int i,j;
    1846   for(i=0; i<IDELEMS(I); i++)
    1847   {
    1848     for(j=0; j<IDELEMS(J); j++)
    1849     {
    1850       MATELEM(r,i+1,j+1)=pDiffOp(I->m[i],J->m[j],multiply);
    1851     }
    1852   }
    1853   return r;
    1854 }
    1855 
    1856 int idElem(const ideal F)
    1857 {
    1858   int i=0,j=IDELEMS(F)-1;
    1859 
    1860   while(j>=0)
    1861   {
    1862     if ((F->m)[j]!=NULL) i++;
    1863     j--;
    1864   }
    1865   return i;
    1866 }
    1867 
    1868 /*
    1869 *computes module-weights for liftings of homogeneous modules
    1870 */
    1871 intvec * idMWLift(ideal mod,intvec * weights)
    1872 {
    1873   if (idIs0(mod)) return new intvec(2);
    1874   int i=IDELEMS(mod);
    1875   while ((i>0) && (mod->m[i-1]==NULL)) i--;
    1876   intvec *result = new intvec(i+1);
    1877   while (i>0)
    1878   {
    1879     (*result)[i]=pFDeg(mod->m[i],currRing)+(*weights)[pGetComp(mod->m[i])];
    1880   }
    1881   return result;
    1882 }
    1883 
    1884 /*2
    1885 *sorts the kbase for idCoef* in a special way (lexicographically
    1886 *with x_max,...,x_1)
    1887 */
    1888 ideal idCreateSpecialKbase(ideal kBase,intvec ** convert)
    1889 {
    1890   int i;
    1891   ideal result;
    1892 
    1893   if (idIs0(kBase)) return NULL;
    1894   result = idInit(IDELEMS(kBase),kBase->rank);
    1895   *convert = idSort(kBase,FALSE);
    1896   for (i=0;i<(*convert)->length();i++)
    1897   {
    1898     result->m[i] = pCopy(kBase->m[(**convert)[i]-1]);
    1899   }
    1900   return result;
    1901 }
    1902 
    1903 /*2
    1904 *returns the index of a given monom in the list of the special kbase
    1905 */
    1906 int idIndexOfKBase(poly monom, ideal kbase)
    1907 {
    1908   int j=IDELEMS(kbase);
    1909 
    1910   while ((j>0) && (kbase->m[j-1]==NULL)) j--;
    1911   if (j==0) return -1;
    1912   int i=rVar(r);
    1913   while (i>0)
    1914   {
    1915     loop
    1916     {
    1917       if (pGetExp(monom,i)>pGetExp(kbase->m[j-1],i)) return -1;
    1918       if (pGetExp(monom,i)==pGetExp(kbase->m[j-1],i)) break;
    1919       j--;
    1920       if (j==0) return -1;
    1921     }
    1922     if (i==1)
    1923     {
    1924       while(j>0)
    1925       {
    1926         if (pGetComp(monom)==pGetComp(kbase->m[j-1])) return j-1;
    1927         if (pGetComp(monom)>pGetComp(kbase->m[j-1])) return -1;
    1928         j--;
    1929       }
    1930     }
    1931     i--;
    1932   }
    1933   return -1;
    1934 }
    1935 
    1936 /*2
    1937 *decomposes the monom in a part of coefficients described by the
    1938 *complement of how and a monom in variables occuring in how, the
    1939 *index of which in kbase is returned as integer pos (-1 if it don't
    1940 *exists)
    1941 */
    1942 poly idDecompose(poly monom, poly how, ideal kbase, int * pos)
    1943 {
    1944   int i;
    1945   poly coeff=p_One(r), base=p_One(r);
    1946 
    1947   for (i=1;i<=rVar(r);i++)
    1948   {
    1949     if (pGetExp(how,i)>0)
    1950     {
    1951       pSetExp(base,i,pGetExp(monom,i));
    1952     }
    1953     else
    1954     {
    1955       pSetExp(coeff,i,pGetExp(monom,i));
    1956     }
    1957   }
    1958   pSetComp(base,pGetComp(monom));
    1959   pSetm(base);
    1960   pSetCoeff(coeff,nCopy(pGetCoeff(monom)));
    1961   pSetm(coeff);
    1962   *pos = idIndexOfKBase(base,kbase);
    1963   if (*pos<0)
    1964     pDelete(&coeff);
    1965   pDelete(&base);
    1966   return coeff;
    1967 }
    1968 
    1969 /*2
    1970 *returns a matrix A of coefficients with kbase*A=arg
    1971 *if all monomials in variables of how occur in kbase
    1972 *the other are deleted
    1973 */
    1974 matrix idCoeffOfKBase(ideal arg, ideal kbase, poly how)
    1975 {
    1976   matrix result;
    1977   ideal tempKbase;
    1978   poly p,q;
    1979   intvec * convert;
    1980   int i=IDELEMS(kbase),j=IDELEMS(arg),k,pos;
    1981 #if 0
    1982   while ((i>0) && (kbase->m[i-1]==NULL)) i--;
    1983   if (idIs0(arg))
    1984     return mpNew(i,1);
    1985   while ((j>0) && (arg->m[j-1]==NULL)) j--;
    1986   result = mpNew(i,j);
    1987 #else
    1988   result = mpNew(i, j);
    1989   while ((j>0) && (arg->m[j-1]==NULL)) j--;
    1990 #endif
    1991 
    1992   tempKbase = idCreateSpecialKbase(kbase,&convert);
    1993   for (k=0;k<j;k++)
    1994   {
    1995     p = arg->m[k];
    1996     while (p!=NULL)
    1997     {
    1998       q = idDecompose(p,how,tempKbase,&pos);
    1999       if (pos>=0)
    2000       {
    2001         MATELEM(result,(*convert)[pos],k+1) =
    2002             pAdd(MATELEM(result,(*convert)[pos],k+1),q);
    2003       }
    2004       else
    2005         pDelete(&q);
    2006       pIter(p);
    2007     }
    2008   }
    2009   idDelete(&tempKbase);
    2010   return result;
    20111520}
    20121521
    … …  
    20971606#endif
    20981607
    2099 static void idDeleteComps(ideal arg,int* red_comp,int del)
    2100 // red_comp is an array [0..args->rank]
    2101 {
    2102   int i,j;
    2103   poly p;
    2104 
    2105   for (i=IDELEMS(arg)-1;i>=0;i--)
    2106   {
    2107     p = arg->m[i];
    2108     while (p!=NULL)
    2109     {
    2110       j = pGetComp(p);
    2111       if (red_comp[j]!=j)
    2112       {
    2113         pSetComp(p,red_comp[j]);
    2114         pSetmComp(p);
    2115       }
    2116       pIter(p);
    2117     }
    2118   }
    2119   (arg->rank) -= del;
    2120 }
    2121 
    2122 /*2
    2123 * returns the presentation of an isomorphic, minimally
    2124 * embedded  module (arg represents the quotient!)
    2125 */
    2126 ideal idMinEmbedding(ideal arg,BOOLEAN inPlace, intvec **w, const ring r)
    2127 {
    2128   if (idIs0(arg)) return idInit(1,arg->rank);
    2129   int i,next_gen,next_comp;
    2130   ideal res=arg;
    2131   if (!inPlace) res = id_Copy(arg,r);
    2132   res->rank=si_max(res->rank,id_RankFreeModule(res,r));
    2133   int *red_comp=(int*)omAlloc((res->rank+1)*sizeof(int));
    2134   for (i=res->rank;i>=0;i--) red_comp[i]=i;
    2135 
    2136   int del=0;
    2137   loop
    2138   {
    2139     next_gen = idReadOutPivot(res,&next_comp);
    2140     if (next_gen<0) break;
    2141     del++;
    2142     syGaussForOne(res,next_gen,next_comp,0,IDELEMS(res));
    2143     for(i=next_comp+1;i<=arg->rank;i++) red_comp[i]--;
    2144     if ((w !=NULL)&&(*w!=NULL))
    2145     {
    2146       for(i=next_comp;i<(*w)->length();i++) (**w)[i-1]=(**w)[i];
    2147     }
    2148   }
    2149 
    2150   idDeleteComps(res,red_comp,del);
    2151   idSkipZeroes(res);
    2152   omFree(red_comp);
    2153 
    2154   if ((w !=NULL)&&(*w!=NULL) &&(del>0))
    2155   {
    2156     intvec *wtmp=new intvec((*w)->length()-del);
    2157     for(i=0;i<res->rank;i++) (*wtmp)[i]=(**w)[i];
    2158     delete *w;
    2159     *w=wtmp;
    2160   }
    2161   return res;
    2162 }
    2163 
    2164 intvec * idQHomWeight(ideal id)
     1608intvec * id_QHomWeight(ideal id, const ring r)
    21651609{
    21661610  poly head, tail;
    … …  
    21811625        all++;
    21821626        for (k=1;k<=coldim;k++)
    2183           IMATELEM(*imat,all,k) = pGetExpDiff(head,tail,k);
     1627          IMATELEM(*imat,all,k) = p_GetExpDiff(head,tail,k,r);
    21841628        if (all==rowmax)
    21851629        {
    … …  
    22091653}
    22101654
    2211 BOOLEAN idIsZeroDim(ideal I)
     1655BOOLEAN id_IsZeroDim(ideal I, const ring r)
    22121656{
    22131657  BOOLEAN *UsedAxis=(BOOLEAN *)omAlloc0(rVar(r)*sizeof(BOOLEAN));
    … …  
    22181662  {
    22191663    po=I->m[i];
    2220     if ((po!=NULL) &&((n=pIsPurePower(po))!=0)) UsedAxis[n-1]=TRUE;
     1664    if ((po!=NULL) &&((n=p_IsPurePower(po,r))!=0)) UsedAxis[n-1]=TRUE;
    22211665  }
    22221666  for(i=rVar(r)-1;i>=0;i--)
    … …  
    23231767}
    23241768#endif
    2325 /* currently unsed:
    2326 ideal idChineseRemainder(ideal *xx, intvec *iv)
    2327 {
    2328   int rl=iv->length();
    2329   number *q=(number *)omAlloc(rl*sizeof(number));
    2330   int i;
    2331   for(i=0; i<rl; i++)
    2332   {
    2333     q[i]=nInit((*iv)[i]);
    2334   }
    2335   return idChineseRemainder(xx,q,rl);
    2336 }
    2337 */
    2338 /*
    2339  * lift ideal with coeffs over Z (mod N) to Q via Farey
    2340  */
    2341 ideal id_Farey(ideal x, number N, const ring r)
    2342 {
    2343   int cnt=IDELEMS(x)*x->nrows;
    2344   ideal result=idInit(cnt,x->rank);
    2345   result->nrows=x->nrows; // for lifting matrices
    2346   result->ncols=x->ncols; // for lifting matrices
    2347 
    2348   int i;
    2349   for(i=cnt-1;i>=0;i--)
    2350   {
    2351     poly h=p_Copy(x->m[i],r);
    2352     result->m[i]=h;
    2353     while(h!=NULL)
    2354     {
    2355       number c=pGetCoeff(h);
    2356       pSetCoeff0(h,nlFarey(c,N));
    2357       n_Delete(&c,r->cf);
    2358       pIter(h);
    2359     }
    2360     while((result->m[i]!=NULL)&&(n_IsZero(pGetCoeff(result->m[i]),r->cf)))
    2361     {
    2362       p_LmDelete(&(result->m[i]),r);
    2363     }
    2364     h=result->m[i];
    2365     while((h!=NULL) && (pNext(h)!=NULL))
    2366     {
    2367       if(n_IsZero(pGetCoeff(pNext(h)),r->cf))
    2368       {
    2369         p_LmDelete(&pNext(h),r);
    2370       }
    2371       else pIter(h);
    2372     }
    2373   }
    2374   return result;
    2375 }
    23761769
    23771770/*2
Note: See TracChangeset for help on using the changeset viewer.