Changeset 4dcfc0f in git


Ignore:
Timestamp:
Feb 5, 2008, 11:44:30 PM (16 years ago)
Author:
Viktor Levandovskyy <levandov@…>
Branches:
(u'spielwiese', '17f1d200f27c5bd38f5dfc6e8a0879242279d1d8')
Children:
a539ad11b73ebc8b93d14f33279ca7d7f3a5412c
Parents:
e671591d1d038c1f88b1a7d2d0e2ff5ad5775d32
Message:
*levandov: updates and fixes for libraries


git-svn-id: file:///usr/local/Singular/svn/trunk@10561 2c84dea3-7e68-4137-9b89-c4e89433aadc
Location:
Singular/LIB
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/dmod.lib

    re67159 r4dcfc0f  
    11//////////////////////////////////////////////////////////////////////////////
    2 version="$Id: dmod.lib,v 1.24 2008-01-17 21:05:31 levandov Exp $";
     2version="$Id: dmod.lib,v 1.25 2008-02-05 22:44:28 levandov Exp $";
    33category="Noncommutative";
    44info="
     
    27832783  dbprint(ppl,"// -1-3- u,v are eliminated");
    27842784  dbprint(ppl-1, K);  // K is without u,v
    2785 
    2786 
    27872785  setring save;
    27882786  // ------------ new ring @R2 ------------------
     
    31463144}
    31473145
     3146// use a finer ordering
     3147
    31483148proc SannfsLOT(poly F, list #)
     3149"USAGE:  SannfsLOT(f [,eng]);  f a poly, eng an optional int
     3150RETURN:  ring
     3151PURPOSE: compute the D-module structure of basering[1/f]*f^s, according to the Levandovskyy's modification of the algorithm by Oaku and Takayama in the ring D[s], where D is the Weyl algebra
     3152NOTE:    activate this ring with the @code{setring} command.
     3153@*       In the ring D[s], the ideal LD (which is NOT a Groebner basis) is the needed D-module structure.
     3154@*       If eng <>0, @code{std} is used for Groebner basis computations,
     3155@*       otherwise, and by default @code{slimgb} is used.
     3156@*       If printlevel=1, progress debug messages will be printed,
     3157@*       if printlevel>=2, all the debug messages will be printed.
     3158EXAMPLE: example SannfsLOT; shows examples
     3159"
     3160{
     3161  int eng = 0;
     3162  if ( size(#)>0 )
     3163  {
     3164    if ( typeof(#[1]) == "int" )
     3165    {
     3166      eng = int(#[1]);
     3167    }
     3168  }
     3169  // returns a list with a ring and an ideal LD in it
     3170  int ppl = printlevel-voice+2;
     3171  //  printf("plevel :%s, voice: %s",printlevel,voice);
     3172  def save = basering;
     3173  int N = nvars(basering);
     3174  //  int Nnew = 2*(N+2);
     3175  int Nnew = 2*(N+1)+1; //removed u,v; added s
     3176  int i,j;
     3177  string s;
     3178  list RL = ringlist(basering);
     3179  list L, Lord;
     3180  list tmp;
     3181  intvec iv;
     3182  L[1] = RL[1]; // char
     3183  L[4] = RL[4]; // char, minpoly
     3184  // check whether vars have admissible names
     3185  list Name  = RL[2];
     3186  list RName;
     3187//   RName[1] = "u";
     3188//   RName[2] = "v";
     3189  RName[1] = "t";
     3190  RName[2] = "Dt";
     3191  for(i=1;i<=N;i++)
     3192  {
     3193    for(j=1; j<=size(RName);j++)
     3194    {
     3195      if (Name[i] == RName[j])
     3196      {
     3197        ERROR("Variable names should not include t,Dt");
     3198      }
     3199    }
     3200  }
     3201  // now, create the names for new vars
     3202//   tmp[1]     = "u";
     3203//   tmp[2]     = "v";
     3204//   list UName = tmp;
     3205  list DName;
     3206  for(i=1;i<=N;i++)
     3207  {
     3208    DName[i] = "D"+Name[i]; // concat
     3209  }
     3210  tmp    =  0;
     3211  tmp[1] = "t";
     3212  tmp[2] = "Dt";
     3213  list SName ; SName[1] = "s";
     3214  //  list NName = tmp + Name + DName + SName;
     3215  list NName = tmp + SName + Name + DName;
     3216  L[2]   = NName;
     3217  tmp    = 0;
     3218  // Name, Dname will be used further
     3219  //  kill UName;
     3220  kill NName;
     3221  // block ord (a(1,1),dp);
     3222  tmp[1]  = "a"; // string
     3223  iv      = 1,1;
     3224  tmp[2]  = iv; //intvec
     3225  Lord[1] = tmp;
     3226  // continue with a(0,0,1)
     3227  tmp[1]  = "a"; // string
     3228  iv      = 0,0,1;
     3229  tmp[2]  = iv; //intvec
     3230  Lord[2] = tmp;
     3231  // continue with dp 1,1,1,1...
     3232  tmp[1]  = "dp"; // string
     3233  s       = "iv=";
     3234  for(i=1;i<=Nnew;i++)
     3235  {
     3236    s = s+"1,";
     3237  }
     3238  s[size(s)]= ";";
     3239  execute(s);
     3240  tmp[2]    = iv;
     3241  Lord[3]   = tmp;
     3242  tmp[1]    = "C";
     3243  iv        = 0;
     3244  tmp[2]    = iv;
     3245  Lord[4]   = tmp;
     3246  tmp       = 0;
     3247  L[3]      = Lord;
     3248  // we are done with the list
     3249  def @R@ = ring(L);
     3250  setring @R@;
     3251  matrix @D[Nnew][Nnew];
     3252  @D[1,2]=1;
     3253  for(i=1; i<=N; i++)
     3254  {
     3255    @D[3+i,N+3+i]=1;
     3256  }
     3257  // ADD [s,t]=-t, [s,Dt]=Dt
     3258  @D[1,3] = -var(1);
     3259  @D[2,3] = var(2);
     3260  //  @D[N+3,2*(N+2)]=1; old t,Dt stuff
     3261  //  L[5] = matrix(UpOneMatrix(Nnew));
     3262  //  L[6] = @D;
     3263  def @R = nc_algebra(1,@D);
     3264  setring @R;
     3265  kill @R@;
     3266  dbprint(ppl,"// -1-1- the ring @R(t,Dt,s,_x,_Dx) is ready");
     3267  dbprint(ppl-1, @R);
     3268  // create the ideal I
     3269  poly  F = imap(save,F);
     3270  //  ideal I = u*F-t,u*v-1;
     3271  ideal I = F-t;
     3272  poly p;
     3273  for(i=1; i<=N; i++)
     3274  {
     3275    //    p = u*Dt; // u*Dt
     3276    p = Dt;
     3277    p = diff(F,var(3+i))*p;
     3278    I = I, var(N+3+i) + p;
     3279  }
     3280  //  I = I, var(1)*var(2) + var(Nnew) +1; // reduce it with t-f!!!
     3281  // t*Dt + s +1 reduced with t-f gives f*Dt + s
     3282  I = I, F*var(2) + var(3);
     3283  // -------- the ideal I is ready ----------
     3284  dbprint(ppl,"// -1-2- starting the elimination of t,Dt in @R");
     3285  dbprint(ppl-1, I);
     3286  ideal J = engine(I,eng);
     3287  ideal K = nselect(J,1,2);
     3288  dbprint(ppl,"// -1-3- t,Dt are eliminated");
     3289  dbprint(ppl-1, K);  // K is without t, Dt
     3290  K = engine(K,eng);  // std does the job too
     3291  // now, we must change the ordering
     3292  // and create a ring without t, Dt
     3293  setring save;
     3294  // ----------- the ring @R3 ------------
     3295  // _x, _Dx,s;  elim.ord for _x,_Dx.
     3296  // keep: N, i,j,s, tmp, RL
     3297  Nnew = 2*N+1;
     3298  kill Lord, tmp, iv, RName;
     3299  list Lord, tmp;
     3300  intvec iv;
     3301  L[1] = RL[1];
     3302  L[4] = RL[4];  // char, minpoly
     3303  // check whether vars hava admissible names -> done earlier
     3304  // now, create the names for new var
     3305  tmp[1] = "s";
     3306  // DName is defined earlier
     3307  list NName = Name + DName + tmp;
     3308  L[2] = NName;
     3309  tmp = 0;
     3310  // block ord (dp(N),dp);
     3311  // string s is already defined
     3312  s = "iv=";
     3313  for (i=1; i<=Nnew-1; i++)
     3314  {
     3315    s = s+"1,";
     3316  }
     3317  s[size(s)]=";";
     3318  execute(s);
     3319  tmp[1] = "dp";  // string
     3320  tmp[2] = iv;   // intvec
     3321  Lord[1] = tmp;
     3322  // continue with dp 1,1,1,1...
     3323  tmp[1] = "dp";  // string
     3324  s[size(s)] = ",";
     3325  s = s+"1;";
     3326  execute(s);
     3327  kill s;
     3328  kill NName;
     3329  tmp[2]      = iv;
     3330  Lord[2]     = tmp;
     3331  tmp[1]      = "C";  iv  = 0;  tmp[2]=iv;
     3332  Lord[3]     = tmp;  tmp = 0;
     3333  L[3]        = Lord;
     3334  // we are done with the list. Now add a Plural part
     3335  def @R2@ = ring(L);
     3336  setring @R2@;
     3337  matrix @D[Nnew][Nnew];
     3338  for (i=1; i<=N; i++)
     3339  {
     3340    @D[i,N+i]=1;
     3341  }
     3342  def @R2 = nc_algebra(1,@D);
     3343  setring @R2;
     3344  kill @R2@;
     3345  dbprint(ppl,"//  -2-1- the ring @R2(_x,_Dx,s) is ready");
     3346  dbprint(ppl-1, @R2);
     3347  ideal MM = maxideal(1);
     3348  //  MM = 0,s,MM;
     3349  MM = 0,0,s,MM[1..size(MM)-1];
     3350  map R01 = @R, MM;
     3351  ideal K = R01(K);
     3352  // total cleanup
     3353  ideal LD = K;
     3354  // make leadcoeffs positive
     3355  for (i=1; i<= ncols(LD); i++)
     3356  {
     3357    if (leadcoef(LD[i]) <0 )
     3358    {
     3359      LD[i] = -LD[i];
     3360    }
     3361  }
     3362  export LD;
     3363  kill @R;
     3364  return(@R2);
     3365}
     3366example
     3367{
     3368  "EXAMPLE:"; echo = 2;
     3369  ring r = 0,(x,y,z),Dp;
     3370  poly F = x^3+y^3+z^3;
     3371  printlevel = 0;
     3372  def A  = SannfsLOT(F);
     3373  setring A;
     3374  LD;
     3375}
     3376
     3377
     3378proc SannfsLOTold(poly F, list #)
    31493379"USAGE:  SannfsLOT(f [,eng]);  f a poly, eng an optional int
    31503380RETURN:  ring
     
    33643594  poly F = x^3+y^3+z^3;
    33653595  printlevel = 0;
    3366   def A  = SannfsLOT(F);
     3596  def A  = SannfsLOTold(F);
    33673597  setring A;
    33683598  LD;
  • Singular/LIB/dmodapp.lib

    re67159 r4dcfc0f  
    11//////////////////////////////////////////////////////////////////////////////
    2 version="$Id: dmodapp.lib,v 1.1 2007-12-11 12:00:01 levandov Exp $";
     2version="$Id: dmodapp.lib,v 1.2 2008-02-05 22:44:29 levandov Exp $";
    33category="Noncommutative";
    44info="
     
    3838"USAGE:  charVariety(I);  I an ideal
    3939RETURN:  ring
    40 PURPOSE: compute the D-module structure of basering[1/f]*f^s with the algorithm given in S and with the Groebner basis engine given in 'eng'
     40PURPOSE: compute the characteristic variety of a D-module D/I
    4141ASSUME: the ground ring is the Weyl algebra with x's before d's
    4242NOTE:    activate the output ring with the @code{setring} command.
     
    7373  }
    7474  L[3] = NEW;
     75  list ncr =ncRelations(save);
     76  matrix @C = ncr[1];
     77  matrix @D = ncr[2];
    7578  def @U = ring(L);
    7679  // 2. create the commutative ring 
     
    8588  // comm ring is ready
    8689  setring @U;
     90  // make @U noncommutative
     91  //  matrix @C = imap(save,@C);
     92  matrix @D = imap(save,@D);
     93  def @@U = nc_algebra(@C,@D);
     94  setring @@U; kill @U;
    8795  // 2. compute Groebner basis
    8896  ideal I = imap(save,I);
    8997  //  I = groebner(I);
    90   I = slimgb(I);
     98  I = slimgb(I); // a bug?
    9199  setring @CU;
    92   ideal CV = imap(@U,I);
     100  ideal CV = imap(@@U,I);
    93101  //  CV = groebner(CV); // cosmetics
    94102  CV = slimgb(CV);
  • Singular/LIB/ratgb.lib

    re67159 r4dcfc0f  
    11//////////////////////////////////////////////////////////////////////////////
    2 version="$Id: ratgb.lib,v 1.9 2008-01-17 21:05:04 levandov Exp $";
     2version="$Id: ratgb.lib,v 1.10 2008-02-05 22:44:30 levandov Exp $";
    33category="Noncommutative";
    44info="
     
    1414LIB "poly.lib";
    1515
    16 static proc rm_content_id(ideal J)
     16//static
     17proc rm_content_id(ideal J)
    1718"USAGE:  rm_content_id(I);  I an ideal
    1819RETURN:  ideal
     
    170171  ideal CG = imap(save,G);
    171172  CG = rm_content_id(CG);
    172   CG = simplify(CG,2+32);
     173  CG = simplify(CG,2);
     174
     175  // warning: a bugfarm! in this ring, the ordering might change!!! (see appelF4)
     176  // so, simplify(32) should take place in the orig. ring! and NOT here
     177  //  CG = simplify(CG,2+32);
    173178
    174179  // 4b. create L(G) with X's as coeffs (for minimization)
    175   int sG  = ncols(CG);
     180  setring save;
     181  G = imap(@RAT,CG);
     182  int sG  = ncols(G);
    176183  ideal LG;
    177184  for (i=1; i<= sG; i++)
    178185  {
    179     LG[i] = lead(CG[i]);
    180   }
    181 
    182 //   ideal SLG = simplify(LG,8+32); //contains zeros
    183 //   intvec islg;
    184 //   if (SLG[1] == 0)
    185 //   {  islg = 0;  }
    186 //   else
    187 //   {    islg = 1;  }
    188 //   for (i=2; i<= ncols(SLG); i++)
    189 //   {
    190 //     if (SLG[i] == 0)
    191 //     {
    192 //       islg = islg, 0;
    193 //     }
    194 //     else
    195 //     {
    196 //       islg = islg, 1;
    197 //     }
    198 //   }
    199 
    200   // compute the D-dimension of the ideal
    201 
    202   LG = groebner(LG); // cosmetics
    203   int d = dim(LG);
     186    LG[i] = lead(G[i]);
     187  }
     188  // compute the D-dimension of the ideal in the ring @RAT
     189  setring @RAT;
     190  ideal LG = imap(save,LG);
     191  ideal LGG = groebner(LG); // cosmetics
     192  int d = dim(LGG);
    204193  int Ddim = d;
    205194  printf("the D-dimension is %s",d);
    206195  if (d==0)
    207196  {
    208     d = vdim(LG);
     197    d = vdim(LGG);
    209198    int Dvdim = d;
    210199    printf("the K-dimension is %s",d);
    211200  }
    212   setring save;
    213 //   for (i=1; i<= ncols(LG); i++)
    214 //   {
    215 //     if (islg[i] == 0)
    216 //     {
    217 //       G[i] = 0;
    218 //     }
    219 //   }
    220 //   G = simplify(G,2); // ready!
    221 
    222   G = imap(@RAT,CG);
     201  ideal SLG = simplify(LG,8+32); //contains zeros
     202  setring save;
     203  ideal SLG = imap(@RAT,SLG);
     204  // simplify(LG,8+32); //contains zeros
     205  intvec islg;
     206  if (SLG[1] == 0)
     207  {  islg = 0;  }
     208  else
     209  {    islg = 1;  }
     210  for (i=2; i<= ncols(SLG); i++)
     211  {
     212    if (SLG[i] == 0)
     213    {
     214      islg = islg, 0;
     215    }
     216    else
     217    {
     218      islg = islg, 1;
     219    }
     220  }
     221  for (i=1; i<= ncols(LG); i++)
     222  {
     223    if (islg[i] == 0)
     224    {
     225      G[i] = 0;
     226    }
     227  }
     228  G = simplify(G,2); // ready!
     229  //  G = imap(@RAT,CG);
    223230  // return the result
    224231  ideal pGBid = G;
     
    227234  //  export Dvdim;
    228235  setring @RAT;
    229   ideal rGBid = CG;
    230   //imap(save,G);
     236  ideal rGBid = imap(save,G);
     237  // CG;
    231238  export rGBid;
    232239  setring save;
Note: See TracChangeset for help on using the changeset viewer.