[3d124a7] | 1 | /////////////////////////////////////////////////////////////////////////////// |
---|
[341696] | 2 | version="$Id$"; |
---|
[49998f] | 3 | category="General purpose"; |
---|
[5480da] | 4 | info=" |
---|
[803c5a1] | 5 | LIBRARY: general.lib Elementary Computations of General Type |
---|
[3d124a7] | 6 | |
---|
[f34c37c] | 7 | PROCEDURES: |
---|
[0b59f5] | 8 | A_Z(\"a\",n); string a,b,... of n comma separated letters |
---|
[63be42] | 9 | ASCII([n,m]); string of printable ASCII characters (number n to m) |
---|
[ca41246] | 10 | absValue(c); absolute value of c |
---|
[fc62b67] | 11 | binomial(n,m[,../..]); n choose m (type int), [type bigint] |
---|
[ebbe4a] | 12 | deleteSublist(iv,l); delete entries given by iv from list l |
---|
[fc62b67] | 13 | factorial(n[,../..]); n factorial (=n!) (type int), [type bigint] |
---|
[e67c6f] | 14 | fibonacci(n); nth Fibonacci number |
---|
[dd2aa36] | 15 | kmemory([n[,v]]); active [allocated] memory in kilobyte |
---|
[3d124a7] | 16 | killall(); kill all user-defined variables |
---|
| 17 | number_e(n); compute exp(1) up to n decimal digits |
---|
| 18 | number_pi(n); compute pi (area of unit circle) up to n digits |
---|
| 19 | primes(n,m); intvec of primes p, n<=p<=m |
---|
| 20 | product(../..[,v]); multiply components of vector/ideal/...[indices v] |
---|
| 21 | sort(ideal/module); sort generators according to monomial ordering |
---|
| 22 | sum(vector/id/..[,v]); add components of vector/ideal/...[with indices v] |
---|
[ebbe4a] | 23 | watchdog(i,cmd); only wait for result of command cmd for i seconds |
---|
[298d0a] | 24 | primecoeffs(J[,q]); primefactors <= min(p,32003) of coeffs of J |
---|
[90d772] | 25 | timeStd(i,d) std(i) if the standard basis computation finished |
---|
| 26 | after d-1 seconds and i otherwhise |
---|
| 27 | timeFactorize(p,d) factorize(p) if the factorization finished after d-1 |
---|
| 28 | seconds otherwhise f is considered to be irreducible |
---|
| 29 | factorH(p) changes variables to become the last variable the |
---|
| 30 | principal one in the multivariate factorization and |
---|
| 31 | factorizes then the polynomial |
---|
| 32 | |
---|
[194f5e5] | 33 | (parameters in square brackets [] are optional) |
---|
[5480da] | 34 | "; |
---|
| 35 | |
---|
[3d124a7] | 36 | LIB "inout.lib"; |
---|
[8b87364] | 37 | LIB "poly.lib"; |
---|
| 38 | LIB "matrix.lib"; |
---|
[3d124a7] | 39 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 40 | |
---|
| 41 | proc A_Z (string s,int n) |
---|
[d2b2a7] | 42 | "USAGE: A_Z(\"a\",n); a any letter, n integer (-26<= n <=26, !=0) |
---|
[3d124a7] | 43 | RETURN: string of n small (if a is small) or capital (if a is capital) |
---|
[0b59f5] | 44 | letters, comma separated, beginning with a, in alphabetical |
---|
[a7a00b] | 45 | order (or reverse alphabetical order if n<0) |
---|
[3d124a7] | 46 | EXAMPLE: example A_Z; shows an example |
---|
[d2b2a7] | 47 | " |
---|
[3d124a7] | 48 | { |
---|
| 49 | if ( n>=-26 and n<=26 and n!=0 ) |
---|
| 50 | { |
---|
| 51 | string alpha = |
---|
| 52 | "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,"+ |
---|
| 53 | "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,"+ |
---|
| 54 | "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,"+ |
---|
| 55 | "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"; |
---|
| 56 | int ii; int aa; |
---|
| 57 | for(ii=1; ii<=51; ii=ii+2) |
---|
| 58 | { |
---|
| 59 | if( alpha[ii]==s ) { aa=ii; } |
---|
| 60 | } |
---|
| 61 | if ( aa==0) |
---|
| 62 | { |
---|
| 63 | for(ii=105; ii<=155; ii=ii+2) |
---|
| 64 | { |
---|
| 65 | if( alpha[ii]==s ) { aa=ii; } |
---|
| 66 | } |
---|
| 67 | } |
---|
| 68 | if( aa!=0 ) |
---|
| 69 | { |
---|
| 70 | string out; |
---|
| 71 | if (n > 0) { out = alpha[aa,2*(n)-1]; return (out); } |
---|
| 72 | if (n < 0) |
---|
| 73 | { |
---|
| 74 | string beta = |
---|
| 75 | "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,"+ |
---|
| 76 | "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,"+ |
---|
| 77 | "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,"+ |
---|
| 78 | "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"; |
---|
| 79 | if ( aa < 52 ) { aa=52-aa; } |
---|
| 80 | if ( aa > 104 ) { aa=260-aa; } |
---|
| 81 | out = beta[aa,2*(-n)-1]; return (out); |
---|
| 82 | } |
---|
| 83 | } |
---|
| 84 | } |
---|
| 85 | } |
---|
| 86 | example |
---|
| 87 | { "EXAMPLE:"; echo = 2; |
---|
| 88 | A_Z("c",5); |
---|
| 89 | A_Z("Z",-5); |
---|
| 90 | string sR = "ring R = (0,"+A_Z("A",6)+"),("+A_Z("a",10)+"),dp;"; |
---|
| 91 | sR; |
---|
[034ce1] | 92 | execute(sR); |
---|
[3d124a7] | 93 | R; |
---|
| 94 | } |
---|
| 95 | /////////////////////////////////////////////////////////////////////////////// |
---|
[63be42] | 96 | proc ASCII (list #) |
---|
| 97 | "USAGE: ASCII([n,m]); n,m= integers (32 <= n <= m <= 126) |
---|
[a7a00b] | 98 | RETURN: string of printable ASCII characters (no native language support)@* |
---|
| 99 | ASCII(): string of all ASCII characters with its numbers,@* |
---|
| 100 | ASCII(n): n-th ASCII character@* |
---|
[b42ab6] | 101 | ASCII(n,m): n-th up to m-th ASCII character (inclusive) |
---|
[63be42] | 102 | EXAMPLE: example ASCII; shows an example |
---|
| 103 | " |
---|
| 104 | { |
---|
| 105 | string s1 = |
---|
[008846] | 106 | " ! \" # $ % & ' ( ) * + , - . |
---|
[63be42] | 107 | 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 |
---|
| 108 | / 0 1 2 3 4 5 6 7 8 9 : ; < = |
---|
| 109 | 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 |
---|
| 110 | > ? @ A B C D E F G H I J K L |
---|
| 111 | 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 |
---|
| 112 | M N O P Q R S T U V W X Y Z [ |
---|
| 113 | 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 |
---|
| 114 | \\ ] ^ _ ` a b c d e f g h i j |
---|
| 115 | 92 93 94 95 96 97 98 99 100 101 102 103 104 105 10 |
---|
| 116 | k l m n o p q r s t u v w x y |
---|
| 117 | 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 |
---|
| 118 | z { | } ~ |
---|
| 119 | 122 123 124 125 126 "; |
---|
| 120 | |
---|
| 121 | string s2 = |
---|
| 122 | " !\"#$%&'()*+,-./0123456789:;<=>?@ABCDEFGHIJKLMNOPQRSTUVWXYZ[\\]^_`abcdefghijklmnopqrstuvwxyz{|}~"; |
---|
| 123 | |
---|
| 124 | if ( size(#) == 0 ) |
---|
| 125 | { |
---|
| 126 | return(s1); |
---|
| 127 | } |
---|
| 128 | if ( size(#) == 1 ) |
---|
| 129 | { |
---|
| 130 | return( s2[#[1]-31] ); |
---|
| 131 | } |
---|
| 132 | if ( size(#) == 2 ) |
---|
| 133 | { |
---|
| 134 | return( s2[#[1]-31,#[2]-#[1]+1] ); |
---|
| 135 | } |
---|
| 136 | } |
---|
| 137 | example |
---|
| 138 | { "EXAMPLE:"; echo = 2; |
---|
| 139 | ASCII();""; |
---|
| 140 | ASCII(42); |
---|
| 141 | ASCII(32,126); |
---|
| 142 | } |
---|
| 143 | /////////////////////////////////////////////////////////////////////////////// |
---|
[3d124a7] | 144 | |
---|
[ca41246] | 145 | proc absValue(def c) |
---|
[298d0a] | 146 | "USAGE: absValue(c); c int, number or poly |
---|
[ca41246] | 147 | RETURN: absValue(c); the absolute value of c |
---|
| 148 | NOTE: absValue(c)=c if c>=0; absValue=-c if c<=0. |
---|
| 149 | @* So the function can be applied to any type, for which comparison |
---|
[298d0a] | 150 | @* operators are defined. |
---|
[ca41246] | 151 | EXAMPLE: example absValue; shows an example |
---|
| 152 | " |
---|
| 153 | { |
---|
| 154 | if (c>=0) { return(c); } |
---|
| 155 | else { return(-c); } |
---|
| 156 | } |
---|
| 157 | example |
---|
| 158 | { "EXAMPLE:"; echo = 2; |
---|
| 159 | ring r1 = 0,x,dp; |
---|
| 160 | absValue(-2002); |
---|
| 161 | |
---|
| 162 | poly f=-4; |
---|
| 163 | absValue(f); |
---|
| 164 | } |
---|
| 165 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 166 | |
---|
[fc62b67] | 167 | proc binomial (int n, int k) |
---|
| 168 | "USAGE: binomial(n,k); n,k integers |
---|
[b42ab6] | 169 | RETURN: binomial(n,k); binomial coefficient n choose k |
---|
[3d124a7] | 170 | EXAMPLE: example binomial; shows an example |
---|
[d2b2a7] | 171 | " |
---|
[3d124a7] | 172 | { |
---|
[fc62b67] | 173 | bigint l; |
---|
| 174 | bigint r=1; |
---|
| 175 | bigint kk=k; |
---|
| 176 | bigint nn=n; |
---|
| 177 | if ( k > n-k ) |
---|
| 178 | { k = n-k; } |
---|
| 179 | if ( k<=0 or k>n ) //trivial cases |
---|
| 180 | { r = (k==0)*r; } |
---|
| 181 | for (l=1; l<=kk; l=l+1 ) |
---|
[c860e9] | 182 | { |
---|
[fc62b67] | 183 | r=r*(nn+1-l)/l; |
---|
[c860e9] | 184 | } |
---|
[fc62b67] | 185 | return(r); |
---|
| 186 | } |
---|
[3d124a7] | 187 | example |
---|
| 188 | { "EXAMPLE:"; echo = 2; |
---|
[883957] | 189 | binomial(200,100);""; //type bigint |
---|
[c860e9] | 190 | int n,k = 200,100; |
---|
[7f3ad4] | 191 | bigint b1 = binomial(n,k); |
---|
[c860e9] | 192 | ring r = 0,x,dp; |
---|
[883957] | 193 | poly b2 = coeffs((x+1)^n,x)[k+1,1]; //coefficient of x^k in (x+1)^n |
---|
| 194 | b1-b2; //b1 and b2 should coincide |
---|
[c860e9] | 195 | } |
---|
| 196 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 197 | |
---|
[3d124a7] | 198 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 199 | |
---|
[fc62b67] | 200 | proc factorial (int n) |
---|
| 201 | "USAGE: factorial(n); n integer |
---|
| 202 | RETURN: factorial(n): n! of type bigint. |
---|
[b42ab6] | 203 | SEE ALSO: prime |
---|
| 204 | EXAMPLE: example factorial; shows an example |
---|
[d2b2a7] | 205 | " |
---|
[fc62b67] | 206 | { |
---|
| 207 | bigint r=1; |
---|
[c860e9] | 208 | //------------------------------ computation -------------------------------- |
---|
[a292211] | 209 | for (int l=2; l<=n; l++) |
---|
[3d124a7] | 210 | { |
---|
[f937e2] | 211 | r=r*l; |
---|
[3d124a7] | 212 | } |
---|
[fc62b67] | 213 | return(r); |
---|
[3d124a7] | 214 | } |
---|
| 215 | example |
---|
| 216 | { "EXAMPLE:"; echo = 2; |
---|
[e67c6f] | 217 | factorial(37); |
---|
[3d124a7] | 218 | } |
---|
| 219 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 220 | |
---|
[fc62b67] | 221 | proc fibonacci (int n) |
---|
| 222 | "USAGE: fibonacci(n); n,p integers |
---|
[b42ab6] | 223 | RETURN: fibonacci(n): nth Fibonacci number, f(0)=f(1)=1, f(i+1)=f(i-1)+f(i) |
---|
[fc62b67] | 224 | @* - computed in characteristic 0, of type bigint |
---|
[b42ab6] | 225 | SEE ALSO: prime |
---|
[3d124a7] | 226 | EXAMPLE: example fibonacci; shows an example |
---|
[d2b2a7] | 227 | " |
---|
[fc62b67] | 228 | { |
---|
| 229 | bigint f,g,h=1,1,1; |
---|
[c860e9] | 230 | //------------------------------ computation -------------------------------- |
---|
[e67c6f] | 231 | for (int ii=3; ii<=n; ii++) |
---|
[3d124a7] | 232 | { |
---|
[f937e2] | 233 | h = f+g; f = g; g = h; |
---|
[fc62b67] | 234 | } |
---|
| 235 | return(h); |
---|
[3d124a7] | 236 | } |
---|
| 237 | example |
---|
| 238 | { "EXAMPLE:"; echo = 2; |
---|
[e67c6f] | 239 | fibonacci(42); |
---|
[3d124a7] | 240 | } |
---|
| 241 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 242 | |
---|
[d694de] | 243 | proc kmemory (list #) |
---|
[b42ab6] | 244 | "USAGE: kmemory([n,[v]]); n,v integers |
---|
[a7a00b] | 245 | RETURN: memory in kilobyte of type bigint |
---|
[b42ab6] | 246 | @* n=0: memory used by active variables (same as no parameters) |
---|
| 247 | @* n=1: total memory allocated by Singular |
---|
[dd2aa36] | 248 | DISPLAY: detailed information about allocated and used memory if v!=0 |
---|
[d694de] | 249 | NOTE: kmemory uses internal function 'memory' to compute kilobyte, and |
---|
| 250 | is the same as 'memory' for n!=0,1,2 |
---|
[3d124a7] | 251 | EXAMPLE: example kmemory; shows an example |
---|
[d2b2a7] | 252 | " |
---|
[917fb5] | 253 | { |
---|
[d694de] | 254 | int n; |
---|
[dd2aa36] | 255 | int verb; |
---|
[d694de] | 256 | if (size(#) != 0) |
---|
| 257 | { |
---|
| 258 | n=#[1]; |
---|
[dd2aa36] | 259 | if (size(#) >1) |
---|
| 260 | { verb=#[2]; } |
---|
[d694de] | 261 | } |
---|
[917fb5] | 262 | |
---|
[dd2aa36] | 263 | if ( verb != 0) |
---|
| 264 | { |
---|
| 265 | if ( n==0) |
---|
| 266 | { dbprint(printlevel-voice+3, |
---|
| 267 | "// memory used, at the moment, by active variables (kilobyte):"); } |
---|
| 268 | if ( n==1 ) |
---|
| 269 | { dbprint(printlevel-voice+3, |
---|
| 270 | "// total memory allocated, at the moment, by SINGULAR (kilobyte):"); } |
---|
[c60d60] | 271 | } |
---|
| 272 | return ((memory(n)+1023) div 1024); |
---|
[3d124a7] | 273 | } |
---|
| 274 | example |
---|
| 275 | { "EXAMPLE:"; echo = 2; |
---|
| 276 | kmemory(); |
---|
[dd2aa36] | 277 | kmemory(1,1); |
---|
[3d124a7] | 278 | } |
---|
| 279 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 280 | |
---|
| 281 | proc killall |
---|
[d2b2a7] | 282 | "USAGE: killall(); (no parameter) |
---|
| 283 | killall(\"type_name\"); |
---|
| 284 | killall(\"not\", \"type_name\"); |
---|
[b42ab6] | 285 | RETURN: killall(); kills all user-defined variables except loaded procedures, |
---|
| 286 | no return value. |
---|
| 287 | @* - killall(\"type_name\"); kills all user-defined variables, |
---|
| 288 | of type \"type_name\" |
---|
| 289 | @* - killall(\"not\", \"type_name\"); kills all user-defined variables, |
---|
| 290 | except those of type \"type_name\" and except loaded procedures |
---|
[65546eb] | 291 | @* - killall(\"not\", \"name_1\", \"name_2\", ...); |
---|
| 292 | kills all user-defined variables, except those of name \"name_i\" |
---|
[b42ab6] | 293 | and except loaded procedures |
---|
[3d124a7] | 294 | NOTE: killall should never be used inside a procedure |
---|
| 295 | EXAMPLE: example killall; shows an example AND KILLS ALL YOUR VARIABLES |
---|
[d2b2a7] | 296 | " |
---|
[3d124a7] | 297 | { |
---|
[01cda73] | 298 | list @marie=names(Top); |
---|
[09f420] | 299 | int j, no_kill, @joni; |
---|
[b42ab6] | 300 | for ( @joni=1; @joni<=size(#); @joni++) |
---|
[5c187b] | 301 | { |
---|
[b42ab6] | 302 | if (typeof(#[@joni]) != "string") |
---|
[5c187b] | 303 | { |
---|
[b42ab6] | 304 | ERROR("Need string as " + string(@joni) + "th argument"); |
---|
[5c187b] | 305 | } |
---|
| 306 | } |
---|
[65546eb] | 307 | |
---|
[5c187b] | 308 | // kills all user-defined variables but not loaded procedures |
---|
| 309 | if( size(#)==0 ) |
---|
| 310 | { |
---|
[b42ab6] | 311 | for ( @joni=size(@marie); @joni>0; @joni-- ) |
---|
[5c187b] | 312 | { |
---|
[48c165a] | 313 | if( @marie[@joni]!="LIB" and typeof(`@marie[@joni]`)!="proc" |
---|
| 314 | and typeof(`@marie[@joni]`)!="package") |
---|
[b42ab6] | 315 | { kill `@marie[@joni]`; } |
---|
[5c187b] | 316 | } |
---|
| 317 | } |
---|
| 318 | else |
---|
| 319 | { |
---|
| 320 | // kills all user-defined variables |
---|
| 321 | if( size(#)==1 ) |
---|
| 322 | { |
---|
| 323 | // of type proc |
---|
| 324 | if( #[1] == "proc" ) |
---|
[3d124a7] | 325 | { |
---|
[b42ab6] | 326 | for ( @joni=size(@marie); @joni>0; @joni-- ) |
---|
[5c187b] | 327 | { |
---|
[298d0a] | 328 | if( (@marie[@joni]!="General") |
---|
[48c165a] | 329 | and (@marie[@joni]!="Top") |
---|
| 330 | and (@marie[@joni]!="killall") |
---|
[298d0a] | 331 | and (@marie[@joni]!="LIB") and |
---|
| 332 | ((typeof(`@marie[@joni]`)=="package") or |
---|
| 333 | (typeof(`@marie[@joni]`)=="proc"))) |
---|
| 334 | { |
---|
| 335 | if (defined(`@marie[@joni]`)) {kill `@marie[@joni]`;} |
---|
| 336 | } |
---|
| 337 | if (!defined(@joni)) break; |
---|
[48c165a] | 338 | } |
---|
[01cda73] | 339 | if (defined(General)) |
---|
[48c165a] | 340 | { |
---|
| 341 | @marie=names(General); |
---|
| 342 | for ( @joni=size(@marie); @joni>0; @joni-- ) |
---|
| 343 | { |
---|
| 344 | if(@marie[@joni]!="killall" |
---|
| 345 | and typeof(`@marie[@joni]`)=="proc") |
---|
| 346 | { kill General::`@marie[@joni]`; } |
---|
| 347 | } |
---|
| 348 | kill General::killall; |
---|
[5c187b] | 349 | } |
---|
[3d124a7] | 350 | } |
---|
[5c187b] | 351 | else |
---|
[65546eb] | 352 | { |
---|
[5c187b] | 353 | // other types |
---|
[b42ab6] | 354 | for ( @joni=size(@marie); @joni>2; @joni-- ) |
---|
[5c187b] | 355 | { |
---|
[65546eb] | 356 | if(typeof(`@marie[@joni]`)==#[1] and @marie[@joni]!="LIB" |
---|
| 357 | and typeof(`@marie[@joni]`)!="proc") |
---|
[b42ab6] | 358 | { kill `@marie[@joni]`; } |
---|
[5c187b] | 359 | } |
---|
| 360 | } |
---|
| 361 | } |
---|
| 362 | else |
---|
| 363 | { |
---|
[65546eb] | 364 | // kills all user-defined variables whose name or type is not #i |
---|
[b42ab6] | 365 | for ( @joni=size(@marie); @joni>2; @joni-- ) |
---|
[5c187b] | 366 | { |
---|
[a1b1dd] | 367 | if ( @marie[@joni] != "LIB" && @marie[@joni] != "Top" |
---|
| 368 | && typeof(`@marie[@joni]`) != "proc") |
---|
[5c187b] | 369 | { |
---|
| 370 | no_kill = 0; |
---|
[09f420] | 371 | for (j=2; j<= size(#); j++) |
---|
[6f2edc] | 372 | { |
---|
[b42ab6] | 373 | if (typeof(`@marie[@joni]`)==#[j] or @marie[@joni] == #[j]) |
---|
[5c187b] | 374 | { |
---|
| 375 | no_kill = 1; |
---|
| 376 | break; |
---|
| 377 | } |
---|
[6f2edc] | 378 | } |
---|
[5c187b] | 379 | if (! no_kill) |
---|
[6f2edc] | 380 | { |
---|
[b42ab6] | 381 | kill `@marie[@joni]`; |
---|
[6f2edc] | 382 | } |
---|
| 383 | } |
---|
[298d0a] | 384 | if (!defined(@joni)) break; |
---|
[5c187b] | 385 | } |
---|
| 386 | } |
---|
[6f2edc] | 387 | } |
---|
[3d124a7] | 388 | } |
---|
| 389 | example |
---|
| 390 | { "EXAMPLE:"; echo = 2; |
---|
[c860e9] | 391 | ring rtest; ideal i=x,y,z; string str="hi"; int j = 3; |
---|
| 392 | export rtest,i,str,j; //this makes the local variables global |
---|
| 393 | listvar(); |
---|
| 394 | killall("ring"); // kills all rings |
---|
| 395 | listvar(); |
---|
| 396 | killall("not", "int"); // kills all variables except int's (and procs) |
---|
| 397 | listvar(); |
---|
| 398 | killall(); // kills all vars except loaded procs |
---|
| 399 | listvar(); |
---|
[3d124a7] | 400 | } |
---|
| 401 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 402 | |
---|
| 403 | proc number_e (int n) |
---|
[d2b2a7] | 404 | "USAGE: number_e(n); n integer |
---|
[b42ab6] | 405 | RETURN: Euler number e=exp(1) up to n decimal digits (no rounding) |
---|
| 406 | @* - of type string if no basering of char 0 is defined |
---|
| 407 | @* - of type number if a basering of char 0 is defined |
---|
| 408 | DISPLAY: decimal format of e if printlevel > 0 (default:printlevel=0 ) |
---|
| 409 | NOTE: procedure uses algorithm of A.H.J. Sale |
---|
[3d124a7] | 410 | EXAMPLE: example number_e; shows an example |
---|
[d2b2a7] | 411 | " |
---|
[3d124a7] | 412 | { |
---|
| 413 | int i,m,s,t; |
---|
| 414 | intvec u,e; |
---|
| 415 | u[n+2]=0; e[n+1]=0; e=e+1; |
---|
| 416 | if( defined(basering) ) |
---|
| 417 | { |
---|
| 418 | if( char(basering)==0 ) { number r=2; t=1; } |
---|
| 419 | } |
---|
| 420 | string result = "2."; |
---|
| 421 | for( i=1; i<=n+1; i=i+1 ) |
---|
| 422 | { |
---|
| 423 | e = e*10; |
---|
| 424 | for( m=n+1; m>=1; m=m-1 ) |
---|
| 425 | { |
---|
| 426 | s = e[m]+u[m+1]; |
---|
[18dd47] | 427 | u[m] = s div (m+1); |
---|
[3d124a7] | 428 | e[m] = s%(m+1); |
---|
| 429 | } |
---|
| 430 | result = result+string(u[1]); |
---|
| 431 | if( t==1 ) { r = r+number(u[1])/number(10)^i; } |
---|
| 432 | } |
---|
[c860e9] | 433 | if( t==1 ) |
---|
| 434 | { dbprint(printlevel-voice+2,"// "+result[1,n+1]); |
---|
| 435 | return(r); |
---|
| 436 | } |
---|
[3d124a7] | 437 | return(result[1,n+1]); |
---|
| 438 | } |
---|
| 439 | example |
---|
| 440 | { "EXAMPLE:"; echo = 2; |
---|
[c860e9] | 441 | number_e(30);""; |
---|
[3d124a7] | 442 | ring R = 0,t,lp; |
---|
[c860e9] | 443 | number e = number_e(30); |
---|
[3d124a7] | 444 | e; |
---|
| 445 | } |
---|
| 446 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 447 | |
---|
| 448 | proc number_pi (int n) |
---|
[d2b2a7] | 449 | "USAGE: number_pi(n); n positive integer |
---|
[b42ab6] | 450 | RETURN: pi (area of unit circle) up to n decimal digits (no rounding) |
---|
| 451 | @* - of type string if no basering of char 0 is defined, |
---|
| 452 | @* - of type number, if a basering of char 0 is defined |
---|
| 453 | DISPLAY: decimal format of pi if printlevel > 0 (default:printlevel=0 ) |
---|
| 454 | NOTE: procedure uses algorithm of S. Rabinowitz |
---|
[3d124a7] | 455 | EXAMPLE: example number_pi; shows an example |
---|
[d2b2a7] | 456 | " |
---|
[3d124a7] | 457 | { |
---|
| 458 | int i,m,t,e,q,N; |
---|
| 459 | intvec r,p,B,Prelim; |
---|
| 460 | string result,prelim; |
---|
[18dd47] | 461 | N = (10*n) div 3 + 2; |
---|
[3d124a7] | 462 | p[N+1]=0; p=p+2; r=p; |
---|
| 463 | for( i=1; i<=N+1; i=i+1 ) { B[i]=2*i-1; } |
---|
| 464 | if( defined(basering) ) |
---|
| 465 | { |
---|
| 466 | if( char(basering)==0 ) { number pi; number pri; t=1; } |
---|
| 467 | } |
---|
| 468 | for( i=0; i<=n; i=i+1 ) |
---|
| 469 | { |
---|
| 470 | p = r*10; |
---|
| 471 | e = p[N+1]; |
---|
| 472 | for( m=N+1; m>=2; m=m-1 ) |
---|
| 473 | { |
---|
| 474 | r[m] = e%B[m]; |
---|
[18dd47] | 475 | q = e div B[m]; |
---|
[3d124a7] | 476 | e = q*(m-1)+p[m-1]; |
---|
| 477 | } |
---|
| 478 | r[1] = e%10; |
---|
[18dd47] | 479 | q = e div 10; |
---|
[3d124a7] | 480 | if( q!=10 and q!=9 ) |
---|
| 481 | { |
---|
| 482 | result = result+prelim; |
---|
| 483 | Prelim = q; |
---|
| 484 | prelim = string(q); |
---|
| 485 | } |
---|
| 486 | if( q==9 ) |
---|
| 487 | { |
---|
| 488 | Prelim = Prelim,9; |
---|
| 489 | prelim = prelim+"9"; |
---|
| 490 | } |
---|
| 491 | if( q==10 ) |
---|
| 492 | { |
---|
[18dd47] | 493 | Prelim = (Prelim+1)-((Prelim+1) div 10)*10; |
---|
[3d124a7] | 494 | for( m=size(Prelim); m>0; m=m-1) |
---|
| 495 | { |
---|
| 496 | prelim[m] = string(Prelim[m]); |
---|
| 497 | } |
---|
| 498 | result = result+prelim; |
---|
| 499 | if( t==1 ) { pi=pi+pri; } |
---|
| 500 | Prelim = 0; |
---|
| 501 | prelim = "0"; |
---|
| 502 | } |
---|
| 503 | if( t==1 ) { pi=pi+number(q)/number(10)^i; } |
---|
| 504 | } |
---|
| 505 | result = result,prelim[1]; |
---|
| 506 | result = "3."+result[2,n-1]; |
---|
[c860e9] | 507 | if( t==1 ) |
---|
| 508 | { dbprint(printlevel-voice+2,"// "+result); |
---|
| 509 | return(pi); |
---|
| 510 | } |
---|
[3d124a7] | 511 | return(result); |
---|
| 512 | } |
---|
| 513 | example |
---|
| 514 | { "EXAMPLE:"; echo = 2; |
---|
[c860e9] | 515 | number_pi(11);""; |
---|
| 516 | ring r = (real,10),t,dp; |
---|
| 517 | number pi = number_pi(11); pi; |
---|
[3d124a7] | 518 | } |
---|
| 519 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 520 | |
---|
| 521 | proc primes (int n, int m) |
---|
[d2b2a7] | 522 | "USAGE: primes(n,m); n,m integers |
---|
[3d124a7] | 523 | RETURN: intvec, consisting of all primes p, prime(n)<=p<=m, in increasing |
---|
[b42ab6] | 524 | order if n<=m, resp. prime(m)<=p<=n, in decreasing order if m<n. |
---|
| 525 | NOTE: prime(n); returns the biggest prime number <= min(n,32003) |
---|
| 526 | if n>=2, else 2 |
---|
[3d124a7] | 527 | EXAMPLE: example primes; shows an example |
---|
[d2b2a7] | 528 | " |
---|
[3d124a7] | 529 | { int change; |
---|
| 530 | if ( n>m ) { change=n; n=m ; m=change; change=1; } |
---|
| 531 | int q,p = prime(m),prime(n); intvec v = q; q = q-1; |
---|
| 532 | while ( q>=p ) { q = prime(q); v = q,v; q = q-1; } |
---|
| 533 | if ( change==1 ) { v = v[size(v)..1]; } |
---|
| 534 | return(v); |
---|
| 535 | } |
---|
| 536 | example |
---|
| 537 | { "EXAMPLE:"; echo = 2; |
---|
[c860e9] | 538 | primes(50,100);""; |
---|
| 539 | intvec v = primes(37,1); v; |
---|
[3d124a7] | 540 | } |
---|
| 541 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 542 | |
---|
| 543 | proc product (id, list #) |
---|
[c860e9] | 544 | "USAGE: product(id[,v]); id ideal/vector/module/matrix/intvec/intmat/list, |
---|
[b42ab6] | 545 | v intvec (default: v=1..number of entries of id) |
---|
| 546 | ASSUME: list members can be multiplied. |
---|
[65546eb] | 547 | RETURN: The product of all entries of id [with index given by v] of type |
---|
[b42ab6] | 548 | depending on the entries of id. |
---|
| 549 | NOTE: If id is not a list, id is treated as a list of polys resp. integers. |
---|
| 550 | A module m is identified with the corresponding matrix M (columns |
---|
| 551 | of M generate m). |
---|
[7708934] | 552 | @* If v is outside the range of id, we have the empty product and the |
---|
| 553 | result will be 1 (of type int). |
---|
[3d124a7] | 554 | EXAMPLE: example product; shows an example |
---|
[d2b2a7] | 555 | " |
---|
[65546eb] | 556 | { |
---|
[7708934] | 557 | //-------------------- initialization and special feature --------------------- |
---|
[c860e9] | 558 | int n,j,tt; |
---|
[65546eb] | 559 | string ty; //will become type of id |
---|
[c860e9] | 560 | list l; |
---|
[7708934] | 561 | |
---|
| 562 | // We wish to allow something like product(x(1..10)) if x(1),...,x(10) are |
---|
[65546eb] | 563 | // variables. x(1..10) is a list of polys and enters the procedure with |
---|
[7708934] | 564 | // id=x(1) and # a list with 9 polys, #[1]= x(2),...,#[9]= x(10). Hence, in |
---|
| 565 | // this case # is never empty. If an additional intvec v is given, |
---|
| 566 | // it is added to #, so we have to separate it first and make |
---|
| 567 | // the rest a list which has to be multiplied. |
---|
| 568 | |
---|
[c860e9] | 569 | int s = size(#); |
---|
| 570 | if( s!=0 ) |
---|
[65546eb] | 571 | { if ( typeof(#[s])=="intvec" or typeof(#[s])=="int") |
---|
| 572 | { |
---|
[7708934] | 573 | intvec v = #[s]; |
---|
[65546eb] | 574 | tt=1; |
---|
[7708934] | 575 | s=s-1; |
---|
[c860e9] | 576 | if ( s>0 ) { # = #[1..s]; } |
---|
| 577 | } |
---|
| 578 | } |
---|
| 579 | if ( s>0 ) |
---|
| 580 | { |
---|
[7708934] | 581 | l = list(id)+#; |
---|
| 582 | kill id; |
---|
| 583 | list id = l; //case: id = list |
---|
| 584 | ty = "list"; |
---|
| 585 | n = size(id); |
---|
[c860e9] | 586 | } |
---|
| 587 | else |
---|
[65546eb] | 588 | { |
---|
[7708934] | 589 | ty = typeof(id); |
---|
[65546eb] | 590 | if( ty == "list" ) |
---|
[7708934] | 591 | { n = size(id); } |
---|
[c860e9] | 592 | } |
---|
[7708934] | 593 | //------------------------------ reduce to 3 cases --------------------------- |
---|
[c860e9] | 594 | if( ty=="poly" or ty=="ideal" or ty=="vector" |
---|
| 595 | or ty=="module" or ty=="matrix" ) |
---|
[3d124a7] | 596 | { |
---|
| 597 | ideal i = ideal(matrix(id)); |
---|
[c860e9] | 598 | kill id; |
---|
[7708934] | 599 | ideal id = i; //case: id = ideal |
---|
| 600 | n = ncols(id); |
---|
[3d124a7] | 601 | } |
---|
[c860e9] | 602 | if( ty=="int" or ty=="intvec" or ty=="intmat" ) |
---|
[3d124a7] | 603 | { |
---|
[c860e9] | 604 | if ( ty == "int" ) { intmat S =id; } |
---|
[3d124a7] | 605 | else { intmat S = intmat(id); } |
---|
| 606 | intvec i = S[1..nrows(S),1..ncols(S)]; |
---|
[c860e9] | 607 | kill id; |
---|
[7708934] | 608 | intvec id = i; //case: id = intvec |
---|
[65546eb] | 609 | n = size(id); |
---|
[7708934] | 610 | } |
---|
| 611 | //--------------- consider intvec v and empty product ----------------------- |
---|
[65546eb] | 612 | if( tt!=0 ) |
---|
[7708934] | 613 | { |
---|
| 614 | for (j=1; j<=size(v); j++) |
---|
| 615 | { |
---|
| 616 | if ( v[j] <= 0 or v[j] > n ) //v outside range of id |
---|
[65546eb] | 617 | { |
---|
[7708934] | 618 | return(1); //empty product is 1 |
---|
[65546eb] | 619 | } |
---|
[7708934] | 620 | } |
---|
| 621 | id = id[v]; //consider part of id |
---|
| 622 | } //corresponding to v |
---|
| 623 | //--------------------- special case: one factor is zero --------------------- |
---|
| 624 | if ( typeof(id) == "ideal") |
---|
| 625 | { |
---|
| 626 | if( size(id) < ncols(id) ) |
---|
| 627 | { |
---|
| 628 | poly f; return(f); |
---|
| 629 | } |
---|
[3d124a7] | 630 | } |
---|
[7708934] | 631 | //-------------------------- finally, multiply objects ----------------------- |
---|
[65546eb] | 632 | n = size(id); |
---|
[7708934] | 633 | def f(1) = id[1]; |
---|
[c860e9] | 634 | for( j=2; j<=n; j=j+1 ) { def f(j)=f(j-1)*id[j]; } |
---|
| 635 | return(f(n)); |
---|
[3d124a7] | 636 | } |
---|
| 637 | example |
---|
| 638 | { "EXAMPLE:"; echo = 2; |
---|
| 639 | ring r= 0,(x,y,z),dp; |
---|
| 640 | ideal m = maxideal(1); |
---|
| 641 | product(m); |
---|
[c860e9] | 642 | product(m[2..3]); |
---|
[3d124a7] | 643 | matrix M[2][3] = 1,x,2,y,3,z; |
---|
| 644 | product(M); |
---|
| 645 | intvec v=2,4,6; |
---|
| 646 | product(M,v); |
---|
| 647 | intvec iv = 1,2,3,4,5,6,7,8,9; |
---|
| 648 | v=1..5,7,9; |
---|
| 649 | product(iv,v); |
---|
| 650 | intmat A[2][3] = 1,1,1,2,2,2; |
---|
| 651 | product(A,3..5); |
---|
| 652 | } |
---|
| 653 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 654 | |
---|
| 655 | proc sort (id, list #) |
---|
[29c4ee] | 656 | "USAGE: sort(id[,v,o,n]); id = ideal/module/intvec/list |
---|
[b42ab6] | 657 | @* sort may be called with 1, 2 or 3 arguments in the following way: |
---|
[6e52cf] | 658 | @* sort(id[,v,n]); v=intvec of positive integers, n=integer, |
---|
| 659 | @* sort(id[,o,n]); o=string (any allowed ordstr of a ring), n=integer |
---|
[b42ab6] | 660 | RETURN: a list l of two elements: |
---|
| 661 | @format |
---|
| 662 | l[1]: object of same type as input but sorted in the following way: |
---|
[3d124a7] | 663 | - if id=ideal/module: generators of id are sorted w.r.t. intvec v |
---|
| 664 | (id[v[1]] becomes 1-st, id[v[2]] 2-nd element, etc.). If no v is |
---|
| 665 | present, id is sorted w.r.t. ordering o (if o is given) or w.r.t. |
---|
| 666 | actual monomial ordering (if no o is given): |
---|
[b42ab6] | 667 | NOTE: generators with SMALLER(!) leading term come FIRST |
---|
| 668 | (e.g. sort(id); sorts backwards to actual monomial ordering) |
---|
[29c4ee] | 669 | - if id=list or intvec: sorted w.r.t. < (indep. of other arguments) |
---|
[3d124a7] | 670 | - if n!=0 the ordering is inverse, i.e. w.r.t. v(size(v)..1) |
---|
| 671 | default: n=0 |
---|
[b42ab6] | 672 | l[2]: intvec, describing the permutation of the input (hence l[2]=v |
---|
| 673 | if v is given (with positive integers)) |
---|
| 674 | @end format |
---|
[63be42] | 675 | NOTE: If v is given id may be any simply indexed object (e.g. any list or |
---|
| 676 | string); if v[i]<0 and i<=size(id) v[i] is set internally to i; |
---|
[f32177] | 677 | entries of v must be pairwise distinct to get a permutation id. |
---|
[3d124a7] | 678 | Zero generators of ideal/module are deleted |
---|
[d4f44f] | 679 | If 'o' is given, the input is sorted by considering leading terms |
---|
| 680 | w.r.t. the new ring ordering given by 'o' |
---|
[3d124a7] | 681 | EXAMPLE: example sort; shows an example |
---|
[d2b2a7] | 682 | " |
---|
[c860e9] | 683 | { int ii,jj,s,n = 0,0,1,0; |
---|
[3d124a7] | 684 | intvec v; |
---|
| 685 | if ( defined(basering) ) { def P = basering; } |
---|
[b5726c] | 686 | if ( size(#)==0 and (typeof(id)=="ideal" or typeof(id)=="module" |
---|
| 687 | or typeof(id)=="matrix")) |
---|
[3d124a7] | 688 | { |
---|
| 689 | id = simplify(id,2); |
---|
[558209] | 690 | for ( ii=1; ii<ncols(id); ii++ ) |
---|
[3d124a7] | 691 | { |
---|
| 692 | if ( id[ii]!=id[ii+1] ) { break;} |
---|
| 693 | } |
---|
[558209] | 694 | if ( ii != ncols(id) ) { v = sortvec(id); } |
---|
| 695 | else { v = ncols(id)..1; } |
---|
[3d124a7] | 696 | } |
---|
[b5726c] | 697 | if ( size(#)>=1 and (typeof(id)=="ideal" or typeof(id)=="module" |
---|
| 698 | or typeof(id)=="matrix") ) |
---|
[3d124a7] | 699 | { |
---|
| 700 | if ( typeof(#[1])=="string" ) |
---|
| 701 | { |
---|
[034ce1] | 702 | execute("ring r1 =("+charstr(P)+"),("+varstr(P)+"),("+#[1]+");"); |
---|
[d4f44f] | 703 | def i = imap(P,id); |
---|
[3d124a7] | 704 | v = sortvec(i); |
---|
| 705 | setring P; |
---|
| 706 | n=2; |
---|
| 707 | } |
---|
| 708 | } |
---|
[29c4ee] | 709 | if ( typeof(id)=="intvec" or typeof(id)=="list" ) |
---|
[3d124a7] | 710 | { |
---|
[29c4ee] | 711 | int Bn,Bi,Bj,Bb; |
---|
| 712 | intvec pivot; |
---|
| 713 | for(Bi=size(id);Bi>0;Bi--) { pivot[Bi]=Bi; } |
---|
| 714 | while(Bj==0) |
---|
[3d124a7] | 715 | { |
---|
[29c4ee] | 716 | Bi++; |
---|
| 717 | Bj=1; |
---|
| 718 | for(Bn=1;Bn<=size(id)-Bi;Bn++) |
---|
| 719 | { |
---|
| 720 | if(id[pivot[Bn]]>id[pivot[Bn+1]]) |
---|
| 721 | { |
---|
| 722 | Bb=pivot[Bn]; |
---|
| 723 | pivot[Bn]=pivot[Bn+1]; |
---|
| 724 | pivot[Bn+1]=Bb; |
---|
| 725 | Bj=0; |
---|
| 726 | } |
---|
| 727 | } |
---|
[3d124a7] | 728 | } |
---|
[29c4ee] | 729 | def Br=id; |
---|
| 730 | for(Bi=size(id);Bi>0;Bi--) { Br[Bi]=id[pivot[Bi]]; } |
---|
| 731 | return(list(Br,pivot)); |
---|
[3d124a7] | 732 | } |
---|
| 733 | if ( size(#)!=0 and n==0 ) { v = #[1]; } |
---|
| 734 | if( size(#)==2 ) |
---|
| 735 | { |
---|
| 736 | if ( #[2] != 0 ) { v = v[size(v)..1]; } |
---|
| 737 | } |
---|
| 738 | s = size(v); |
---|
[63be42] | 739 | if( size(id) < s ) { s = size(id); } |
---|
[3d124a7] | 740 | def m = id; |
---|
[63be42] | 741 | if ( size(m) != 0 ) |
---|
| 742 | { |
---|
[558209] | 743 | for ( jj=1; jj<=s; jj++) |
---|
[63be42] | 744 | { |
---|
| 745 | if ( v[jj]<=0 ) { v[jj]=jj; } |
---|
| 746 | m[jj] = id[v[jj]]; |
---|
| 747 | } |
---|
| 748 | } |
---|
| 749 | if ( v == 0 ) { v = 1; } |
---|
[3d124a7] | 750 | list L=m,v; |
---|
| 751 | return(L); |
---|
| 752 | } |
---|
| 753 | example |
---|
| 754 | { "EXAMPLE:"; echo = 2; |
---|
[c860e9] | 755 | ring r0 = 0,(x,y,z,t),lp; |
---|
| 756 | ideal i = x3,z3,xyz; |
---|
[584f84d] | 757 | sort(i); //sorts using lex ordering, smaller polys come first |
---|
[65546eb] | 758 | |
---|
[c860e9] | 759 | sort(i,3..1); |
---|
[b42ab6] | 760 | |
---|
[584f84d] | 761 | sort(i,"ls")[1]; //sort w.r.t. negative lex ordering |
---|
[b42ab6] | 762 | |
---|
| 763 | intvec v =1,10..5,2..4;v; |
---|
[584f84d] | 764 | sort(v)[1]; // sort v lexicographically |
---|
[b42ab6] | 765 | |
---|
[584f84d] | 766 | sort(v,"Dp",1)[1]; // sort v w.r.t (total sum, reverse lex) |
---|
[d4f44f] | 767 | |
---|
| 768 | // Note that in general: lead(sort(M)) != sort(lead(M)), e.g: |
---|
| 769 | module M = [0, 1, 1, 0], [1, 0, 0, 1]; M; |
---|
[3f4e52] | 770 | sort(lead(M), "c, dp")[1]; |
---|
| 771 | lead(sort(M, "c, dp")[1]); |
---|
[d4f44f] | 772 | |
---|
| 773 | // In order to sort M wrt a NEW ordering by considering OLD leading |
---|
| 774 | // terms use one of the following equivalent commands: |
---|
[3f4e52] | 775 | module( M[ sort(lead(M), "c,dp")[2] ] ); |
---|
| 776 | sort( M, sort(lead(M), "c,dp")[2] )[1]; |
---|
[4de1db] | 777 | |
---|
[29c4ee] | 778 | // BUG: Please, don't use this sort for integer vectors or lists |
---|
[4de1db] | 779 | // with them if there can be negative integers! |
---|
| 780 | // TODO: for some HiWi |
---|
[29c4ee] | 781 | sort(3..-3)[1]; |
---|
| 782 | sort(list(-v, v))[1]; |
---|
[4de1db] | 783 | |
---|
[3d124a7] | 784 | } |
---|
| 785 | /////////////////////////////////////////////////////////////////////////////// |
---|
[15e59a] | 786 | |
---|
[84375a] | 787 | static proc lsum (int n, list l) |
---|
[15e59a] | 788 | { if (n>10) |
---|
[4173c7] | 789 | { return( lsum(n div 2,list(l[1..(n div 2)])) + lsum(n-n div 2, list(l[(n div 2+1)..n])) ); |
---|
[15e59a] | 790 | } |
---|
| 791 | else |
---|
| 792 | { def Summe=l[1]; |
---|
| 793 | for (int i=2;i<=n;i++) |
---|
| 794 | { Summe=Summe+l[i]; |
---|
| 795 | } |
---|
| 796 | return(Summe); |
---|
| 797 | } |
---|
| 798 | } |
---|
| 799 | |
---|
| 800 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 801 | |
---|
[3d124a7] | 802 | proc sum (id, list #) |
---|
[b42ab6] | 803 | "USAGE: sum(id[,v]); id ideal/vector/module/matrix/intvec/intmat/list, |
---|
| 804 | v intvec (default: v=1..number of entries of id) |
---|
| 805 | ASSUME: list members can be added. |
---|
[65546eb] | 806 | RETURN: The sum of all entries of id [with index given by v] of type |
---|
[b42ab6] | 807 | depending on the entries of id. |
---|
[7708934] | 808 | NOTE: If id is not a list, id is treated as a list of polys resp. integers. |
---|
[b42ab6] | 809 | A module m is identified with the corresponding matrix M (columns |
---|
| 810 | of M generate m). |
---|
[7708934] | 811 | @* If v is outside the range of id, we have the empty sum and the |
---|
| 812 | result will be 0 (of type int). |
---|
[3d124a7] | 813 | EXAMPLE: example sum; shows an example |
---|
[d2b2a7] | 814 | " |
---|
[3d124a7] | 815 | { |
---|
[7708934] | 816 | //-------------------- initialization and special feature --------------------- |
---|
[b42ab6] | 817 | int n,j,tt; |
---|
[7708934] | 818 | string ty; // will become type of id |
---|
[b42ab6] | 819 | list l; |
---|
[7708934] | 820 | |
---|
| 821 | // We wish to allow something like sum(x(1..10)) if x(1),...,x(10) are |
---|
[65546eb] | 822 | // variables. x(1..10) is a list of polys and enters the procedure with |
---|
[7708934] | 823 | // id=x(1) and # a list with 9 polys, #[1]= x(2),...,#[9]= x(10). Hence, in |
---|
| 824 | // this case # is never empty. If an additional intvec v is given, |
---|
| 825 | // it is added to #, so we have to separate it first and make |
---|
| 826 | // the rest a list which has to be added. |
---|
| 827 | |
---|
[b42ab6] | 828 | int s = size(#); |
---|
| 829 | if( s!=0 ) |
---|
[7708934] | 830 | { if ( typeof(#[s])=="intvec" or typeof(#[s])=="int") |
---|
[b42ab6] | 831 | { intvec v = #[s]; |
---|
[65546eb] | 832 | tt=1; |
---|
[7708934] | 833 | s=s-1; |
---|
[b42ab6] | 834 | if ( s>0 ) { # = #[1..s]; } |
---|
| 835 | } |
---|
| 836 | } |
---|
| 837 | if ( s>0 ) |
---|
| 838 | { |
---|
[7708934] | 839 | l = list(id)+#; |
---|
| 840 | kill id; |
---|
| 841 | list id = l; //case: id = list |
---|
| 842 | ty = "list"; |
---|
[b42ab6] | 843 | } |
---|
| 844 | else |
---|
[65546eb] | 845 | { |
---|
[7708934] | 846 | ty = typeof(id); |
---|
[b42ab6] | 847 | } |
---|
[7708934] | 848 | //------------------------------ reduce to 3 cases --------------------------- |
---|
[b42ab6] | 849 | if( ty=="poly" or ty=="ideal" or ty=="vector" |
---|
| 850 | or ty=="module" or ty=="matrix" ) |
---|
[7708934] | 851 | { //case: id = ideal |
---|
[3d124a7] | 852 | ideal i = ideal(matrix(id)); |
---|
[b42ab6] | 853 | kill id; |
---|
[7708934] | 854 | ideal id = simplify(i,2); //delete 0 entries |
---|
[3d124a7] | 855 | } |
---|
[b42ab6] | 856 | if( ty=="int" or ty=="intvec" or ty=="intmat" ) |
---|
[15e59a] | 857 | { //case: id = intvec |
---|
[b42ab6] | 858 | if ( ty == "int" ) { intmat S =id; } |
---|
[3d124a7] | 859 | else { intmat S = intmat(id); } |
---|
| 860 | intvec i = S[1..nrows(S),1..ncols(S)]; |
---|
[b42ab6] | 861 | kill id; |
---|
[15e59a] | 862 | intvec id = i; |
---|
[3d124a7] | 863 | } |
---|
[7708934] | 864 | //------------------- consider intvec v and empty sum ----------------------- |
---|
[65546eb] | 865 | if( tt!=0 ) |
---|
[7708934] | 866 | { |
---|
| 867 | for (j=1; j<=size(v); j++) |
---|
| 868 | { |
---|
| 869 | if ( v[j] <= 0 or v[j] > size(id) ) //v outside range of id |
---|
[65546eb] | 870 | { |
---|
[7708934] | 871 | return(0); //empty sum is 0 |
---|
[65546eb] | 872 | } |
---|
[7708934] | 873 | } |
---|
| 874 | id = id[v]; //consider part of id |
---|
| 875 | } //corresponding to v |
---|
| 876 | |
---|
| 877 | //-------------------------- finally, add objects --------------------------- |
---|
[65546eb] | 878 | n = size(id); |
---|
[15e59a] | 879 | if (n>10) |
---|
[4173c7] | 880 | { return( lsum(n div 2,list(id[1..(n div 2)])) + lsum(n-n div 2, list(id[(n div 2+1)..n])) ); |
---|
[15e59a] | 881 | } |
---|
| 882 | else |
---|
| 883 | { def Summe=id[1]; |
---|
| 884 | for (int lauf=2;lauf<=n;lauf++) |
---|
| 885 | { Summe=Summe+id[lauf]; |
---|
| 886 | } |
---|
| 887 | return(Summe); |
---|
| 888 | } |
---|
[b42ab6] | 889 | } |
---|
[3d124a7] | 890 | example |
---|
| 891 | { "EXAMPLE:"; echo = 2; |
---|
[15e59a] | 892 | ring r1 = 0,(x,y,z),dp; |
---|
[3d124a7] | 893 | vector pv = [xy,xz,yz,x2,y2,z2]; |
---|
| 894 | sum(pv); |
---|
[c860e9] | 895 | sum(pv,2..5); |
---|
| 896 | matrix M[2][3] = 1,x,2,y,3,z; |
---|
| 897 | intvec w=2,4,6; |
---|
| 898 | sum(M,w); |
---|
| 899 | intvec iv = 1,2,3,4,5,6,7,8,9; |
---|
| 900 | sum(iv,2..4); |
---|
[15e59a] | 901 | iv = intvec(1..100); |
---|
| 902 | sum(iv); |
---|
| 903 | ring r2 = 0,(x(1..10)),dp; |
---|
| 904 | sum(x(3..7),intvec(1,3,5)); |
---|
[3d124a7] | 905 | } |
---|
| 906 | /////////////////////////////////////////////////////////////////////////////// |
---|
[6f2edc] | 907 | |
---|
[15e59a] | 908 | |
---|
| 909 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 910 | |
---|
[ebbe4a] | 911 | proc watchdog(int i, string cmd) |
---|
[a7a00b] | 912 | "USAGE: watchdog(i,cmd); i integer, cmd string |
---|
[b42ab6] | 913 | RETURN: Result of cmd, if the result can be computed in i seconds. |
---|
| 914 | Otherwise the computation is interrupted after i seconds, |
---|
| 915 | the string "Killed" is returned and the global variable |
---|
| 916 | 'watchdog_interrupt' is defined. |
---|
[65546eb] | 917 | NOTE: * the MP package must be enabled |
---|
| 918 | * the current basering should not be watchdog_rneu, since |
---|
[b42ab6] | 919 | watchdog_rneu will be killed |
---|
[ebbe4a] | 920 | * if there are variable names of the structure x(i) all |
---|
| 921 | polynomials have to be put into eval(...) in order to be |
---|
| 922 | interpreted correctly |
---|
[a6606e6] | 923 | * a second Singular process is started by this procedure |
---|
[ebbe4a] | 924 | EXAMPLE: example watchdog; shows an example |
---|
| 925 | " |
---|
| 926 | { |
---|
| 927 | string rname=nameof(basering); |
---|
[17ee4f] | 928 | def rsave=basering; |
---|
[ebbe4a] | 929 | if (defined(watchdog_rneu)) |
---|
| 930 | { |
---|
| 931 | kill watchdog_rneu; |
---|
| 932 | } |
---|
[1b2216] | 933 | if ( i > 0 ) |
---|
[ebbe4a] | 934 | { |
---|
[1b2216] | 935 | int j=10; |
---|
| 936 | int k=999999; |
---|
[65546eb] | 937 | // fork, get the pid of the child and send it the command |
---|
[1b2216] | 938 | link l_fork="ssi:fork"; |
---|
| 939 | open(l_fork); |
---|
| 940 | execute("write(l_fork,quote(" + cmd + "));"); |
---|
[ebbe4a] | 941 | |
---|
| 942 | |
---|
| 943 | // sleep in small, but growing intervals for appr. 1 second |
---|
[1b2216] | 944 | while(j < k) |
---|
| 945 | { |
---|
| 946 | if (status(l_fork, "read", "ready", j)) {break;} |
---|
| 947 | j = j + j; |
---|
| 948 | } |
---|
[ebbe4a] | 949 | |
---|
| 950 | // sleep in intervals of one second |
---|
[1b2216] | 951 | j = 1; |
---|
| 952 | if (!status(l_fork,"read","ready")) |
---|
| 953 | { |
---|
| 954 | while (j < i) |
---|
[ebbe4a] | 955 | { |
---|
[1b2216] | 956 | if (status(l_fork, "read", "ready", k)) {break;} |
---|
| 957 | j = j + 1; |
---|
[ebbe4a] | 958 | } |
---|
[1b2216] | 959 | } |
---|
[ebbe4a] | 960 | // check, whether we have a result, and return it |
---|
[1b2216] | 961 | if (status(l_fork, "read", "ready")) |
---|
| 962 | { |
---|
| 963 | def result = read(l_fork); |
---|
| 964 | if (nameof(basering)!=rname) |
---|
[ebbe4a] | 965 | { |
---|
[1b2216] | 966 | def watchdog_rneu=basering; |
---|
| 967 | setring rsave; |
---|
| 968 | if (!defined(result)) |
---|
[ebbe4a] | 969 | { |
---|
[1b2216] | 970 | def result=fetch(watchdog_rneu,result); |
---|
[ebbe4a] | 971 | } |
---|
| 972 | } |
---|
[1b2216] | 973 | if(defined(watchdog_interrupt)) |
---|
[ebbe4a] | 974 | { |
---|
[1b2216] | 975 | kill watchdog_interrupt; |
---|
[ebbe4a] | 976 | } |
---|
[1b2216] | 977 | close(l_fork); |
---|
[ebbe4a] | 978 | } |
---|
| 979 | else |
---|
| 980 | { |
---|
[1b2216] | 981 | string result="Killed"; |
---|
| 982 | if(!defined(watchdog_interrupt)) |
---|
| 983 | { |
---|
| 984 | int watchdog_interrupt=1; |
---|
| 985 | export watchdog_interrupt; |
---|
| 986 | } |
---|
| 987 | close(l_fork); |
---|
[ebbe4a] | 988 | } |
---|
[1b2216] | 989 | return(result); |
---|
[50cbdc] | 990 | } |
---|
| 991 | else |
---|
| 992 | { |
---|
[1b2216] | 993 | ERROR("First argument of watchdog has to be a positive integer."); |
---|
[65546eb] | 994 | } |
---|
[ebbe4a] | 995 | } |
---|
| 996 | example |
---|
| 997 | { "EXAMPLE:"; echo=2; |
---|
| 998 | ring r=0,(x,y,z),dp; |
---|
| 999 | poly f=x^30+y^30; |
---|
| 1000 | watchdog(1,"factorize(eval("+string(f)+"))"); |
---|
| 1001 | watchdog(100,"factorize(eval("+string(f)+"))"); |
---|
| 1002 | } |
---|
| 1003 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 1004 | |
---|
| 1005 | proc deleteSublist(intvec v,list l) |
---|
[803c5a1] | 1006 | "USAGE: deleteSublist(v,l); intvec v; list l |
---|
[ebbe4a] | 1007 | where the entries of the integer vector v correspond to the |
---|
| 1008 | positions of the elements to be deleted |
---|
| 1009 | RETURN: list without the deleted elements |
---|
| 1010 | EXAMPLE: example deleteSublist; shows an example" |
---|
| 1011 | { |
---|
| 1012 | list k; |
---|
| 1013 | int i,j,skip; |
---|
| 1014 | j=1; |
---|
| 1015 | skip=0; |
---|
| 1016 | intvec vs=sort(v)[1]; |
---|
| 1017 | for ( i=1 ; i <=size(vs) ; i++) |
---|
| 1018 | { |
---|
| 1019 | while ((j+skip) < vs[i]) |
---|
| 1020 | { |
---|
| 1021 | k[j] = l[j+skip]; |
---|
| 1022 | j++; |
---|
| 1023 | } |
---|
| 1024 | skip++; |
---|
| 1025 | } |
---|
| 1026 | if(vs[size(vs)]<size(l)) |
---|
| 1027 | { |
---|
| 1028 | k=k+list(l[(vs[size(vs)]+1)..size(l)]); |
---|
| 1029 | } |
---|
| 1030 | return(k); |
---|
| 1031 | } |
---|
| 1032 | example |
---|
| 1033 | { "EXAMPLE:"; echo=2; |
---|
| 1034 | list l=1,2,3,4,5; |
---|
| 1035 | intvec v=1,3,4; |
---|
| 1036 | l=deleteSublist(v,l); |
---|
| 1037 | l; |
---|
| 1038 | } |
---|
[8b87364] | 1039 | |
---|
| 1040 | /////////////////////////////////////////////////////////////////////////////// |
---|
| 1041 | proc primecoeffs(J, list #) |
---|
[a7a00b] | 1042 | "USAGE: primecoeffs(J[,p]); J any type which can be converted to a matrix |
---|
[8b87364] | 1043 | e.g. ideal, matrix, vector, module, int, intvec |
---|
[a7a00b] | 1044 | p = integer |
---|
[c94e60f] | 1045 | COMPUTE: primefactors <= p of coeffs of J (default p = 32003) |
---|
[a7a00b] | 1046 | RETURN: a list, say l, of two intvectors:@* |
---|
| 1047 | l[1] : the different primefactors of all coefficients of J@* |
---|
[8b87364] | 1048 | l[2] : the different remaining factors |
---|
| 1049 | EXAMPLE: example primecoeffs; shows an example |
---|
| 1050 | " |
---|
| 1051 | { |
---|
| 1052 | int q,ii,n,mark;; |
---|
[298d0a] | 1053 | if (size(#) == 0) |
---|
[8b87364] | 1054 | { |
---|
[298d0a] | 1055 | q=32003; |
---|
[8b87364] | 1056 | } |
---|
[298d0a] | 1057 | else |
---|
[8b87364] | 1058 | { |
---|
| 1059 | if( typeof(#[1]) != "int") |
---|
| 1060 | { |
---|
| 1061 | ERROR("2nd parameter must be of type int"+newline); |
---|
| 1062 | } |
---|
| 1063 | q=#[1]; |
---|
[be8dbb] | 1064 | if (q > 32003) { q = 32003; } |
---|
[8b87364] | 1065 | } |
---|
| 1066 | |
---|
| 1067 | if (defined(basering) == 0) |
---|
| 1068 | { |
---|
| 1069 | mark=1; |
---|
| 1070 | ring r = 0,x,dp; |
---|
| 1071 | } |
---|
| 1072 | def I = ideal(matrix(J)); |
---|
| 1073 | poly p = product(maxideal(1)); |
---|
[298d0a] | 1074 | matrix Coef=coef(I[1],p); |
---|
[8b87364] | 1075 | ideal id, jd, rest; |
---|
| 1076 | intvec v,re; |
---|
| 1077 | list result,l; |
---|
| 1078 | for(ii=2; ii<=ncols(I); ii++) |
---|
| 1079 | { |
---|
| 1080 | Coef=concat(Coef,coef(I[ii],p)); |
---|
| 1081 | } |
---|
| 1082 | id = Coef[2,1..ncols(Coef)]; |
---|
| 1083 | id = simplify(id,6); |
---|
[298d0a] | 1084 | for (ii=1; ii<=size(id); ii++) |
---|
| 1085 | { |
---|
| 1086 | l = primefactors(number(id[ii]),q); |
---|
[8acd267] | 1087 | list jdl=l[1]; |
---|
| 1088 | jd = jd,jdl[1..size(jdl)]; |
---|
[c94e60f] | 1089 | kill jdl; |
---|
[8acd267] | 1090 | rest=rest, l[3]; |
---|
[298d0a] | 1091 | } |
---|
[8b87364] | 1092 | jd = simplify(jd,6); |
---|
[298d0a] | 1093 | for (ii=1; ii<=size(jd); ii++) |
---|
| 1094 | { |
---|
[8b87364] | 1095 | v[ii]=int(jd[ii]); |
---|
| 1096 | } |
---|
| 1097 | v = sort(v)[1]; |
---|
| 1098 | rest = simplify(rest,6); |
---|
| 1099 | id = sort(id)[1]; |
---|
| 1100 | if (mark) |
---|
| 1101 | { |
---|
| 1102 | for (ii=1; ii<=size(rest); ii++) |
---|
| 1103 | { |
---|
| 1104 | re[ii] = int(rest[ii]); |
---|
| 1105 | } |
---|
| 1106 | result = v,re; |
---|
| 1107 | } |
---|
| 1108 | else |
---|
| 1109 | { |
---|
[298d0a] | 1110 | result = v,rest; |
---|
[8b87364] | 1111 | } |
---|
| 1112 | return(result); |
---|
| 1113 | } |
---|
| 1114 | example |
---|
| 1115 | { "EXAMPLE:"; echo = 2; |
---|
| 1116 | primecoeffs(intvec(7*8*121,7*8));""; |
---|
| 1117 | ring r = 0,(b,c,t),dp; |
---|
| 1118 | ideal I = -13b6c3t+4b5c4t,-10b4c2t-5b4ct2; |
---|
| 1119 | primecoeffs(I); |
---|
[298d0a] | 1120 | } |
---|
[8b87364] | 1121 | /////////////////////////////////////////////////////////////////////////////// |
---|
[90d772] | 1122 | proc timeFactorize(poly i,list #) |
---|
[a7a00b] | 1123 | "USAGE: timeFactorize(p,d); poly p , integer d |
---|
[3c4dcc] | 1124 | RETURN: factorize(p) if the factorization finished after d-1 |
---|
[90d772] | 1125 | seconds otherwhise f is considered to be irreducible |
---|
| 1126 | EXAMPLE: example timeFactorize; shows an example |
---|
| 1127 | " |
---|
| 1128 | { |
---|
| 1129 | def P=basering; |
---|
| 1130 | if (size(#) > 0) |
---|
| 1131 | { |
---|
[1b2216] | 1132 | if ((typeof(#[1]) == "int")&&(#[1])) |
---|
[90d772] | 1133 | { |
---|
[1b2216] | 1134 | int wait = #[1]; |
---|
| 1135 | int j = 10; |
---|
[90d772] | 1136 | |
---|
[1b2216] | 1137 | string bs = nameof(basering); |
---|
| 1138 | link l_fork = "ssi:fork"; |
---|
| 1139 | open(l_fork); |
---|
| 1140 | write(l_fork, quote(timeFactorize(eval(i)))); |
---|
[90d772] | 1141 | |
---|
[1b2216] | 1142 | // sleep in small intervalls for appr. one second |
---|
| 1143 | if (wait > 0) |
---|
| 1144 | { |
---|
| 1145 | while(j < 1000000) |
---|
[90d772] | 1146 | { |
---|
[1b2216] | 1147 | if (status(l_fork, "read", "ready", j)) {break;} |
---|
| 1148 | j = j + j; |
---|
[90d772] | 1149 | } |
---|
[1b2216] | 1150 | } |
---|
[90d772] | 1151 | |
---|
[1b2216] | 1152 | // sleep in intervalls of one second from now on |
---|
| 1153 | j = 1; |
---|
| 1154 | while (j < wait) |
---|
| 1155 | { |
---|
| 1156 | if (status(l_fork, "read", "ready", 1000000)) {break;} |
---|
| 1157 | j = j + 1; |
---|
| 1158 | } |
---|
[90d772] | 1159 | |
---|
[1b2216] | 1160 | if (status(l_fork, "read", "ready")) |
---|
| 1161 | { |
---|
| 1162 | def result = read(l_fork); |
---|
| 1163 | if (bs != nameof(basering)) |
---|
[90d772] | 1164 | { |
---|
[1b2216] | 1165 | def PP = basering; |
---|
| 1166 | setring P; |
---|
| 1167 | def result = imap(PP, result); |
---|
| 1168 | kill PP; |
---|
[90d772] | 1169 | } |
---|
[1b2216] | 1170 | close(l_fork); |
---|
| 1171 | } |
---|
| 1172 | else |
---|
| 1173 | { |
---|
| 1174 | list result; |
---|
| 1175 | intvec v=1,1; |
---|
| 1176 | result[1]=list(1,i); |
---|
| 1177 | result[2]=v; |
---|
| 1178 | close(l_fork); |
---|
[90d772] | 1179 | } |
---|
[1b2216] | 1180 | return (result); |
---|
[90d772] | 1181 | } |
---|
| 1182 | } |
---|
| 1183 | return(factorH(i)); |
---|
| 1184 | } |
---|
| 1185 | example |
---|
| 1186 | { "EXAMPLE:"; echo = 2; |
---|
| 1187 | ring r=0,(x,y),dp; |
---|
| 1188 | poly p=((x2+y3)^2+xy6)*((x3+y2)^2+x10y); |
---|
| 1189 | p=p^2; |
---|
| 1190 | //timeFactorize(p,2); |
---|
| 1191 | //timeFactorize(p,20); |
---|
| 1192 | } |
---|
| 1193 | |
---|
| 1194 | proc timeStd(ideal i,list #) |
---|
| 1195 | "USAGE: timeStd(i,d), i ideal, d integer |
---|
| 1196 | RETURN: std(i) if the standard basis computation finished after |
---|
| 1197 | d-1 seconds and i otherwhise |
---|
| 1198 | EXAMPLE: example timeStd; shows an example |
---|
| 1199 | " |
---|
| 1200 | { |
---|
| 1201 | def P=basering; |
---|
| 1202 | if (size(#) > 0) |
---|
| 1203 | { |
---|
[1b2216] | 1204 | if ((typeof(#[1]) == "int")&&(#[1])) |
---|
[90d772] | 1205 | { |
---|
[1b2216] | 1206 | int wait = #[1]; |
---|
| 1207 | int j = 10; |
---|
[90d772] | 1208 | |
---|
[1b2216] | 1209 | string bs = nameof(basering); |
---|
| 1210 | link l_fork = "ssi:fork"; |
---|
| 1211 | open(l_fork); |
---|
| 1212 | write(l_fork, quote(system("pid"))); |
---|
| 1213 | int pid = read(l_fork); |
---|
| 1214 | write(l_fork, quote(timeStd(eval(i)))); |
---|
[90d772] | 1215 | |
---|
[1b2216] | 1216 | // sleep in small intervalls for appr. one second |
---|
| 1217 | if (wait > 0) |
---|
| 1218 | { |
---|
| 1219 | while(j < 1000000) |
---|
[90d772] | 1220 | { |
---|
[1b2216] | 1221 | if (status(l_fork, "read", "ready", j)) {break;} |
---|
| 1222 | j = j + j; |
---|
[90d772] | 1223 | } |
---|
[1b2216] | 1224 | } |
---|
| 1225 | j = 1; |
---|
| 1226 | while (j < wait) |
---|
| 1227 | { |
---|
| 1228 | if (status(l_fork, "read", "ready", 1000000)) {break;} |
---|
| 1229 | j = j + 1; |
---|
| 1230 | } |
---|
| 1231 | if (status(l_fork, "read", "ready")) |
---|
| 1232 | { |
---|
| 1233 | def result = read(l_fork); |
---|
| 1234 | if (bs != nameof(basering)) |
---|
[90d772] | 1235 | { |
---|
[1b2216] | 1236 | def PP = basering; |
---|
| 1237 | setring P; |
---|
| 1238 | def result = imap(PP, result); |
---|
| 1239 | kill PP; |
---|
[90d772] | 1240 | } |
---|
[1b2216] | 1241 | close(l_fork); |
---|
| 1242 | } |
---|
| 1243 | else |
---|
| 1244 | { |
---|
| 1245 | ideal result=i; |
---|
| 1246 | close(l_fork); |
---|
[90d772] | 1247 | } |
---|
[1b2216] | 1248 | return (result); |
---|
[90d772] | 1249 | } |
---|
| 1250 | } |
---|
| 1251 | return(std(i)); |
---|
| 1252 | } |
---|
| 1253 | example |
---|
| 1254 | { "EXAMPLE:"; echo = 2; |
---|
| 1255 | ring r=32003,(a,b,c,d,e),dp; |
---|
[cb8103a] | 1256 | int n=7; |
---|
[90d772] | 1257 | ideal i= |
---|
| 1258 | a^n-b^n, |
---|
| 1259 | b^n-c^n, |
---|
| 1260 | c^n-d^n, |
---|
| 1261 | d^n-e^n, |
---|
| 1262 | a^(n-1)*b+b^(n-1)*c+c^(n-1)*d+d^(n-1)*e+e^(n-1)*a; |
---|
[cb8103a] | 1263 | def i1=timeStd(i,1); |
---|
| 1264 | def i2=timeStd(i,20); |
---|
| 1265 | listvar(); |
---|
[90d772] | 1266 | } |
---|
| 1267 | |
---|
| 1268 | proc factorH(poly p) |
---|
| 1269 | "USAGE: factorH(p) p poly |
---|
| 1270 | RETURN: factorize(p) |
---|
[a7a00b] | 1271 | NOTE: changes variables to make the last variable the principal |
---|
[90d772] | 1272 | one in the multivariate factorization and factorizes then |
---|
| 1273 | the polynomial |
---|
| 1274 | EXAMPLE: example factorH; shows an example |
---|
| 1275 | " |
---|
| 1276 | { |
---|
| 1277 | def R=basering; |
---|
| 1278 | int i,j; |
---|
| 1279 | int n=1; |
---|
| 1280 | int d=nrows(coeffs(p,var(1))); |
---|
| 1281 | for(i=1;i<=nvars(R);i++) |
---|
| 1282 | { |
---|
| 1283 | j=nrows(coeffs(p,var(i))); |
---|
| 1284 | if(d>j) |
---|
| 1285 | { |
---|
| 1286 | n=i; |
---|
| 1287 | d=j; |
---|
| 1288 | } |
---|
| 1289 | } |
---|
| 1290 | ideal ma=maxideal(1); //die letzte Variable ist die Hauptvariable |
---|
| 1291 | ma[nvars(R)]=var(n); |
---|
| 1292 | ma[n]=var(nvars(R)); |
---|
| 1293 | map phi=R,ma; |
---|
| 1294 | list fac=factorize(phi(p)); |
---|
| 1295 | list re=phi(fac); |
---|
| 1296 | return(re); |
---|
| 1297 | } |
---|
| 1298 | example |
---|
| 1299 | { "EXAMPLE:"; echo = 2; |
---|
| 1300 | system("random",992851144); |
---|
| 1301 | ring r=32003,(x,y,z,w,t),lp; |
---|
| 1302 | poly p=y2w9+yz7t-yz5w4-z2w4t4-w8t3; |
---|
| 1303 | factorize(p); //fast |
---|
| 1304 | system("random",992851262); |
---|
| 1305 | //factorize(p); //slow |
---|
| 1306 | system("random",992851262); |
---|
| 1307 | factorH(p); |
---|
| 1308 | } |
---|