source:git/Singular/LIB/general.lib@803c5a1

jengelh-datetimespielwiese
Last change on this file since 803c5a1 was 803c5a1, checked in by Gert-Martin Greuel <greuel@…>, 22 years ago
* GMG: Kosmetik git-svn-id: file:///usr/local/Singular/svn/trunk@4981 2c84dea3-7e68-4137-9b89-c4e89433aadc
• Property mode set to `100644`
File size: 34.0 KB
Line
2//anne, added deleteSublist and watchdog 12.12.2000
3///////////////////////////////////////////////////////////////////////////////
4version="\$Id: general.lib,v 1.30 2000-12-22 13:52:31 greuel Exp \$";
5category="General purpose";
6info="
7LIBRARY:  general.lib   Elementary Computations of General Type
8
9PROCEDURES:
10 A_Z(\"a\",n);          string a,b,... of n comma separated letters
11 ASCII([n,m]);          string of printable ASCII characters (number n to m)
12 binomial(n,m[,../..]); n choose m (type int), [type string/type number]
13 deleteSublist(iv,l);   delete entries given by iv from list l
14 factorial(n[,../..]);  n factorial (=n!) (type int), [type string/number]
15 fibonacci(n[,p]);      nth Fibonacci number [char p]
16 kmemory([n[,v]]);      active [allocated] memory in kilobyte
17 killall();             kill all user-defined variables
18 number_e(n);           compute exp(1) up to n decimal digits
19 number_pi(n);          compute pi (area of unit circle) up to n digits
20 primes(n,m);           intvec of primes p, n<=p<=m
21 product(../..[,v]);    multiply components of vector/ideal/...[indices v]
22 ringweights(r);        intvec of weights of ring variables of ring r
23 sort(ideal/module);    sort generators according to monomial ordering
24 sum(vector/id/..[,v]); add components of vector/ideal/...[with indices v]
25 watchdog(i,cmd);       only wait for result of command cmd for i seconds
26 which(command);        search for command and return absolute path, if found
27           (parameters in square brackets [] are optional)
28";
29
30LIB "inout.lib";
31///////////////////////////////////////////////////////////////////////////////
32
33proc A_Z (string s,int n)
34"USAGE:   A_Z(\"a\",n);  a any letter, n integer (-26<= n <=26, !=0)
35RETURN:  string of n small (if a is small) or capital (if a is capital)
36         letters, comma separated, beginning with a, in alphabetical
37         order (or revers alphabetical order if n<0)
38EXAMPLE: example A_Z; shows an example
39"
40{
41  if ( n>=-26 and n<=26 and n!=0 )
42  {
43    string alpha =
44    "a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x,y,z,"+
45    "a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x,y,z,"+
46    "A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T,U,V,W,X,Y,Z,"+
47    "A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T,U,V,W,X,Y,Z";
48    int ii; int aa;
49    for(ii=1; ii<=51; ii=ii+2)
50    {
51       if( alpha[ii]==s ) { aa=ii; }
52    }
53    if ( aa==0)
54    {
55      for(ii=105; ii<=155; ii=ii+2)
56      {
57        if( alpha[ii]==s ) { aa=ii; }
58      }
59    }
60    if( aa!=0 )
61    {
62      string out;
63      if (n > 0) { out = alpha[aa,2*(n)-1];  return (out); }
64      if (n < 0)
65      {
66        string beta =
67        "z,y,x,w,v,u,t,s,r,q,p,o,n,m,l,k,j,i,h,g,f,e,d,c,b,a,"+
68        "z,y,x,w,v,u,t,s,r,q,p,o,n,m,l,k,j,i,h,g,f,e,d,c,b,a,"+
69        "Z,Y,X,W,V,U,T,S,R,Q,P,O,N,M,L,K,J,I,H,G,F,E,D,C,B,A,"+
70        "Z,Y,X,W,V,U,T,S,R,Q,P,O,N,M,L,K,J,I,H,G,F,E,D,C,B,A";
71        if ( aa < 52 ) { aa=52-aa; }
72        if ( aa > 104 ) { aa=260-aa; }
73        out = beta[aa,2*(-n)-1]; return (out);
74      }
75    }
76  }
77}
78example
79{ "EXAMPLE:"; echo = 2;
80   A_Z("c",5);
81   A_Z("Z",-5);
82   string sR = "ring R = (0,"+A_Z("A",6)+"),("+A_Z("a",10)+"),dp;";
83   sR;
84   execute(sR);
85   R;
86}
87///////////////////////////////////////////////////////////////////////////////
88proc ASCII (list #)
89"USAGE:   ASCII([n,m]); n,m= integers (32 <= n <= m <= 126)
90RETURN:   printable ASCII characters (no native language support)
91          ASCII():    string of  all ASCII characters with its numbers,
92                      no return value
93          ASCII(n):   string, n-th ASCII character
94          ASCII(n,m): list, n-th up to m-th ASCII character (inclusive)
95EXAMPLE: example ASCII; shows an example
96"
97{
98  string s1 =
99 "     !    \"    #    \$    %    &    '    (    )    *    +    ,    -    .
10032   33   34   35   36   37   38   39   40   41   42   43   44   45   46
101/    0    1    2    3    4    5    6    7    8    9    :    ;    <    =
10247   48   49   50   51   52   53   54   55   56   57   58   59   60   61
103>    ?    @    A    B    C    D    E    F    G    H    I    J    K    L
10462   63   64   65   66   67   68   69   70   71   72   73   74   75   76
105M    N    O    P    Q    R    S    T    U    V    W    X    Y    Z    [
10677   78   79   80   81   82   83   84   85   86   87   88   89   90   91
107\\    ]    ^    _    `    a    b    c    d    e    f    g    h    i    j
10892   93   94   95   96   97   98   99  100  101  102  103  104  105  10
109k    l    m    n    o    p    q    r    s    t    u    v    w    x    y
110107  108  109  110  111  112  113  114  115  116  117  118  119  120  121
111z    {    |    }    ~
112122  123  124  125  126 ";
113
114   string s2 =
115 " !\"#\$%&'()*+,-./0123456789:;<=>?@ABCDEFGHIJKLMNOPQRSTUVWXYZ[\\]^_`abcdefghijklmnopqrstuvwxyz{|}~";
116
117   if ( size(#) == 0 )
118   {
119      return(s1);
120   }
121   if ( size(#) == 1 )
122   {
123      return( s2[#[1]-31] );
124   }
125   if ( size(#) == 2 )
126   {
127      return( s2[#[1]-31,#[2]-#[1]+1] );
128   }
129}
130example
131{ "EXAMPLE:"; echo = 2;
132   ASCII();"";
133   ASCII(42);
134   ASCII(32,126);
135}
136///////////////////////////////////////////////////////////////////////////////
137
138proc binomial (int n, int k, list #)
139"USAGE:   binomial(n,k[,p]); n,k,p integers
140RETURN:  binomial(n,k); binomial coefficient n choose k,
141@*         - of type string (computed in characteristic 0)
142         binomial(n,k,p); n choose k, computed in characteristic prime(p)
143@*         - of type number if a basering is present and prime(p)=char(basering)
144@*         - of type string else
145NOTE:    In any characteristic, binomial(n,k) = coefficient of x^k in (1+x)^n
146EXAMPLE: example binomial; shows an example
147"
148{
149   int str,p;
150//---------------------------- initialization -------------------------------
151   if ( size(#) == 0 )
152   {  str = 1;
153      ring bin = 0,x,dp;
154      number r=1;
155   }
156   if ( size(#) > 0 )
157   {
158      p = (#[1]!=0)*prime(#[1]);
159      if ( defined(basering) )
160      {
161         if ( p == char(basering) )
162         {  number r=1;
163         }
164         else
165         {  str = 1;
166            ring bin = p,x,dp;
167            number r=1;
168         }
169      }
170      else
171      {  str = 1;
172         ring bin = p,x,dp;
173         number r=1;
174      }
175   }
176//-------------------------------- char 0 -----------------------------------
177   if ( p==0 )
178   {
179      r = binom0(n,k);
180   }
181//-------------------------------- char p -----------------------------------
182   else
183   {
184      r = binomp(n,k,p);
185   }
186//-------------------------------- return -----------------------------------
187   if ( str==1 ) { return(string(r)); }
188   else { return(r); }
189 }
190example
191{ "EXAMPLE:"; echo = 2;
192   binomial(200,100);"";                   //type string, computed in char 0
193   binomial(200,100,3);"";                 //type string, computed in char 3
194   int n,k = 200,100;
195   ring r = 0,x,dp;
196   number b1 = binomial(n,k,0);            //type number, computed in ring r
197   poly b2 = coeffs((x+1)^n,x)[k+1,1];     //coefficient of x^k in (x+1)^n
198   b1-b2;                                  //b1 and b2 should coincide
199}
200///////////////////////////////////////////////////////////////////////////////
201
202static proc binom0 (int n, int k)
203 //computes binomial coefficient n choose k in basering, assume 0<k<=n
204 //and char(basering) = 0 or n < char(basering)
205{
206   int l;
207   number r=1;
208   if ( k > n-k )
209   { k = n-k;
210   }
211   if ( k<=0 or k>n )               //trivial cases
212   { r = (k==0)*r;
213   }
214   for (l=1; l<=k; l++ )
215   {
216      r=r*(n+1-l)/l;
217   }
218   return(r);
219}
220///////////////////////////////////////////////////////////////////////////////
221
222static proc binomp (int n, int k, int p)
223 //computes binomial coefficient n choose k in basering of char p > 0
224 //binomial(n,k) = coefficient of x^k in (1+x)^n.
225 //Let n=q*p^j, gcd(q,p)=1, then (1+x)^n = (1 + x^(p^j))^q. We have
226 //binomial(n,k)=0 if k!=l*p^j and binomial(n,l*p^j) = binomial(q,l).
227 //Do this reduction first. Then, in denominator and numerator
228 //of defining formula for binomial coefficient, reduce those factors
229 //mod p which are not divisible by p and cancel common factors p. Hence,
230 //if n = h*p+r, k=l*p+s, r,s<p, binomial(n,k) = binomial(r,s)*binomial(h,l)
231{
232   int l,q,i= 1,n,1;
233   number r=1;
234   if ( k > n-k )
235   { k = n-k;
236   }
237   if ( k<=0 or k>n)               //trivial cases
238   { r = (k==0)*r;
239   }
240   else
241   {
242      while ( q mod p == 0 )
243      {  l = l*p;
244         q = q div p;
245      }                            //we have now n=q*l, l=p^j, gcd(q,p)=1;
246      if (k mod l != 0 )
247      { r = 0;
248      }
249      else
250      {  l = k div l;
251         n = q mod p;
252         k = l mod p;              //now 0<= k,n <p, use binom0 for n,k
253         q = q div p;              //recursion for q,l
254         l = l div p;              //use binomp for q,l
255         r = binom0(n,k)*binomp(q,l,p);
256      }
257   }
258   return(r);
259}
260///////////////////////////////////////////////////////////////////////////////
261
262proc factorial (int n, list #)
263"USAGE:   factorial(n[,p]);  n,p integers
264RETURN:  factorial(n):   n! (computed in characteristic 0), of type string
265         factorial(n,p): n! computed in characteristic prime(p)
266         - of type number if a basering is present and prime(p)=char(basering)
267         - of type string else
268EXAMPLE: example factorial; shows an example
269"
270{   int str,l,p;
271//---------------------------- initialization -------------------------------
272   if ( size(#) == 0 )
273   {  str = 1;
274      ring bin = 0,x,dp;
275      number r=1;
276   }
277   if ( size(#) > 0 )
278   {
279      p = (#[1]!=0)*prime(#[1]);
280      if ( defined(basering) )
281      {
282         if ( p == char(basering) )
283         {  number r=1;
284         }
285         else
286         {  str = 1;
287            ring bin = p,x,dp;
288            number r=1;
289         }
290      }
291      else
292      {  str = 1;
293         ring bin = p,x,dp;
294         number r=1;
295      }
296   }
297//------------------------------ computation --------------------------------
298   for (l=2; l<=n; l++)
299   {
300      r=r*l;
301   }
302   if ( str==1 ) { return(string(r)); }
303   else { return(r); }
304}
305example
306{ "EXAMPLE:"; echo = 2;
307   factorial(37);"";                   //37! of type string (as long integer)
308   ring r1 = 0,x,dp;
309   number p = factorial(37,0);         //37! of type number, computed in r
310   p;
311}
312///////////////////////////////////////////////////////////////////////////////
313
314proc fibonacci (int n, list #)
315"USAGE:   fibonacci(n);  n,p integers
316RETURN:  fibonacci(n): nth Fibonacci number, f(0)=f(1)=1, f(i+1)=f(i-1)+f(i)
317         - computed in characteristic 0, of type string
318         of type number computed in char(basering) if n is of type number
319         fibonacci(n,p): f(n) computed in characteristic prime(p)
320         - of type number if a basering is present and prime(p)=char(basering)
321         - of type string else
322EXAMPLE: example fibonacci; shows an example
323"
324{   int str,ii,p;
325//---------------------------- initialization -------------------------------
326   if ( size(#) == 0 )
327   {  str = 1;
328      ring bin = 0,x,dp;
329      number f,g,h=1,1,1;
330   }
331   if ( size(#) > 0 )
332   {
333      p = (#[1]!=0)*prime(#[1]);
334      if ( defined(basering) )
335      {
336         if ( p == char(basering) )
337         {  number f,g,h=1,1,1;
338         }
339         else
340         {  str = 1;
341            ring bin = p,x,dp;
342            number f,g,h=1,1,1;
343         }
344      }
345      else
346      {  str = 1;
347         ring bin = p,x,dp;
348         number f,g,h=1,1,1;
349      }
350   }
351//------------------------------ computation --------------------------------
352   for (ii=3; ii<=n; ii=ii+1)
353   {
354      h = f+g; f = g; g = h;
355    }
356   if ( str==1 ) { return(string(h)); }
357   else { return(h); }
358}
359example
360{ "EXAMPLE:"; echo = 2;
361   fibonacci(333); "";              //f(333) of type string (as long integer)
362   ring r = 17,x,dp;
363   number b = fibonacci(333,17);    //f(333) of type number, computed in r
364   b;
365}
366///////////////////////////////////////////////////////////////////////////////
367
368proc kmemory (list #)
369"USAGE:   kmemory([n,[v]]); n = int
370RETURN:  memory in kilobyte of type int
371         n=0: memory used by active variables (same as no parameters)
372         n=1: total memory allocated by Singular
373         n=2: difference between top and init memory adress (sbrk memory)
374         n!=0,1,2: 0
375DISPLAY: detailed information about allocated and used memory if v!=0
376NOTE:    kmemory uses internal function 'memory' to compute kilobyte, and
377         is the same as 'memory' for n!=0,1,2
378EXAMPLE: example kmemory; shows an example
379"
380{
381   int n;
382   int verb;
383   if (size(#) != 0)
384   {
385     n=#[1];
386     if (size(#) >1)
387     { verb=#[2]; }
388   }
389
390  if ( verb != 0)
391  {
392    if ( n==0)
393    { dbprint(printlevel-voice+3,
394      "// memory used, at the moment, by active variables (kilobyte):"); }
395    if ( n==1 )
396    { dbprint(printlevel-voice+3,
397      "// total memory allocated, at the moment, by SINGULAR (kilobyte):"); }
398   }
399   return ((memory(n)+1023)/1024);
400}
401example
402{ "EXAMPLE:"; echo = 2;
403   kmemory();
404   kmemory(1,1);
405}
406///////////////////////////////////////////////////////////////////////////////
407
408proc killall
409"USAGE:   killall(); (no parameter)
410         killall(\"type_name\");
411         killall(\"not\", \"type_name\");
412COMPUTE: killall(); kills all user-defined variables but not loaded procedures
413         killall(\"type_name\"); kills all user-defined variables,
414         of type \"type_name\"
415         killall(\"not\", \"type_name\"); kills all user-defined variables,
416         except those of type \"type_name\" and except loaded procedures
417         killall(\"not\", \"name_1\", \"name_2\", ...);
418         kills all user-defined variables, except those of name \"name_i\"
420RETURN:  no return value
421NOTE:    killall should never be used inside a procedure
422EXAMPLE: example killall; shows an example AND KILLS ALL YOUR VARIABLES
423"
424{
425  list L=names(); int joni=size(L);
426  int no_kill, j;
427  for (j=1; j<=size(#); j++)
428  {
429    if (typeof(#[j]) != "string")
430    {
431      ERROR("Need string as " + string(j) + "th argument");
432    }
433  }
434
435  // kills all user-defined variables but not loaded procedures
436  if( size(#)==0 )
437  {
438    for ( ; joni>0; joni-- )
439    {
440      if( L[joni]!="LIB" and typeof(`L[joni]`)!="proc" ) { kill `L[joni]`; }
441    }
442  }
443  else
444  {
445    // kills all user-defined variables
446    if( size(#)==1 )
447    {
448      // of type proc
449      if( #[1] == "proc" )
450      {
451        for ( joni=size(L); joni>0; joni-- )
452        {
453          if((L[joni]!="killall") and (L[joni]=="LIB" or typeof(`L[joni]`)=="proc"))
454          { kill `L[joni]`; }
455        }
456      }
457      else
458      {
459        // other types
460        for ( ; joni>2; joni-- )
461        {
462          if(typeof(`L[joni]`)==#[1] and L[joni]!="LIB" and typeof(`L[joni]`)!="proc") { kill `L[joni]`; }
463        }
464      }
465    }
466    else
467    {
468      // kills all user-defined variables whose name or type is not #i
469      for ( ; joni>2; joni-- )
470      {
471        if ( L[joni] != "LIB" && typeof(`L[joni]`) != "proc")
472        {
473          no_kill = 0;
474          for (j=2; j<= size(#); j++)
475          {
476            if (typeof(`L[joni]`)==#[j] or L[joni] == #[j])
477            {
478              no_kill = 1;
479              break;
480            }
481          }
482          if (! no_kill)
483          {
484            kill `L[joni]`;
485          }
486        }
487      }
488    }
489  }
490}
491example
492{ "EXAMPLE:"; echo = 2;
493   ring rtest; ideal i=x,y,z; string str="hi"; int j = 3;
494   export rtest,i,str,j;       //this makes the local variables global
495   listvar();
496   killall("ring");            // kills all rings
497   listvar();
498   killall("not", "int");      // kills all variables except int's (and procs)
499   listvar();
500   killall();                  // kills all vars except loaded procs
501   listvar();
502}
503///////////////////////////////////////////////////////////////////////////////
504
505proc number_e (int n)
506"USAGE:   number_e(n);  n integer
507COMPUTE: Euler number e=exp(1) up to n decimal digits (no rounding)
508         by A.H.J. Sale's algorithm
509RETURN:  - string of exp(1) if no basering of char 0 is defined;
510         - exp(1), type number, if a basering of char 0 is defined, display its
511         decimal format if printlevel >= voice (default:printlevel=voice-1 )
512EXAMPLE: example number_e; shows an example
513"
514{
515   int i,m,s,t;
516   intvec u,e;
517   u[n+2]=0; e[n+1]=0; e=e+1;
518   if( defined(basering) )
519   {
520      if( char(basering)==0 ) { number r=2; t=1; }
521   }
522   string result = "2.";
523   for( i=1; i<=n+1; i=i+1 )
524   {
525      e = e*10;
526      for( m=n+1; m>=1; m=m-1 )
527      {
528         s    = e[m]+u[m+1];
529         u[m] = s div (m+1);
530         e[m] = s%(m+1);
531      }
532      result = result+string(u[1]);
533      if( t==1 ) { r = r+number(u[1])/number(10)^i; }
534   }
535   if( t==1 )
536   { dbprint(printlevel-voice+2,"// "+result[1,n+1]);
537     return(r);
538   }
539   return(result[1,n+1]);
540}
541example
542{ "EXAMPLE:"; echo = 2;
543   number_e(30);"";
544   ring R = 0,t,lp;
545   number e = number_e(30);
546   e;
547}
548///////////////////////////////////////////////////////////////////////////////
549
550proc number_pi (int n)
551"USAGE:   number_pi(n);  n positive integer
552COMPUTE: pi (area of unit circle) up to n decimal digits (no rounding)
553         by algorithm of S. Rabinowitz
554RETURN:  - string of pi if no basering of char 0 is defined,
555         - pi, of type number, if a basering of char 0 is defined, display its
556         decimal format if printlevel >= voice (default:printlevel=voice-1 )
557EXAMPLE: example number_pi; shows an example
558"
559{
560   int i,m,t,e,q,N;
561   intvec r,p,B,Prelim;
562   string result,prelim;
563   N = (10*n) div 3 + 2;
564   p[N+1]=0; p=p+2; r=p;
565   for( i=1; i<=N+1; i=i+1 ) { B[i]=2*i-1; }
566   if( defined(basering) )
567   {
568      if( char(basering)==0 ) { number pi; number pri; t=1; }
569   }
570   for( i=0; i<=n; i=i+1 )
571   {
572      p = r*10;
573      e = p[N+1];
574      for( m=N+1; m>=2; m=m-1 )
575      {
576         r[m] = e%B[m];
577         q    = e div B[m];
578         e    = q*(m-1)+p[m-1];
579      }
580      r[1] = e%10;
581      q    = e div 10;
582      if( q!=10 and q!=9 )
583      {
584         result = result+prelim;
585         Prelim = q;
586         prelim = string(q);
587      }
588      if( q==9 )
589      {
590         Prelim = Prelim,9;
591         prelim = prelim+"9";
592      }
593      if( q==10 )
594      {
595         Prelim = (Prelim+1)-((Prelim+1) div 10)*10;
596         for( m=size(Prelim); m>0; m=m-1)
597         {
598            prelim[m] = string(Prelim[m]);
599         }
600         result = result+prelim;
601         if( t==1 ) { pi=pi+pri; }
602         Prelim = 0;
603         prelim = "0";
604      }
605      if( t==1 ) { pi=pi+number(q)/number(10)^i; }
606   }
607   result = result,prelim[1];
608   result = "3."+result[2,n-1];
609   if( t==1 )
610   { dbprint(printlevel-voice+2,"// "+result);
611     return(pi);
612   }
613   return(result);
614}
615example
616{ "EXAMPLE:"; echo = 2;
617   number_pi(11);"";
618   ring r = (real,10),t,dp;
619   number pi = number_pi(11); pi;
620}
621///////////////////////////////////////////////////////////////////////////////
622
623proc primes (int n, int m)
624"USAGE:   primes(n,m);  n,m integers
625RETURN:  intvec, consisting of all primes p, prime(n)<=p<=m, in increasing
626         order if n<=m, resp. prime(m)<=p<=n, in decreasing order if m<n
627NOTE:    prime(n); returns the biggest prime number <= n (if n>=2, else 2)
628EXAMPLE: example primes; shows an example
629"
630{  int change;
631   if ( n>m ) { change=n; n=m ; m=change; change=1; }
632   int q,p = prime(m),prime(n); intvec v = q; q = q-1;
633   while ( q>=p ) { q = prime(q); v = q,v; q = q-1; }
634   if ( change==1 ) { v = v[size(v)..1]; }
635   return(v);
636}
637example
638{  "EXAMPLE:"; echo = 2;
639    primes(50,100);"";
640    intvec v = primes(37,1); v;
641}
642///////////////////////////////////////////////////////////////////////////////
643
644proc product (id, list #)
645"USAGE:    product(id[,v]); id ideal/vector/module/matrix/intvec/intmat/list,
646          v intvec  (default: v=1.. number of entries of id)
647RETURN:   - if id is not a list: poly resp. int, the product of all entries of
648          id with index given by v.
649          id is treated as a list of polys resp. integers. A module m is
650          identified with corresponding matrix M (columns of M generate m)
651          - if id is a list: product of list entries, with index given by v.
652          Assume that list members can be multiplied
653EXAMPLE:  example product; shows an example
654"
655{
656   int n,j,tt;
657   string ty;
658   list l;
659   int s = size(#);
660   if( s!=0 )
661   {  if ( typeof(#[s])=="intvec" )
662      {  intvec v = #[s];
663         tt=1; s=s-1;
664         if ( s>0 ) { # = #[1..s]; }
665      }
666   }
667   if ( s>0 )
668   {
669     l = list(id)+#;
670     kill id;
671     list id = l;
672     ty = "list";
673   }
674   else
675   { ty = typeof(id);
676   }
677   if( ty=="list" )
678   { n = size(id);
679     def f(1) = id[1];
680     for( j=2; j<=n; j=j+1 ) { def f(j)=f(j-1)*id[j]; }
681     return(f(n));
682   }
683   if( ty=="poly" or ty=="ideal" or ty=="vector"
684       or ty=="module" or ty=="matrix" )
685   {
686      ideal i = ideal(matrix(id));
687      kill id;
688      ideal id = i;
689      if( tt!=0 ) { id = id[v]; }
690      n = ncols(id); poly f(1)=id[1];
691   }
692   if( ty=="int" or ty=="intvec" or ty=="intmat" )
693   {
694      if ( ty == "int" ) { intmat S =id; }
695      else { intmat S = intmat(id); }
696      intvec i = S[1..nrows(S),1..ncols(S)];
697      kill id;
698      intvec id = i;
699      if( tt!=0 ) { id = id[v]; }
700      n = size(id); int f(1)=id[1];
701   }
702   for( j=2; j<=n; j=j+1 ) { def f(j)=f(j-1)*id[j]; }
703   return(f(n));
704}
705example
706{  "EXAMPLE:"; echo = 2;
707   ring r= 0,(x,y,z),dp;
708   ideal m = maxideal(1);
709   product(m);
710   product(m[2..3]);
711   matrix M[2][3] = 1,x,2,y,3,z;
712   product(M);
713   intvec v=2,4,6;
714   product(M,v);
715   intvec iv = 1,2,3,4,5,6,7,8,9;
716   v=1..5,7,9;
717   product(iv,v);
718   intmat A[2][3] = 1,1,1,2,2,2;
719   product(A,3..5);
720}
721///////////////////////////////////////////////////////////////////////////////
722proc ringweights (list # )
723"USAGE:   ringweights (P); P=name of an existing ring (true name, not a string)
724RETURN:  intvec, size=nvars(P), consisting of the weights of the variables of P
725NOTE:    This is useful when enlarging P but keeping the weights of the old
726         variables
727EXAMPLE: example ringweights;  shows an example
728"
729{
730   int ii,q,fi,fo,fia;
731   intvec rw,nw;
732   string os;
733   def P = #[1];
734   string osP = ordstr(P);
735   fo  = 1;
736//------------------------- find weights in ordstr(P) -------------------------
737   fi  = find(osP,"(",fo);
738   fia = find(osP,"a",fo)+find(osP,"w",fo)+find(osP,"W",fo);
739   while ( fia )
740   {
741      os = osP[fi+1,find(osP,")",fi)-fi-1];
742      if( find(os,",") )
743      {
744         execute("nw = "+os+";");
745         if( size(nw) > ii )
746         {
747             rw = rw,nw[ii+1..size(nw)];
748         }
749         else  {  ii = ii - size(nw); }
750
751         if( find(osP[1,fi],"a",fo) ) { ii = size(nw); }
752      }
753      else
754      {
755         execute("q = "+os+";");
756         if( q > ii )
757         {
758            nw = 0; nw[q-ii] = 0;
759            nw = nw + 1;          //creates an intvec 1,...,1 of length q-ii
760            rw = rw,nw;
761         }
762         else { ii = ii - q; }
763      }
764      fo  = fi+1;
765      fi  = find(osP,"(",fo);
766      fia = find(osP,"a",fo)+find(osP,"w",fo)+find(osP,"W",fo);
767   }
768//-------------- adjust weight vector to length = nvars(P)  -------------------
769   if( fo > 1 )
770   {                                            // case when weights were found
771      rw = rw[2..size(rw)];
772      if( size(rw) > nvars(P) )
773      {
774         rw = rw[1..nvars(P)];
775      }
776      if( size(rw) < nvars(P) )
777      {
778         nw=0; nw[nvars(P)-size(rw)]=0; nw=nw+1; rw=rw,nw;
779      }
780   }
781   else
782   {                                         // case when no weights were found
783      rw[nvars(P)]= 0; rw=rw+1;
784   }
785   return(rw);
786}
787example
788{"EXAMPLE:";  echo = 2;
789  ring r0 = 0,(x,y,z),dp;
790  ringweights(r0);
791  ring r1 = 0,x(1..5),(ds(3),wp(2,3));
792  ringweights(r1);
793  ring r2 = 0,x(1..5),(a(1,2,3,0),dp);
794  ringweights(r2);
795  ring r3 = 0,x(1..10),(a(1..5),dp(5),a(10..13),Wp(5..9));
796  ringweights(r3);
797  // an example for enlarging the ring:
798  intvec v = 6,2,3,4,5;
799  ring R = 0,x(1..10),(a(ringweights(r1),v),dp);
800  ordstr(R);
801}
802
803///////////////////////////////////////////////////////////////////////////////
804
805proc sort (id, list #)
806"USAGE:   sort(id[v,o,n]); id=ideal/module/intvec/list (of intvec's or int's)
807         sort may be called with 1, 2 or 3 arguments in the following way:
808         sort(id[v,n]); v=intvec of positive integers, n=integer,
809         sort(id[o,n]); o=string (any allowed ordstr of a ring), n=integer
810RETURN:  a list of two elements:
811         [1]: object of same type as input but sorted in the following manner:
812           - if id=ideal/module: generators of id are sorted w.r.t. intvec v
813             (id[v[1]] becomes 1-st, id[v[2]]  2-nd element, etc.). If no v is
814             present, id is sorted w.r.t. ordering o (if o is given) or w.r.t.
815             actual monomial ordering (if no o is given):
816                    generators with smaller leading term come first
817             (e.g. sort(id); sorts w.r.t actual monomial ordering)
818           - if id=list of intvec's or int's: consider a list element, say
819             id[1]=3,2,5, as exponent vector of the monomial x^3*y^2*z^5;
820             the corresponding monomials are ordered w.r.t. intvec v (s.a.).
821             If no v is present, the monomials are sorted w.r.t. ordering o
822             (if o is given) or w.r.t. lexicographical ordering (if no o is
823             given). The corresponding ordered list of exponent vectors is
824             returned.
825             (e.g. sort(id); sorts lexicographically, smaller int's come first)
826             WARNING: Since negative exponents create the 0 polynomial in
827             Singular, id should not contain negative integers: the result
828             might not be as expected
829           - if id=intvec: id is treated as list of integers
830           - if n!=0 the ordering is inverse, i.e. w.r.t. v(size(v)..1)
831             default: n=0
832         [2]: intvec, describing the permutation of the input (hence [2]=v if
833             v is given (with positive integers)
834NOTE:    If v is given id may be any simply indexed object (e.g. any list or
835         string); if v[i]<0 and i<=size(id) v[i] is set internally to i;
836         entries of v must be pairwise distinct to get a permutation if id.
837         Zero generators of ideal/module are deleted
838EXAMPLE: example sort; shows an example
839"
840{  int ii,jj,s,n = 0,0,1,0;
841   intvec v;
842   if ( defined(basering) ) { def P = basering; }
843   if ( size(#)==0 and (typeof(id)=="ideal" or typeof(id)=="module") )
844   {
845      id = simplify(id,2);
846      for ( ii=1; ii<size(id); ii++ )
847      {
848         if ( id[ii]!=id[ii+1] ) { break;}
849      }
850      if ( ii != size(id) ) { v = sortvec(id); }
851      else  { v = size(id)..1; }
852   }
853   if ( size(#)>=1 and (typeof(id)=="ideal" or typeof(id)=="module") )
854   {
855      if ( typeof(#[1])=="string" )
856      {
857         execute("ring r1 =("+charstr(P)+"),("+varstr(P)+"),("+#[1]+");");
858         def i = imap(P,id);
859         v = sortvec(i);
860         setring P;
861         n=2;
862      }
863   }
864   if ( typeof(id)=="intvec" or typeof(id)=="list" and n==0 )
865   {
866      string o;
867      if ( size(#)==0 ) { o = "lp"; n=1; }
868      if ( size(#)>=1 )
869      {
870         if ( typeof(#[1])=="string" ) { o = #[1]; n=1; }
871      }
872   }
873   if ( typeof(id)=="intvec" or typeof(id)=="list" and n==1 )
874   {
875      if ( typeof(id)=="list" )
876      {
877         for (ii=1; ii<=size(id); ii++)
878         {
879            if (typeof(id[ii]) != "intvec" and typeof(id[ii]) != "int")
880               { "// list elements must be intvec/int"; return(); }
881            else
882               { s=size(id[ii])*(s < size(id[ii])) + s*(s >= size(id[ii])); }
883         }
884      }
885      execute("ring r=0,x(1..s),("+o+");");
886      ideal i;
887      poly f;
888      for (ii=1; ii<=size(id); ii++)
889      {
890         f=1;
891         for (jj=1; jj<=size(id[ii]); jj++)
892         {
893            f=f*x(jj)^(id[ii])[jj];
894         }
895         i[ii]=f;
896      }
897      v = sort(i)[2];
898   }
899   if ( size(#)!=0 and n==0 ) { v = #[1]; }
900   if( size(#)==2 )
901   {
902      if ( #[2] != 0 ) { v = v[size(v)..1]; }
903   }
904   s = size(v);
905   if( size(id) < s ) { s = size(id); }
906   def m = id;
907   if ( size(m) != 0 )
908   {
909      for ( jj=1; jj<=s; jj=jj+1)
910      {
911         if ( v[jj]<=0 ) { v[jj]=jj; }
912         m[jj] = id[v[jj]];
913      }
914   }
915   if ( v == 0 ) { v = 1; }
916   list L=m,v;
917   return(L);
918}
919example
920{  "EXAMPLE:"; echo = 2;
921   ring r0 = 0,(x,y,z,t),lp;
922   ideal i = x3,z3,xyz;
923   sort(i);                // sort w.r.t. lex ordering
924   sort(i,3..1);
925   sort(i,"ls")[1];        // sort w.r.t. negative lex ordering
926   list L =1,8..5,3..10;
927   sort(L)[1];             // sort L lexicographically
928   sort(L,"Dp",1)[1];      // sort L w.r.t (total sum, reverse lex)
929}
930///////////////////////////////////////////////////////////////////////////////
931
932proc sum (id, list #)
933"USAGE:    sum(id[,v]); id=ideal/vector/module/matrix resp. id=intvec/intmat,
934                       v=intvec (e.g. v=1..n, n=integer)
935RETURN:   poly resp. int which is the sum of all entries of id, with index
936          given by v (default: v=1..number of entries of id)
937NOTE:     id is treated as a list of polys resp. integers. A module m is
938          identified with corresponding matrix M (columns of M generate m)
939EXAMPLE:  example sum; shows an example
940"
941{
942   if( typeof(id)=="poly" or typeof(id)=="ideal" or typeof(id)=="vector"
943       or typeof(id)=="module" or typeof(id)=="matrix" )
944   {
945      ideal i = ideal(matrix(id));
946      if( size(#)!=0 ) { i = i[#[1]]; }
947      matrix Z = matrix(i);
948   }
949   if( typeof(id)=="int" or typeof(id)=="intvec" or typeof(id)=="intmat" )
950   {
951      if ( typeof(id) == "int" ) { intmat S =id; }
952      else { intmat S = intmat(id); }
953      intvec i = S[1..nrows(S),1..ncols(S)];
954      if( size(#)!=0 ) { i = i[#[1]]; }
955      intmat Z=transpose(i);
956   }
957   intvec v; v[ncols(Z)]=0; v=v+1;
958   return((Z*v)[1,1]);
959 }
960example
961{  "EXAMPLE:"; echo = 2;
962   ring r= 0,(x,y,z),dp;
963   vector pv = [xy,xz,yz,x2,y2,z2];
964   sum(pv);
965   sum(pv,2..5);
966   matrix M[2][3] = 1,x,2,y,3,z;
967   intvec w=2,4,6;
968   sum(M,w);
969   intvec iv = 1,2,3,4,5,6,7,8,9;
970   sum(iv,2..4);
971}
972///////////////////////////////////////////////////////////////////////////////
973
974proc which (command)
975"USAGE:    which(command); command = string expression
976RETURN:   Absolute pathname of command, if found in search path.
977          Empty string, otherwise.
978NOTE:     Based on the Unix command 'which'.
979EXAMPLE:  example which; shows an example
980"
981{
982   int rs;
983   int i;
984   string fn = "which_" + string(system("pid"));
985   string pn;
986   string cmd;
987   if( typeof(command) != "string")
988   {
989     return (pn);
990   }
991   if (system("uname") != "ix86-Win")
992   {
993     cmd = "which ";
994   }
995   else
996   {
997     // unfortunately, it does not take -path
998     cmd = "type ";
999   }
1000   i = system("sh", cmd + command + " > " + fn);
1002   if (system("uname") != "ix86-Win")
1003   {
1004     // TBC: Hmm... should parse output to get rid of 'command is '
1005     pn[size(pn)] = "";
1006     i = 1;
1007     while ((pn[i] != " ") and (pn[i] != ""))
1008     {
1009       i = i+1;
1010     }
1011     if (pn[i] == " ") {pn[i] = "";}
1012     rs = system("sh", "ls " + pn + " > " + fn + " 2>&1 ");
1013   }
1014   else
1015   {
1016     rs = 0;
1017   }
1018   i = system("sh", "rm " + fn);
1019   if (rs == 0) {return (pn);}
1020   else
1021   {
1023     return ("");
1024   }
1025}
1026example
1027{  "EXAMPLE:"; echo = 2;
1028    which("sh");
1029}
1030///////////////////////////////////////////////////////////////////////////////
1031
1032proc watchdog(int i, string cmd)
1033"USAGE  : watchdog(i,cmd); i integer; cmd string
1034RETURNS: Result of cmd, if the result can be computed in
1035         i seconds. Otherwise the computation is interrupted after
1036         i seconds, the string "Killed" is returned and the global
1037         variable 'watchdog_interrupt' is defined.
1038NOTE:    * the MP package must be enabled
1039         * the current basering should not be watchdog_rneu
1040         * if there are variable names of the structure x(i) all
1041           polynomials have to be put into eval(...) in order to be
1042           interpreted correctly
1043         * a second Singular process is started by this procedure
1044EXAMPLE: example watchdog; shows an example
1045"
1046{
1047  string rname=nameof(basering);
1048  if (defined(watchdog_rneu))
1049  {
1050    kill watchdog_rneu;
1051  }
1052// If we do not have MP-links, watchdog cannot be used
1053  if (system("with","MP"))
1054  {
1055    if ( i > 0 )
1056    {
1057      int j=10;
1058      int k=999999;
1059// fork, get the pid of the child and send it the command
1061      open(l_fork);
1062      write(l_fork,quote(system("pid")));
1064      execute("write(l_fork,quote(" + cmd + "));");
1065
1066
1067// sleep in small, but growing intervals for appr. 1 second
1068      while(j < k)
1069      {
1071        j = j + j;
1072      }
1073
1074// sleep in intervals of one second
1075      j = 1;
1077      {
1078        while (j < i)
1079        {
1081          j = j + 1;
1082        }
1083      }
1084// check, whether we have a result, and return it
1086      {
1088        if (nameof(basering)!=rname)
1089        {
1090          def watchdog_rneu=basering;
1091        }
1092        if(defined(watchdog_interrupt))
1093        {
1094          kill (watchdog_interrupt);
1095        }
1096        close(l_fork);
1097      }
1098      else
1099      {
1100        string result="Killed";
1101        if(!defined(watchdog_interrupt))
1102        {
1103          int watchdog_interrupt=1;
1104          export watchdog_interrupt;
1105        }
1106        close(l_fork);
1107        j = system("sh","kill " + string(pid));
1108      }
1109      if (defined(watchdog_rneu))
1110      {
1111        keepring watchdog_rneu;
1112      }
1113      return(result);
1114    }
1115    else
1116    {
1117      ERROR("First argument of watchdog has to be a positive integer.");
1118    }
1119    ERROR("MP-support is not enabled in this version of Singular.");
1120  }
1121}
1122example
1123{ "EXAMPLE:"; echo=2;
1124  ring r=0,(x,y,z),dp;
1125  poly f=x^30+y^30;
1126  watchdog(1,"factorize(eval("+string(f)+"))");
1127  watchdog(100,"factorize(eval("+string(f)+"))");
1128}
1129///////////////////////////////////////////////////////////////////////////////
1130
1131proc deleteSublist(intvec v,list l)
1132"USAGE:   deleteSublist(v,l); intvec v; list l
1133         where the entries of the integer vector v correspond to the
1134         positions of the elements to be deleted
1135RETURN:  list without the deleted elements
1136EXAMPLE: example deleteSublist; shows an example"
1137{
1138  list k;
1139  int i,j,skip;
1140  j=1;
1141  skip=0;
1142  intvec vs=sort(v)[1];
1143  for ( i=1 ; i <=size(vs) ; i++)
1144  {
1145    while ((j+skip) < vs[i])
1146    {
1147      k[j] = l[j+skip];
1148      j++;
1149    }
1150    skip++;
1151  }
1152  if(vs[size(vs)]<size(l))
1153  {
1154    k=k+list(l[(vs[size(vs)]+1)..size(l)]);
1155  }
1156  return(k);
1157}
1158example
1159{ "EXAMPLE:"; echo=2;
1160   list l=1,2,3,4,5;
1161   intvec v=1,3,4;
1162   l=deleteSublist(v,l);
1163   l;
1164}
1165///////////////////////////////////////////////////////////////////////////////
Note: See TracBrowser for help on using the repository browser.