Changeset ebdaa1 in git for libpolys/polys/sparsmat.cc


Ignore:
Timestamp:
Nov 9, 2011, 4:24:50 PM (12 years ago)
Author:
Mohamed Barakat <mohamed.barakat@…>
Branches:
(u'spielwiese', '8e0ad00ce244dfd0756200662572aef8402f13d5')
Children:
b2929e71ee342c64f5189c222795ac906cb21ee9
Parents:
8f40312428677d88ab37a31824f2ec53f144a18a
Message:
minor fixes after last rebase

make works
make check works
File:
1 edited

Legend:

Unmodified
Added
Removed

  • libpolys/polys/sparsmat.cc

    r8f4031 rebdaa1  
    6060
    6161
     62/* declare internal 'C' stuff */
     63static void sm_ExactPolyDiv(poly, poly, const ring);
     64static BOOLEAN sm_IsNegQuot(poly, const poly, const poly, const ring);
     65static void sm_ExpMultDiv(poly, const poly, const poly, const ring);
     66static void sm_PolyDivN(poly, const number, const ring);
     67static BOOLEAN smSmaller(poly, poly);
     68static void sm_CombineChain(poly *, poly, const ring);
     69static void sm_FindRef(poly *, poly *, poly, const ring);
     70
     71static void sm_ElemDelete(smpoly *, const ring);
     72static smpoly smElemCopy(smpoly);
     73static float sm_PolyWeight(smpoly, const ring);
     74static smpoly sm_Poly2Smpoly(poly, const ring);
     75static poly sm_Smpoly2Poly(smpoly, const ring);
     76static BOOLEAN sm_HaveDenom(poly, const ring);
     77static number sm_Cleardenom(ideal, const ring);
     78
    6279omBin smprec_bin = omGetSpecBin(sizeof(smprec));
    63 
    64 static void smExpMultDiv(poly t, const poly b, const poly c)
    65 {
    66   int i;
    67   pTest(t);
    68   pLmTest(b);
    69   pLmTest(c);
    70   poly bc = p_New(currRing);
    71 
    72   p_ExpVectorDiff(bc, b, c, currRing);
    73 
    74   while(t!=NULL)
    75   {
    76     pExpVectorAdd(t, bc);
    77     pIter(t);
    78   }
    79   p_LmFree(bc, currRing);
    80 }
    8180
    8281static poly pp_Mult_Coeff_mm_DivSelect_MultDiv(poly p, int &lp, poly m,
    … …  
    353352
    354353/* ----------------- basics (used from 'C') ------------------ */
    355 static BOOLEAN smHaveDenom(poly a)
    356 {
    357   BOOLEAN sw;
    358   number x;
    359 
    360   while (a != NULL)
    361   {
    362     x = nGetDenom(pGetCoeff(a));
    363     sw = nIsOne(x);
    364     nDelete(&x);
    365     if (!sw)
    366     {
    367       return TRUE;
    368     }
    369     pIter(a);
    370   }
    371   return FALSE;
    372 }
    373 
    374 static number smCleardenom(ideal id)
    375 {
    376   poly a;
    377   number x,y,res=nInit(1);
    378   BOOLEAN sw=FALSE;
    379 
    380   for (int i=0; i<IDELEMS(id); i++)
    381   {
    382     a = id->m[i];
    383     sw = smHaveDenom(a);
    384     if (sw) break;
    385   }
    386   if (!sw) return res;
    387   for (int i=0; i<IDELEMS(id); i++)
    388   {
    389     a = id->m[i];
    390     if (a!=NULL)
    391     {
    392       x = nCopy(pGetCoeff(a));
    393       p_Cleardenom(a, currRing);
    394       y = nDiv(x,pGetCoeff(a));
    395       nDelete(&x);
    396       x = nMult(res,y);
    397       nNormalize(x);
    398       nDelete(&res);
    399       res = x;
    400     }
    401   }
    402   return res;
    403 }
    404 
    405354/*2
    406355*returns the determinant of the module I;
    … …  
    486435  sm_KillModifiedRing(tmpR);
    487436  M=II;
    488 }
    489 
    490 /*
    491 * from poly to smpoly
    492 * do not destroy p
    493 */
    494 static smpoly smPoly2Smpoly(poly q)
    495 {
    496   poly pp;
    497   smpoly res, a;
    498   long x;
    499 
    500   if (q == NULL)
    501     return NULL;
    502   a = res = (smpoly)omAllocBin(smprec_bin);
    503   a->pos = x = pGetComp(q);
    504   a->m = q;
    505   a->e = 0;
    506   loop
    507   {
    508     pSetComp(q,0);
    509     pp = q;
    510     pIter(q);
    511     if (q == NULL)
    512     {
    513       a->n = NULL;
    514       return res;
    515     }
    516     if (pGetComp(q) != x)
    517     {
    518       a = a->n = (smpoly)omAllocBin(smprec_bin);
    519       pNext(pp) = NULL;
    520       a->pos = x = pGetComp(q);
    521       a->m = q;
    522       a->e = 0;
    523     }
    524   }
    525437}
    526438
    … …  
    587499
    588500/*
    589 * from smpoly to poly
    590 * destroy a
    591 */
    592 static poly smSmpoly2Poly(smpoly a)
    593 {
    594   smpoly b;
    595   poly res, pp, q;
    596   long x;
    597 
    598   if (a == NULL)
    599     return NULL;
    600   x = a->pos;
    601   q = res = a->m;
    602   loop
    603   {
    604     pSetComp(q,x);
    605     pp = q;
    606     pIter(q);
    607     if (q == NULL)
    608       break;
    609   }
    610   loop
    611   {
    612     b = a;
    613     a = a->n;
    614     omFreeBin((void *)b,  smprec_bin);
    615     if (a == NULL)
    616       return res;
    617     x = a->pos;
    618     q = pNext(pp) = a->m;
    619     loop
    620     {
    621       pSetComp(q,x);
    622       pp = q;
    623       pIter(q);
    624       if (q == NULL)
    625         break;
    626     }
    627   }
    628 }
    629 
    630 /*
    631501* transform the result to a module
    632502*/
    … …  
    789659
    790660/*
    791 * weigth of a polynomial, for pivot strategy
    792 */
    793 static float smPolyWeight(smpoly a)
    794 {
    795   poly p = a->m;
    796   int i;
    797   float res = (float)nSize(pGetCoeff(p));
    798 
    799   if (pNext(p) == NULL)
    800   {
    801     for(i=pVariables; i>0; i--)
    802     {
    803       if (pGetExp(p,i) != 0) return res+1.0;
    804     }
    805     return res;
    806   }
    807   else
    808   {
    809     i = 0;
    810     res = 0.0;
    811     do
    812     {
    813       i++;
    814       res += (float)nSize(pGetCoeff(p));
    815       pIter(p);
    816     }
    817     while (p);
    818     return res+(float)i;
    819   }
    820 }
    821 
    822 /*
    823661* prepare smPivot, compute weights for rows and columns
    824662* and the weight for all points
    … …  
    1011849
    1012850/* ----------------- elimination ------------------ */
    1013 static void smElemDelete(smpoly *r)
    1014 {
    1015   smpoly a = *r, b = a->n;
    1016 
    1017   pDelete(&a->m);
    1018   omFreeBin((void *)a,  smprec_bin);
    1019   *r = b;
    1020 }
    1021 
    1022 static smpoly smElemCopy(smpoly a)
    1023 {
    1024   smpoly r = (smpoly)omAllocBin(smprec_bin);
    1025   memcpy(r, a, sizeof(smprec));
    1026 /*  r->m = pCopy(r->m); */
    1027   return r;
    1028 }
    1029851
    1030852/* first step of elimination */
    … …  
    18431665
    18441666
    1845 static BOOLEAN smIsNegQuot(poly a, const poly b, const poly c)
    1846 {
    1847   if (pLmDivisibleByNoComp(c, b))
    1848   {
    1849     pExpVectorDiff(a, b, c);
    1850     // Hmm: here used to be a pSetm(a): but it is unnecessary,
    1851     // if b and c are correct
    1852     return FALSE;
    1853   }
    1854   else
    1855   {
    1856     int i;
    1857     for (i=pVariables; i>0; i--)
    1858     {
    1859       if(pGetExp(c,i) > pGetExp(b,i))
    1860         pSetExp(a,i,pGetExp(c,i)-pGetExp(b,i));
    1861       else
    1862         pSetExp(a,i,0);
    1863     }
    1864     // here we actually might need a pSetm, if a is to be used in
    1865     // comparisons
    1866     return TRUE;
    1867   }
    1868 }
    1869 
    1870 static BOOLEAN smSmaller(poly a, poly b)
    1871 {
    1872   loop
    1873   {
    1874     pIter(b);
    1875     if (b == NULL) return TRUE;
    1876     pIter(a);
    1877     if (a == NULL) return FALSE;
    1878   }
    1879 }
    1880 
    1881 static void smCombineChain(poly *px, poly r)
    1882 {
    1883   poly pa = *px, pb;
    1884   number x;
    1885   int i;
    1886 
    1887   loop
    1888   {
    1889     pb = pNext(pa);
    1890     if (pb == NULL)
    1891     {
    1892       pa = pNext(pa) = r;
    1893       break;
    1894     }
    1895     i = pLmCmp(pb, r);
    1896     if (i > 0)
    1897       pa = pb;
    1898     else
    1899     {
    1900       if (i == 0)
    1901       {
    1902         x = nAdd(pGetCoeff(pb), pGetCoeff(r));
    1903         pLmDelete(&r);
    1904         if (nIsZero(x))
    1905         {
    1906           pLmDelete(&pb);
    1907           pNext(pa) = pAdd(pb,r);
    1908         }
    1909         else
    1910         {
    1911           pa = pb;
    1912           pSetCoeff(pa,x);
    1913           pNext(pa) = pAdd(pNext(pa), r);
    1914         }
    1915       }
    1916       else
    1917       {
    1918         pa = pNext(pa) = r;
    1919         pNext(pa) = pAdd(pb, pNext(pa));
    1920       }
    1921       break;
    1922     }
    1923   }
    1924   *px = pa;
    1925 }
    1926 
    1927 static void smFindRef(poly *ref, poly *px, poly r)
    1928 {
    1929   number x;
    1930   int i;
    1931   poly pa = *px, pp = NULL;
    1932 
    1933   loop
    1934   {
    1935     i = pLmCmp(pa, r);
    1936     if (i > 0)
    1937     {
    1938       pp = pa;
    1939       pIter(pa);
    1940       if (pa==NULL)
    1941       {
    1942         pNext(pp) = r;
    1943         break;
    1944       }
    1945     }
    1946     else
    1947     {
    1948       if (i == 0)
    1949       {
    1950         x = nAdd(pGetCoeff(pa), pGetCoeff(r));
    1951         pLmDelete(&r);
    1952         if (nIsZero(x))
    1953         {
    1954           pLmDelete(&pa);
    1955           if (pp!=NULL)
    1956             pNext(pp) = pAdd(pa,r);
    1957           else
    1958             *px = pAdd(pa,r);
    1959         }
    1960         else
    1961         {
    1962           pp = pa;
    1963           pSetCoeff(pp,x);
    1964           pNext(pp) = pAdd(pNext(pp), r);
    1965         }
    1966       }
    1967       else
    1968       {
    1969         if (pp!=NULL)
    1970           pp = pNext(pp) = r;
    1971         else
    1972           *px = pp = r;
    1973         pNext(pp) = pAdd(pa, pNext(r));
    1974       }
    1975       break;
    1976     }
    1977   }
    1978   *ref = pp;
    1979 }
    1980 
    19811667//disable that, as it fails with coef buckets
    19821668//#define X_MAS
    … …  
    19941680  if ((c == NULL) || pLmIsConstantComp(c))
    19951681  {
    1996     return ppMult_qq(a, b);
     1682    return pp_Mult_qq(a, b);
    19971683  }
    19981684
    … …  
    26182304    }
    26192305  }
    2620   smExactPolyDiv(a, b);
     2306  return res;
    26212307}
    26222308
Note: See TracChangeset for help on using the changeset viewer.