source: git/Singular/LIB/general.lib @ 371848

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