Changeset 6942859 in git


Ignore:
Timestamp:
Oct 2, 2013, 3:24:29 PM (9 years ago)
Author:
Andreas Steenpass <steenpass@…>
Branches:
(u'jengelh-datetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', 'a800fe4b3e9d37a38c5a10cc0ae9dfa0c15a4ee6')
Children:
2a13ece6182151684381e63ed79c0ab959d69b85
Parents:
c5a262ad53df1f22857ebbf33ab1df008ab05831
git-author:
Andreas Steenpass <steenpass@mathematik.uni-kl.de>2013-10-02 15:24:29+02:00
git-committer:
Andreas Steenpass <steenpass@mathematik.uni-kl.de>2013-10-04 16:55:05+02:00
Message:
chg: update realclassify.lib
(cherry picked from commit 9eb2f474dffc7af6a2fe0ffdd8d65a80f37027c8)

Signed-off-by: Andreas Steenpass <steenpass@mathematik.uni-kl.de>
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • Singular/LIB/realclassify.lib

    rc5a262 r6942859  
    2828";
    2929
     30LIB "elim.lib";
     31LIB "primdec.lib";
    3032LIB "classify.lib";
    3133LIB "rootsur.lib";
     34LIB "rootsmr.lib";
    3235LIB "atkins.lib";
    3336LIB "solve.lib";
    … …  
    3538proc realclassify(poly f, list #)
    3639"
    37 USAGE:    realclassify(f[, printcomments]); f poly, printcomments int
     40USAGE:    realclassify(f[, format]); f poly, format string
    3841RETURN:   A list containing (in this order)
    3942          @* - the type of the singularity as a string,
    … …  
    4649             the normal form is given as a string with an \"a\" as the
    4750             parameter. Otherwise, it is given as a polynomial.
    48           @* An optional integer printcomments can be provided. If its value
    49              is non-zero, a string will be added at the end of the returned
    50              list, containing the result in more readable form and in some
    51              cases also more comments on how to interpret the result. The
    52              default is zero.
     51          @* An optional string @code{format} can be provided. Its default
     52             value is \"short\" in which case the return value is the list
     53             described above. If set to \"nice\", a string is added at the end
     54             of this list, containing the result in a more readable form.
    5355NOTE:     The classification is done over the real numbers, so in contrast to
    5456          classify.lib, the signs of coefficients of monomials where even
    … …  
    7476  if(size(#) > 0)
    7577  {
    76     if(size(#) > 1 || typeof(#[1]) != "int")
     78    if(size(#) > 1 || typeof(#[1]) != "string")
    7779    {
    7880      ERROR("Wrong optional parameters.");
    7981    }
    80     printcomments = #[1];
     82    if(#[1] != "short" && #[1] != "nice")
     83    {
     84      ERROR("Wrong optional parameters.");
     85    }
     86    if(#[1] == "nice")
     87    {
     88      printcomments = 1;
     89    }
    8190  }
    8291
    … …  
    128137
    129138  /* determine the type */
    130   string typeofsing, typeofsing_alternative;
    131   poly nf, nf_alternative;
     139  string typeofsing;
     140  poly nf;
    132141  poly monparam;   // the monomial whose coefficient is the parameter
    133142                   // in the modality 1 cases, 0 otherwise
    134143  string morecomments = newline;
    135   map phi;
     144
    136145  if(cr == 0)   // case A[1]
    137146  {
    138     if(lambda < n)
    139     {
    140       typeofsing = "A[1]+";
    141       typeofsing_alternative = "A[1]-";
    142     }
    143     else
    144     {
    145       typeofsing = "A[1]-";
    146       typeofsing_alternative = "A[1]+";
    147     }
     147    typeofsing, nf = caseA1(rf, lambda, n);
    148148  }
    149149  if(cr == 1)   // case A[k], k > 1
    150150  {
    151     int k = deg(lead(rf), 1:n)-1;
    152     if(k%2 == 0)
    153     {
    154       nf = var(1)^(k+1);
    155       nf_alternative = nf;
    156       typeofsing = "A["+string(k)+"]";
    157       typeofsing_alternative = typeofsing;
    158     }
    159     else
    160     {
    161       if(leadcoef(rf) > 0)
    162       {
    163         nf = var(1)^(k+1);
    164         typeofsing = "A["+string(k)+"]+";
    165         typeofsing_alternative = "A["+string(k)+"]-";
    166       }
    167       else
    168       {
    169         nf = -var(1)^(k+1);
    170         typeofsing = "A["+string(k)+"]-";
    171         typeofsing_alternative = "A["+string(k)+"]+";
    172       }
    173       nf_alternative = -nf;
    174     }
     151    typeofsing, nf = caseAk(rf, n);
    175152  }
    176153  if(cr == 2)
    … …  
    180157      if(mu == 4)   // case D[4]
    181158      {
    182         rf = jet(rf, 3);
    183         number s1 = number(rf/(var(1)^3));
    184         number s2 = number(rf/(var(2)^3));
    185         if(s2 == 0 && s1 != 0)
    186         {
    187           phi = br, var(2), var(1);
    188           rf = phi(rf);
    189         }
    190         if(s1 == 0 && s2 == 0)
    191         {
    192           number t1 = number(rf/(var(1)^2*var(2)));
    193           number t2 = number(rf/(var(2)^2*var(1)));
    194           if(t1+t2 == 0)
    195           {
    196             phi = br, var(1), 2*var(2);
    197             rf = phi(rf);
    198           }
    199           phi = br, var(1)+var(2), var(2);
    200           rf = phi(rf);
    201         }
    202         ring R = 0, y, dp;
    203         map phi = br, 1, y;
    204         poly rf = phi(rf);
    205         int k = nrroots(rf);
    206         setring(br);
    207         if(k == 3)
    208         {
    209           nf = var(1)^2*var(2)-var(2)^3;
    210           typeofsing = "D[4]-";
    211         }
    212         else
    213         {
    214           nf = var(1)^2*var(2)+var(2)^3;
    215           typeofsing = "D[4]+";
    216         }
    217         typeofsing_alternative = typeofsing;
    218         nf_alternative = nf;
     159        typeofsing, nf = caseD4(rf);
    219160      }
    220161      else   // case D[k], k > 4
    221162      {
    222         rf = jet(rf, d);
    223         list factorization = factorize(jet(rf, 3));
    224         list factors = factorization[1][2];
    225         if(factorization[2][2] == 2)
    226         {
    227           factors = insert(factors, factorization[1][3], 1);
    228         }
    229         else
    230         {
    231           factors = insert(factors, factorization[1][3]);
    232         }
    233         factors[2] = factorization[1][1]*factors[2];
    234         matrix T[2][2] = factors[1]/var(1), factors[1]/var(2),
    235                          factors[2]/var(1), factors[2]/var(2);
    236         phi = br, luinverse(T)[2]*matrix(ideal(var(1), var(2)), 2, 1);
    237         rf = phi(rf);
    238         rf = jet(rf, d);
    239         poly g;
    240         for(i = 4; i < mu; i++)
    241         {
    242           g = jet(rf, i) - var(1)^2*var(2);
    243           if(g != 0)
    244           {
    245             phi = br, var(1)-(g/(var(1)*var(2)))/2,
    246                       var(2)-(g/var(1)^i)*var(1)^(i-2);
    247             rf = phi(rf);
    248             rf = jet(rf, d);
    249           }
    250         }
    251         number a = number(rf/var(2)^(mu-1));
    252         if(a > 0)
    253         {
    254           typeofsing = "D["+string(mu)+"]+";
    255           nf = var(1)^2*var(2)+var(2)^(mu-1);
    256           if(mu%2 == 0)
    257           {
    258             nf_alternative = nf;
    259             typeofsing_alternative = typeofsing;
    260           }
    261           else
    262           {
    263             nf_alternative = var(1)^2*var(2)-var(2)^(mu-1);
    264             typeofsing_alternative = "D["+string(mu)+"]-";
    265           }
    266         }
    267         else
    268         {
    269           typeofsing = "D["+string(mu)+"]-";
    270           nf = var(1)^2*var(2)-var(2)^(mu-1);
    271           if(mu%2 == 0)
    272           {
    273             nf_alternative = nf;
    274             typeofsing_alternative = typeofsing;
    275           }
    276           else
    277           {
    278             nf_alternative = var(1)^2*var(2)+var(2)^(mu-1);
    279             typeofsing_alternative = "D["+string(mu)+"]+";
    280           }
    281         }
     163        typeofsing, nf = caseDk(rf, mu);
    282164      }
    283165    }
    284     if(complextype == "E[6]")  // case E[6] ;
    285     {
    286       poly g = jet(rf,3);
    287       number s = number(g/(var(1)^3));
    288       if(s == 0)
    289       {
    290         phi = br, var(2), var(1);
    291         rf = phi(rf);
    292         g = jet(rf,3);
    293       }
    294       list Factors = factorize(g);
    295       poly g1 = Factors[1][2];
    296       phi = br, (var(1)-(g1/var(2))*var(2))/(g1/var(1)), var(2);
    297       rf = phi(rf);
    298       rf = jet(rf,4);
    299       number w = number(rf/(var(2)^4));
    300       if(w > 0)
    301       {
    302         typeofsing = "E[6]+";
    303         nf = var(1)^3+var(2)^4;
    304         typeofsing_alternative = "E[6]-";
    305         nf_alternative = var(1)^3-var(2)^4;
    306       }
    307       else
    308       {
    309         typeofsing = "E[6]-";
    310         nf = var(1)^3-var(2)^4;
    311         typeofsing_alternative = "E[6]+";
    312         nf_alternative = var(1)^3+var(2)^4;
    313       }
     166    if(complextype == "E[6]")   // case E[6]
     167    {
     168      typeofsing, nf = caseE6(rf);
    314169    }
    315170    if(complextype == "E[7]")   // case E[7]
    316171    {
    317       typeofsing = "E[7]";
    318       nf = var(1)^3+var(1)*var(2)^3;
    319       typeofsing_alternative = typeofsing;
    320       nf_alternative = nf;
     172      typeofsing, nf = caseE7();
    321173    }
    322174    if(complextype == "E[8]")   // case E[8]
    323175    {
    324       typeofsing = "E[8]";
    325       nf = var(1)^3+var(2)^5;
    326       typeofsing_alternative = typeofsing;
    327       nf_alternative = nf;
    328     }
    329     if(complextype == "J[2,0]")  // case J[10]
    330     {
    331       int signforJ10;
    332       poly g = jet(rf,3);
    333       number s = number(g/(var(1)^3));
    334       if (s == 0)
    335       {
    336         phi = br, var(2), var(1);
    337         rf = phi(rf);
    338         g = jet(rf,3);
    339       }
    340       list Factors = factorize(g);
    341       poly g1 = Factors[1][2];
    342       phi = br, (var(1)-(g1/var(2))*var(2))/(g1/var(1)), var(2);
    343       rf = phi(rf);
    344       poly rf3 = jet(rf,3);
    345       number w0 = number(rf3/(var(1)^3));
    346       if(w0 < 0)
    347       {
    348         phi = br, -var(1), var(2);
    349         rf = phi(rf);
    350       }
    351       rf3 = jet(rf,3);
    352       poly rf4 = jet(rf,4);
    353       poly rf5 = jet(rf,5);
    354       poly rf6 = jet(rf,6);
    355       poly b1 = 3*(rf3/(var(1)^3))*var(1)^2+2*(rf4/(var(1)^2*var(2)^2))
    356                 *var(1)+(rf5/(var(1)*var(2)^4));
    357       ring R = 0,x,dp;
    358       map phi = br, x, 1;
    359       poly b11 = phi(b1);
    360       int r = nrroots(b11);
    361       if( r == 0 || r == 1)
    362       {
    363         setring(br);
    364         signforJ10 = 1;
    365       }
    366       else
    367       {
    368         setring(br);
    369         poly c1 = (rf3/(var(1)^3))*var(1)^3+(rf4/(var(1)^2*var(2)^2))*
    370                   var(1)^2+(rf5/(var(1)*var(2)^4))*
    371                   var(1)+rf6/(var(2)^6);
    372         list Factors = factorize(c1);
    373         if( size(Factors[1])>2)
    374         {
    375           if( deg(Factors[1][2]) == 1)
    376           {
    377             poly g1 = Factors[1][2];
    378           }
    379           else
    380           {
    381             poly g1 = Factors[1][3];
    382           }
    383           phi = br, ((g1/var(1))*var(1)-g1)/(g1/var(1)),
    384                   var(2);
    385           number b = number(phi(b1));
    386           if(b > 0)
    387           {
    388             signforJ10 = 1;
    389           }
    390           else
    391           {
    392             signforJ10 = -1;
    393           }
    394         }
    395         else
    396         {
    397           ring R = (complex,40,40),x,lp;
    398           phi = br, x, 1;
    399           poly c1 = phi(c1);
    400           poly b1 = phi(b1);
    401           list L = laguerre_solve(c1,30);
    402           list LL;
    403           for(i = 1; i <= size(L); i++)
    404           {
    405             if(impart(L[i]) == 0)
    406             {
    407               LL = insert(LL,L[i]);
    408             }
    409           }
    410           number r1 = LL[1];
    411           for(j = 1; j <= size(LL); j++)
    412           {
    413             r1 = round(r1*10000000000)/10000000000;
    414             number c1min = number(subst(c1,x,r1-0.0000000001));
    415             number c1plus = number(subst(c1,x,r1+0.0000000001));
    416             number b1min = number(subst(b1,x,r1-0.00000000001));
    417             number b1plus = number(subst(b1,x,r1+0.00000000001));
    418             if(c1min*c1plus<0)
    419             {
    420               int c = -1;
    421             }
    422             if(c1min*c1plus>0)
    423             {
    424               int c = 1;
    425             }
    426             if(b1min>0 && b1plus>0)
    427             {
    428               int b = 1;
    429             }
    430             if(b1min<0 && b1plus<0)
    431             {
    432               int b = -1;
    433             }
    434             if(b1min*b1plus<=0)
    435             {
    436               int b = 0;
    437             }
    438             if( c < 0 && b != 0)
    439             {
    440               r1 = LL[j];
    441               break;
    442             }
    443           }
    444           setring(br);
    445           if (c == -1 && b == 1)
    446           {
    447             signforJ10 = 1;
    448           }
    449           if (c == -1 && b == -1)
    450           {
    451             signforJ10 = -1;
    452           }
    453           if (c == 1 || b == 0)
    454           {
    455             ERROR("Ask Arnold the normal form.)");
    456           }
    457         }
    458       }
    459       if(signforJ10 == 1)
    460       {
    461         typeofsing = "J[10]+";
    462       }
    463       else
    464       {
    465         typeofsing = "J[10]-";
    466       }
    467       nf = var(1)^3+var(1)^2*var(2)^2+signforJ10*var(1)*var(2)^4;
    468       monparam = var(1)^2*var(2)^2;
    469       typeofsing_alternative = typeofsing;
    470       nf_alternative = nf;
    471     }
    472     if(complextype[1,3] == "X[1") //case X[1,k]
    473     {
    474       if(mu > 9)
    475       {
    476         rf = jet(rf,4);
    477 
    478         number s1 = number(rf/(var(1)^4));
    479         number s2 = number(rf/(var(2)^4));
    480         if(s2 != 0 && s1 == 0)
    481         {
    482           phi = br, var(2), var(1);
    483           rf = phi(rf);
    484         }
    485         if(s2 == 0 && s1 == 0)
    486         {
    487           number t1 = number(rf/(var(1)^3*var(2)));
    488           number t2 = number(rf/(var(1)^2*var(2)^2));
    489           number t3 = number(rf/(var(1)*var(2)^3));
    490           if(t1+t2+t3 == 0)
    491           {
    492             if(2*t1+4*t2+8*t3 != 0)
    493             {
    494               phi = br, var(1), 2*var(2);
    495               rf = phi(rf);
    496             }
    497             else
    498             {
    499               phi = br, var(1), 3*var(2);
    500               rf = phi(rf);
    501             }
    502           }
    503           phi = br, var(1), var(1)+var(2);
    504           rf = phi(rf);
    505         }
    506         ring R = 0, x, dp;
    507         map phi = br, var(1), 1;
    508         int k = nrroots(phi(rf));
    509         setring(br);
    510         if(k == 1)
    511         {
    512           number w = number(rf/(var(1)^4));
    513           if(w > 0)
    514           {
    515             typeofsing = "X["+string(mu)+"]++";
    516             nf = var(1)^4+var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    517             typeofsing_alternative = "X["+string(mu)+"]--";
    518             nf_alternative = -var(1)^4-var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    519           }
    520           else
    521           {
    522             typeofsing = "X["+string(mu)+"]--";
    523             nf = -var(1)^4-var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    524             typeofsing_alternative = "X["+string(mu)+"]++";
    525             nf_alternative = var(1)^4+var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    526           }
    527         }
    528         if(k == 3)
    529         {
    530           list Factors = factorize(rf);
    531           for(i = 2; i <= size(Factors[1]); i++)
    532           {
    533             if(Factors[2][i] == 2)
    534             {
    535               poly g1 = Factors[1][i];
    536               break;
    537             }
    538           }
    539           map phi = br, (var(1)-(g1/(var(2))*var(2)))/(g1/var(1)), var(2);
    540           rf = phi(rf);
    541           number w = number(rf/(var(1)^2*var(2)^2));
    542           if(w > 0)
    543           {
    544             typeofsing = "X["+string(mu)+"]-+";
    545             nf = -var(1)^4+var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    546             typeofsing_alternative = "X["+string(mu)+"]+-";
    547             nf_alternative = var(1)^4-var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    548           }
    549           else
    550           {
    551             typeofsing = "X["+string(mu)+"]+-";
    552             nf = var(1)^4-var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    553             typeofsing_alternative = "X["+string(mu)+"]-+";
    554             nf_alternative = -var(1)^4+var(1)^2*var(2)^2+var(2)^(4+(mu-9));
    555           }
    556         }
    557         monparam = var(2)^(4+(mu-9));
    558       }
    559     }
    560     if(complextype == "E[12]")   // case E[12]
    561     {
    562       typeofsing = "E[12]";
    563       nf = var(1)^3+var(2)^7+var(1)*var(2)^5;
    564       monparam = var(1)*var(2)^5;
    565       typeofsing_alternative = typeofsing;
    566       nf_alternative = nf;
    567     }
    568     if(complextype == "E[13]")   // case E[13]
    569     {
    570       typeofsing = "E[13]";
    571       nf = var(1)^3+var(1)*var(2)^5+var(2)^8;
    572       monparam = var(2)^8;
    573       typeofsing_alternative = typeofsing;
    574       nf_alternative = nf;
    575     }
    576     if(complextype == "E[14]") //case E[14]
    577     {
    578       poly g = jet(rf,3);
    579       number s = number(g/(var(1)^3));
    580       if(s == 0)
    581       {
    582         phi = br, var(2), var(1);
    583         rf = phi(rf);
    584         g = jet(rf,3);
    585         s = number(g/(var(1)^3));
    586       }
    587       rf = rf/s;
    588       list Factors = factorize(g);
    589       poly g1 = Factors[1][2];
    590       phi = br, (var(1)-(g1/var(2))*var(2))/(g1/var(1)), var(2);
    591       rf = phi(rf);
    592       g = jet(rf,3);
    593       number w0 = number(g/(var(1)^3));
    594       phi = br, var(1)-((jet(rf,4)-(w0*var(1)^3))/(3*var(1)^2)), var(2);
    595       rf = phi(rf);
    596       phi = br, var(1)-((jet(rf,5)-(w0*var(1)^3))/(3*var(1)^2)), var(2);
    597       rf = phi(rf);
    598       rf = s*rf;
    599       rf = jet(rf,8);
    600       number w = number(rf/(var(2)^8));
    601       if(w > 0)
    602       {
    603         typeofsing = "E[14]+";
    604         nf = var(1)^3+var(2)^8+var(1)*var(2)^6;
    605         typeofsing_alternative = "E[14]-";
    606         nf_alternative = var(1)^3-var(2)^8+var(1)*var(2)^6;
    607       }
    608       if(w < 0)
    609       {
    610         typeofsing = "E[14]-";
    611         nf = var(1)^3-var(2)^8+var(1)*var(2)^6;
    612         typeofsing_alternative = "E[14]+";
    613         nf_alternative = var(1)^3+var(2)^8+var(1)*var(2)^6;
    614       }
    615       monparam = var(1)*var(2)^6;
    616     }
    617     if(complextype == "Z[11]")   // case Z[11]
    618     {
    619       typeofsing = "Z[11]";
    620       nf = var(1)^3*var(2)+var(2)^5+var(1)*var(2)^4;
    621       monparam = var(1)*var(2)^4;
    622       typeofsing_alternative = typeofsing;
    623       nf_alternative = nf;
    624     }
    625     if(complextype == "Z[12]")   // case Z[12]
    626     {
    627       typeofsing = "Z[12]";
    628       nf = var(1)^3*var(2)+var(1)*var(2)^4+var(1)^2*var(2)^3;
    629       monparam = var(1)^2*var(2)^3;
    630       typeofsing_alternative = typeofsing;
    631       nf_alternative = nf;
    632     }
    633     if(complextype == "Z[13]")
    634     {
    635       poly g = jet(rf,4);
    636       number s = number(g/var(1)^3*var(2));
    637       if(s == 0)
    638       {
    639         phi = br, var(2), var(1);
    640         rf = phi(rf);
    641         g = jet(rf,4);
    642       }
    643       list Factors = factorize(g);
    644       if(Factors[2][2] == 3)
    645       {
    646         poly g1 = Factors[1][2];
    647       }
    648       else
    649       {
    650         poly g1 = Factors[1][3];
    651       }
    652       phi = br, var(1)-(g1/var(2))*var(2), var(2);
    653       rf = phi(rf);
    654       rf = jet(rf,6);
    655       number w = number(rf/var(2)^6);
    656       if(w > 0)
    657       {
    658         typeofsing = "Z[13]+";
    659         nf = var(1)^3*var(2)+var(2)^6+var(1)*var(2)^5;
    660         typeofsing_alternative = "Z[13]-";
    661         nf_alternative = var(1)^3*var(2)-var(2)^6+var(1)*var(2)^5;
    662       }
    663       else
    664       {
    665         typeofsing = "Z[13]-";
    666         nf = var(1)^3*var(2)-var(2)^6+var(1)*var(2)^5;
    667         typeofsing_alternative = "Z[13]+";
    668         nf_alternative = var(1)^3*var(2)+var(2)^6+var(1)*var(2)^5;
    669       }
    670       monparam = var(1)*var(2)^5;
    671     }
    672     if(complextype == "W[12]") //case W[12]
    673     {
    674       poly g = jet(rf, 4);
    675       number s = number(g/(var(1)^4));
    676       if(s == 0)
    677       {
    678         s = number(g/(var(2)^4));
    679         phi = br, var(2), var(1);   // maybe we'll need this transformation
    680         rf = phi(rf);               // later
    681       }
    682       if(s > 0)
    683       {
    684         typeofsing = "W[12]+";
    685         nf = var(1)^4+var(2)^5+var(1)^2*var(2)^3;
    686         typeofsing_alternative = "W[12]-";
    687         nf_alternative = -var(1)^4+var(2)^5+var(1)^2*var(2)^3;
    688       }
    689       else
    690       {
    691         typeofsing = "W[12]-";
    692         nf = -var(1)^4+var(2)^5+var(1)^2*var(2)^3;
    693         typeofsing_alternative = "W[12]+";
    694         nf_alternative = var(1)^4+var(2)^5+var(1)^2*var(2)^3;
    695       }
    696       monparam = var(1)^2*var(2)^3;
    697     }
    698     if(complextype == "W[13]")   //case W[13]
    699     {
    700       poly g = jet(rf, 4);
    701       number s = number(g/(var(1)^4));
    702       if(s == 0)
    703       {
    704         s = number(g/(var(2)^4));
    705         phi = br, var(2), var(1);   // maybe we'll need this transformation
    706         rf = phi(rf);               // later
    707       }
    708       if(s > 0)
    709       {
    710         typeofsing = "W[13]+";
    711         nf = var(1)^4+var(1)*var(2)^4+var(2)^6;
    712         typeofsing_alternative = "W[13]-";
    713         nf_alternative = -var(1)^4+var(1)*var(2)^4+var(2)^6;
    714       }
    715       else
    716       {
    717         typeofsing = "W[13]-";
    718         nf = -var(1)^4+var(1)*var(2)^4+var(2)^6;
    719         typeofsing_alternative = "W[13]+";
    720         nf_alternative = var(1)^4+var(1)*var(2)^4+var(2)^6;
    721       }
    722       monparam = var(2)^6;
     176      typeofsing, nf = caseE8();
    723177    }
    724178    if(typeofsing == "")
    … …  
    734188  /* add the non-corank variables to the normal forms */
    735189  nf = addnondegeneratevariables(nf, lambda, cr);
    736   nf_alternative = addnondegeneratevariables(nf_alternative, n-cr-lambda, cr);
    737190
    738191  /* write normal form as a string in the cases with modality greater than 0 */
    … …  
    740193  {
    741194    poly nf_tmp = nf;
    742     poly nf_alternative_tmp = nf_alternative;
     195    kill nf;
    743196    def nf = modality1NF(nf_tmp, monparam);
    744     def nf_alternative = modality1NF(nf_alternative_tmp, monparam);
    745197  }
    746198
    … …  
    755207                       +"Inertia index:       "   +string(lambda)+newline
    756208                       +"Determinacy:         <= "+string(d)     +newline;
    757     if(typeofsing != typeofsing_alternative || nf != nf_alternative
    758        || lambda != n-cr-lambda)
    759     {
    760       comments = comments+newline
    761         +"Note: By multiplying the input polynomial with -1,"+newline
    762         +"      it can also be regarded as of the following case:"+newline
    763         +"Type of singularity: "+typeofsing_alternative+newline
    764         +"Normal form:         "+string(nf_alternative)+newline
    765         +"Inertia index:       "+string(n-cr-lambda)   +newline;
    766     }
    767209    if(morecomments != newline)
    768210    {
    … …  
    787229  ring r = 0, (x,y,z), ds;
    788230  poly f = (x2+3y-2z)^2+xyz-(x-y3+x2z3)^3;
    789   realclassify(f, 1);
     231  realclassify(f, "nice");
     232}
     233
     234///////////////////////////////////////////////////////////////////////////////
     235static proc caseA1(poly rf, int lambda, int n)
     236{
     237  string typeofsing = "A[1]";
     238  poly nf = 0;
     239  return(typeofsing, nf);
     240}
     241
     242///////////////////////////////////////////////////////////////////////////////
     243static proc caseAk(poly rf, int n)
     244{
     245  /* preliminaries */
     246  string typeofsing;
     247  poly nf;
     248
     249  int k = deg(lead(rf), 1:n)-1;
     250  if(k%2 == 0)
     251  {
     252    nf = var(1)^(k+1);
     253    typeofsing = "A["+string(k)+"]";
     254  }
     255  else
     256  {
     257    if(leadcoef(rf) > 0)
     258    {
     259      nf = var(1)^(k+1);
     260      typeofsing = "A["+string(k)+"]+";
     261    }
     262    else
     263    {
     264      nf = -var(1)^(k+1);
     265      typeofsing = "A["+string(k)+"]-";
     266    }
     267  }
     268  return(typeofsing, nf);
     269}
     270
     271///////////////////////////////////////////////////////////////////////////////
     272static proc caseD4(poly rf)
     273{
     274  /* preliminaries */
     275  string typeofsing;
     276  poly nf;
     277  def br = basering;
     278  map phi;
     279
     280  rf = jet(rf, 3);
     281  number s1 = number(rf/(var(1)^3));
     282  number s2 = number(rf/(var(2)^3));
     283  if(s2 == 0 && s1 != 0)
     284  {
     285    phi = br, var(2), var(1);
     286    rf = phi(rf);
     287  }
     288  if(s1 == 0 && s2 == 0)
     289  {
     290    number t1 = number(rf/(var(1)^2*var(2)));
     291    number t2 = number(rf/(var(2)^2*var(1)));
     292    if(t1+t2 == 0)
     293    {
     294      phi = br, var(1)+2*var(2), var(2);
     295      rf = phi(rf);
     296    }
     297    else
     298    {
     299      phi = br, var(1)+var(2), var(2);
     300      rf = phi(rf);
     301    }
     302  }
     303  ring R = 0, y, dp;
     304  map phi = br, 1, y;
     305  poly rf = phi(rf);
     306  int k = nrroots(rf);
     307  setring(br);
     308  if(k == 3)
     309  {
     310    nf = var(1)^2*var(2)-var(2)^3;
     311    typeofsing = "D[4]-";
     312  }
     313  else
     314  {
     315    nf = var(1)^2*var(2)+var(2)^3;
     316    typeofsing = "D[4]+";
     317  }
     318  return(typeofsing, nf);
     319}
     320
     321///////////////////////////////////////////////////////////////////////////////
     322static proc caseDk(poly rf, int mu)
     323{
     324  /* preliminaries */
     325  string typeofsing;
     326  poly nf;
     327  def br = basering;
     328  map phi;
     329
     330  rf = jet(rf, mu-1);
     331  list factorization = factorize(jet(rf, 3));
     332  list factors = factorization[1][2];
     333  if(factorization[2][2] == 2)
     334  {
     335    factors = insert(factors, factorization[1][3], 1);
     336  }
     337  else
     338  {
     339    factors = insert(factors, factorization[1][3]);
     340  }
     341  factors[2] = factorization[1][1]*factors[2];
     342  matrix T[2][2] = factors[1]/var(1), factors[1]/var(2),
     343         factors[2]/var(1), factors[2]/var(2);
     344  phi = br, luinverse(T)[2]*matrix(ideal(var(1), var(2)), 2, 1);
     345  rf = phi(rf);
     346  rf = jet(rf, mu-1);
     347  poly g;
     348  int i;
     349  for(i = 4; i < mu; i++)
     350  {
     351    g = jet(rf, i) - var(1)^2*var(2);
     352    if(g != 0)
     353    {
     354      phi = br, var(1)-(g/(var(1)*var(2)))/2,
     355          var(2)-(g/var(1)^i)*var(1)^(i-2);
     356      rf = phi(rf);
     357      rf = jet(rf, mu-1);
     358    }
     359  }
     360  number a = number(rf/var(2)^(mu-1));
     361  if(a > 0)
     362  {
     363    typeofsing = "D["+string(mu)+"]+";
     364    nf = var(1)^2*var(2)+var(2)^(mu-1);
     365  }
     366  else
     367  {
     368    typeofsing = "D["+string(mu)+"]-";
     369    nf = var(1)^2*var(2)-var(2)^(mu-1);
     370  }
     371  return(typeofsing, nf);
     372}
     373
     374///////////////////////////////////////////////////////////////////////////////
     375static proc caseE6(poly rf)
     376{
     377  /* preliminaries */
     378  string typeofsing;
     379  poly nf;
     380  def br = basering;
     381  map phi;
     382
     383  poly g = jet(rf,3);
     384  number s = number(g/(var(1)^3));
     385  if(s == 0)
     386  {
     387    phi = br, var(2), var(1);
     388    rf = phi(rf);
     389    g = jet(rf,3);
     390  }
     391  list Factors = factorize(g);
     392  poly g1 = Factors[1][2];
     393  phi = br, (var(1)-(g1/var(2))*var(2))/(g1/var(1)), var(2);
     394  rf = phi(rf);
     395  rf = jet(rf,4);
     396  number w = number(rf/(var(2)^4));
     397  if(w > 0)
     398  {
     399    typeofsing = "E[6]+";
     400    nf = var(1)^3+var(2)^4;
     401  }
     402  else
     403  {
     404    typeofsing = "E[6]-";
     405    nf = var(1)^3-var(2)^4;
     406  }
     407  return(typeofsing, nf);
     408}
     409
     410///////////////////////////////////////////////////////////////////////////////
     411static proc caseE7()
     412{
     413  string typeofsing = "E[7]";
     414  poly nf = var(1)^3+var(1)*var(2)^3;
     415  return(typeofsing, nf);
     416}
     417
     418///////////////////////////////////////////////////////////////////////////////
     419static proc caseE8()
     420{
     421  string typeofsing = "E[8]";
     422  poly nf = var(1)^3+var(2)^5;
     423  return(typeofsing, nf);
    790424}
    791425
    … …  
    1153787  determinacy(f);
    1154788}
     789

  • Tst/Short/realclassify_s.res.gz.uu

    rc5a262 r6942859  
    1 begin 640 realclassify_s.res.gz
    2 M'XL("`<"!T\``W)E86QC;&%S<VEF>5]S+G)E<P"54\%NHS`0O?LK1M4>B`)L
    3 M;$B:IAL.W5PJ;7MI;RB-G,2DUAJ(;$?!^?H=(`&ZVLMR@,%OYGEFWLS;^^KY
    4 M%0!H`K^>G^#.&ALJN;U[)&]7A"6`AQM92.N-'DG]A20!+;C:*6Z,S-S&A(4X
    5 MA\9RVX5%"71V'+;DPYB_;ID.W&<A:%D<0,,2)CYXE>_\R\B'O>G][Q,XELI!
    6 MACY>Q<:1"]AE],'&E;L$7A6X:%RQ2S3ZB/J8^=>LO<P'BA6E=+T@"*_2>!V0
    7 ME+5_%Z1B#GE(&K4GC*!#8\4DG;;6A*2SUD*_^]8B[^XHH,S`8!$GQ;6T;G%E
    8 M?RUUSA5D^%G`[>FO^EEJ7OSND>9AY$6JHM10G/*MT#<T)L^%T%9RD,5>5+?C
    9 M"5D)*W0N"[YS/=./)6:(UUNQ@"<'^4E9>52N;K/]%$AQ/-FFHT692\SP+.TG
    10 M!-0G;3AJON,%<&5*V`ILXX'KO=@#-W6A-4-6*E6>:[X=-V+QWTT(^B[\LRY*
    11 M.AT?^EFAD[#1-"^U$>:H<$:S@:2LDW/2R1AU,@95-&;?H\JQ3L\-1BZKUF+K
    12 M):)N?`DH.C%\TWGE6BQ:+UV7`QWDPZZS:_K9_3JY-!I,+B;@XAZ*$\@;J5NE
    13 MFU)F'3I8$8HKTI`<6I+*]6-.<37V_0AX!R29=N!\0/(0-IM=K^W)>+@*W\@?
    14 (S]NXW`\$````
     1begin 600 realclassify_s.res.gz
     2M'XL("&H<3%(``W)E86QC;&%S<VEF>5]S+G)E<P!=4\%NHS`0O?LK1M$>0`$*
     3M-DG39.&PFTND;2_M#:4132"R"B:RB=;FZ]?%B4W7!SSXS3S//!ZO;]O="P`D
     4M.?S9_8)9+_JHH1^S#7J](3@'?7B@C/:>OT%?.^0Y\*ILCDTI!*W5042L^AN)
     5MONQM&<G!QFEDR*<U_]VRF*0O(^"4G8%#!G$`G@Q4,/@!G(3+?\SATC4*:IWC
     6M23PG*L2#_X[G4@VA)T-%YA(/Q'\GKF;UO6NO#F#&Z+&:Z;&*9+]&.F=;I/L0
     7M%=B\#9H/*TV&"F).,-()8Y2B8F&B&!5+$^F\1Q.A-W6IH*M!Z$FN3<EIK]8W
     8M]I>.MV4#M=[6<%_NJM\=+]FG0\:%T3-M6,>!7=N/BM_1%.U8Q7M:`F6G2MZ/
     9M8[2M^HJWE)5'Y9A^9KI#*\>3DSR)HU&:MN.B$I=&?^IZ(@JV@L16"&*%""69
     10MXP<B%;:*''1E)DV$]YE&U7P($YV$]3-9264PLL^4[2&9](-O%A#.`M\-D)")
     11M`70#*G50FD,[BF6T&D=96G3BM$0[;20Y&Q*IG%L2[;"3$]$[:Y*%!5<3DJ=H
     124_$&^W'\57N)O?J!_US0V458#````
    1513`
    1614end

  • Tst/Short/realclassify_s.stat

    rc5a262 r6942859  
    1 1 >> tst_memory_0 :: 1325859335:3130- 14004 :3-1-3:ix86-Linux:mamawutz:320108
    2 1 >> tst_memory_1 :: 1325859335:3130- 14004 :3-1-3:ix86-Linux:mamawutz:2355200
    3 1 >> tst_memory_2 :: 1325859335:3130- 14004 :3-1-3:ix86-Linux:mamawutz:2433352
    4 1 >> tst_timer_1 :: 1325859335:3130- 14004 :3-1-3:ix86-Linux:mamawutz:16
     11 >> tst_memory_0 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:323576
     21 >> tst_memory_1 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:2322432
     31 >> tst_memory_2 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:2401344
     41 >> tst_timer_1 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:16

  • Tst/Short/realclassify_s.tst

    rc5a262 r6942859  
    66ring r = 0, (x,y,z), ds;
    77poly f = (x2+3y-2z)^2+xyz-(x-y3+x2z3)^3;
    8 realclassify(f, 1);
     8realclassify(f, "nice");
    99
    1010realmorsesplit(f);
Note: See TracChangeset for help on using the changeset viewer.