Annotation of OpenXM_contrib/PHC/Ada/Schubert/ts_org_pieri.adb, Revision 1.1
1.1 ! maekawa 1: with text_io,integer_io; use text_io,integer_io;
! 2: with Communications_with_User; use Communications_with_User;
! 3: with Timing_Package; use Timing_Package;
! 4: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
! 5: with Standard_Natural_Vectors; use Standard_Natural_Vectors;
! 6: with Standard_Natural_Vectors_io; use Standard_Natural_Vectors_io;
! 7: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
! 8: with Standard_Complex_Vectors;
! 9: with Standard_Complex_Vectors_io; use Standard_Complex_Vectors_io;
! 10: with Standard_Complex_Matrices;
! 11: with Standard_Complex_Matrices_io; use Standard_Complex_Matrices_io;
! 12: with Standard_Complex_VecMats; use Standard_Complex_VecMats;
! 13: with Standard_Complex_VecMats_io; use Standard_Complex_VecMats_io;
! 14: with Standard_Random_Matrices; use Standard_Random_Matrices;
! 15: with Standard_Complex_Poly_Matrices;
! 16: with Standard_Complex_Poly_Matrices_io; use Standard_Complex_Poly_Matrices_io;
! 17: with Symbol_Table; use Symbol_Table;
! 18: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
! 19: with Standard_Complex_Polynomials_io; use Standard_Complex_Polynomials_io;
! 20: with Standard_Complex_Poly_Systems; use Standard_Complex_Poly_Systems;
! 21: with Standard_Complex_Poly_Systems_io; use Standard_Complex_Poly_Systems_io;
! 22: with Standard_Complex_Poly_SysFun; use Standard_Complex_Poly_SysFun;
! 23: with Matrix_Indeterminates; use Matrix_Indeterminates;
! 24: with Drivers_for_Poly_Continuation; use Drivers_for_Poly_Continuation;
! 25: with Brackets,Brackets_io; use Brackets,Brackets_io;
! 26: with Bracket_Monomials; use Bracket_Monomials;
! 27: with Bracket_Monomials_io; use Bracket_Monomials_io;
! 28: with Bracket_Expansions; use Bracket_Expansions;
! 29: with Bracket_Polynomials; use Bracket_Polynomials;
! 30: with Bracket_Polynomials_io; use Bracket_Polynomials_io;
! 31: with Bracket_Systems; use Bracket_Systems;
! 32: with Pieri_Trees,Pieri_Trees_io; use Pieri_Trees,Pieri_Trees_io;
! 33: with Pieri_Root_Counts; use Pieri_Root_Counts;
! 34: with Symbolic_Minor_Equations; use Symbolic_Minor_Equations;
! 35: with Numeric_Minor_Equations; use Numeric_Minor_Equations;
! 36: with Solve_Pieri_Leaves;
! 37: with Specialization_of_Planes; use Specialization_of_Planes;
! 38: with Pieri_Deformations; use Pieri_Deformations;
! 39:
! 40: procedure ts_org_pieri is
! 41:
! 42: -- DESCRIPTION :
! 43: -- This program executes the original implementation of the Pieri homotopy
! 44: -- algorithm, implemented strictly along the lines of the description in
! 45: -- the paper "Numerical Schubert Calculus" of Birkett Huber, Frank Sottile
! 46: -- and Bernd Sturmfels.
! 47:
! 48: function Maximum ( n1,n2 : natural ) return natural is
! 49: begin
! 50: if n1 >= n2
! 51: then return n1;
! 52: else return n2;
! 53: end if;
! 54: end Maximum;
! 55:
! 56: procedure Add_t_Symbol is
! 57:
! 58: -- DESCRIPTION :
! 59: -- Adds the symbol for the continuation parameter t to the symbol table.
! 60:
! 61: tsb : Symbol;
! 62:
! 63: begin
! 64: Symbol_Table.Enlarge(1);
! 65: tsb(1) := 't';
! 66: for i in 2..tsb'last loop
! 67: tsb(i) := ' ';
! 68: end loop;
! 69: Symbol_Table.Add(tsb);
! 70: end Add_t_Symbol;
! 71:
! 72: -- DISPLAYING THE REPRESENTATIONS OF THE PLANES :
! 73:
! 74: procedure Display_Polynomial_Pattern
! 75: ( file : in file_type; n : in natural; b1,b2 : in Bracket ) is
! 76:
! 77: -- DESCRIPTION :
! 78: -- Displays the pattern by a polynomial matrix.
! 79:
! 80: pm : Standard_Complex_Poly_Matrices.Matrix(1..n,b1'range)
! 81: := Schubert_Pattern(n,b1,b2);
! 82:
! 83: begin
! 84: put(file,pm);
! 85: Standard_Complex_Poly_Matrices.Clear(pm);
! 86: end Display_Polynomial_Pattern;
! 87:
! 88: -- AUXILIARIES TO SET UP EQUATIONS :
! 89:
! 90: procedure Expand_Minors ( file : in file_type;
! 91: mat : in Standard_Complex_Poly_Matrices.Matrix;
! 92: bm : in Bracket_Monomial ) is
! 93:
! 94: -- DESCRIPTION :
! 95: -- Expands the quadratic bracket monomial.
! 96:
! 97: first : boolean := true;
! 98: lb : Link_to_Bracket;
! 99:
! 100: procedure Visit_Bracket ( b : in Bracket; continue : out boolean ) is
! 101:
! 102: p : Poly := Null_Poly;
! 103:
! 104: begin
! 105: if first
! 106: then lb := new Bracket'(b);
! 107: first := false;
! 108: else p := Expanded_Minor(mat,b);
! 109: if p /= Null_Poly
! 110: then put(file,lb.all);
! 111: put(file,"*("); put(file,p); put(file,")");
! 112: end if;
! 113: Clear(lb); Clear(p);
! 114: end if;
! 115: continue := true;
! 116: end Visit_Bracket;
! 117: procedure Visit_Brackets is new Enumerate_Brackets(Visit_Bracket);
! 118:
! 119: begin
! 120: Visit_Brackets(bm);
! 121: end Expand_Minors;
! 122:
! 123: procedure Expand_Minors ( file : in file_type;
! 124: mat : in Standard_Complex_Poly_Matrices.Matrix;
! 125: bp : Bracket_Polynomial ) is
! 126:
! 127: -- DESCRIPTION :
! 128: -- Writes the expansion of the matrix, using the bracket polynomial which
! 129: -- is a list of quadratic monomials that represent the Laplace expansion.
! 130:
! 131: procedure Visit_Term ( t : in Bracket_Term; continue : out boolean ) is
! 132: begin
! 133: if REAL_PART(t.coeff) > 0.0
! 134: then put("+");
! 135: else put("-");
! 136: end if;
! 137: Expand_Minors(file,mat,t.monom);
! 138: continue := true;
! 139: end Visit_Term;
! 140: procedure Visit_Terms is new Enumerate_Terms(Visit_Term);
! 141:
! 142: begin
! 143: Visit_Terms(bp);
! 144: end Expand_Minors;
! 145:
! 146: procedure Pieri_Equations
! 147: ( file : in file_type; n,d : in natural; bs : in Bracket_System;
! 148: b1,b2 : in Bracket ) is
! 149:
! 150: -- DESCRIPTION :
! 151: -- Writes the Pieri equations corresponding to the pair of brackets.
! 152:
! 153: cffmat : Standard_Complex_Matrices.Matrix(1..n,1..(n-d))
! 154: := Random_Matrix(n,n-d);
! 155: polmat : Standard_Complex_Poly_Matrices.Matrix(1..n,1..d)
! 156: := Schubert_Pattern(n,b1,b2);
! 157: sys : Poly_Sys(bs'first+1..bs'last);
! 158: sol : Standard_Complex_Matrices.Matrix(1..n,1..d);
! 159:
! 160: begin
! 161: put(file,"Plane X for (");
! 162: put(file,b1); put(file,","); put(file,b2); put_line(file,") :");
! 163: put(file,polmat);
! 164: put_line(file,"The system with expanded minors of X : ");
! 165: for i in 1..bs'last loop
! 166: Expand_Minors(file,polmat,bs(i)); new_line(file);
! 167: end loop;
! 168: -- put("Give a "); put(n,1); put("x"); put(n-d,1);
! 169: -- put_line("-matrix of complex numbers : ");
! 170: -- get(cffmat);
! 171: -- put("Your "); put(n-d,1); put_line("-plane : "); put(cffmat);
! 172: -- put("A random "); put(n-d,1); put_line("-plane : "); put(cffmat);
! 173: put("The minors evaluated at a random ");
! 174: put(n-d,1); put_line("-plane :");
! 175: sys := Expanded_Minors(cffmat,polmat,bs);
! 176: put_line(sys);
! 177: Standard_Complex_Poly_Matrices.Clear(polmat);
! 178: if Pieri_Condition(n,b1,b2)
! 179: then put("The "); put(d,1); put_line("-plane at the leaves :");
! 180: sol := Solve_Pieri_Leaves(Standard_Output,b1,b2,cffmat); put(sol);
! 181: put_line("The solution evaluated at the system : ");
! 182: put(Evaluate(sys,sol)); new_line;
! 183: else put_line("Pair of leaves does not satisfy Pieri's condition.");
! 184: end if;
! 185: end Pieri_Equations;
! 186:
! 187: procedure Expand_Minors ( n,d : in natural; bs : in Bracket_System ) is
! 188:
! 189: -- DESCRIPTION :
! 190: -- Expands the minors to obtain a symbolic formulation of the equations.
! 191:
! 192: b1,b2 : Bracket(1..d);
! 193:
! 194: begin
! 195: put("Give 1st bracket : "); get(b1);
! 196: put("Give 2nd bracket : "); get(b2);
! 197: Pieri_Equations(Standard_Output,n,d,bs,b1,b2);
! 198: end Expand_Minors;
! 199:
! 200: procedure Pieri_Equations_for_Paired_Chains
! 201: ( file : in file_type; n,d : in natural; bs : in Bracket_System;
! 202: b1,b2 : in bracket_Array ) is
! 203:
! 204: -- DESCRIPTION :
! 205: -- Writes the equations for each node along a pair of chains.
! 206:
! 207: maxlen : constant natural := Maximum(b1'last,b2'last);
! 208: lb1,lb2 : Link_to_Bracket;
! 209:
! 210: begin
! 211: for i in b1'first..maxlen loop
! 212: if i <= b1'last
! 213: then lb1 := b1(i);
! 214: else lb1 := b1(b1'last);
! 215: end if;
! 216: if i <= b2'last
! 217: then lb2 := b2(i);
! 218: else lb2 := b2(b2'last);
! 219: end if;
! 220: Pieri_Equations(file,n,d,bs,lb1.all,lb2.all);
! 221: end loop;
! 222: end Pieri_Equations_for_Paired_Chains;
! 223:
! 224: procedure Write_Pieri_Equations
! 225: ( file : in file_type; n,d : in natural; t1,t2 : in Pieri_Tree;
! 226: bs : in Bracket_System ) is
! 227:
! 228: -- DESCRIPTION :
! 229: -- Writes the Pieri Equations for all pairs of chains.
! 230:
! 231: procedure Visit_Chain ( b1,b2 : in Bracket_Array; cont : out boolean ) is
! 232: begin
! 233: Pieri_Equations_for_Paired_Chains(file,n,d,bs,b1,b2);
! 234: cont := true;
! 235: end Visit_Chain;
! 236: procedure Visit_Chains is new Enumerate_Paired_Chains(Visit_Chain);
! 237:
! 238: begin
! 239: Visit_Chains(t1,t2);
! 240: end Write_Pieri_Equations;
! 241:
! 242: procedure Write_Pieri_Equations
! 243: ( n,d : in natural; t1,t2 : in Pieri_Tree ) is
! 244:
! 245: k,kd : natural;
! 246: bm : Bracket_Monomial;
! 247: file : file_type;
! 248:
! 249: begin
! 250: new_line; skip_line;
! 251: put_line("Reading the name of the output file.");
! 252: Read_Name_and_Create_File(file);
! 253: put("Give k to determine (m-k+1)-plane : "); get(k);
! 254: kd := n-k+1;
! 255: bm := Maximal_Minors(n,kd); -- because n = m+p
! 256: put(file,"All maximal minors : "); put(file,bm); new_line(file);
! 257: declare
! 258: bs : constant Bracket_System := Minor_Equations(kd,kd-d,bm);
! 259: begin
! 260: put_line(file,"The generic equation in the Laplace expansion : ");
! 261: put(file,bs(0));
! 262: put_line(file,"The specific equations in the system : ");
! 263: for i in 1..bs'last loop
! 264: put(file,bs(i));
! 265: end loop;
! 266: Write_Pieri_Equations(file,n,d,t1,t2,bs);
! 267: end;
! 268: Close(file);
! 269: Clear(bm);
! 270: end Write_Pieri_Equations;
! 271:
! 272: -- AUXILIARIES TO REPRESENT PIERI TREES :
! 273:
! 274: procedure Write_a_Chain ( file : in file_type; b : in Bracket_Array ) is
! 275: begin
! 276: put(b(b'first).all);
! 277: for i in b'first+1..b'last loop
! 278: put(" < "); put(b(i).all);
! 279: end loop;
! 280: new_line;
! 281: end Write_a_Chain;
! 282:
! 283: procedure Write_Chains ( t : in Pieri_Tree ) is
! 284:
! 285: procedure Write_Chain ( b : in Bracket_Array; cont : out boolean ) is
! 286: begin
! 287: Write_a_Chain(Standard_Output,b);
! 288: cont := true;
! 289: end Write_Chain;
! 290: procedure Write is new Enumerate_Chains(Write_Chain);
! 291:
! 292: begin
! 293: Write(t);
! 294: end Write_Chains;
! 295:
! 296: -- AUXILIARIES TO COUNT THE ROOTS :
! 297:
! 298: procedure Write_Polynomial_Patterns
! 299: ( file : in file_type;
! 300: n : in natural; b1,b2 : in Bracket_Array ) is
! 301:
! 302: -- DESCRIPTION :
! 303: -- Writes the two chains in a paired fashion. If they have unequal
! 304: -- length, then the last element of the shortest chain appears repeated.
! 305:
! 306: maxlen : constant natural := Maximum(b1'last,b2'last);
! 307: lb1,lb2 : Link_to_Bracket;
! 308:
! 309: begin
! 310: for i in b1'first..maxlen loop
! 311: put(file,"(");
! 312: if i <= b1'last
! 313: then lb1 := b1(i);
! 314: else lb1 := b1(b1'last);
! 315: end if;
! 316: put(file,lb1.all); put(file,",");
! 317: if i <= b2'last
! 318: then lb2 := b2(i);
! 319: else lb2 := b2(b2'last);
! 320: end if;
! 321: put(file,lb2.all); put_line(file,") has pattern : ");
! 322: Display_Polynomial_Pattern(file,n,lb1.all,lb2.all);
! 323: end loop;
! 324: end Write_Polynomial_Patterns;
! 325:
! 326: procedure Write_Paired_Chain
! 327: ( file : in file_type;
! 328: n : in natural; b1,b2 : in Bracket_Array ) is
! 329:
! 330: -- DESCRIPTION :
! 331: -- Writes the two chains in a paired fashion. If they have unequal
! 332: -- length, then the last element of the shortest chain appears repeated.
! 333:
! 334: maxlen : constant natural := Maximum(b1'last,b2'last);
! 335:
! 336: begin
! 337: put(file,"("); put(file,b1(b1'first).all); put(file,",");
! 338: put(file,b2(b2'first).all); put(file,")");
! 339: for i in b1'first+1..maxlen loop
! 340: put(file," < "); put(file,"(");
! 341: if i <= b1'last
! 342: then put(file,b1(i).all);
! 343: else put(file,b1(b1'last).all);
! 344: end if;
! 345: put(file,",");
! 346: if i <= b2'last
! 347: then put(file,b2(i).all);
! 348: else put(file,b2(b2'last).all);
! 349: end if;
! 350: put(file,")");
! 351: end loop;
! 352: new_line(file);
! 353: Write_Polynomial_Patterns(Standard_Output,n,b1,b2);
! 354: if Pieri_Condition(n,b1(b1'last).all,b2(b2'last).all)
! 355: then put_line("Leaves satisfy Pieri's condition.");
! 356: else put_line("Leaves do not satisfy Pieri's condition.");
! 357: end if;
! 358: end Write_Paired_Chain;
! 359:
! 360: procedure Write_Pieri_Chains ( n : in natural; t1,t2 : in Pieri_Tree ) is
! 361:
! 362: -- DESCRIPTION :
! 363: -- Enumerates all pairs of chains and checks Pieri's condition
! 364: -- at the leaves.
! 365:
! 366: procedure Visit_Chain ( b1,b2 : in Bracket_Array; cont : out boolean ) is
! 367: begin
! 368: Write_Paired_Chain(Standard_Output,n,b1,b2);
! 369: cont := true;
! 370: end Visit_Chain;
! 371: procedure Visit_Chains is new Enumerate_Paired_Chains(Visit_Chain);
! 372:
! 373: begin
! 374: Visit_Chains(t1,t2);
! 375: end Write_Pieri_Chains;
! 376:
! 377: function First_Standard_Plane
! 378: ( n,m,r : natural ) return Standard_Complex_Matrices.Matrix is
! 379:
! 380: -- DESCRIPTION :
! 381: -- Returns the plane spanned by the first m+1-r standard basis vectors.
! 382:
! 383: res : Standard_Complex_Matrices.Matrix := Random_Matrix(n,m+1-r);
! 384:
! 385: begin
! 386: for i in res'range(1) loop
! 387: for j in res'range(2) loop
! 388: if i = j
! 389: then res(i,j) := Create(1.0);
! 390: else res(i,j) := Create(0.0);
! 391: end if;
! 392: end loop;
! 393: end loop;
! 394: return res;
! 395: end First_Standard_Plane;
! 396:
! 397: function Last_Standard_Plane
! 398: ( n,m,r : natural ) return Standard_Complex_Matrices.Matrix is
! 399:
! 400: -- DESCRIPTION :
! 401: -- Returns the plane spanned by the first m+1-r standard basis vectors.
! 402:
! 403: res : Standard_Complex_Matrices.Matrix := Random_Matrix(n,m+1-r);
! 404:
! 405: begin
! 406: for i in res'range(1) loop
! 407: for j in res'range(2) loop
! 408: if i = res'last(2) + 1 - j
! 409: then res(i,j) := Create(1.0);
! 410: else res(i,j) := Create(0.0);
! 411: end if;
! 412: end loop;
! 413: end loop;
! 414: return res;
! 415: end Last_Standard_Plane;
! 416:
! 417: function First_Random_Input_Sequence
! 418: ( n,m,a : natural; kp : Vector ) return VecMat is
! 419:
! 420: -- DESCRIPTION :
! 421: -- Returns the first sequence of random input planes. The first plane
! 422: -- in the sequence is spanned by the first standard basis vectors.
! 423: -- The dimensions of the planes are m+1-kp(i), for the appropriate i.
! 424:
! 425: res : VecMat(0..a-1);
! 426:
! 427: begin
! 428: res(0) := new Standard_Complex_Matrices.Matrix'
! 429: (First_Standard_Plane(n,m,kp(1)));
! 430: for i in 1..res'last loop
! 431: res(i) := new Standard_Complex_Matrices.Matrix'
! 432: (Random_Matrix(n,m+1-kp(i+1)));
! 433: end loop;
! 434: return res;
! 435: end First_Random_Input_Sequence;
! 436:
! 437: function Second_Random_Input_Sequence
! 438: ( n,m,a : natural; kp : Vector ) return VecMat is
! 439:
! 440: -- DESCRIPTION :
! 441: -- Returns the second sequence of random input planes. The first plane
! 442: -- in the sequence is spanned by the last standard basis vectors.
! 443:
! 444: res : VecMat(0..kp'last-a-2);
! 445:
! 446: begin
! 447: res(0) := new Standard_Complex_Matrices.Matrix'
! 448: (Last_Standard_Plane(n,m,kp(1)));
! 449: for i in 1..res'last loop
! 450: res(i) := new Standard_Complex_Matrices.Matrix'
! 451: (Random_Matrix(n,m+1-kp(a+1+i)));
! 452: end loop;
! 453: return res;
! 454: end Second_Random_Input_Sequence;
! 455:
! 456: procedure Reallify ( c : in out Complex_Number ) is
! 457:
! 458: -- DESCRIPTION :
! 459: -- Sets the imaginary part of c to zero.
! 460:
! 461: begin
! 462: c := Create(REAL_PART(c),0.0);
! 463: end Reallify;
! 464:
! 465: procedure Reallify ( m : in out Standard_Complex_Matrices.Matrix ) is
! 466:
! 467: -- DESCRIPTION :
! 468: -- Sets the imaginary part of every entry in the matrix to zero.
! 469:
! 470: begin
! 471: for i in m'range(1) loop
! 472: for j in m'range(2) loop
! 473: Reallify(m(i,j));
! 474: end loop;
! 475: end loop;
! 476: end Reallify;
! 477:
! 478: procedure Reallify ( v : in out VecMat ) is
! 479:
! 480: -- DESCRIPTION :
! 481: -- Sets the imaginary part of every entry of every matrix in v to zero
! 482: -- and makes the matrices orthonormal.
! 483:
! 484: begin
! 485: for i in v'range loop
! 486: Reallify(v(i).all);
! 487: v(i).all := Orthogonalize(v(i).all);
! 488: end loop;
! 489: end Reallify;
! 490:
! 491: procedure Orthogonalize ( v : in out VecMat ) is
! 492:
! 493: -- DESCRIPTION :
! 494: -- Orthonormalizes every matrix in the array.
! 495:
! 496: begin
! 497: for i in v'range loop
! 498: v(i).all := Orthogonalize(v(i).all);
! 499: end loop;
! 500: end Orthogonalize;
! 501:
! 502: procedure Set_Parameters ( file : in file_type; report : out boolean ) is
! 503:
! 504: -- DESCRIPTION :
! 505: -- Interactive determination of the continuation and output parameters.
! 506:
! 507: oc : natural;
! 508:
! 509: begin
! 510: new_line;
! 511: Driver_for_Continuation_Parameters(file);
! 512: new_line;
! 513: Driver_for_Process_io(file,oc);
! 514: report := not (oc = 0);
! 515: end Set_Parameters;
! 516:
! 517: function Select_Pairs ( lps : List_of_Paired_Nodes )
! 518: return List_of_Paired_Nodes is
! 519:
! 520: -- DESCRIPTION :
! 521: -- Returns a selection of the list of pairs.
! 522:
! 523: res,res_last : List_of_Paired_Nodes;
! 524: k : natural;
! 525: tmp : List_of_Paired_Nodes := lps;
! 526:
! 527: begin
! 528: put("Give the number of pairs : "); get(k);
! 529: declare
! 530: sel : Standard_Natural_Vectors.Vector(1..k);
! 531: ind : natural := 1;
! 532: begin
! 533: put("Give an increasing sequence of "); put(k,1); put(" numbers : ");
! 534: get(sel);
! 535: for i in 1..Length_Of(lps) loop
! 536: if i = sel(ind)
! 537: then ind := ind+1;
! 538: Append(res,res_last,Head_Of(tmp));
! 539: end if;
! 540: tmp := Tail_Of(tmp);
! 541: end loop;
! 542: end;
! 543: skip_line;
! 544: return res;
! 545: end Select_Pairs;
! 546:
! 547: procedure put ( file : in file_type; lp : in List_of_Paired_Nodes ) is
! 548:
! 549: tmp : List_of_Paired_Nodes := lp;
! 550: pnd : Paired_Nodes;
! 551:
! 552: begin
! 553: put_line(file,"The pairs of leaves that satisfy Pieri's condition :");
! 554: for i in 1..Length_Of(lp) loop
! 555: pnd := Head_Of(tmp);
! 556: put(file,"Leaf "); put(file,i,3); put(file," : ");
! 557: put(file,"("); put(file,pnd.left.node); put(file,",");
! 558: put(file,pnd.right.node); put_line(file,")");
! 559: tmp := Tail_Of(tmp);
! 560: end loop;
! 561: end put;
! 562:
! 563: procedure Root_Count ( n,d,a : in natural; kp : in Vector ) is
! 564:
! 565: -- DESCRIPTION :
! 566: -- Set up of Pieri trees from a partition of the planes.
! 567:
! 568: -- ON ENTRY :
! 569: -- n the dimension of the space we are working in;
! 570: -- d dimension of the brackets, of the output planes;
! 571: -- a number of planes in the first set of the partition;
! 572: -- kp the dimensions of the input planes are m+1-kp(i),
! 573: -- where m = n-d.
! 574:
! 575: file : file_type;
! 576: timer : Timing_Widget;
! 577: ans : character;
! 578: m : constant natural := n-d;
! 579: v1 : Vector(0..a-1) := kp(1..a);
! 580: t1 : Pieri_Tree(d,a-1) := Create(n,d,v1);
! 581: a2 : constant natural := kp'last-a-1;
! 582: v2 : Vector(0..a2-1) := kp(a+1..kp'last-1);
! 583: t2 : Pieri_Tree(d,a2) := Create(n,d,v2);
! 584: lp : List_of_Paired_Nodes := Create(n,d,t1,t2);
! 585: np : Nodal_Pair(d) := Create(d,lp);
! 586: sel_lp : List_of_Paired_Nodes;
! 587: rc : constant natural := Length_Of(lp);
! 588: nb : constant natural := Number_of_Paths(np);
! 589: l1 : VecMat(0..a-1) := First_Random_Input_Sequence(n,m,a,kp);
! 590: l2 : VecMat(0..kp'last-a-2) := Second_Random_Input_Sequence(n,m,a,kp);
! 591: ln : Standard_Complex_Matrices.Matrix(1..n,1..m+1-kp(kp'last))
! 592: := Random_Matrix(n,m+1-kp(kp'last));
! 593: report,outlog : boolean;
! 594: sols : VecMat(1..rc);
! 595:
! 596: begin
! 597: skip_line;
! 598: new_line;
! 599: put_line("Reading the name of the output file.");
! 600: Read_Name_and_Create_File(file);
! 601: new_line;
! 602: put_line("See the output file for results.");
! 603: new_line;
! 604: put(file,"(m = "); put(file,m,1); put(file,",p = "); put(file,d,1);
! 605: put(file,")");
! 606: put(file," k = {");
! 607: for i in 1..a-1 loop
! 608: put(file,kp(i),1); put(file,",");
! 609: end loop;
! 610: put(file,kp(a),1);
! 611: put(file,"}");
! 612: put(file,"{");
! 613: for i in a+1..kp'last-2 loop
! 614: put(file,kp(i),1); put(file,",");
! 615: end loop;
! 616: put(file,kp(kp'last-1),1);
! 617: put(file,"} ");
! 618: put(file,kp(kp'last),1);
! 619: new_line(file);
! 620: put(file,"The first Pieri tree has height "); put(file,Height(t1),1);
! 621: put_line(file," :"); Write_Tree(file,t1); -- put(file,t1);
! 622: put(file,"The second Pieri tree has height "); put(file,Height(t2),1);
! 623: put_line(file," :"); Write_Tree(file,t2); -- put(file,t2);
! 624: put(file,"The root count equals : "); put(file,rc,1); new_line(file);
! 625: put("The root count equals : "); put(rc,1); new_line;
! 626: put(file,lp);
! 627: put_line(file,"The tree of paired nodes :");
! 628: Write(file,np);
! 629: put("The number of paths : "); put(nb,1); new_line;
! 630: put(file,"The number of paths : "); put(file,nb,1); new_line(file);
! 631: new_line;
! 632: put("Do you want real random input planes ? (y/n) "); get(ans);
! 633: if ans = 'y'
! 634: then Reallify(l1); Reallify(l2); Reallify(ln);
! 635: else Orthogonalize(l1); Orthogonalize(l2);
! 636: ln := Orthogonalize(ln);
! 637: end if;
! 638: put("Do you want moving cycles/polynomial systems on file ? (y/n) ");
! 639: Ask_Yes_or_No(ans);
! 640: outlog := (ans = 'y');
! 641: put("Do you want to select only certain pairs ? (y/n) ? ");
! 642: Ask_Yes_or_No(ans);
! 643: if ans = 'y'
! 644: then sel_lp := Select_Pairs(lp);
! 645: else sel_lp := lp;
! 646: end if;
! 647: Set_Parameters(file,report);
! 648: new_line(file);
! 649: put_line(file,"The first sequence of random input planes :");
! 650: put(file,l1,2);
! 651: put_line(file,"The second sequence of random input planes :");
! 652: put(file,l2,2);
! 653: put_line(file,"The last random input plane :"); put(file,ln,2);
! 654: put_line(file,"Starting the deformations at the paired leaves.");
! 655: Matrix_Indeterminates.Initialize_Symbols(n,d);
! 656: Add_t_Symbol;
! 657: tstart(timer);
! 658: Deform_Pairs(file,n,d,sel_lp,l1,l2,ln,report,outlog,sols);
! 659: tstop(timer);
! 660: new_line(file);
! 661: print_times(file,timer,"Pieri Deformations");
! 662: Matrix_Indeterminates.Clear_Symbols;
! 663: Clear(t1); Clear(t2); Clear(lp);
! 664: Clear(l1); Clear(l2); Clear(sols);
! 665: Close(file);
! 666: end Root_Count;
! 667:
! 668: -- FIVE MAJOR TEST PROGRAMS :
! 669:
! 670: procedure Test_Pieri_Condition ( n,d : in natural ) is
! 671:
! 672: -- DESCRIPTION :
! 673: -- Reads two brackets and tests Pieri's condition.
! 674:
! 675: b1,b2 : Bracket(1..d);
! 676: ans : character;
! 677: cnt : natural := 0;
! 678:
! 679: begin
! 680: new_line;
! 681: put_line("Interactive test on Pieri's condition for pair of brackets.");
! 682: Matrix_Indeterminates.Initialize_Symbols(n,d);
! 683: loop
! 684: new_line;
! 685: put("Give first bracket : "); get(b1);
! 686: put("Give second bracket : "); get(b2);
! 687: put("("); put(b1); put(","); put(b2); put(")");
! 688: if Pieri_Condition(n,b1,b2)
! 689: then put_line(" satisfies Pieri's condition.");
! 690: else put_line(" does not satisfy Pieri's condition.");
! 691: end if;
! 692: put("Do you want more tests ? (y/n) "); get(ans);
! 693: exit when ans /= 'y';
! 694: end loop;
! 695: Matrix_Indeterminates.Clear_Symbols;
! 696: end Test_Pieri_Condition;
! 697:
! 698: procedure Test_Pieri_Tree ( n,d : in natural ) is
! 699:
! 700: -- DESCRIPTION : Constructs T(r_0,r_1,..,r_a).
! 701:
! 702: a : natural;
! 703: ans : character;
! 704:
! 705: begin
! 706: new_line;
! 707: put_line("Interactive test on the construction of Pieri trees");
! 708: loop
! 709: new_line;
! 710: put("Give number of jumped-branching levels : "); get(a);
! 711: declare
! 712: v : Vector(0..a);
! 713: t : Pieri_Tree(d,a);
! 714: begin
! 715: put("Give "); put(a+1,1); put(" natural numbers : "); get(v);
! 716: t := Create(n,d,v);
! 717: put_line("The Pieri tree : "); Write_Tree(t); --put(t);
! 718: put_line("The chains in the Pieri tree :"); Write_Chains(t);
! 719: Clear(t);
! 720: end;
! 721: put("Do you want more tests ? (y/n) "); get(ans);
! 722: exit when ans /= 'y';
! 723: end loop;
! 724: end Test_Pieri_Tree;
! 725:
! 726: procedure Test_Root_Count ( n,d : in natural ) is
! 727:
! 728: -- DESCRIPTION :
! 729: -- Reads the input planes and sets up a pair of trees to perform
! 730: -- the combinatorial root count.
! 731:
! 732: m : constant natural := n-d;
! 733: np,sum,a : natural;
! 734: ans : character;
! 735:
! 736: begin
! 737: new_line;
! 738: put_line("Counting roots combinatorially by Pieri trees");
! 739: loop
! 740: new_line;
! 741: put("Give the number of planes to intersect : "); get(np);
! 742: declare
! 743: kp : Vector(1..np);
! 744: begin
! 745: put("Give "); put(np,1); put(" co-dimensions of the planes : ");
! 746: get(kp);
! 747: put("m = "); put(m,1); put(" p = "); put(d,1);
! 748: put(" Sum : "); put(kp(kp'first),1);
! 749: sum := kp(kp'first);
! 750: for i in kp'first+1..kp'last loop
! 751: put(" + "); put(kp(i),1);
! 752: sum := sum + kp(i);
! 753: end loop;
! 754: put(" = "); put(sum,1);
! 755: if sum = m*d
! 756: then put_line(" = m*p");
! 757: loop
! 758: put("Give number of elements in first set : "); get(a);
! 759: Root_Count(n,d,a,kp);
! 760: new_line;
! 761: put("Do you want to test other partitions ? (y/n) "); get(ans);
! 762: exit when ans /= 'y';
! 763: end loop;
! 764: else put_line(" /= m*p");
! 765: end if;
! 766: end;
! 767: put("Do you want more tests ? (y/n) "); get(ans);
! 768: exit when ans /= 'y';
! 769: end loop;
! 770: end Test_Root_Count;
! 771:
! 772: procedure Test_Pieri_Equations ( n,d : in natural ) is
! 773:
! 774: -- DESCRIPTION :
! 775: -- Set up of the expansions of the maximal minors.
! 776:
! 777: k,kd,nb : natural;
! 778: bm : Bracket_Monomial;
! 779: ans : character;
! 780: b1,b2 : Bracket(1..d);
! 781:
! 782: begin
! 783: new_line;
! 784: put_line("Set up equations for certain Schubert conditions.");
! 785: Matrix_Indeterminates.Initialize_Symbols(n,d);
! 786: loop
! 787: new_line;
! 788: put("Give k to determine (m-k+1)-plane : "); get(k);
! 789: kd := n-k+1;
! 790: bm := Maximal_Minors(n,kd); -- because n = m+p
! 791: put("All maximal minors : "); put(bm); new_line;
! 792: declare
! 793: bs : constant Bracket_System := Minor_Equations(kd,kd-d,bm);
! 794: begin
! 795: put_line("The generic equation in the Laplace expansion : ");
! 796: put(bs(0));
! 797: put_line("The specific equations in the system : ");
! 798: for i in 1..bs'last loop
! 799: put(bs(i));
! 800: end loop;
! 801: Expand_Minors(n,d,bs);
! 802: end;
! 803: put("Do you want more tests ? (y/n) "); get(ans);
! 804: Clear(bm);
! 805: exit when (ans /= 'y');
! 806: end loop;
! 807: Matrix_Indeterminates.Clear_Symbols;
! 808: new_line;
! 809: put("Do you want to test memory consumption ? (y/n) "); get(ans);
! 810: if ans = 'y'
! 811: then put("Give number of times : "); get(nb);
! 812: k := 1;
! 813: for i in 1..d loop
! 814: b1(i) := i;
! 815: b2(i) := i;
! 816: end loop;
! 817: for i in 1..nb loop
! 818: kd := n-k+1;
! 819: bm := Maximal_Minors(n,kd);
! 820: declare
! 821: nva : constant natural := n*d+1;
! 822: bs : Bracket_System(0..Number_of_Brackets(bm))
! 823: := Minor_Equations(kd,kd-d,bm);
! 824: cffmat : Standard_Complex_Matrices.Matrix(1..n,1..(n-d))
! 825: := Random_Matrix(n,n-d);
! 826: stamat : Standard_Complex_Matrices.Matrix(1..n,1..(n-d))
! 827: := Random_Matrix(n,n-d);
! 828: movmat : Standard_Complex_Poly_Matrices.Matrix
! 829: (cffmat'range(1),cffmat'range(2))
! 830: := Moving_U_Matrix(nva,stamat,cffmat);
! 831: polmat : Standard_Complex_Poly_Matrices.Matrix(1..n,1..d)
! 832: := Schubert_Pattern(n,b1,b2);
! 833: sys : Poly_Sys(1..bs'last)
! 834: := Expanded_Minors(cffmat,polmat,bs);
! 835: movcyc : Poly_Sys(1..bs'last)
! 836: := Lifted_Expanded_Minors(movmat,polmat,bs);
! 837: begin
! 838: Clear(bm); Clear(bs);
! 839: Standard_Complex_Poly_Matrices.Clear(movmat);
! 840: Standard_Complex_Poly_Matrices.Clear(polmat);
! 841: Clear(sys); Clear(movcyc);
! 842: end;
! 843: end loop;
! 844: end if;
! 845: end Test_Pieri_Equations;
! 846:
! 847: procedure Test_Pieri_Homotopies ( n,d : in natural ) is
! 848:
! 849: -- DESCRIPTION :
! 850: -- Set up of the parametric families in the Pieri homotopy algorithm.
! 851:
! 852: b1,b2 : Bracket(1..d);
! 853: m : constant natural := n-d;
! 854: nva : constant natural := n*d+1;
! 855: k,kd : natural;
! 856: bm : Bracket_Monomial;
! 857: F1 : constant Standard_Complex_Matrices.Matrix
! 858: := Random_Upper_Triangular(n);
! 859: F2 : constant Standard_Complex_Matrices.Matrix
! 860: := Random_Lower_Triangular(n);
! 861:
! 862: begin
! 863: new_line;
! 864: put_line("Set up the moving cycles in the Pieri homotopy algorithm.");
! 865: new_line;
! 866: Matrix_Indeterminates.Initialize_Symbols(n,d);
! 867: Add_t_Symbol;
! 868: put_line("The general upper triangular matrix F : "); put(F1,2);
! 869: put_line("The general lower triangular matrix F' : "); put(F2,2);
! 870: put("Give first "); put(d,1); put("-bracket : "); get(b1);
! 871: put("Give second "); put(d,1); put("-bracket : "); get(b2);
! 872: put("The pair of brackets : ");
! 873: put("("); put(b1); put(","); put(b2); put_line(")");
! 874: put("Give k to determine (m-k+1)-plane : "); get(k);
! 875: kd := n-k+1;
! 876: bm := Maximal_Minors(n,kd);
! 877: declare
! 878: L : constant Standard_Complex_Matrices.Matrix := Random_Matrix(n,m+1-k);
! 879: X : Standard_Complex_Poly_Matrices.Matrix(1..n,b1'range);
! 880: U1 : constant Standard_Complex_Matrices.Matrix := U_Matrix(F1,b1);
! 881: movU1 : Standard_Complex_Poly_Matrices.Matrix(L'range(1),L'range(2));
! 882: U2 : constant Standard_Complex_Matrices.Matrix := U_Matrix(F2,b2);
! 883: movU2 : Standard_Complex_Poly_Matrices.Matrix(L'range(1),L'range(2));
! 884: bs : constant Bracket_System := Minor_Equations(kd,kd-d,bm);
! 885: movcyc1,movcyc2 : Poly_Sys(1..bs'last);
! 886: begin
! 887: put("The U-matrix for F and "); put(b1); put_line(" :"); put(U1,2);
! 888: put("The U-matrix for F' and "); put(b2); put_line(" :"); put(U2,2);
! 889: put("A random "); put(m+1-k,1); put_line("-plane :"); put(L,2);
! 890: movU1 := Moving_U_Matrix(nva,U1,L);
! 891: put_line("The first moving U-matrix :"); put(movU1);
! 892: movU2 := Moving_U_Matrix(nva,U2,L);
! 893: put_line("The second moving U-matrix :"); put(movU2);
! 894: X := Schubert_Pattern(n,b1,b2);
! 895: put("The unknown "); put(d,1); put_line("-plane X :"); put(X);
! 896: movcyc1 := Lifted_Expanded_Minors(movU1,X,bs);
! 897: put_line("The polynomial system for the first moving U-matrix :");
! 898: put_line(movcyc1);
! 899: movcyc2 := Lifted_Expanded_Minors(movU2,X,bs);
! 900: put_line("The polynomial system for the second moving U-matrix :");
! 901: put_line(movcyc2);
! 902: put_line("Target system at t = 1 :");
! 903: put_line(Eval(movcyc2,Create(1.0),Number_of_Unknowns(movcyc2(1))));
! 904: put_line("System that must be solved :");
! 905: put_line(Expanded_Minors(L,X,bs));
! 906: end;
! 907: Matrix_Indeterminates.Clear_Symbols;
! 908: end Test_Pieri_Homotopies;
! 909:
! 910: procedure Main is
! 911:
! 912: m,p,n : natural;
! 913: ans : character;
! 914:
! 915: begin
! 916: new_line;
! 917: put_line("Interactive test for setting up of Pieri homotopies.");
! 918: new_line;
! 919: put("Give m, m+p = n, dimension of space : "); get(m);
! 920: put("Give p, number of entries in brackets : "); get(p);
! 921: n := m+p;
! 922: new_line;
! 923: put_line("MENU for testing the Pieri Homotopy Algorithm : ");
! 924: put_line(" 1. Construction of a Pieri tree. ");
! 925: put_line(" 2. Interactively test Pieri's condition on pairs of brackets.");
! 926: put_line(" 3. Set up equations for certain Schubert conditions.");
! 927: put_line(" 4. Set up the moving cycles in the Pieri homotopy algorithm.");
! 928: put_line(" 5. Test combinatorial root count and start deforming.");
! 929: put("Make your choice (1, 2, 3, 4 or 5) : "); get(ans);
! 930: case ans is
! 931: when '1' => Test_Pieri_Tree(n,p);
! 932: when '2' => Test_Pieri_Condition(n,p);
! 933: when '3' => Test_Pieri_Equations(n,p);
! 934: when '4' => Test_Pieri_Homotopies(n,p);
! 935: when '5' => Test_Root_Count(n,p);
! 936: when others => null;
! 937: end case;
! 938: end Main;
! 939:
! 940: begin
! 941: Main;
! 942: end ts_org_pieri;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to ts_org_pieri.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Schubert