Changeset d5f47f in git
 Timestamp:
 May 5, 2006, 4:36:12 PM (17 years ago)
 Branches:
 (u'spielwiese', '8d54773d6c9e2f1d2593a28bc68b7eeab54ed529')
 Children:
 999b3a5588fa00975695577b7240befd77a2ca12
 Parents:
 933165be31db0f88135bcfa7c8151c9624bd34c5
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

Singular/LIB/center.lib
r933165b rd5f47f 1 1 /////////////////////////////////////////////////////////////////////////////// 2 version="$Id: center.lib,v 1.1 8 20060403 13:16:53motsak Exp $"2 version="$Id: center.lib,v 1.19 20060505 14:36:12 motsak Exp $" 3 3 category="Noncommutative" 4 4 info=" 5 5 LIBRARY: center.lib computation of central elements of GRalgebras 6 6 AUTHOR: Oleksandr Motsak, motsak@mathematik.unikl.de. 7 OVERVIEW: A library for computing elements of the center and centralizers of elements in various noncommutative algebras. 8 9 SUPPORT: Forschungsschwerpunkt 'Mathematik und Praxis', University of Kaiserslautern 7 OVERVIEW: Computing elements of center and centralizers of sets of elements 8 in noncommutative algebras. 10 9 11 10 MAIN PROCEDURES: 12 @format 13 centralizeSet(F, V): computes a V.S. basis of the centralizer of set F within V, 14 centralizerVS(F, d): computes a V.S. basis of the centralizer of set F, 15 centralizerRed(F, D[, N]): computes reduced elements of the centralizer of set F, 16 centerVS(D): computes a V.S. basis of the center up to degree D, 17 centerRed(D[, k]): computes reduced elements of the center up to degree D, 18 19 center(D[, k]): computes reduced elements of the center, 20 centralizer(F, D[, k]): computes reduced elements of the centralizer of set F, 21 22 sa_reduce(V): 'subalgebra reduction' of a set of pairwise commuting polynomials V, 23 sa_poly_reduce(p, V): 'subalgebra reduction' of a polynomial p wrt a set of pairwise commuting polynomials V. 24 25 inCenter(T): checks the centrality of polynomials of a list/ideal/poly T, 26 inCentralizer(T, S): checks whether polynomials of list/ideal/poly T commute with polynomials of S, 27 isCartan(p): checks whether polynomial p is a Cartan element, 28 @end format 11 centralizeSet(F, V): v.s. basis of the centralizer of F within V 12 centralizerVS(F, D): v.s. basis of the centralizer of F 13 centralizerRed(F, D[, N]): reduced basis of the centralizer of F 14 centerVS(D): v.s. basis of the center 15 centerRed(D[, k]): reduced basis of the center 16 17 center(D[, k]): reduced basis of the center 18 centralizer(F, D[, k]): reduced bais of the centralizer of F 19 20 sa_reduce(V): 's.a. reduction' of pairwise commuting elements 21 sa_poly_reduce(p, V): 's.a. reduction' of p by pairwise commuting elements 22 23 inCenter(T): checks the centrality of list/ideal/poly T 24 inCentralizer(T, S): checks whether list/ideal/poly T commute with S 25 isCartan(p): checks whether polynomial p is a Cartan element 29 26 30 27 AUXILIARY PROCEDURES: 31 @format 32 applyAdF(Basis, f): computes images of basis elements under the linear map Ad_{f}, 33 linearMapKernel(Images): computes the kernel of a linear map given by image of basis, 34 linearCombinations(Basis, C): computes linear combinations of Basis vectors with the coefficients from C, 35 36 variablesStandard(): computes the set of algebra generators in their natural order, 37 variablesSorted(): computes the sorted set of algebra generators, 38 39 PBW_eqDeg(Deg): computes PBW monomials of a given degree Deg, 40 PBW_maxDeg(MaxDeg): computes PBW monomials up to a given degree MaxDeg, 41 PBW_maxMonom(MaxMonom): computes PBW monomials up to a given maximal monomial MaxMonom. 42 @end format 28 applyAdF(Basis, f): images of elements under the klinear map Ad_{f} 29 linearMapKernel(Images): kernel of a linear map given by images 30 linearCombinations(Basis, C): klinear combinations of elements 31 32 variablesStandard(): set of algebra generators in their natural order 33 variablesSorted(): heuristically sorted set of algebra generators 34 35 PBW_eqDeg(Deg): PBW monomials of given degree 36 PBW_maxDeg(MaxDeg): PBW monomials up to given degree 37 PBW_maxMonom(MaxMonom): PBW monomials up to given maximal monomial 43 38 44 39 KEYWORDS: center; centralizer; cartan; reduce; centralize; PBW … … 56 51 /******************************************************/ 57 52 static proc DefaultValue ( def s, list # ) // Process general variable parameters list 58 "53 " 59 54 RETURN: s or (typeof(s))(#) 60 55 " 61 56 { 62 63 64 65 66 67 68 69 70 57 def @p = s; 58 if ( size(#) > 0 ) 59 { 60 if ( typeof(#[1]) == typeof(s) ) 61 { 62 @p = #[1]; 63 } 64 } 65 return( @p ); 71 66 } 72 67 73 68 /******************************************************/ 74 69 static proc DefaultInt( list # ) // Process variable parameters list with 'int' default value 75 "70 " 76 71 RETURN: 0 or int(#) 77 72 " 78 73 { 79 80 74 int @p = 0; 75 return( DefaultValue( @p, # ) ); 81 76 } 82 77 83 78 /******************************************************/ 84 79 static proc DefaultIdeal ( list # ) // Process variable parameters list with 'ideal' default value 85 "80 " 86 81 RETURN: 0 or ideal(#) 87 82 " 88 83 { 89 90 84 ideal @p = 0; 85 return( DefaultValue( @p, # ) ); 91 86 } 92 87 … … 100 95 /******************************************************/ 101 96 static proc toprint( int pl ) // To print or not to print?!? The answer is here! 102 "97 " 103 98 RETURN: 1 means to print, otherwise 0. 104 99 " 105 100 { 106 101 return( printlevel >= ( 3  pl) ); // voice + ? 107 102 } 108 103 109 104 /******************************************************/ 110 105 static proc DBPrint( int pl, list # ) // My 'dbprint' which uses toprint(i). 111 "106 " 112 107 USAGE: 113 108 " 114 109 { 115 116 { 117 110 if( toprint(pl) ) 111 { 112 dbprint(1, #); 118 113 } 119 114 } … … 121 116 /******************************************************/ 122 117 static proc BCall( string Name, list # ) // This function must be called at the beginning of every 'mathematical' function. 123 "118 " 124 119 USAGE: Name is a name of a mathematical function to trace. # means parameters into the function. 125 120 " 126 121 { 127 128 { 129 130 131 122 if( toprint(0) ) 123 { 124 "CALL: ", Name, "( "; 125 dbprint(1, #); 126 " )"; 132 127 } 133 128 } … … 135 130 /******************************************************/ 136 131 static proc ECall(string Name, list #) // This function must be called at the end of every 'mathematical' function. 137 "132 " 138 133 USAGE: Name is a name of a mathematical function to trace. # means result of the function. 139 134 " 140 135 { 141 142 { 143 144 145 136 if( toprint(0) ) 137 { 138 "RET : ", Name, " => Result: {"; 139 dbprint(1, #); 140 " }"; 146 141 } 147 142 } … … 155 150 /******************************************************/ 156 151 static proc makeNice( def l ) 157 "152 " 158 153 RETURN: the same as input 159 154 PURPOSE: make 'nice' polynomials, kill scalars 160 155 " 161 156 { 162 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "makeNice", l ); }; /*4DEBUG*/ 163 164 poly p; 165 166 if( typeof(l) == "poly" ) 167 { 168 // "normal" polynomial form == no denominators, gcd of coeffs is a unit 169 p = cleardenom( l ); 170 l = cleardenom( p / leadcoef(p) ); 157 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "makeNice", l ); }; /*4DEBUG*/ 158 159 poly p; 160 161 if( typeof(l) == "poly" ) 162 { 163 // "normal" polynomial form == no denominators, gcd of coeffs is a unit 164 l = cleardenom( l ); 165 if ( deg(p) > 0 ) 166 { 167 l = cleardenom( l / leadcoef(l) ); 168 } 171 169 } else 172 {170 { 173 171 if( typeof(l) == "ideal" ) 174 172 { 175 173 for( int i = 1; i <= size(l); i++ ) 176 174 { 177 175 p = l[i]; 176 p = cleardenom( p ); 178 177 179 178 // Now make polynomials look nice: 180 179 if ( deg(p) > 0 ) // throw away scalars! 181 { 182 // "normal" polynomial form == no denominators, gcd of coeffs is a unit 183 p = cleardenom( p ); 184 l[i] = cleardenom( p / leadcoef(p) ); 185 } else 186 { 187 l[i] = 0; 188 } 189 } 180 { 181 // "normal" polynomial form == no denominators, gcd of coeffs is a unit 182 p = cleardenom( p / leadcoef(p) ); 183 } else 184 { 185 p = 0; // no constants 186 } 187 l[i] = p; 188 189 } 190 190 191 191 l = simplify(l, 2 + 8); 192 193 }194 195 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "makeNice", l ); }; /*4DEBUG*/196 192 } 193 } 194 195 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "makeNice", l ); }; /*4DEBUG*/ 196 return( l ); 197 197 } 198 198 … … 201 201 /******************************************************/ 202 202 static proc monomialForm( def p ) 203 "203 " 204 204 : p is either poly or ideal 205 205 RETURN: polynomial with all monomials from p but without coefficients. 206 206 " 207 207 { 208 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "monomialForm", p ); }; /*4DEBUG*/209 210 211 212 213 { 214 215 { 216 217 218 208 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "monomialForm", p ); }; /*4DEBUG*/ 209 210 poly result = 0; int k, j; poly m; 211 212 if( typeof(p) == "ideal" ) // 213 { 214 if( ncols(p) > 0 ) 215 { 216 result = uni_poly( p[1] ); 217 218 for ( k = ncols(p); k > 1; k  ) 219 219 { 220 220 for( j = size(p[k]); j > 0; j ) 221 221 { 222 222 m = leadmonom( p[k][j] ); 223 223 224 224 if( size(result + m) > size(result) ) // trick! 225 225 { 226 226 result = result + m; 227 227 } 228 228 } … … 231 231 } 232 232 } 233 234 { 235 236 { 237 233 else 234 { 235 if( typeof(p) == "poly" ) // 236 { 237 result = uni_poly(p); 238 238 } else 239 239 { 240 240 ERROR( "ERROR: Wrong input! Expected polynomial or ideal!" ); 241 242 } 243 244 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "monomialForm", result ); }; /*4DEBUG*/245 241 } 242 } 243 244 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "monomialForm", result ); }; /*4DEBUG*/ 245 return( result ); 246 246 } 247 247 248 248 /******************************************************/ 249 249 static proc uni_poly( poly p ) 250 "250 " 251 251 returns polynomial with the same monomials but without coefficients. 252 252 " 253 253 { 254 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "uni_poly", p ); }; /*4DEBUG*/255 256 257 258 259 { 260 261 } 262 263 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "uni_poly", result ); }; /*4DEBUG*/264 254 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "uni_poly", p ); }; /*4DEBUG*/ 255 256 poly result = 0; 257 258 for ( int k = size(p); k > 0; k  ) 259 { 260 result = result + leadmonom( p[k] ); 261 } 262 263 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "uni_poly", result ); }; /*4DEBUG*/ 264 return( result ); 265 265 } 266 266 … … 271 271 /******************************************************/ 272 272 static proc smoothQideal( ideal I, list #) 273 "273 " 274 274 PURPOSE: smooths the ideal in a current QUOTIENT(!) ring. 275 275 " 276 276 { 277 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "smoothQideal", I ); }; /*4DEBUG*/278 279 ideal A = NF( I, twostd(DefaultIdeal(#)), 1 );280 281 ideal D = I  A;282 283 if( size(D) > 0 ) // Were there any changes???284 { 285 ideal T = ideal();int j = 1;286 287 288 { 289 if( size(D[i]) == 0 ) // keep only unchanged elements277 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "smoothQideal", I ); }; /*4DEBUG*/ 278 279 ideal A = I  NF( I, twostd(DefaultIdeal(#)), 1 ); // component wise 280 281 if( size(A) > 0 ) // Were there any changes (any nonzero component)? 282 { 283 ideal T; 284 285 int j = 1; 286 287 for( int i = 1; i <= ncols(I); i++ ) 288 { 289 if( size(A[i]) == 0 ) // keep only unchanged elements 290 290 { 291 291 T[ j ] = I[i]; j++; 292 292 } 293 293 } 294 295 kill A; 296 297 ideal A = T; 298 299 kill T; 300 } 301 302 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "smoothQideal", A ); }; /*4DEBUG*/ 303 304 return( A ); 294 I = T; 295 } 296 297 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "smoothQideal", I ); }; /*4DEBUG*/ 298 299 return( I ); 305 300 } 306 301 … … 318 313 proc PBW_maxDeg( int MaxDeg ) 319 314 "USAGE: PBW_maxDeg(MaxDeg); int MaxDeg 320 PURPOSE: Compute the PBW basis (up to a given maximal degree) of a current algebra.321 RETURN: ideal consisting of PBWelements.315 PURPOSE: Compute PBW elements up to a given maximal degree. 316 RETURN: ideal consisting of found elements. 322 317 NOTE: unit is omitted. Weights are ignored! 323 318 " 324 319 { 325 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "PBW_maxDeg", MaxDeg ); }; /*4DEBUG*/326 327 ideal Basis = ideal();328 329 330 { 331 ideal T = Basis + maxideal(k); kill Basis; ideal Basis = T; kill T; // ?332 } 333 334 ideal T = smoothQideal( Basis ); kill Basis;335 336 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "PBW_maxDeg", T); }; /*4DEBUG*/337 return( T);320 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "PBW_maxDeg", MaxDeg ); }; /*4DEBUG*/ 321 322 ideal Basis = 0; 323 324 for (int k = 1; k <= MaxDeg; k ++ ) 325 { 326 Basis = Basis + maxideal(k); 327 } 328 329 Basis = smoothQideal( Basis ); 330 331 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "PBW_maxDeg", Basis ); }; /*4DEBUG*/ 332 return( Basis ); 338 333 } 339 334 example 340 335 { 341 "EXAMPLE:"; echo = 2;342 ring A = 0,(e,f,h),dp;343 matrix D[3][3]=0;344 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f;345 ncalgebra(1,D); // this algebra is U(sl_2)346 347 // PBW Basis of A_2  monomials of degree <= 2, without unit:348 PBW_maxDeg( 2 );336 "EXAMPLE:"; echo = 2; 337 ring A = 0,(e,f,h),dp; 338 matrix D[3][3]=0; 339 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f; 340 ncalgebra(1,D); // this algebra is U(sl_2) 341 342 // PBW Basis of A_2  monomials of degree <= 2, without unit: 343 PBW_maxDeg( 2 ); 349 344 } 350 345 … … 353 348 proc PBW_eqDeg( int Deg ) 354 349 "USAGE: PBW_eqDeg(Deg); int Deg 355 PURPOSE: Compute the PBW basis (of a given degree) of a current algebra.356 RETURN: ideal consisting of PBWelements.350 PURPOSE: Compute PBW elements of a given degree. 351 RETURN: ideal consisting of found elements. 357 352 NOTE: Unit is omitted. Weights are ignored! 358 353 " 359 354 { 360 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "PBW_eqDeg", Deg ); }; /*4DEBUG*/361 362 363 364 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "PBW_eqDeg", Basis ); }; /*4DEBUG*/365 355 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "PBW_eqDeg", Deg ); }; /*4DEBUG*/ 356 357 ideal Basis = smoothQideal( maxideal( Deg ) ); 358 359 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "PBW_eqDeg", Basis ); }; /*4DEBUG*/ 360 return( Basis ); 366 361 } 367 362 example 368 363 { 369 "EXAMPLE:"; echo = 2;370 ring A = 0,(e,f,h),dp;371 matrix D[3][3]=0;372 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f;373 ncalgebra(1,D); // this algebra is U(sl_2)374 375 // PBW Basis of A_2 \ A_1  monomials of degree == 2:376 PBW_eqDeg( 2 );364 "EXAMPLE:"; echo = 2; 365 ring A = 0,(e,f,h),dp; 366 matrix D[3][3]=0; 367 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f; 368 ncalgebra(1,D); // this algebra is U(sl_2) 369 370 // PBW Basis of A_2 \ A_1  monomials of degree == 2: 371 PBW_eqDeg( 2 ); 377 372 } 378 373 … … 381 376 proc PBW_maxMonom( poly MaxMonom ) 382 377 "USAGE: PBW_maxMonom(m); poly m 383 PURPOSE: Compute the PBW basis, up to a given maximal exponent, of a current algebra. 384 INPUT: Maximal exponent is given by the corresponding monomial. 385 RETURN: ideal consisting of PBW elements. 378 PURPOSE: Compute PBW elements up to a given maximal one. 379 RETURN: ideal consisting of found elements. 386 380 NOTE: Unit is omitted. Weights are ignored! 387 381 " 388 382 { 389 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "PBW_maxMonom", MaxMonom ); }; /*4DEBUG*/ 390 391 ideal K = ideal(); 392 393 intvec exp = leadexp( MaxMonom ); 394 395 for ( int k = 1; k <= size(exp); k ++ ) 396 { 397 K[ 1 + size(K) ] = var(k)^( 1 + exp[k] ); 398 } 399 400 attrib(K, "isSB", 1); 401 402 ideal Basis = kbase( K ); 403 404 kill K; 405 406 ideal K = Basis[ (ncols(Basis)1)..1]; // reverse, kill last 1 407 408 kill Basis; 409 410 ideal T = smoothQideal( K ); kill K; 411 412 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "PBW_maxMonom", T ); }; /*4DEBUG*/ 413 414 return( T ); 383 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "PBW_maxMonom", MaxMonom ); }; /*4DEBUG*/ 384 385 ideal K = 0; 386 387 intvec exp = leadexp( MaxMonom ); 388 389 for ( int k = 1; k <= size(exp); k ++ ) 390 { 391 K[ 1 + size(K) ] = var(k)^( 1 + exp[k] ); 392 } 393 394 attrib(K, "isSB", 1); 395 396 K = kbase(K); 397 398 K = K[ (ncols(K)1)..1]; // reverse, kill last 1 399 400 K = smoothQideal( K ); 401 402 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "PBW_maxMonom", K ); }; /*4DEBUG*/ 403 404 return( K ); 415 405 } 416 406 example 417 407 { 418 "EXAMPLE:"; echo = 2;419 ring A = 0,(e,f,h),dp;420 matrix D[3][3]=0;421 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f;422 ncalgebra(1,D); // this algebra is U(sl_2)423 424 // At most 1st degree in e, h and at most 2nd degree in f, unit is omitted:425 PBW_maxMonom( e*(f^2)* h );408 "EXAMPLE:"; echo = 2; 409 ring A = 0,(e,f,h),dp; 410 matrix D[3][3]=0; 411 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f; 412 ncalgebra(1,D); // this algebra is U(sl_2) 413 414 // At most 1st degree in e, h and at most 2nd degree in f, unit is omitted: 415 PBW_maxMonom( e*(f^2)* h ); 426 416 } 427 417 … … 437 427 /******************************************************/ 438 428 proc applyAdF( ideal I, poly p ) 439 "429 " 440 430 USAGE: applyAdF( Basis, f); ideal Basis, poly f 441 431 PURPOSE: Apply Ad_{f} to every element of Basis … … 444 434 " 445 435 { 446 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "applyAdF", I, p ); }; /*4DEBUG*/ 447 448 poly t; 449 450 ideal II = ideal(); 451 452 for ( int k = ncols(I); k > 0; k ) 453 { 454 t = I[k]; 455 II[k] = p * t  t * p; // we have to reduce smooth images in Qrings... 456 } 457 458 ideal J = NF( II, twostd(0) ); // ... now! 459 460 kill II; 461 462 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "applyAdF", J ); }; /*4DEBUG*/ 463 return( J ); 436 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "applyAdF", I, p ); }; /*4DEBUG*/ 437 438 poly t; ideal II = 0; 439 440 for ( int k = ncols(I); k > 0; k ) 441 { 442 t = I[k]; 443 II[k] = p * t  t * p; // we have to reduce smooth images in Qrings... 444 } 445 446 II = NF( II, twostd(0) ); // ... now! 447 448 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "applyAdF", II ); }; /*4DEBUG*/ 449 return( II ); 464 450 } 465 451 example 466 452 { 467 "EXAMPLE:"; echo = 2;468 ring A = 0,(e,f,h),dp;469 matrix D[3][3]=0;470 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f;471 ncalgebra(1,D); // this algebra is U(sl_2)472 473 // Let us consider the linear map Ad_{e} from A_2 into A.474 475 // Compute the PBW basis of A_2:476 ideal Basis = PBW_maxDeg( 2 ); Basis;477 478 // Compute images of basis elements under the linear map Ad_e:479 ideal Image = applyAdF( Basis, e ); Image;480 481 // Now we have a linear map given by: Basis_i > Image_i482 // Let's compute its kernel:483 module C = linearMapKernel( Image ); C;484 485 // Now we can compute the kernel of Ad_e by means of basis vectors:486 ideal K = linearCombinations(Basis, C); K;487 488 // Let's check that Ad_e(K) is zero:489 applyAdF( K, e );453 "EXAMPLE:"; echo = 2; 454 ring A = 0,(e,f,h),dp; 455 matrix D[3][3]=0; 456 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f; 457 ncalgebra(1,D); // this algebra is U(sl_2) 458 459 // Let us consider the linear map Ad_{e} from A_2 into A. 460 461 // Compute the PBW basis of A_2: 462 ideal Basis = PBW_maxDeg( 2 ); Basis; 463 464 // Compute images of basis elements under the linear map Ad_e: 465 ideal Image = applyAdF( Basis, e ); Image; 466 467 // Now we have a linear map given by: Basis_i > Image_i 468 // Let's compute its kernel: 469 module C = linearMapKernel( Image ); C; 470 471 // Now we can compute the kernel of Ad_e by means of basis vectors: 472 ideal K = linearCombinations(Basis, C); K; 473 474 // Let's check that Ad_e(K) is zero: 475 applyAdF( K, e ); 490 476 } 491 477 … … 495 481 proc linearMapKernel( ideal Images ) 496 482 "USAGE: linearMapKernel( Images ); ideal Images 497 PURPOSE: Computes the kernel of a linear map given by its images on certain basis vectors 483 PURPOSE: Computes the kernel of a linear map: e_i \mapsto Images[i], 484 @* where e_i are certain basis vectors 498 485 RETURN: syzygy module, or 0 if all images are zeroes 499 486 SEE ALSO: applyAdF; linearMapKernel 500 487 " 501 488 { 502 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "linearMapKernel", Images ); }; /*4DEBUG*/ 503 504 // This must be a list of monomials in a form of polynomial (sum with coeffs == 1) 505 poly Monomials = monomialForm( Images ); 506 507 int N = size( Monomials ); // number of different monomials 508 509 if ( N == 0 ) // & ncols( Images ) > 0 => all Images == 0 510 { 511 int result = 0; 512 513 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "linearMapKernel", result ); }; /*4DEBUG*/ 514 return( result ); 515 } 516 517 // Compute matrix MD 518 module MD; // zero 519 520 int x, y; 521 522 vector w; 523 524 poly p, m; 525 526 int V = ncols(Images); // must be equal to ncols(Basis) and size(Basis)! 527 528 // We take monomials as vector space basis of <Image>_k... 529 530 // TODO: Is there any other way to compute a basis of it and represent images as 531 // linear combination of them??? 532 533 // Maybe some 'free resolution' stuff??? But we need linear maps only!!! 534 535 for ( x = 1; x <= V; x++ ) // For every Image vector 536 { 537 w = 0; 538 539 540 p = Images[x]; 541 542 y = 1; // from 1st monomial in Monomials... 543 544 while( size(p) > 0 ) 545 { 546 m = leadmonom(p); 547 548 // y < N! 549 while( Monomials[y] != m ) 489 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "linearMapKernel", Images ); }; /*4DEBUG*/ 490 491 // This must be a list of monomials in a form of polynomial (sum with coeffs == 1) 492 poly Monomials = monomialForm( Images ); 493 494 int N = size( Monomials ); // number of different monomials 495 496 if ( N == 0 ) // & ncols( Images ) > 0 => all Images == 0 497 { 498 int result = 0; 499 500 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "linearMapKernel", result ); }; /*4DEBUG*/ 501 return( result ); 502 } 503 504 // Compute matrix MD 505 module MD; // zero 506 507 int x, y; 508 509 vector w; 510 511 poly p, m; 512 513 int V = ncols(Images); // must be equal to ncols(Basis) and size(Basis)! 514 515 // We take monomials as vector space basis of <Image>_k... 516 517 // TODO: Is there any other way to compute a basis of it and represent images as 518 // linear combination of them??? 519 520 // Maybe some 'free resolution' stuff??? But we need linear maps only!!! 521 522 for ( x = 1; x <= V; x++ ) // For every Image vector 523 { 524 w = 0; 525 526 p = Images[x]; 527 528 y = 1; // from 1st monomial in Monomials... 529 530 while( size(p) > 0 ) 531 { 532 m = leadmonom(p); 533 534 // y < N! 535 while( Monomials[y] != m ) 550 536 // There MUST be this monomial because of the construction of Monomials polynomial! 551 537 { 552 538 y++; // to size() 553 539 } 554 540 555 556 557 558 541 // found monomial m at position y. 542 543 w = w + leadcoef(p) * gen(y); // leadcoef(p) MUST be nonzero! 544 p = p  lead(p); // kill lead term 559 545 } 560 546 561 MD[x] = w; 562 } 563 564 /*******************************************/ 565 566 // save options 567 intvec v = option( get ); 568 569 // set right options 570 option( redSB ); 571 option( redTail ); 572 573 // compute everything in a right form 574 module linearMapKernel = simplify( std( syz(MD) ), 1 + 2 + 8 ); 575 // note that MD is a matrix of numbers  no polynomials... 576 577 // restore options 578 option( set, v ); 579 580 // kill used structure 581 kill MD; 582 583 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "linearMapKernel", linearMapKernel ); }; /*4DEBUG*/ 584 585 return( linearMapKernel ); 547 MD[x] = w; 548 } 549 550 /*******************************************/ 551 552 // save options 553 intvec v = option( get ); 554 555 // set right options 556 option( redSB ); 557 option( redTail ); 558 559 // compute everything in a right form 560 MD = simplify( std( syz(MD) ), 1 + 2 + 8 ); 561 // note that MD is a matrix of numbers  no polynomials... 562 563 // restore options 564 option( set, v ); 565 566 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "linearMapKernel", MD ); }; /*4DEBUG*/ 567 568 return( MD ); 586 569 } 587 570 example 588 571 { 589 "EXAMPLE:"; echo = 2;590 ring A = 0,(e,f,h),dp;591 matrix D[3][3]=0;592 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f;593 ncalgebra(1,D); // this algebra is U(sl_2)594 595 // Let us consider the linear map Ad_{e} from A_2 into A.596 597 // Compute the PBW basis of A_2:598 ideal Basis = PBW_maxDeg( 2 ); Basis;599 600 // Compute images of basis elements under the linear map Ad_e:601 ideal Image = applyAdF( Basis, e ); Image;602 603 // Now we have a linear map given by: Basis_i > Image_i604 // Let's compute its kernel:605 module C = linearMapKernel( Image ); C;606 607 // Now we can compute the kernel of Ad_e by means of basis vectors:608 ideal K = linearCombinations(Basis, C); K;609 610 // Let's check that Ad_e(K) is zero:611 ideal Z = applyAdF( K, e ); Z;612 613 // Now linearMapKernel will return a single integer 0:614 def CC = linearMapKernel(Z); typeof(CC); CC;572 "EXAMPLE:"; echo = 2; 573 ring A = 0,(e,f,h),dp; 574 matrix D[3][3]=0; 575 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f; 576 ncalgebra(1,D); // this algebra is U(sl_2) 577 578 // Let us consider the linear map Ad_{e} from A_2 into A. 579 580 // Compute the PBW basis of A_2: 581 ideal Basis = PBW_maxDeg( 2 ); Basis; 582 583 // Compute images of basis elements under the linear map Ad_e: 584 ideal Image = applyAdF( Basis, e ); Image; 585 586 // Now we have a linear map given by: Basis_i > Image_i 587 // Let's compute its kernel: 588 module C = linearMapKernel( Image ); C; 589 590 // Now we can compute the kernel of Ad_e by means of basis vectors: 591 ideal K = linearCombinations(Basis, C); K; 592 593 // Let's check that Ad_e(K) is zero: 594 ideal Z = applyAdF( K, e ); Z; 595 596 // Now linearMapKernel will return a single integer 0: 597 def CC = linearMapKernel(Z); typeof(CC); CC; 615 598 } 616 599 … … 618 601 /******************************************************/ 619 602 proc linearCombinations( ideal Basis, module KER ) 620 "603 " 621 604 USAGE: linearCombinations( Basis, C ); ideal Basis, module C 622 PURPOSE: computes linear combinations of Basis vectors with the coefficients from C.623 RETURN: ideal of linear combinations of Basis vectors with the coefficients from C.605 PURPOSE: computes linear combinations of Basis vectors with coefficients from C 606 RETURN: ideal generated by computed linear combinations 624 607 SEE ALSO: linearMapKernel; applyAdF 625 608 " 626 609 { 627 610 628 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "linearCombinations", Basis, KER ); }; /*4DEBUG*/629 630 631 632 633 634 635 636 637 638 639 ideal result = ideal();640 641 642 643 { 644 645 646 647 648 { 649 650 651 611 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "linearCombinations", Basis, KER ); }; /*4DEBUG*/ 612 613 614 number c; 615 616 int x, y; 617 618 vector w; 619 620 poly p; 621 622 ideal result = 0; 623 624 // Kernel' basis translation 625 for ( x = 1; x <= ncols(KER); x++ ) 626 { 627 p = 0; 628 w = KER[x]; 629 630 for ( y = 1; y <= nrows(w); y++ ) 631 { 632 c = leadcoef( w[y] ); 633 634 if ( c != 0 ) 652 635 { 653 636 p = p + c * Basis[y]; // linear combination of base vectors { Basis[y] } 654 637 } 655 638 } 656 657 result[ 1 + size(result) ] = p;658 }659 660 661 662 663 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "linearCombinations", result ); }; /*4DEBUG*/664 665 639 result[ x ] = p; 640 } 641 642 643 // no reduction in quotient algebras is needed. No multiplications were done! 644 645 646 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "linearCombinations", result ); }; /*4DEBUG*/ 647 648 return( result ); 666 649 } 667 650 example 668 651 { 669 "EXAMPLE:"; echo = 2;670 ring A = 0,(e,f,h),dp;671 matrix D[3][3]=0;672 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f;673 ncalgebra(1,D); // this algebra is U(sl_2)674 675 // Let us consider the linear map Ad_{e} from A_2 into A.676 677 // Compute the PBW basis of A_2:678 ideal Basis = PBW_maxDeg( 2 ); Basis;679 680 // Compute images of basis elements under the linear map Ad_e:681 ideal Image = applyAdF( Basis, e ); Image;682 683 // Now we have a linear map given by: Basis_i > Image_i684 // Let's compute its kernel:685 module C = linearMapKernel( Image ); C;686 687 // Now we can compute the kernel of Ad_e by means of basis vectors:688 ideal K = linearCombinations(Basis, C); K;689 690 // Let's check that Ad_e(K) is zero:691 applyAdF( K, e );652 "EXAMPLE:"; echo = 2; 653 ring A = 0,(e,f,h),dp; 654 matrix D[3][3]=0; 655 D[1,2]=h; D[1,3]=2*e; D[2,3]=2*f; 656 ncalgebra(1,D); // this algebra is U(sl_2) 657 658 // Let us consider the linear map Ad_{e} from A_2 into A. 659 660 // Compute the PBW basis of A_2: 661 ideal Basis = PBW_maxDeg( 2 ); Basis; 662 663 // Compute images of basis elements under the linear map Ad_e: 664 ideal Image = applyAdF( Basis, e ); Image; 665 666 // Now we have a linear map given by: Basis_i > Image_i 667 // Let's compute its kernel: 668 module C = linearMapKernel( Image ); C; 669 670 // Now we can compute the kernel of Ad_e by means of basis vectors: 671 ideal K = linearCombinations(Basis, C); K; 672 673 // Let's check that Ad_e(K) is zero: 674 applyAdF( K, e ); 692 675 } 693 676 … … 696 679 /******************************************************/ 697 680 static proc LINEAR_MAP_KERNEL(ideal Basis, ideal Images ) // Ker of the linear map given by its values on basis vectors 698 "681 " 699 682 PURPOSE: Computation of the kernel basis of the linear map given by the list @given 700 683 " 701 684 { 702 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "LINEAR_MAP_KERNEL", Basis, Images ); }; /*4DEBUG*/703 704 ideal result = ideal();705 706 707 { 708 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "LINEAR_MAP_KERNEL", result ); }; /*4DEBUG*/709 710 } 711 712 713 714 715 716 717 718 { 719 720 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "LINEAR_MAP_KERNEL", Basis ); }; /*4DEBUG*/721 685 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "LINEAR_MAP_KERNEL", Basis, Images ); }; /*4DEBUG*/ 686 687 ideal result = 0; 688 689 if ( size( Basis ) == 0 ) 690 { 691 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "LINEAR_MAP_KERNEL", result ); }; /*4DEBUG*/ 692 return( result ); 693 } 694 695 // compute fundamental solutions system 696 def T = linearMapKernel( Images ); 697 698 699 // check result of linearMapKernel 700 if( (typeof(T) == "int") and (T == 0) ) 701 { 702 // All zeroes! Return Basis: 703 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "LINEAR_MAP_KERNEL", Basis ); }; /*4DEBUG*/ 704 return( Basis ); 722 705 } 723 724 { 725 726 { 727 706 else 707 { 708 if( typeof(T) != "module" ) 709 { 710 ERROR( "Wrong output from the 'linearMapKernel' function!" ); 728 711 } 729 712 } 730 713 731 result = linearCombinations( Basis, T ); 732 733 kill T; 734 735 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "LINEAR_MAP_KERNEL", result ); }; /*4DEBUG*/ 736 return( result ); 714 result = linearCombinations( Basis, T ); 715 716 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "LINEAR_MAP_KERNEL", result ); }; /*4DEBUG*/ 717 return( result ); 737 718 } 738 719 … … 747 728 " 748 729 { 749 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "ZeroKer", Basis, Images ); }; /*4DEBUG*/750 751 ideal result = ideal();752 753 754 { 755 756 { 757 730 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "ZeroKer", Basis, Images ); }; /*4DEBUG*/ 731 732 ideal result = 0; 733 734 for( int i = 1; i <= ncols( Basis ); i++ ) 735 { 736 if( size( Images[i] ) == 0 ) // zero image? 737 { 738 result[ 1 + size(result) ] = Basis[i]; // take this basis vector! 758 739 } 759 740 } 760 741 761 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "ZeroKer", result ); }; /*4DEBUG*/762 742 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "ZeroKer", result ); }; /*4DEBUG*/ 743 return( result ); 763 744 } 764 745 … … 776 757 "USAGE: variablesStandard(); 777 758 RETURN: ideal, generated by algebra variables 778 PURPOSE: computes the ideal generated byalgebra variables taken in their natural order759 PURPOSE: computes the set of algebra variables taken in their natural order 779 760 SEE ALSO: variablesSorted 780 761 EXAMPLE: example variablesStandard; shows an example 781 762 " 782 763 { 783 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "variablesStandard" ); }; /*4DEBUG*/784 785 786 787 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "variablesStandard", result ); }; /*4DEBUG*/788 764 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "variablesStandard" ); }; /*4DEBUG*/ 765 766 ideal result = maxideal(1); 767 768 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "variablesStandard", result ); }; /*4DEBUG*/ 769 return( result ); 789 770 } 790 771 example 791 772 { 792 "EXAMPLE:"; echo = 2;793 ring A = 0,(x,y,z),dp;794 matrix D[3][3]=0;795 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;796 ncalgebra(1,D); // this algebra is U(sl_2)797 // Variables in their natural order:798 variablesStandard();773 "EXAMPLE:"; echo = 2; 774 ring A = 0,(x,y,z),dp; 775 matrix D[3][3]=0; 776 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 777 ncalgebra(1,D); // this algebra is U(sl_2) 778 // Variables in their natural order: 779 variablesStandard(); 799 780 } 800 781 … … 803 784 "USAGE: variablesSorted(); 804 785 RETURN: ideal, generated by sorted algebra variables 805 PURPOSE: computes the ideal generated by algebra variables sorted so that Cartan variables are first and all other variables are behind. 806 NOTE: This is a heuristics for the computation of center: it is better to compute centralizers of Cartan variables first since we can omit solving the system of equations. 786 PURPOSE: computes the set of algebra variables sorted so that 787 @* Cartan variables go first 788 NOTE: This is a heuristics for the computation of center: 789 @* it is better to compute centralizers of Cartan variables first since in this 790 @* case we can omit solving the system of equations. 807 791 SEE ALSO: variablesStandard 808 792 EXAMPLE: example variablesSorted; shows an example 809 793 "{ 810 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "variablesSorted" ); }; /*4DEBUG*/811 812 813 814 815 ideal result;816 817 818 819 820 821 822 823 { 824 825 826 827 { 828 829 794 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "variablesSorted" ); }; /*4DEBUG*/ 795 796 ideal V = variablesStandard(); 797 int N = size( V ); // == nvars( basering ) 798 799 ideal result = 0; 800 801 int r_begin = 1; 802 int r_end = N; 803 804 poly v; 805 806 for( int k = 1; k <= N; k++ ) 807 { 808 v = V[k]; 809 810 if( isCartan(v) == 1 ) // Cartan elements go 1st 811 { 812 result[r_begin] = v; 813 r_begin++; 830 814 } else // Other  in the end... 831 815 { 832 816 result[r_end] = v; 833 817 r_end; 834 835 } 836 837 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "variablesSorted", result ); }; /*4DEBUG*/838 818 } 819 } 820 821 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "variablesSorted", result ); }; /*4DEBUG*/ 822 return( result ); 839 823 } 840 824 example 841 825 { 842 "EXAMPLE:"; echo = 2;843 ring A = 0,(x,y,z),dp;844 matrix D[3][3]=0;845 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;846 ncalgebra(1,D); // this algebra is U(sl_2)847 // There is only one Cartan variable  z in U(sl_2),848 // it must go 1st:849 variablesSorted();826 "EXAMPLE:"; echo = 2; 827 ring A = 0,(x,y,z),dp; 828 matrix D[3][3]=0; 829 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 830 ncalgebra(1,D); // this algebra is U(sl_2) 831 // There is only one Cartan variable  z in U(sl_2), 832 // it must go 1st: 833 variablesSorted(); 850 834 } 851 835 … … 867 851 proc centralizeSet( ideal F, ideal V ) // HL 'core' function 868 852 "USAGE: centralizeSet( F, V ); ideal F, ideal V 869 INPUT: 853 INPUT: a finite set of elements F, vector space basis V 870 854 RETURN: ideal, generated by base elements 871 PURPOSE: computes the vector space basis of the centralizer of F in the vector space spanned by V, that is, Cen(F[N],Cen(F[N1],...,Cen(F[1],V)...)) 855 PURPOSE: computes a v.s. basis of the centralizer of F in v.s. V 856 NOTE: Cen(F,V) is computed by the formula Cen(F[N],..,Cen(F[1],V)..) 872 857 SEE ALSO: centralizerVS; centralizer; inCentralizer 873 858 EXAMPLE: example centralizeSet; shows an example 874 859 " 875 860 { 876 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centralizeSet", F, V ); }; /*4DEBUG*/877 878 879 880 881 { 882 883 } 884 885 886 887 for( int v = 1; (v <= N) and (size(V) > 0); v++ )888 {889 DBPrint(1, "Centralizing " + string(F[v]) );890 891 //892 ideal Images = applyAdF( V, F[v] );893 894 ideal K;895 896 { 897 K= ZeroKer( V, Images );861 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centralizeSet", F, V ); }; /*4DEBUG*/ 862 863 int N = size(F); 864 865 if( N == 0) 866 { 867 ERROR( "F MUST be non empty!!!" ); 868 } 869 870 DBPrint(1, "BasisSize: " + string(size(V)) ); 871 872 ideal Images; 873 874 for( int v = 1; (v <= N) and (size(V) > 0); v++ ) 875 { 876 DBPrint(1, "Centralizing " + string(F[v]) ); 877 878 Images = applyAdF( V, F[v] ); 879 880 if( (isCartan(F[v]) == 1) or (size(V) == 1) ) 881 { 882 V = ZeroKer( V, Images ); 898 883 } else 899 { 900 K = LINEAR_MAP_KERNEL( V, Images ); 901 } 902 903 kill Images; 904 kill V; ideal V = K; kill K; 905 906 // Printing... 907 DBPrint(1, "Progress: [ " + string(v) + " / " + string(N) + " ]"+ 908 " => BasisSize: " + string(size(V)) ); 909 } 910 911 ideal result = makeNice(V); kill V; 912 913 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centralizeSet", result ); }; /*4DEBUG*/ 914 915 return( result ); 884 { 885 V = LINEAR_MAP_KERNEL( V, Images ); 886 } 887 888 // Printing... 889 DBPrint(1, "Progress: [ " + string(v) + " / " + string(N) + " ]"+ 890 " => BasisSize: " + string(size(V)) ); 891 } 892 893 V = makeNice(V); 894 895 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centralizeSet", V ); }; /*4DEBUG*/ 896 897 return( V ); 916 898 } 917 899 example 918 900 { 919 "EXAMPLE:"; echo = 2;920 ring A_4_1 = 0,(e(1..4)),dp;921 matrix D[4][4]=0;922 D[2,4] = e(1);923 D[3,4] = e(2);924 // This is $A_{41}$  the first real Lie algebra of dimension $4$.925 ncalgebra(1,D);901 "EXAMPLE:"; echo = 2; 902 ring A_4_1 = 0,(e(1..4)),dp; 903 matrix D[4][4]=0; 904 D[2,4] = e(1); 905 D[3,4] = e(2); 906 // This is $A_{41}$  the first real Lie algebra of dimension $4$. 907 ncalgebra(1,D); 926 908 927 ideal F = variablesSorted(); F;909 ideal F = variablesSorted(); F; 928 910 929 // the center of $A_{41}$ is generated by930 // $e(1)$ and $1/2*e(2)^2+e(1)*e(3)$931 // therefore one may consider computing it in the following way:911 // the center of $A_{41}$ is generated by 912 // $e(1)$ and $1/2*e(2)^2+e(1)*e(3)$ 913 // therefore one may consider computing it in the following way: 932 914 933 // 1. Compute PBW basis consisting of934 // monomials of exponent <= (1,2,1,0)935 ideal V = PBW_maxMonom( e(1) * e(2)^ 2 * e(3) );915 // 1. Compute PBW basis consisting of 916 // monomials of exponent <= (1,2,1,0) 917 ideal V = PBW_maxMonom( e(1) * e(2)^ 2 * e(3) ); 936 918 937 // 2. Compute the centralizer of F within vector space938 // spanned by these monomials:939 ideal C = centralizeSet( F, V ); C;919 // 2. Compute the centralizer of F within vector space 920 // spanned by these monomials: 921 ideal C = centralizeSet( F, V ); C; 940 922 941 inCenter(C);923 inCenter(C); 942 924 } 943 925 … … 946 928 /******************************************************/ 947 929 proc centralizerVS( ideal F, int d ) 948 "USAGE: centralizerVS( F, D ); ideal F, int D930 "USAGE: centralizerVS( F, D ); ideal F, int D 949 931 RETURN: ideal, generated by elements of degree <= D 950 PURPOSE: computes a v ector spacebasis of the centralizer of F up to degree D.932 PURPOSE: computes a v.s. basis of the centralizer of F up to degree D. 951 933 NOTE: D must be nonnegative 952 934 SEE ALSO: centerVS; centralizer; inCentralizer … … 954 936 " 955 937 { 956 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centralizerVS", F, d ); }; /*4DEBUG*/957 958 959 { 960 938 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centralizerVS", F, d ); }; /*4DEBUG*/ 939 940 if( size(F) == 0) 941 { 942 ERROR( "F MUST be nonempty!!!" ); 961 943 } 962 944 963 ideal V = PBW_maxDeg( d ); // PBW basis 964 965 ideal result = centralizeSet( F, V ); // basis of the Centralizer of S in VS <V> 966 967 kill V; 968 969 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centralizerVS", result ); }; /*4DEBUG*/ 970 971 return( result ); 945 ideal V = centralizeSet( F, PBW_maxDeg( d ) ); // basis of the Centralizer of S in PBW basis 946 947 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centralizerVS", V ); }; /*4DEBUG*/ 948 949 return( V ); 972 950 } 973 951 example 974 952 { 975 "EXAMPLE:"; echo = 2;976 ring A = 0,(x,y,z),dp;977 matrix D[3][3]=0;978 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;979 ncalgebra(1,D); // this algebra is U(sl_2)980 ideal F = x, y;981 // find all elements commuting with x and y of degree <= 4:982 ideal C = centralizerVS(F, 4);983 C;984 inCentralizer(C, F);953 "EXAMPLE:"; echo = 2; 954 ring A = 0,(x,y,z),dp; 955 matrix D[3][3]=0; 956 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 957 ncalgebra(1,D); // this algebra is U(sl_2) 958 ideal F = x, y; 959 // find all elements commuting with x and y of degree <= 4: 960 ideal C = centralizerVS(F, 4); 961 C; 962 inCentralizer(C, F); 985 963 } 986 964 … … 999 977 "USAGE: centerVS( D ); int D 1000 978 RETURN: ideal, generated by elements of degree <= D 1001 PURPOSE: computes a v ector space basis of the center of the current algebraup to degree D.979 PURPOSE: computes a v.s. basis of the center up to degree D. 1002 980 NOTE: D must be nonnegative 1003 981 SEE ALSO: centralizerVS; center; inCenter … … 1005 983 " 1006 984 { 1007 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centerVS", D ); }; /*4DEBUG*/1008 1009 1010 1011 { 1012 ERROR( "No current ring!" );1013 } 1014 1015 1016 { 1017 1018 } 1019 1020 985 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centerVS", D ); }; /*4DEBUG*/ 986 987 988 if( nameof( basering ) == "basering" ) 989 { 990 // ERROR( "No current ring!" ); 991 } 992 993 if( D < 0 ) 994 { 995 ERROR( "Degree D must be nonnegative!" ); 996 } 997 998 ideal result = centralizerVS( variablesSorted(), D ); 1021 999 1022 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centerVS", result ); }; /*4DEBUG*/1023 1024 1000 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centerVS", result ); }; /*4DEBUG*/ 1001 1002 return( result ); 1025 1003 } 1026 1004 example 1027 1005 { 1028 "EXAMPLE:"; echo = 2;1029 ring A = 0,(x,y,z),dp;1030 matrix D[3][3]=0;1031 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;1032 ncalgebra(1,D); // this algebra is U(sl_2)1033 // find all central elements of degree <= 41034 ideal Z = centerVS(4);1035 Z;1036 // note that the second element is the square of the first1037 // plus the multiple of the first:1038 Z[2]  Z[1]^2 + 8*Z[1];1039 inCenter(Z);1006 "EXAMPLE:"; echo = 2; 1007 ring A = 0,(x,y,z),dp; 1008 matrix D[3][3]=0; 1009 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 1010 ncalgebra(1,D); // this algebra is U(sl_2) 1011 // find all central elements of degree <= 4 1012 ideal Z = centerVS(4); 1013 Z; 1014 // note that the second element is the square of the first 1015 // plus the multiple of the first: 1016 Z[2]  Z[1]^2 + 8*Z[1]; 1017 inCenter(Z); 1040 1018 } 1041 1019 … … 1045 1023 "USAGE: centralizerRed( F, D[, N] ); ideal F, int D[, int N] 1046 1024 RETURN: ideal, generated by computed generators 1047 PURPOSE: if N is absent and D >= 0 computes a subalgebra generators of the centralizer of F up to degree D, otherwise if N is present computes N(at least) first generators of the centralizer, if moreover D > 0 it will be used as the first maximal degree estimation. 1048 NOTE: Current ordering must be a degree compatible wellordering. 1025 PURPOSE: computes a subalgebra generators of centralizer(F) up to degree D. 1026 NOTE: In general, one cannot compute the whole centralizer(F). 1027 @* Hence, one has to specify a termination condition via arguments D and/or N. 1028 @* If D is positive, only centralizing elements up to degree D will be found. 1029 @* If D is negative, the termination is determined by N only. 1030 @* If N is given, the computation stops if at least N elements has been found. 1031 @* Warning: if N is given and bigger than the actual number of generators, 1032 @* the procedure may not terminate. 1033 @* Current ordering must be a degree compatible wellordering. 1049 1034 SEE ALSO: centralizerVS; centerRed; centralizer; inCentralizer 1050 1035 EXAMPLE: example centralizerRed; shows an example 1051 1036 " 1052 1037 { 1053 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centralizerRed", F, D, # ); }; /*4DEBUG*/1054 1055 1056 { 1057 ERROR( "No current ring!" );1038 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centralizerRed", F, D, # ); }; /*4DEBUG*/ 1039 1040 if( nameof( basering ) == "basering" ) 1041 { 1042 // ERROR( "No current ring!" ); 1058 1043 } 1059 1044 1060 1061 { 1062 1045 if( size(F) == 0) 1046 { 1047 ERROR( "F MUST be nonempty!!!" ); 1063 1048 } 1064 1065 def NCRING = basering; // Noncommutative ring 1066 list L = ringlist( NCRING ); 1067 def L1, L2, L3, L4 = L[1..4]; // General components 1068 def COMMRING = ring( list( L1, L2, L3, L4 ) ); // Underlying commutative ring 1069 kill L1, L2, L3, L4, L; 1070 1071 1072 int k = DefaultInt(#); 1073 1074 1075 int m = (k > 0); 1076 1077 int MinDeg = 6; // starting guess for Maximal Bounding Degree, 6 1078 int Delta = 4; // increment of it, 4 1079 1080 if( m and (D <= 0) ) 1081 { 1082 // minimal guess 1083 D = MinDeg; 1084 } 1085 1086 if( !m and D < 0) 1087 { 1088 ERROR("Wrong bounding condition!"); 1089 } 1090 1091 ideal result = ideal(); 1092 1093 int i, j, l, d; 1094 1095 // we keep the list of found leading monomials in the commutative ring! 1096 setring COMMRING; 1097 1098 // Init 1099 list FOUND_LEADING_MONOMIALS = list(); 1100 1101 for( i = 1; i <= D; i++ ) 1102 { 1103 FOUND_LEADING_MONOMIALS[i] = ideal(); 1104 } 1105 1106 setring NCRING; 1107 1108 // Main loop: 1109 i = 1; 1110 1111 ideal PBW = ideal(); 1112 ideal NEW; 1113 1114 while( i <= D ) 1115 { 1116 DBPrint( 1, "Current degree is " + string(i) ); 1117 1118 // Compute current "reduced" PBW basis... 1119 1120 // Prepare current found leading monomials 1121 setring COMMRING; 1122 1123 if( defined(FLM) ) 1124 { 1125 kill FLM; 1049 1050 ///////////////////////////////////////////////////////////////////////////// 1051 1052 int i, j, l, d; 1053 1054 ///////////////////////////////////////////////////////////////////////////// 1055 1056 int k = DefaultInt(#); 1057 1058 int m = (k > 0); 1059 1060 int MinDeg = 6; // starting guess for Maximal Bounding Degree, 6 1061 int Delta = 4; // increment of it, 4 1062 1063 if( m and (D <= 0) ) 1064 { 1065 // minimal guess 1066 D = MinDeg; 1067 } 1068 1069 if( !m and D < 0) 1070 { 1071 ERROR("Wrong bounding condition!"); 1072 } 1073 1074 ///////////////////////////////////////////////////////////////////////////// 1075 1076 def NCRING = basering; // Noncommutative ring 1077 list L = ringlist( NCRING ); 1078 def L1, L2, L3, L4 = L[1..4]; // General components 1079 1080 def COMMRING = ring( list( L1, L2, L3, L4 ) ); // Underlying commutative ring 1081 1082 1083 ///////////////////////////////////////////////////////////////////////////// 1084 1085 // we keep the list of found leading monomials in the commutative ring! 1086 setring COMMRING; 1087 1088 // Init 1089 list FOUND_LEADING_MONOMIALS = list(); 1090 1091 for( i = 1; i <= D; i++ ) 1092 { 1093 FOUND_LEADING_MONOMIALS[i] = ideal(); 1094 } 1095 1096 ideal FLM, NEW, T; // in COMMRING 1097 1098 ///////////////////////////////////////////////////////////////////////////// 1099 1100 setring NCRING; 1101 1102 ideal result, FLM, PBW, NEW, T, P; // in NCRING 1103 1104 // Main loop: 1105 i = 1; 1106 1107 while( i <= D ) 1108 { 1109 DBPrint( 1, "Current degree is " + string(i) ); 1110 1111 ///////////////////////////////////////////////////////////////////////////// 1112 1113 // Compute current "reduced" PBW basis... 1114 1115 // Prepare current found leading monomials 1116 setring COMMRING; 1117 FLM = FOUND_LEADING_MONOMIALS[i]; 1118 1119 // And back to NCRing 1120 setring NCRING; 1121 1122 FLM = imap(COMMRING, FLM); // We cannot write imap(COMMRING, FOUND_LEADING_MONOMIALS[i]) :((( 1123 1124 attrib(FLM, "isSB", 1); // just to avoid "no standard basis" warning. 1125 1126 // degrees should not change, 1127 // no monomials should be multiplied here 1128 T = reduce( PBW_eqDeg( i ), FLM, 1 ); 1129 1130 // we simply kill in T monomials occurring in FOUND_LEADING_MONOMIALS[i] 1131 P = PBW + T; // + simplifies 1132 1133 // Compute current centralizer 1134 NEW = centralizeSet( F, P ); 1135 1136 if( size(NEW) > 0 ) 1137 { 1138 // In order to speedup multiplications we are going into a commutative ring: 1139 setring COMMRING; 1140 1141 // we can perform commutative interreduction 1142 // since no monomials should be multiplied! 1143 // degrees should not change 1144 NEW = interred( imap( NCRING, NEW ) ); 1145 1146 // Go back! 1147 setring NCRING; 1148 1149 NEW = imap( COMMRING, NEW ); 1150 1151 DBPrint( 1, "Found: ", NEW ); 1152 1153 // Add them to result... 1154 result = result + NEW; 1126 1155 } 1127 1156 1128 ideal FLM = FOUND_LEADING_MONOMIALS[i]; 1129 1130 // And back to NCRing 1131 setring NCRING; 1132 1133 1134 1135 ideal FLM = imap(COMMRING, FLM); 1136 attrib(FLM, "isSB", 1); // just to avoid "no standard basis" warning. 1137 1138 // degrees should not change, 1139 // no monomials should be multiplied here 1140 ideal T = reduce( PBW_eqDeg( i ), FLM, 1 ); 1141 1142 kill FLM; 1143 1144 // we simply kill in T monomials occurring in FOUND_LEADING_MONOMIALS[i] 1145 ideal P = PBW + T; // here we simplify T 1146 1147 // Compute current centralizer 1148 NEW = centralizeSet( F, P ); 1149 1150 if( size(NEW) > 0 ) 1151 { 1152 // In order to speedup multiplications we are going into a commutative ring: 1153 setring COMMRING; 1154 1155 if( defined(NEW) ) 1157 // Did we find needed number of generators? Or reached the bound? 1158 if( (m and (size(result) >= k)) or (!m and (i == D)) ) 1159 { 1160 break; // Get out of here!!! 1161 } 1162 1163 // otherwise we must update FOUND_LEADING_MONOMIALS 1164 if( size(NEW) > 0 ) 1165 { 1166 setring COMMRING; 1167 1168 FLM = 0; 1169 1170 // We must update FOUND_LEADING_MONOMIALS!!! 1171 for( j = 1; j <= size(NEW); j++ ) 1156 1172 { 1157 kill NEW; 1173 FLM[j] = leadmonom( NEW[j] ); // we are interested in leading monomials only! 1158 1174 } 1159 1175 1160 // we can perform commutative interreduction 1161 // since no monomials should be multiplied! 1162 // degrees should not change 1163 ideal NEW = interred( imap( NCRING, NEW ) ); 1164 1165 // Go back! 1166 setring NCRING; 1167 1168 kill NEW; ideal NEW = imap( COMMRING, NEW ); 1169 1170 DBPrint( 1, "Found: ", NEW ); 1171 1172 // Add them to result... 1173 result = result + NEW; 1174 } 1175 1176 // Did we find needed number of generators? Or reached the bound? 1177 if( (m and (size(result) >= k)) or (!m and (i == D)) ) 1178 { 1179 break; // Get out of here!!! 1180 } 1181 1182 // otherwise we must update FOUND_LEADING_MONOMIALS 1183 if( size(NEW) > 0 ) 1184 { 1185 setring COMMRING; 1186 1187 if( defined(FLM) ) 1188 { 1189 kill FLM; 1190 } 1191 1192 ideal FLM = ideal(); 1193 1194 // We must update FOUND_LEADING_MONOMIALS!!! 1195 for( j = 1; j <= size(NEW); j++ ) 1196 { 1197 FLM[j] = leadmonom( NEW[j] ); // we are interested in leading monomials only! 1198 } 1199 1200 kill NEW; 1201 1202 FOUND_LEADING_MONOMIALS[i] = FOUND_LEADING_MONOMIALS[i] + FLM; 1203 1204 for( j = 1; j <= D; j = j + i ) // For every degree (j*i) of LNEW, do: 1176 FOUND_LEADING_MONOMIALS[i] = FOUND_LEADING_MONOMIALS[i] + FLM; 1177 1178 for( j = 1; j <= D; j = j + i ) // For every degree (j*i) of LNEW, do: 1205 1179 { 1206 1180 for( l = j; (l+i) <= D; l++ ) 1207 1181 { 1208 1209 1182 FOUND_LEADING_MONOMIALS[l+i] = 1183 FOUND_LEADING_MONOMIALS[l+i] + FOUND_LEADING_MONOMIALS[l] * FLM; 1210 1184 } 1211 1185 } 1212 1186 1213 1214 1187 // Return to NCRING 1188 setring NCRING; 1215 1189 1216 // And refine T one more: 1217 1218 1219 ideal FLM = imap(COMMRING, FLM); 1220 attrib(FLM, "isSB", 1);// just to avoid "no standard basis" warning. 1221 1222 // we simply kill in T monomials occurring in FOUND_LEADING_MONOMIALS[i] 1223 kill P; ideal P = PBW + reduce( T, FLM, 1 );1224 kill FLM; 1225 } 1226 1227 kill T, PBW; ideal PBW = P; kill P;1228 1229 1190 FLM = imap(COMMRING, FLM); 1191 attrib(FLM, "isSB", 1);// just to avoid "no standard basis" warning. 1192 1193 // we simply kill in T monomials occurring in FOUND_LEADING_MONOMIALS[i] 1194 T = reduce( T, FLM, 1 ); 1195 1196 PBW = PBW + T; 1197 } else 1198 { 1199 PBW = P; 1200 } 1201 1202 1203 if( m and (i == D) ) // Was the previous estimation too small??? 1230 1204 { 1231 // We must update FOUND_LEADING_MONOMIALS in theirworld:1232 1233 1234 1235 1205 // We must update FOUND_LEADING_MONOMIALS in their Commutative world: 1206 setring COMMRING; 1207 1208 // Init new grades: 1209 for( j = D + 1; j <= (D + Delta); j++ ) 1236 1210 { 1237 1211 FOUND_LEADING_MONOMIALS[j] = ideal(); 1238 1212 } 1239 1213 1240 if( defined(FLM) ) 1214 FLM = 0; 1215 1216 // All previously computed elements in their order! 1217 NEW = imap( NCRING, result ); 1218 1219 for( j = 1; j <= size(NEW); j++ ) 1241 1220 { 1242 kill FLM; 1243 } 1244 1245 ideal FLM = ideal(); 1246 1247 // All previously computed elements in their order! 1248 ideal NEW = imap( NCRING, result ); 1249 1250 for( j = 1; j <= size(NEW); j++ ) 1221 FLM[j] = leadmonom( NEW[j] ); // we are interested in leading monomials only! 1222 } 1223 1224 while( size(FLM) > 0 ) 1251 1225 { 1252 FLM[j] = leadmonom( NEW[j] ); // we are interested in leading monomials only! 1253 } 1254 1255 kill NEW; 1256 1257 while( size(FLM) > 0 ) 1258 { 1259 // minimal degree: 1260 d = mindeg(FLM); 1226 // minimal degree: 1227 d = mindeg(FLM); 1261 1228 1262 1263 idealT = jet( FLM, d );1229 // take all of minimal degree: 1230 T = jet( FLM, d ); 1264 1231 1265 1232 // there are size(T) elements of smallest degree (deg(FLM[1])) in FLM! 1266 1233 1267 1268 1234 // Add them in the same way: 1235 for( j = 1; j <= (D + Delta); j = j + d ) // For every degree (j*d) of T, do: 1269 1236 { 1270 1237 for( l = j; (l + d) <= (D + Delta); l++ ) 1271 1238 { 1272 1239 if( (l + d) > D ) // Only new should be updated! 1273 1240 { 1274 1275 1241 FOUND_LEADING_MONOMIALS[l+d] = 1242 FOUND_LEADING_MONOMIALS[l+d] + FOUND_LEADING_MONOMIALS[l] * T; 1276 1243 } 1277 1244 } 1278 1245 } 1279 1246 1280 1281 1247 // Kill them from FLM: 1248 if( size(T) < size(FLM) ) 1282 1249 { 1283 1250 FLM = FLM[ (size(T)+1) .. size(FLM) ]; 1284 1251 } else 1285 { 1286 FLM = ideal(0); // break; 1287 } 1288 1289 kill T; 1252 { 1253 FLM = 0; // break; 1254 } 1255 1290 1256 } 1291 1257 1292 1293 1294 1295 1296 1258 // Go back... 1259 setring NCRING; 1260 1261 // And set new Bound 1262 D = D + Delta; 1297 1263 } 1298 1264 1299 i++; 1300 } 1301 1302 kill COMMRING; 1303 1304 result = makeNice(result); 1305 1306 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centralizerRed", result ); }; /*4DEBUG*/ 1307 1308 return( result ); 1265 i++; 1266 } 1267 1268 result = makeNice(result); 1269 1270 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centralizerRed", result ); }; /*4DEBUG*/ 1271 1272 return( result ); 1309 1273 } 1310 1274 example 1311 1275 { 1312 "EXAMPLE:"; echo = 2;1313 ring A = 0,(x,y,z),dp;1314 matrix D[3][3]=0;1315 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;1316 ncalgebra(1,D); // this algebra is U(sl_2)1317 ideal F = x, y;1318 // find subalgebra generators degree <= 4 of an algebra of1319 // all elements commuting with x and y:1320 ideal C = centralizerRed(F, 4);1321 C;1322 inCentralizer(C, F);1276 "EXAMPLE:"; echo = 2; 1277 ring A = 0,(x,y,z),dp; 1278 matrix D[3][3]=0; 1279 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 1280 ncalgebra(1,D); // this algebra is U(sl_2) 1281 ideal F = x, y; 1282 // find subalgebra generators degree <= 4 of an algebra of 1283 // all elements commuting with x and y: 1284 ideal C = centralizerRed(F, 4); 1285 C; 1286 inCentralizer(C, F); 1323 1287 } 1324 1288 … … 1326 1290 /******************************************************/ 1327 1291 proc centerRed( int D, list # ) 1328 "USAGE: centerRed( D[, k] ); int D[, int k]1292 "USAGE: centerRed( D[, N] ); int D[, int N] 1329 1293 RETURN: ideal, generated by computed generators 1330 PURPOSE: if N is absent and D >= 0 computes a subalgebra generators of the center up to degree D, otherwise if N is present computes N(at least) first generators of the center, if moreover D > 0 it will be used as the first maximal degree estimation. 1331 NOTE: Current ordering must be a degree compatible wellordering. 1294 PURPOSE: computes a subalgebra generators of the center up to degree D. 1295 NOTE: In general, one cannot compute the whole center. 1296 @* Hence, one has to specify a termination condition via arguments D and/or N. 1297 @* If D is positive, only central elements up to degree D will be found. 1298 @* If D is negative, the termination is determined by N only. 1299 @* If N is given, the computation stops if at least N elements has been found. 1300 @* Warning: if N is given and bigger than the actual number of generators, 1301 @* the procedure may not terminate. 1302 @* Current ordering must be a degree compatible wellordering. 1332 1303 SEE ALSO: centralizerRed; centerVS; center; inCenter 1333 1304 EXAMPLE: example centerRed; shows an example 1334 1305 " 1335 1306 { 1336 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centerRed", D); }; /*4DEBUG*/1337 1338 1339 { 1340 ERROR( "No current ring!" );1341 } 1342 1343 1307 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "centerRed", D, # ); }; /*4DEBUG*/ 1308 1309 if( nameof( basering ) == "basering" ) 1310 { 1311 // ERROR( "No current ring!" ); 1312 } 1313 1314 ideal result = centralizerRed( variablesSorted(), D, # ); 1344 1315 1345 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centerRed", result ); }; /*4DEBUG*/1346 1347 1316 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "centerRed", result ); }; /*4DEBUG*/ 1317 1318 return( result ); 1348 1319 } 1349 1320 example 1350 1321 { 1351 "EXAMPLE:"; echo = 2;1352 ring A = 0,(x,y,z),dp;1353 matrix D[3][3]=0;1354 D[1,2]=z;1355 ncalgebra(1,D); // it is a Heisenberg algebra1356 // find vector space basis of center of degree <= 31357 ideal VSZ = centerVS(3);1358 // There should be 3 degrees of z.1359 VSZ;1360 inCenter(VSZ);1361 // find "minimal" central elements of degree <= 31362 ideal SAZ = centerRed(3);1363 // Only 'z' must be computed1364 SAZ;1365 inCenter(SAZ);1322 "EXAMPLE:"; echo = 2; 1323 ring A = 0,(x,y,z),dp; 1324 matrix D[3][3]=0; 1325 D[1,2]=z; 1326 ncalgebra(1,D); // it is a Heisenberg algebra 1327 // find vector space basis of center of degree <= 3 1328 ideal VSZ = centerVS(3); 1329 // There should be 3 degrees of z. 1330 VSZ; 1331 inCenter(VSZ); 1332 // find "minimal" central elements of degree <= 3 1333 ideal SAZ = centerRed(3); 1334 // Only 'z' must be computed 1335 SAZ; 1336 inCenter(SAZ); 1366 1337 } 1367 1338 … … 1375 1346 /******************************************************/ 1376 1347 static proc INTERRED( ideal S ) 1377 "USAGE: INTERRED( S ); ideal S1348 "USAGE: INTERRED( S ); ideal S 1378 1349 RETURN: ideal, interreduced S 1379 PURPOSE: interreduction without monomial multiplication, just make every leading monomial occur in a single polynomial 1380 " 1381 { 1382 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "INTERRED", S ); }; /*4DEBUG*/ 1383 1384 ideal result = S; 1385 1386 int flag = 1; 1387 1388 int i, j, N; 1389 1390 poly p, lm; 1391 1392 while( flag == 1 ) 1350 PURPOSE: interreduction without monomial multiplication, 1351 just make every leading monomial occur in a single polynomial 1352 " 1353 { 1354 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "INTERRED", S ); }; /*4DEBUG*/ 1355 1356 ideal result = S; 1357 1358 int flag = 1; 1359 1360 int i, j, N; 1361 1362 poly p, lm; 1363 1364 while( flag == 1 ) 1393 1365 { 1394 1395 1396 1397 1398 1399 1400 1401 1402 1403 1404 1405 1406 { 1407 1408 1409 1410 1411 1366 flag = 0; 1367 1368 result = sort( simplify( result, 1 + 2 + 8) )[1]; 1369 // sorting w.r.t. actual monomial ordering 1370 // generators with SMALLER(!) leading term come FIRST 1371 1372 N = size(result); 1373 1374 // kill leading monomials: 1375 1376 i = 1; 1377 while( i < N ) 1378 { 1379 p = result[i]; 1380 lm = leadmonom(p); 1381 1382 j = i + 1; 1383 while( leadmonom(result[j]) == lm ) 1412 1384 { 1413 1414 1415 1385 result[j] = result[j]  p; // leadcoefs are 1 because of simplify. 1386 flag = 1; // we have changed something => we do still need to care about it... 1387 j++; 1416 1388 1417 1389 if( j > N ) 1418 1390 { 1419 1391 break; 1420 1392 } 1421 1393 } 1422 1394 1423 1395 i = j; 1424 1396 } 1425 1397 } 1426 1398 1427 1428 1429 1430 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "INTERRED", result ); }; /*4DEBUG*/1431 1432 1399 // We are done! No common leading monomials! 1400 // The result is sorted 1401 1402 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "INTERRED", result ); }; /*4DEBUG*/ 1403 1404 return( result ); 1433 1405 } 1434 1406 … … 1436 1408 /******************************************************/ 1437 1409 static proc SANF( poly p, list FOUND_LEADING_MONOMIALS ) 1438 "1410 " 1439 1411 reduce p wrt found multiples without ANY polynomial multiplications! 1440 1412 " 1441 1413 { 1442 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "SANF", p, FOUND_LEADING_MONOMIALS); }; /*4DEBUG*/1443 1444 1445 1446 1447 1448 1449 1450 { 1451 1452 1453 1454 { 1455 1414 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "SANF", p, FOUND_LEADING_MONOMIALS); }; /*4DEBUG*/ 1415 1416 poly q = p; 1417 poly head = 0; 1418 1419 int d; int N = size(FOUND_LEADING_MONOMIALS); 1420 1421 while( size(q) > 0 ) 1422 { 1423 d = maxdeg(p); 1424 1425 if( (0 < d) and (d <= N) ) 1426 { 1427 if( size(FOUND_LEADING_MONOMIALS[d]) > 0 ) 1456 1428 { 1457 1458 1429 attrib( FOUND_LEADING_MONOMIALS[d], "isSB", 1); 1430 q = reduce( p, FOUND_LEADING_MONOMIALS[d] ); 1459 1431 } 1460 1432 1461 1433 DBPrint(1, string(p) + " > " + string(q) ); 1462 1434 } 1463 1435 1464 1465 { 1466 1467 1468 1436 if( q == p ) 1437 { 1438 p = lead(q); 1439 1440 if( d > 0 ) 1469 1441 { 1470 1471 1442 // No scalars! 1443 head = head + p; 1472 1444 } 1473 1445 1474 1446 q = q  p; 1475 1447 } 1476 1448 1477 1478 } 1479 1480 1481 1482 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "SANF", head ); }; /*4DEBUG*/1483 1484 1449 p = q; 1450 } 1451 1452 1453 1454 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "SANF", head ); }; /*4DEBUG*/ 1455 1456 return( head ); 1485 1457 } 1486 1458 … … 1489 1461 static proc maxdegInt( ideal I ) 1490 1462 { 1491 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "maxdegInt", I ); }; /*4DEBUG*/1492 1493 1494 1495 1496 1497 1498 { 1499 1500 1501 1502 { 1503 1463 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "maxdegInt", I ); }; /*4DEBUG*/ 1464 1465 intmat D = maxdeg(I); 1466 1467 int max = D[1, 1]; int m; 1468 1469 for( int c = 2; c <= ncols(D); c++ ) 1470 { 1471 m = D[1, c]; 1472 1473 if( m > max ) 1474 { 1475 max = m; 1504 1476 } 1505 1477 } 1506 1478 1507 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "maxdegInt", max ); }; /*4DEBUG*/1508 1509 1479 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "maxdegInt", max ); }; /*4DEBUG*/ 1480 1481 return( max ); 1510 1482 } 1511 1483 … … 1514 1486 "USAGE: sa_reduce(V); ideal V 1515 1487 RETURN: ideal, generated by found elements 1516 PURPOSE: compute a s ubalgebra basis of an algebra generated by polynomial fromV1488 PURPOSE: compute a s.a. basis of an algebra generated by V 1517 1489 NOTE: May produce wrong result in quotient algebras. 1518 1490 SEE ALSO: sa_poly_reduce … … 1520 1492 " 1521 1493 { 1522 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "sa_reduce", V ); }; /*4DEBUG*/1523 1524 1525 1526 1527 1528 1529 1530 1531 1532 1533 1534 1535 1536 1537 { 1538 1494 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "sa_reduce", V ); }; /*4DEBUG*/ 1495 1496 ideal result = ideal(); 1497 1498 ideal FLM = INTERRED( V ); // The output is sorted "[1]<[2]<[3]<..." 1499 1500 // We are bounded by maximal degree!!! 1501 int D = maxdegInt( FLM ); 1502 1503 // Init 1504 list FOUND_LEADING_MONOMIALS = list(); 1505 1506 int i; 1507 1508 for( i = 1; i <= D; i++ ) 1509 { 1510 FOUND_LEADING_MONOMIALS[i] = ideal(); 1539 1511 } 1540 1512 1541 1542 1543 1544 1545 1546 1547 1548 1549 { 1550 1551 1552 1553 { 1554 1555 1556 1513 int d, j, l; 1514 1515 poly p, q; ideal T; 1516 1517 1518 int c = 1; // polynomials in FLM commute pairwise 1519 1520 for( j = 1; (j < size(FLM)) and (c == 1); j++ ) 1521 { 1522 p = FLM[j]; 1523 1524 for( l = j+1; (l <= size(FLM)) and (c == 1); l++ ) 1525 { 1526 q = FLM[l]; 1527 1528 if( NF(p*q  q*p, twostd(0)) != 0 ) 1557 1529 { 1558 1530 c = 0; // There exists noncommuting pair 1559 1531 } 1560 1532 } 1561 1533 } 1562 1534 1563 1564 { 1565 1566 1567 1568 1569 1570 1571 1572 1573 1574 1575 1535 while( size(FLM) > 0 ) 1536 { 1537 // Take the 1st element of FLM... 1538 p = FLM[1]; // SANF( FLM[1], FOUND_LEADING_MONOMIALS ); 1539 1540 FLM[1] = 0; // ...and kill it from FLM 1541 1542 d = maxdeg( p ); 1543 T = ideal(p); 1544 1545 FOUND_LEADING_MONOMIALS[d] = FOUND_LEADING_MONOMIALS[d] + T; 1546 1547 for( j = 1; j <= D; j = j + d ) // For every degree (j*d) of T, do: 1576 1548 { 1577 1549 for( l = j; (l + d) <= D; l++ ) 1578 1550 { 1579 1580 1551 FOUND_LEADING_MONOMIALS[l+d] = 1552 FOUND_LEADING_MONOMIALS[l+d] + FOUND_LEADING_MONOMIALS[l] * T; 1581 1553 1582 1554 if( c != 1 ) 1583 1555 { 1584 1585 1556 FOUND_LEADING_MONOMIALS[l+d] = 1557 FOUND_LEADING_MONOMIALS[l+d] + T * FOUND_LEADING_MONOMIALS[l]; 1586 1558 } 1587 1559 } 1588 1560 } 1589 1561 1590 1591 { 1592 1562 if( size(FLM) > 0 ) 1563 { 1564 for( i = 2; i <= ncols(FLM); i++ ) 1593 1565 { 1594 1566 FLM[i] = SANF( FLM[i], FOUND_LEADING_MONOMIALS ); 1595 1567 } 1596 1568 FLM = INTERRED( FLM ); 1597 1569 } 1598 1570 1599 1600 1601 1602 1603 } 1604 1605 1606 1607 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "sa_reduce", result ); }; /*4DEBUG*/1608 1609 1571 DBPrint(1, "Found: " + string(T) ); 1572 1573 result = result + T; 1574 1575 } 1576 1577 result = makeNice(result); 1578 1579 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "sa_reduce", result ); }; /*4DEBUG*/ 1580 1581 return( result ); 1610 1582 } 1611 1583 example … … 1626 1598 "USAGE: sa_poly_reduce(p, V); poly p, ideal V 1627 1599 RETURN: polynomial, a reduction of p wrt V 1628 PURPOSE: computes a reduction of a given polynomial p wrt a set of polynomialsV1600 PURPOSE: computes a reduction of polynomial p wrt V 1629 1601 NOTE: May produce wrong result in quotient algebras. 1630 1602 SEE ALSO: sa_reduce … … 1632 1604 " 1633 1605 { 1634 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "sa_poly_reduce", p, V ); }; /*4DEBUG*/1635 1636 1637 1638 1639 1640 1641 1642 1643 1644 1645 1646 1647 1648 { 1649 1606 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "sa_poly_reduce", p, V ); }; /*4DEBUG*/ 1607 // As previous... 1608 1609 ideal FLM = INTERRED( V ); // The output is sorted "[1]<[2]<[3]<..." 1610 1611 // We are bounded by maximal degree!!! 1612 int D = maxdegInt( FLM + ideal(p) ); 1613 1614 // Init 1615 list FOUND_LEADING_MONOMIALS = list(); 1616 1617 int i; 1618 1619 for( i = 1; i <= D; i++ ) 1620 { 1621 FOUND_LEADING_MONOMIALS[i] = ideal(); 1650 1622 } 1651 1623 1652 1653 1654 1655 1656 1657 1658 1659 1660 { 1661 1662 1663 1664 { 1665 1666 1667 1624 int d, j, l; 1625 1626 poly f, q; ideal T; 1627 1628 1629 int c = 1; // polynomials in FLM commute pairwise 1630 1631 for( j = 1; (j < size(FLM)) and (c == 1); j++ ) 1632 { 1633 f = FLM[j]; 1634 1635 for( l = j+1; (l <= size(FLM)) and (c == 1); l++ ) 1636 { 1637 q = FLM[l]; 1638 1639 if( NF(f*q  q*f, twostd(0)) != 0 ) 1668 1640 { 1669 1641 c = 0; 1670 1642 } 1671 1643 } … … 1673 1645 1674 1646 1675 1676 { 1677 1678 1679 1680 1681 1682 1683 1684 1685 1686 1687 1647 while( size(FLM) > 0 ) 1648 { 1649 // Take the 1st element of FLM... 1650 q = SANF( FLM[1], FOUND_LEADING_MONOMIALS ); 1651 1652 FLM[1] = 0; // ...and kill it from FLM 1653 1654 d = maxdeg(q); 1655 T = ideal(q); 1656 1657 FOUND_LEADING_MONOMIALS[d] = FOUND_LEADING_MONOMIALS[d] + T; 1658 1659 for( j = 1; j <= D; j = j + d ) // For every degree (j*d) of T, do: 1688 1660 { 1689 1661 for( l = j; (l + d) <= D; l++ ) 1690 1662 { 1691 1692 1663 FOUND_LEADING_MONOMIALS[l+d] = 1664 FOUND_LEADING_MONOMIALS[l+d] + FOUND_LEADING_MONOMIALS[l] * T; 1693 1665 1694 1666 if( c != 1 ) 1695 1667 { 1696 1697 1668 FOUND_LEADING_MONOMIALS[l+d] = 1669 FOUND_LEADING_MONOMIALS[l+d] + T * FOUND_LEADING_MONOMIALS[l]; 1698 1670 } 1699 1671 } 1700 1672 } 1701 1673 1702 1703 { 1704 1674 if( size(FLM) > 0 ) 1675 { 1676 for( i = 2; i <= ncols(FLM); i++ ) 1705 1677 { 1706 1678 FLM[i] = SANF( FLM[i], FOUND_LEADING_MONOMIALS ); 1707 1679 } 1708 1680 FLM = INTERRED( FLM ); 1709 1681 } 1710 1682 } 1711 1683 1712 1713 1714 1715 1716 1717 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "sa_poly_reduce", result ); }; /*4DEBUG*/1718 1719 1684 poly result = SANF(p, FOUND_LEADING_MONOMIALS); 1685 1686 result = makeNice( result ); 1687 1688 1689 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "sa_poly_reduce", result ); }; /*4DEBUG*/ 1690 1691 return( result ); 1720 1692 } 1721 1693 example … … 1742 1714 /******************************************************/ 1743 1715 static proc inCentralizer_poly( poly p, ideal S ) 1744 "1716 " 1745 1717 if p in centralizer(S) => return 1, otherwise return 0 1746 1718 " 1747 1719 { 1748 1749 1750 1751 { 1752 1753 1754 1755 { 1756 1757 1758 1720 poly f; 1721 1722 for( int k = 1; k <= size(S); k++ ) 1723 { 1724 f = S[k]; 1725 1726 if( NF( f * p  p * f, twostd(0) ) != 0 ) 1727 { 1728 DBPrint( 1, "POLY: " + string (p) + 1729 " is NOT in the centralizer of poly {" + string(f) + "}" ); 1730 return (0); 1759 1731 } 1760 1732 } 1761 1733 1762 1734 return( 1 ); 1763 1735 } 1764 1736 … … 1766 1738 static proc inCentralizer_list( def l, ideal S ) 1767 1739 { 1768 1769 { 1770 1771 { 1772 1740 for( int @i = 1; @i <= size(l); @i++ ) 1741 { 1742 if( (typeof(l[@i])=="poly") or (typeof(l[@i]) == "int") or (typeof(l[@i]) == "number") ) 1743 { 1744 if(! inCentralizer_poly(l[@i], S) ) 1773 1745 { 1774 1746 return(0); 1775 1747 } 1776 1748 1777 1749 } else 1778 1750 { 1779 1751 if( (typeof(l[@i])=="list") or (typeof(l[@i])=="ideal") ) 1780 1752 { 1781 1753 if(! inCentralizer_list(l[@i], S) ) 1782 1754 { 1783 1755 return(0); 1784 1785 1786 1787 } 1788 1756 } 1757 } 1758 } 1759 } 1760 return(1); 1789 1761 } 1790 1762 … … 1794 1766 "USAGE: inCentralizer(a, S); a poly/list/ideal, S poly/ideal 1795 1767 RETURN: integer, 1 if a in the centralizer(S), 0 otherwise 1796 PURPOSE: check whether a given elementis centralizing with respect to elements of S1768 PURPOSE: check whether a is centralizing with respect to elements of S 1797 1769 EXAMPLE: example inCentralizer; shows examples 1798 1770 " 1799 1771 { 1800 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "inCentralizer", a, S ); }; /*4DEBUG*/1801 1802 1803 { 1804 ERROR( "No current ring!" );1805 } 1806 1807 1808 1809 1810 1811 { 1812 1772 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "inCentralizer", a, S ); }; /*4DEBUG*/ 1773 1774 if( nameof( basering ) == "basering" ) 1775 { 1776 // ERROR( "No current ring!" ); 1777 } 1778 1779 1780 int res; 1781 1782 if( (typeof(a) == "poly") or (typeof(a) == "int") or (typeof(a) == "number") ) 1783 { 1784 res = inCentralizer_poly(a, S); 1813 1785 } else 1814 {1786 { 1815 1787 if( (typeof(a)=="list") or (typeof(a)=="ideal") ) 1816 1788 { 1817 1789 res = inCentralizer_list(a, S); 1818 1819 1820 1821 1822 }1823 1824 1825 { 1826 1827 } 1828 1829 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "inCentralizer", res ); }; /*4DEBUG*/1830 1831 1790 } else 1791 { 1792 res = 1; 1793 } 1794 } 1795 1796 if( res == 1 ) 1797 { 1798 ERROR( "Wrong argument!" ); 1799 } 1800 1801 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "inCentralizer", res ); }; /*4DEBUG*/ 1802 1803 return (res); 1832 1804 } 1833 1805 example 1834 1806 { 1835 "EXAMPLE:";echo=2;1836 ring r=0,(x,y,z),dp;1837 matrix D[3][3]=0;1838 D[1,2]=z;1839 ncalgebra(1,D); // the Heisenberg algebra1840 poly f = x^2;1841 poly a = z; // we know this element if central 1842 poly b = y^2;1843 inCentralizer(a, f);1844 inCentralizer(b, f);1845 list l = list(1, a);1846 inCentralizer(l, f);1847 ideal I = a, b;1848 inCentralizer(I, f);1849 printlevel = 2;1850 inCentralizer(a, f); // yes1851 inCentralizer(b, f); // no1807 "EXAMPLE:";echo=2; 1808 ring r=0,(x,y,z),dp; 1809 matrix D[3][3]=0; 1810 D[1,2]=z; 1811 ncalgebra(1,D); // the Heisenberg algebra 1812 poly f = x^2; 1813 poly a = z; // 'z' is central => it lies in any centralizer! 1814 poly b = y^2; 1815 inCentralizer(a, f); 1816 inCentralizer(b, f); 1817 list l = list(1, a); 1818 inCentralizer(l, f); 1819 ideal I = a, b; 1820 inCentralizer(I, f); 1821 printlevel = 2; 1822 inCentralizer(a, f); // yes 1823 inCentralizer(b, f); // no 1852 1824 } 1853 1825 1854 1826 /******************************************************/ 1855 1827 proc inCenter( def a ) // Checks the centrality of a 1856 "USAGE: inCenter(a); a poly/list/ideal1828 "USAGE: inCenter(a); a poly/list/ideal 1857 1829 RETURN: integer, 1 if a in the center, 0 otherwise 1858 1830 PURPOSE: check whether a given element is central … … 1860 1832 " 1861 1833 { 1862 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "inCenter", a ); }; /*4DEBUG*/1863 1864 1865 { 1866 ERROR( "No current ring!" );1867 } 1868 1869 1870 1871 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "inCenter", result ); }; /*4DEBUG*/1872 1873 1834 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "inCenter", a ); }; /*4DEBUG*/ 1835 1836 if( nameof( basering ) == "basering" ) 1837 { 1838 // ERROR( "No current ring!" ); 1839 } 1840 1841 int result = inCentralizer( a, variablesStandard() ); 1842 1843 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "inCenter", result ); }; /*4DEBUG*/ 1844 1845 return( result ); 1874 1846 } 1875 1847 example 1876 1848 { 1877 "EXAMPLE:";echo=2;1878 ring r=0,(x,y,z),dp;1879 matrix D[3][3]=0;1880 D[1,2]=z;1881 D[1,3]=2*x;1882 D[2,3]=2*y;1883 ncalgebra(1,D); // this is U(sl_2)1884 poly p=4*x*y+z^22*z;1885 inCenter(p);1886 poly f=4*x*y;1887 inCenter(f);1888 list l= list( 1, p, p^2, p^3);1889 inCenter(l);1890 ideal I= p, f;1891 inCenter(I);1849 "EXAMPLE:";echo=2; 1850 ring r=0,(x,y,z),dp; 1851 matrix D[3][3]=0; 1852 D[1,2]=z; 1853 D[1,3]=2*x; 1854 D[2,3]=2*y; 1855 ncalgebra(1,D); // this is U(sl_2) 1856 poly p=4*x*y+z^22*z; 1857 inCenter(p); 1858 poly f=4*x*y; 1859 inCenter(f); 1860 list l= list( 1, p, p^2, p^3); 1861 inCenter(l); 1862 ideal I= p, f; 1863 inCenter(I); 1892 1864 } 1893 1865 … … 1895 1867 /******************************************************/ 1896 1868 proc isCartan( poly f ) // Checks whether f is a Cartan element. 1897 " 1898 PURPOSE: check whether f is a Cartan element. 1899 RETURN: 1 if f is a Cartan element and 0 otherwise. 1900 NOTE: f is a Cartan element iff for all g in A there exists C in K such that [f, g] = C * g 1901 iff for all variables v_i of A there exist C in K such that [f, v_i] = C * v_i. 1902 " 1903 { 1904 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "isCartan", f ); }; /*4DEBUG*/ 1905 1906 if( nameof( basering ) == "basering" ) 1907 { 1908 ERROR( "No current ring!" ); 1909 } 1910 1911 1912 ideal V = variablesStandard(); 1913 1914 int r = 1; poly v, g; 1915 1916 for( int i = size(V); i > 0; i ) 1917 { 1918 v = leadmonom(V[i]); // V[i] must be just a variable, but... 1919 1920 g = NF( f*v  v*f, twostd(0) ); // [f, V[i]] 1921 1922 if( size(g) > 0 ) 1923 { 1924 if( size(g) > 1 ) // it is not just \alpha * v_i. 1869 "USAGE: isCartan(f); poly f 1870 PURPOSE: check whether f is a Cartan element. 1871 RETURN: integer, 1 if f is a Cartan element and 0 otherwise. 1872 NOTE: f is a Cartan element 1873 @* iff for all g in A there exists C in K such that [f, g] = C * g 1874 @* iff for all variables v_i there exist C in K such that [f, v_i] = C * v_i. 1875 " 1876 { 1877 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ BCall( "isCartan", f ); }; /*4DEBUG*/ 1878 1879 if( nameof( basering ) == "basering" ) 1880 { 1881 // ERROR( "No current ring!" ); 1882 } 1883 1884 1885 ideal V = variablesStandard(); 1886 1887 int r = 1; poly v, g; 1888 1889 for( int i = size(V); i > 0; i ) 1890 { 1891 v = leadmonom(V[i]); // V[i] must be just a variable, but... 1892 1893 g = NF( f*v  v*f, twostd(0) ); // [f, V[i]] 1894 1895 if( size(g) > 0 ) 1896 { 1897 if( size(g) > 1 ) // it is not just \alpha * v_i. 1925 1898 { 1926 1927 1899 r = 0; 1900 break; 1928 1901 } 1929 1902 1930 1903 if( leadmonom(g) != v ) // g = \alpha * v_j, j != i. 1931 1904 { 1932 1933 1905 r = 0; 1906 break; 1934 1907 } 1935 1908 … … 1937 1910 } 1938 1911 1939 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "isCartan", r ); }; /*4DEBUG*/1940 1912 /*4DEBUG*/ if( defined( @@@DEBUG ) ){ ECall( "isCartan", r ); }; /*4DEBUG*/ 1913 return( r ); 1941 1914 } 1942 1915 example 1943 1916 { 1944 "EXAMPLE:";echo=2;1945 ring r=0,(x,y,z),dp;1946 matrix D[3][3]=0;1947 D[1,2]=z;1948 D[1,3]=2*x;1949 D[2,3]=2*y;1950 ncalgebra(1,D); // this is U(sl_2) with cartan  z1951 isCartan(z); // yes!1952 poly p=4*x*y+z^22*z;1953 isCartan(p); // central elements are Cartan elements!1954 poly f=4*x*y;1955 isCartan(f); // no way!1956 isCartan( 10 + p + z ); // scalar + central + cartan1957 } 1958 1959 1960 1961 1962 /******************************************************/ 1963 /******************************************************/ 1964 // ::MainAliases:: The only nonstatic functions, visible to user are here. And they are high level wrappers around basic static functions.1917 "EXAMPLE:";echo=2; 1918 ring r=0,(x,y,z),dp; 1919 matrix D[3][3]=0; 1920 D[1,2]=z; 1921 D[1,3]=2*x; 1922 D[2,3]=2*y; 1923 ncalgebra(1,D); // this is U(sl_2) with cartan  z 1924 isCartan(z); // yes! 1925 poly p=4*x*y+z^22*z; 1926 isCartan(p); // central elements are Cartan elements! 1927 poly f=4*x*y; 1928 isCartan(f); // no way! 1929 isCartan( 10 + p + z ); // scalar + central + cartan 1930 } 1931 1932 1933 1934 1935 /******************************************************/ 1936 /******************************************************/ 1937 // ::MainAliases:: The main nonstatic functions, visible to user are here. They are wrappers around basic functions. 1965 1938 /******************************************************/ 1966 1939 /******************************************************/ … … 1973 1946 "USAGE: center(D[, N]); int D, int N 1974 1947 RETURN: ideal, generated by elements of degree at most D 1975 PURPOSE: computes a minimal set of central elements up to degree D. 1976 NOTE: In general, one cannot predict the number or the highest degree of 1977 central elements. Hence, one has to specify a termination condition via arguments D and/or N. 1978 @* If D is positive, the computation stops after all central elements of degree at most D has been found. 1979 @* If D is negative, the termination is determined by N only. 1980 @* If N is given, the computation stops if at least N central elements has been found. 1981 @* Warning: if N is given and bigger than the actual number of generators, the procedure may not terminate. 1948 PURPOSE: computes a subalgebra generators of the center up to degree D. 1949 NOTE: In general, one cannot compute the whole center. 1950 @* Hence, one has to specify a termination condition via arguments D and/or N. 1951 @* If D is positive, only central elements up to degree D will be found. 1952 @* If D is negative, the termination is determined by N only. 1953 @* If N is given, the computation stops if at least N elements has been found. 1954 @* Warning: if N is given and bigger than the actual number of generators, 1955 @* the procedure may not terminate. 1956 @* Current ordering must be a degree compatible wellordering. 1982 1957 SEE ALSO: centralizer; inCenter 1983 1958 EXAMPLE: example center; shows an example 1984 1959 " 1985 1960 { 1986 1987 { 1988 ERROR( "No current ring!" );1989 } 1990 1991 1992 { 1993 1994 } 1995 1996 1997 { 1998 1999 } 2000 2001 1961 if( nameof( basering ) == "basering" ) 1962 { 1963 // ERROR( "No current ring!" ); 1964 } 1965 1966 if( DefaultInt( # ) > 0 ) 1967 { 1968 return( centerRed( D, # ) ); 1969 } 1970 1971 if( D >= 0 ) 1972 { 1973 return( sa_reduce( centerVS(D) ) ); // Experimental! May be wrong!!! 1974 } 1975 1976 ERROR( "Wrong arguments!" ); 2002 1977 } 2003 1978 example 2004 1979 { 2005 "EXAMPLE:"; echo = 2;2006 ring A = 0,(x,y,z,t),dp;2007 matrix D[4][4]=0;2008 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;2009 ncalgebra(1,D); // this algebra is U(sl_2) tensored with K[t]2010 ideal Z = center(3); // find all central elements of degree <= 32011 Z;2012 inCenter(Z);2013 ideal ZZ = center(1, 1); // find one central element of the lowest degree2014 ZZ;2015 inCenter(ZZ);1980 "EXAMPLE:"; echo = 2; 1981 ring A = 0,(x,y,z,t),dp; 1982 matrix D[4][4]=0; 1983 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 1984 ncalgebra(1,D); // this algebra is U(sl_2) tensored with K[t] 1985 ideal Z = center(3); // find all central elements of degree <= 3 1986 Z; 1987 inCenter(Z); 1988 ideal ZZ = center(1, 1); // find one central element of the lowest degree 1989 ZZ; 1990 inCenter(ZZ); 2016 1991 } 2017 1992 2018 1993 /******************************************************/ 2019 1994 proc centralizer( ideal S, int D, list # ) // Computes the generators of the centralizer of S in a basering 2020 "USAGE: centralizer(S, MaxDeg[, N]); poly/ideal S, int MaxDeg, int N 2021 RETURN: ideal, generated by elements of degree <= MaxDeg 2022 PURPOSE: computes a minimal set of elements centralizer(S) up to degree MaxDeg. 2023 NOTE: In general, one cannot predict the number or the highest degree of 2024 centralizing elements. Hence, one has to specify a termination condition via arguments MaxDeg and/or N. 2025 @* If MaxDeg is positive, the computation stops after all centralizing elements of degree at most MaxDeg has been found. 2026 @* If MaxDeg is negative, the termination is determined by N only. 2027 @* If N is given, the computation stops if at least N centralizing elements has been found. 2028 @* Warning: if N is given and bigger than the actual number of generators, the procedure may not terminate. 1995 "USAGE: centralizer(F, D[, N]); poly/ideal F, int D[, int N] 1996 RETURN: ideal, generated by computed generators 1997 PURPOSE: computes a subalgebra generators of centralizer(F) up to degree D. 1998 NOTE: In general, one cannot compute the whole centralizer(F). 1999 @* Hence, one has to specify a termination condition via arguments D and/or N. 2000 @* If D is positive, only centralizing elements up to degree D will be found. 2001 @* If D is negative, the termination is determined by N only. 2002 @* If N is given, the computation stops if at least N elements has been found. 2003 @* Warning: if N is given and bigger than the actual number of generators, 2004 @* the procedure may not terminate. 2005 @* Current ordering must be a degree compatible wellordering. 2029 2006 SEE ALSO: center; inCentralizer 2030 2007 EXAMPLE: example centralizer; shows an example 2031 2008 " 2032 2009 { 2033 2034 { 2035 ERROR( "No current ring!" );2036 } 2037 2038 2039 { 2040 2041 } 2042 2043 2044 { 2045 2046 } 2047 2048 2010 if( nameof( basering ) == "basering" ) 2011 { 2012 // ERROR( "No current ring!" ); 2013 } 2014 2015 if( DefaultInt( # ) > 0 ) 2016 { 2017 return( centralizerRed( S, D, # ) ); 2018 } 2019 2020 if( D >= 0 ) 2021 { 2022 return( sa_reduce( centralizerVS(S, D) ) ); // Experimental! May be wrong!!! 2023 } 2024 2025 ERROR( "Wrong arguments!" ); 2049 2026 } 2050 2027 example 2051 2028 { 2052 "EXAMPLE:"; echo = 2;2053 ring A = 0,(x,y,z),dp;2054 matrix D[3][3]=0;2055 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y;2056 ncalgebra(1,D); // this algebra is U(sl_2)2057 poly f = 4*x*y+z^22*z; // a central polynomial2058 f;2059 ideal c = centralizer(f, 2); // find all elements of the centralizer of f2060 2061 c; // since f is central, the answer consists of generators of A2062 inCentralizer(c, f);2063 ideal cc = centralizer(f,1,2); // find at least two elements of the centralizer of f2064 cc;2065 inCentralizer(cc, f);2066 poly g = z^22*z; // some noncentral polynomial2067 c = centralizer(g, 2); // find all elements of the centralizer of g2068 2069 c;2070 inCentralizer(c, g);2071 centralizer(g,1,1); // find the element of the lowest degree in the centralizer2072 cc = centralizer(g,1,2); // find at least two elements of the centralizer of g2073 cc;2074 inCentralizer(cc, g);2029 "EXAMPLE:"; echo = 2; 2030 ring A = 0,(x,y,z),dp; 2031 matrix D[3][3]=0; 2032 D[1,2]=z; D[1,3]=2*x; D[2,3]=2*y; 2033 ncalgebra(1,D); // this algebra is U(sl_2) 2034 poly f = 4*x*y+z^22*z; // a central polynomial 2035 f; 2036 ideal c = centralizer(f, 2); // find all elements of the centralizer of f 2037 // of degree <= 2 2038 c; // since f is central, the answer consists of generators of A 2039 inCentralizer(c, f); 2040 ideal cc = centralizer(f,1,2); // find at least two elements of the centralizer of f 2041 cc; 2042 inCentralizer(cc, f); 2043 poly g = z^22*z; // some noncentral polynomial 2044 c = centralizer(g, 2); // find all elements of the centralizer of g 2045 // of degree <= 2 2046 c; 2047 inCentralizer(c, g); 2048 centralizer(g,1,1); // find the element of the lowest degree in the centralizer 2049 cc = centralizer(g,1,2); // find at least two elements of the centralizer of g 2050 cc; 2051 inCentralizer(cc, g); 2075 2052 } 2076 2053 2077 2054 2078 2055 /******************************************************* 2079 // normally one should use this library together with ncalg.lib in the following way:2056 // normally one should use this library together with ncalg.lib in the following way: 2080 2057 2081 2058 LIB "ncalg.lib";
Note: See TracChangeset
for help on using the changeset viewer.