Annotation of OpenXM_contrib/PHC/Ada/Schubert/driver_for_pieri_homotopies.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;
! 6: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
! 7: with Standard_Complex_Vectors;
! 8: with Standard_Complex_Vectors_io; use Standard_Complex_Vectors_io;
! 9: with Standard_Natural_Matrices;
! 10: with Standard_Natural_Matrices_io; use Standard_Natural_Matrices_io;
! 11: with Standard_Floating_Matrices;
! 12: with Standard_Complex_Matrices;
! 13: with Standard_Complex_Matrices_io; use Standard_Complex_Matrices_io;
! 14: with Standard_Complex_Norms_Equals; use Standard_Complex_Norms_Equals;
! 15: with Standard_Random_Numbers; use Standard_Random_Numbers;
! 16: with Standard_Random_Vectors; use Standard_Random_Vectors;
! 17: with Standard_Matrix_Inversion; use Standard_Matrix_Inversion;
! 18: with Standard_Complex_VecMats; use Standard_Complex_VecMats;
! 19: with Symbol_Table; use Symbol_Table;
! 20: with Matrix_Indeterminates;
! 21: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
! 22: with Standard_Complex_Poly_Matrices;
! 23: with Standard_Complex_Poly_Matrices_io; use Standard_Complex_Poly_Matrices_io;
! 24: with Standard_Complex_Poly_Systems; use Standard_Complex_Poly_Systems;
! 25: with Standard_Complex_Poly_SysFun; use Standard_Complex_Poly_SysFun;
! 26: with Standard_Complex_Poly_Systems_io; use Standard_Complex_Poly_Systems_io;
! 27: with Homotopy;
! 28: with Standard_Complex_Solutions; use Standard_Complex_Solutions;
! 29: with Standard_Complex_Solutions_io; use Standard_Complex_Solutions_io;
! 30: with Standard_Root_Refiners; use Standard_Root_Refiners;
! 31: with Increment_and_Fix_Continuation; use Increment_and_Fix_Continuation;
! 32: with Drivers_for_Poly_Continuation; use Drivers_for_Poly_Continuation;
! 33: with Brackets,Brackets_io; use Brackets,Brackets_io;
! 34: with Bracket_Monomials,Bracket_Systems; use Bracket_Monomials,Bracket_Systems;
! 35: with Specialization_of_Planes; use Specialization_of_Planes;
! 36: with Symbolic_Minor_Equations; use Symbolic_Minor_Equations;
! 37: with Numeric_Minor_Equations; use Numeric_Minor_Equations;
! 38: with Determinantal_Systems; use Determinantal_Systems;
! 39: with Plane_Representations; use Plane_Representations;
! 40: with Localization_Posets; use Localization_Posets;
! 41: with Localization_Posets_io; use Localization_Posets_io;
! 42: with Deformation_Posets; use Deformation_Posets;
! 43: with Drivers_for_Input_Planes; use Drivers_for_Input_Planes;
! 44:
! 45: procedure Driver_for_Pieri_Homotopies
! 46: ( file : in file_type; n,d : in natural ) is
! 47:
! 48: -- AUXILIARIES IN THE SOLVERS :
! 49:
! 50: procedure Add_t_Symbol is
! 51:
! 52: -- DESCRIPTION :
! 53: -- Adds the symbol for the continuation parameter t to the symbol table.
! 54:
! 55: tsb : Symbol;
! 56:
! 57: begin
! 58: Symbol_Table.Enlarge(1);
! 59: tsb(1) := 't';
! 60: for i in 2..tsb'last loop
! 61: tsb(i) := ' ';
! 62: end loop;
! 63: Symbol_Table.Add(tsb);
! 64: end Add_t_Symbol;
! 65:
! 66: procedure Set_Parameters ( file : in file_type; report : out boolean ) is
! 67:
! 68: -- DESCRIPTION :
! 69: -- Interactive determination of the continuation and output parameters.
! 70:
! 71: oc : natural;
! 72:
! 73: begin
! 74: new_line;
! 75: Driver_for_Continuation_Parameters(file);
! 76: new_line;
! 77: Driver_for_Process_io(file,oc);
! 78: report := not (oc = 0);
! 79: new_line;
! 80: put_line("No more input expected. See output file for results...");
! 81: new_line;
! 82: new_line(file);
! 83: end Set_Parameters;
! 84:
! 85: function First_Standard_Plane
! 86: ( m,p : natural ) return Standard_Complex_Matrices.Matrix is
! 87:
! 88: -- DESCRIPTION :
! 89: -- Returns the m-plane spanned by the first m standard basis vectors.
! 90:
! 91: res : Standard_Complex_Matrices.Matrix(1..m+p,1..m);
! 92:
! 93: begin
! 94: for j in 1..m loop -- assign j-th column
! 95: for i in 1..m+p loop
! 96: res(i,j) := Create(0.0); -- initialize with zeros
! 97: end loop;
! 98: res(j,j) := Create(1.0); -- j-th column = j-th basis vector
! 99: end loop;
! 100: return res;
! 101: end First_Standard_Plane;
! 102:
! 103: function Last_Standard_Plane
! 104: ( m,p : natural ) return Standard_Complex_Matrices.Matrix is
! 105:
! 106: -- DESCRIPTION :
! 107: -- Returns the m-plane spanned by the last m standard basis vectors.
! 108:
! 109: res : Standard_Complex_Matrices.Matrix(1..m+p,1..m);
! 110:
! 111: begin
! 112: for j in 1..m loop -- assign j-th column
! 113: for i in 1..m+p loop
! 114: res(i,j) := Create(0.0); -- initialize with zeros
! 115: end loop;
! 116: res(p+j,j) := Create(1.0); -- j-th vector = (p+j)-th basis vector
! 117: end loop;
! 118: return res;
! 119: end Last_Standard_Plane;
! 120:
! 121: procedure Basis_Change ( n : in natural; vm : in VecMat; transvm : out VecMat;
! 122: trans : out Standard_Complex_Matrices.Matrix ) is
! 123:
! 124: -- DESCRIPTION :
! 125: -- Changes basis so that the last two planes in vm are spanned by
! 126: -- the first and last standard basis vectors respectively.
! 127:
! 128: wrk : Standard_Complex_Matrices.Matrix(1..n,1..n);
! 129: -- penuplane : constant Standard_Complex_Matrices.Matrix
! 130: -- := vm(vm'last-1).all;
! 131: frstplane : constant Standard_Complex_Matrices.Matrix := vm(vm'first).all;
! 132: lastplane : constant Standard_Complex_Matrices.Matrix := vm(vm'last).all;
! 133: colind : natural := 0;
! 134: use Standard_Complex_Matrices;
! 135: ranflt : double_float;
! 136:
! 137: begin
! 138: for j in frstplane'range(2) loop -- fill up the work matrix
! 139: colind := colind + 1;
! 140: for i in frstplane'range(1) loop
! 141: wrk(i,colind) := frstplane(i,j);
! 142: end loop;
! 143: end loop;
! 144: for j in colind+1..(n-lastplane'length(2)) loop -- random spacing
! 145: for i in 1..n loop
! 146: ranflt := Random;
! 147: wrk(i,j) := Create(ranflt);
! 148: end loop;
! 149: end loop;
! 150: colind := n+1;
! 151: for j in reverse lastplane'range(2) loop
! 152: colind := colind - 1;
! 153: for i in lastplane'range(1) loop
! 154: wrk(i,colind) := lastplane(i,j);
! 155: end loop;
! 156: end loop;
! 157: trans := Inverse(wrk); -- transformation = inverse
! 158: for i in vm'first+1..vm'last-1 loop
! 159: transvm(i-1) := new Standard_Complex_Matrices.Matrix'
! 160: (trans*vm(i).all);
! 161: end loop;
! 162: transvm(transvm'last-1)
! 163: := new Standard_Complex_Matrices.Matrix'(trans*vm(vm'first).all);
! 164: transvm(transvm'last)
! 165: := new Standard_Complex_Matrices.Matrix'(trans*vm(vm'last).all);
! 166: end Basis_Change;
! 167:
! 168: function Solution_Plane ( locmap : Standard_Natural_Matrices.Matrix;
! 169: mat : Standard_Complex_Matrices.Matrix )
! 170: return Solution is
! 171:
! 172: -- DESCRIPTION :
! 173: -- Returns the representation of the solution plane as a vector.
! 174:
! 175: solloc : constant Standard_Complex_Vectors.Vector
! 176: := Vector_Rep(locmap,Localize(locmap,mat));
! 177: sol : Solution(solloc'length);
! 178:
! 179: begin
! 180: sol.m := 1;
! 181: sol.t := Create(0.0);
! 182: sol.res := 0.0;
! 183: sol.err := 0.0;
! 184: sol.rco := 0.0;
! 185: sol.v := solloc;
! 186: return sol;
! 187: end Solution_Plane;
! 188:
! 189: function Solution_Planes ( locmap : Standard_Natural_Matrices.Matrix;
! 190: vm : VecMat ) return Solution_List is
! 191:
! 192: -- DESCRIPTION :
! 193: -- Returns the representation of the vector of planes as a solution list.
! 194:
! 195: res,res_last : Solution_List;
! 196:
! 197: begin
! 198: for i in vm'range loop
! 199: Append(res,res_last,Solution_Plane(locmap,vm(i).all));
! 200: end loop;
! 201: return res;
! 202: end Solution_Planes;
! 203:
! 204: function Create_Hypersurface_System
! 205: ( n : natural; locmap : Standard_Natural_Matrices.Matrix;
! 206: xpm : Standard_Complex_Poly_Matrices.Matrix; planes : VecMat )
! 207: return Poly_Sys is
! 208:
! 209: -- DESCRIPTION :
! 210: -- Returns the polynomial system that collects the hypersurface
! 211: -- intersection conditions for meeting the given m-planes.
! 212: -- The system is localized according to the given localization map.
! 213:
! 214: wrksys : Poly_Sys(1..n) := Polynomial_Equations(planes(1..n),xpm);
! 215: res : Poly_Sys(1..n) := Localize(locmap,wrksys);
! 216:
! 217: begin
! 218: Clear(wrksys);
! 219: return res;
! 220: end Create_Hypersurface_System;
! 221:
! 222: function Create_General_System
! 223: ( locmap : Standard_Natural_Matrices.Matrix;
! 224: xpm : Standard_Complex_Poly_Matrices.Matrix; planes : VecMat )
! 225: return Poly_Sys is
! 226:
! 227: -- DESCRIPTION :
! 228: -- Returns the polynomial system that collects the general
! 229: -- intersection conditions for meeting the given (m+1-k(i))-planes.
! 230: -- The system is localized according to the given localization map.
! 231:
! 232: wrksys : constant Poly_Sys := Polynomial_Equations(planes,xpm);
! 233: res : constant Poly_Sys := Localize(locmap,wrksys);
! 234: tmp : Poly;
! 235:
! 236: begin
! 237: for i in wrksys'range loop
! 238: tmp := wrksys(i);
! 239: Clear(tmp);
! 240: end loop;
! 241: return res;
! 242: end Create_General_System;
! 243:
! 244: procedure Evaluate_Roots ( file : in file_type;
! 245: p : in Poly_sys; sols : in out Solution_List ) is
! 246:
! 247: -- DESCRIPTION :
! 248: -- Evaluates the roots at the system, good for overdetermined systems.
! 249:
! 250: tmp : Solution_List := sols;
! 251:
! 252: begin
! 253: for i in 1..Length_Of(sols) loop
! 254: put(file,"Solution "); put(file,i,1);
! 255: put_line(file," evaluated at the system : ");
! 256: put_line(file,Eval(p,Head_Of(tmp).v));
! 257: tmp := Tail_Of(tmp);
! 258: end loop;
! 259: end Evaluate_Roots;
! 260:
! 261: function Square ( n : natural; p : Poly_Sys ) return Poly_Sys is
! 262:
! 263: -- DESCRIPTION :
! 264: -- Returns a n-by-n system, by adding random linear combinations of
! 265: -- the polynomials p(i), i > n, to the first n polynomials.
! 266:
! 267: res : Poly_Sys(1..n);
! 268: acc : Poly;
! 269:
! 270: begin
! 271: for i in res'range loop
! 272: Copy(p(i),res(i));
! 273: for j in n+1..p'last loop
! 274: acc := Random1*p(j);
! 275: Add(res(i),acc);
! 276: Clear(acc);
! 277: end loop;
! 278: end loop;
! 279: return res;
! 280: end Square;
! 281:
! 282: procedure Refine_Roots ( file : in file_type;
! 283: p : in Poly_Sys; sols : in out Solution_List ) is
! 284:
! 285: -- DESCRIPTION :
! 286: -- Calls the root refiner, for square polynomial systems only.
! 287:
! 288: epsxa : constant double_float := 10.0**(-8);
! 289: epsfa : constant double_float := 10.0**(-8);
! 290: tolsing : constant double_float := 10.0**(-8);
! 291: max : constant natural := 3;
! 292: numit : natural := 0;
! 293:
! 294: begin
! 295: Reporting_Root_Refiner(file,p,sols,epsxa,epsfa,tolsing,numit,max,false);
! 296: end Refine_Roots;
! 297:
! 298: function General_Homotopy
! 299: ( xpm : Standard_Complex_Poly_Matrices.Matrix;
! 300: locmap : Standard_Natural_Matrices.Matrix;
! 301: start,target : VecMat ) return Poly_Sys is
! 302:
! 303: -- DESCRIPTION :
! 304: -- Returns a general homotopy between start and target planes.
! 305: -- The continuation parameter is inside the determinants.
! 306:
! 307: -- REQUIRED : dimensions of start and target matrices must correspond.
! 308:
! 309: res : Link_to_Poly_Sys;
! 310: n : constant natural := xpm'length(1);
! 311: p : constant natural := xpm'length(2);
! 312: m : constant natural := n-p;
! 313: nva : constant natural := n*p + 1;
! 314:
! 315: begin
! 316: for i in start'range loop
! 317: declare
! 318: moving : Standard_Complex_Poly_Matrices.Matrix
! 319: (start(i)'range(1),start(i)'range(2))
! 320: := Moving_U_Matrix(nva,start(i).all,target(i).all);
! 321: kd : constant natural := p + start(i)'length(2);
! 322: bm : Bracket_Monomial := Maximal_Minors(n,kd);
! 323: bs : Bracket_System(0..Number_of_Brackets(bm))
! 324: := Minor_Equations(kd,kd-p,bm);
! 325: sys : Poly_Sys(1..bs'last)
! 326: := Lifted_Expanded_Minors(moving,xpm,bs);
! 327: begin
! 328: Concat(res,sys);
! 329: Standard_Complex_Poly_Matrices.Clear(moving);
! 330: end;
! 331: end loop;
! 332: declare
! 333: locres : constant Poly_Sys := Localize(locmap,res.all);
! 334: begin
! 335: Clear(res);
! 336: return locres;
! 337: end;
! 338: end General_Homotopy;
! 339:
! 340: procedure Continuation ( file : in file_type; sols : in out Solution_List;
! 341: report : in boolean ) is
! 342:
! 343: -- DESCRIPTION :
! 344: -- Calls the path trackers starting at the given solutions.
! 345:
! 346: -- REQUIRED : Homotopy is properly created.
! 347:
! 348: procedure Sil_Cont is
! 349: new Silent_Continue(Max_Norm,
! 350: Homotopy.Eval,Homotopy.Diff,Homotopy.Diff);
! 351: procedure Rep_Cont is
! 352: new Reporting_Continue(Max_Norm,
! 353: Homotopy.Eval,Homotopy.Diff,Homotopy.Diff);
! 354:
! 355: begin
! 356: Set_Continuation_Parameter(sols,Create(0.0));
! 357: if report
! 358: then Rep_Cont(file,sols,false,Create(1.0));
! 359: else Sil_Cont(sols,false,Create(1.0));
! 360: end if;
! 361: end Continuation;
! 362:
! 363: procedure General_Path_Following
! 364: ( file : in file_type; target,homsys : in Poly_Sys;
! 365: sols : in out Solution_List; report : in boolean ) is
! 366:
! 367: -- DESCRIPTION :
! 368: -- Given a homotopy with start solutions, the path trackers will
! 369: -- trace the paths defined by this homotopy.
! 370:
! 371: -- REQUIRED : not Is_Null(sols).
! 372:
! 373: timer : Timing_Widget;
! 374: neq : constant natural := homsys'last; -- number of equations
! 375: dim : constant natural := Head_Of(sols).n; -- actual dimension
! 376: squhom : Poly_Sys(1..dim) := Square(dim,homsys);
! 377:
! 378: begin
! 379: tstart(timer);
! 380: Homotopy.Create(squhom,dim+1);
! 381: Continuation(file,sols,report);
! 382: tstop(timer);
! 383: new_line(file);
! 384: print_times(file,timer,"Determinantal Cheater's homotopy to target system");
! 385: new_line(file);
! 386: Evaluate_Roots(file,target,sols);
! 387: end General_Path_Following;
! 388:
! 389: procedure Path_Following_with_Cheater
! 390: ( file : in file_type; start,target : in Poly_Sys;
! 391: sols : in out Solution_List; report : in boolean ) is
! 392:
! 393: -- DESCRIPTION :
! 394: -- Calls the standard continuation routines to solve a specific
! 395: -- target system, starting at the solutions of the start system.
! 396: -- This is the usual linear cheater between start and target system,
! 397: -- although it will take nonsquare inputs, it should only be used
! 398: -- for hypersurface intersection conditions.
! 399:
! 400: -- REQUIRED : not Is_Null(sols).
! 401:
! 402: timer : Timing_Widget;
! 403: n : constant natural := target'last;
! 404: a : Standard_Complex_Vectors.Vector(1..n) := Random_Vector(1,n);
! 405: b : Standard_Complex_Vectors.Vector(1..n) := Random_Vector(1,n);
! 406: dimsol : constant natural := Head_Of(sols).n;
! 407: squ_target,squ_start : Poly_Sys(1..dimsol);
! 408:
! 409: begin
! 410: tstart(timer);
! 411: if dimsol = n
! 412: then Homotopy.Create(target,start,2,a,b,true); -- linear cheater
! 413: else squ_target := Square(dimsol,target);
! 414: squ_start := Square(dimsol,start);
! 415: Homotopy.Create(squ_target,squ_start,2,a,b,true);
! 416: end if;
! 417: Continuation(file,sols,report);
! 418: tstop(timer);
! 419: new_line(file);
! 420: print_times(file,timer,"Linear Cheater's homotopy to target system");
! 421: if dimsol = n
! 422: then Refine_Roots(file,target,sols);
! 423: else Refine_Roots(file,squ_target,sols);
! 424: Evaluate_Roots(file,target,sols);
! 425: end if;
! 426: end Path_Following_with_Cheater;
! 427:
! 428: procedure Solve_Hypersurface_Target_System
! 429: ( file : in file_type; m,p : in natural;
! 430: start_planes,target_planes : in VecMat;
! 431: index_poset : in Array_of_Array_of_Nodes;
! 432: deform_poset : in Array_of_Array_of_VecMats;
! 433: target_level : in natural; report : in boolean ) is
! 434:
! 435: -- DESCRIPTION :
! 436: -- This procedure tests the output of the deformation poset,
! 437: -- creating polynomial representations of the intersection conditions
! 438: -- and solution lists representing the solution planes.
! 439:
! 440: -- ON ENTRY :
! 441: -- file to write intermediate output on;
! 442: -- m dimension of the input planes;
! 443: -- p dimension of the solution planes;
! 444: -- start_planes input m-planes in general position;
! 445: -- target_planes specific input m-planes;
! 446: -- index_poset indexed localization poset;
! 447: -- deform_poset poset with the solution p-planes;
! 448: -- target_level indicates lowest node in deform_poset that is filled;
! 449: -- report switch for intermediate output during continuation.
! 450:
! 451: dim : constant natural := m*p;
! 452: xpm : Standard_Complex_Poly_Matrices.Matrix(1..m+p,1..p)
! 453: := Localization_Pattern(m+p,index_poset(dim)(1).top,
! 454: index_poset(dim)(1).bottom);
! 455: solplanes : constant VecMat := deform_poset(target_level)(1).all;
! 456: locmap : Standard_Natural_Matrices.Matrix(1..m+p,1..p)
! 457: := Standard_Coordinate_Frame(xpm,solplanes(1).all);
! 458: locsys : Poly_Sys(start_planes'range)
! 459: := Create_Hypersurface_System(dim,locmap,xpm,start_planes);
! 460: target : Poly_Sys(target_planes'range)
! 461: := Create_Hypersurface_System(dim,locmap,xpm,target_planes);
! 462: sols : Solution_List := Solution_Planes(locmap,solplanes);
! 463:
! 464: begin
! 465: Matrix_Indeterminates.Reduce_Symbols(locmap);
! 466: put_line(file,"THE GENERIC SYSTEM : ");
! 467: put_line(file,locsys);
! 468: new_line(file);
! 469: put_line(file,"The localization map :"); put(file,locmap);
! 470: Refine_Roots(file,locsys,sols);
! 471: new_line(file);
! 472: put_line(file,"THE TARGET SYSTEM : ");
! 473: put_line(file,target);
! 474: new_line(file);
! 475: Path_Following_with_Cheater(file,locsys,target,sols,report);
! 476: end Solve_Hypersurface_Target_System;
! 477:
! 478: procedure Solve_General_Target_System
! 479: ( file : in file_type; m,p : in natural; k : in Bracket;
! 480: start_planes,target_planes : in VecMat;
! 481: index_poset : in Array_of_Array_of_Nodes;
! 482: deform_poset : in Array_of_Array_of_VecMats;
! 483: target_level : in natural; report : in boolean ) is
! 484:
! 485: -- DESCRIPTION :
! 486: -- This procedure tests the output of the deformation poset,
! 487: -- creating polynomial representations of the intersection conditions
! 488: -- and solution lists representing the solution planes.
! 489:
! 490: -- ON ENTRY :
! 491: -- file to write intermediate output on;
! 492: -- m dimension of the input planes;
! 493: -- p dimension of the solution planes;
! 494: -- k co-dimension conditions;
! 495: -- start_planes input (m+1-k(i))-planes in general position;
! 496: -- target_planes specific input (m+1-k(i))-planes;
! 497: -- index_poset indexed localization poset;
! 498: -- deform_poset poset with the solution p-planes;
! 499: -- target_level indicates lowest node in deform_poset that is filled;
! 500: -- report switch for intermediate output during continuation.
! 501:
! 502: dim : constant natural := m*p;
! 503: xpm : Standard_Complex_Poly_Matrices.Matrix(1..m+p,1..p)
! 504: := Localization_Pattern(m+p,index_poset(dim)(1).top,
! 505: index_poset(dim)(1).bottom);
! 506: solplanes : constant VecMat := deform_poset(target_level)(1).all;
! 507: locmap : Standard_Natural_Matrices.Matrix(1..m+p,1..p)
! 508: := Standard_Coordinate_Frame(xpm,solplanes(1).all);
! 509: homtpy : constant Poly_Sys
! 510: := General_Homotopy(xpm,locmap,start_planes,target_planes);
! 511: locsys : constant Poly_Sys
! 512: := Create_General_System(locmap,xpm,start_planes);
! 513: target : constant Poly_Sys
! 514: := Create_General_System(locmap,xpm,target_planes);
! 515: sols : Solution_List := Solution_Planes(locmap,solplanes);
! 516:
! 517: begin
! 518: Matrix_Indeterminates.Reduce_Symbols(locmap);
! 519: put_line(file,"THE GENERIC SYSTEM : ");
! 520: put_line(file,locsys);
! 521: new_line(file);
! 522: put_line(file,"The localization map :"); put(file,locmap);
! 523: Evaluate_Roots(file,locsys,sols);
! 524: new_line(file);
! 525: put_line(file,"THE TARGET SYSTEM : ");
! 526: put_line(file,target);
! 527: new_line(file);
! 528: General_Path_Following(file,target,homtpy,sols,report);
! 529: end Solve_General_Target_System;
! 530:
! 531: procedure Write_Poset_Times
! 532: ( file : in file_type; timer : in Timing_Widget;
! 533: npaths : in Standard_Natural_Vectors.Vector;
! 534: timings : in Duration_Array ) is
! 535:
! 536: -- DESCRIPTION :
! 537: -- Writes a overview of #paths and timings spent during deformations
! 538: -- along the poset structure.
! 539:
! 540: begin
! 541: new_line(file);
! 542: put_line(file,"--------------------------------------");
! 543: put_line(file,"| TIMING INFORMATION OVERVIEW |");
! 544: put_line(file,"--------------------------------------");
! 545: put_line(file,"| n | #paths | user cpu time |");
! 546: put_line(file,"--------------------------------------");
! 547: for i in npaths'range loop
! 548: if npaths(i) /= 0
! 549: then put(file,"|"); put(file,i,4); put(file," |");
! 550: put(file,npaths(i),7); put(file," | ");
! 551: print_hms(file,timings(i)); put(file," |"); new_line(file);
! 552: end if;
! 553: end loop;
! 554: put_line(file,"--------------------------------------");
! 555: put(file,"| total |");
! 556: put(file,Standard_Natural_Vectors.Sum(npaths),7);
! 557: put(file," | ");
! 558: print_hms(file,Elapsed_User_Time(timer));
! 559: put(file," |"); new_line(file);
! 560: put_line(file,"--------------------------------------");
! 561: end Write_Poset_Times;
! 562:
! 563: procedure Solve_Hypersurface_Deformation_Poset
! 564: ( file : in file_type; m,p : in natural;
! 565: level_poset : in Array_of_Nodes;
! 566: index_poset : in Array_of_Array_of_Nodes ) is
! 567:
! 568: -- DESCRIPTION :
! 569: -- Creates a deformation poset and applies the Solve operator.
! 570:
! 571: n : constant natural := m+p;
! 572: deform_poset : Array_of_Array_of_VecMats(index_poset'range)
! 573: := Create(index_poset);
! 574: planes : VecMat(1..m*p) := Random_Complex_Planes(m,p);
! 575: target_planes : VecMat(1..m*p);
! 576: ans : character;
! 577: report,outlog : boolean;
! 578: timer : Timing_Widget;
! 579: target_level : natural := m*p;
! 580: nbp : natural := 0;
! 581: npaths : Standard_Natural_Vectors.Vector(1..target_level)
! 582: := (1..target_level => 0);
! 583: timings : Duration_Array(1..target_level) := (1..target_level => 0.0);
! 584:
! 585: begin
! 586: put_line("The size of the deformation poset : ");
! 587: put_line(file,"The size of the deformation poset : ");
! 588: put_roco(index_poset);
! 589: put_roco(file,index_poset);
! 590: new_line;
! 591: put("Give target level <= "); put(target_level,1);
! 592: put(" = root level by default : "); get(target_level);
! 593: for i in 1..target_level loop
! 594: nbp := nbp + Row_Root_Count_Sum(level_poset,i);
! 595: end loop;
! 596: put("The number of paths : "); put(nbp,1); new_line;
! 597: put(file,"The number of paths : "); put(file,nbp,1); new_line(file);
! 598: new_line;
! 599: put("Do you want to have the homotopies on file ? (y/n) ");
! 600: Ask_Yes_or_No(ans);
! 601: outlog := (ans = 'y');
! 602: if target_level < m*p
! 603: then planes(m*p).all := Last_Standard_Plane(m,p);
! 604: if target_level < m*p-1
! 605: then planes(m*p-1).all := First_Standard_Plane(m,p);
! 606: end if;
! 607: end if;
! 608: Driver_for_Input_Planes(file,m,p,target_planes);
! 609: Matrix_Indeterminates.Initialize_Symbols(m+p,p);
! 610: Add_t_Symbol;
! 611: Set_Parameters(file,report);
! 612: tstart(timer);
! 613: for i in index_poset(target_level)'range loop
! 614: declare
! 615: root : Node := index_poset(target_level)(i).all;
! 616: begin
! 617: Solve(file,m+p,deform_poset,root,planes(1..target_level),
! 618: report,outlog,npaths,timings);
! 619: end;
! 620: end loop;
! 621: tstop(timer);
! 622: new_line(file);
! 623: print_times(file,timer,"Solving along the deformation poset");
! 624: Write_Poset_Times(file,timer,
! 625: npaths(1..target_level),timings(1..target_level));
! 626: if deform_poset(target_level)(1) /= null
! 627: then new_line(file);
! 628: Solve_Hypersurface_Target_System
! 629: (file,m,p,planes,target_planes,index_poset,deform_poset,
! 630: target_level,report);
! 631: end if;
! 632: end Solve_Hypersurface_Deformation_Poset;
! 633:
! 634: procedure Solve_General_Deformation_Poset
! 635: ( file : in file_type; m,p : in natural; k : in Bracket;
! 636: index_poset : in Array_of_Array_of_Nodes ) is
! 637:
! 638: -- DESCRIPTION :
! 639: -- Applies the solver to general intersection conditions.
! 640:
! 641: deform_poset : Array_of_Array_of_VecMats(index_poset'range)
! 642: := Create(index_poset);
! 643: planes : VecMat(k'range) := Random_Complex_Planes(m,p,k);
! 644: target_planes : VecMat(k'range);
! 645: ans : character;
! 646: report,outlog : boolean;
! 647: timer : Timing_Widget;
! 648: target_level : natural := m*p;
! 649: npaths : Standard_Natural_Vectors.Vector(1..target_level)
! 650: := (1..target_level => 0);
! 651: timings : Duration_Array(1..target_level) := (1..target_level => 0.0);
! 652:
! 653: begin
! 654: put_line("The size of the deformation poset : ");
! 655: put_line(file,"The size of the deformation poset : ");
! 656: put_roco(index_poset);
! 657: put_roco(file,index_poset);
! 658: new_line;
! 659: put("Give target level <= "); put(target_level,1);
! 660: put(" = root level by default : "); get(target_level);
! 661: new_line;
! 662: put("Do you want to have the homotopies on file ? (y/n) ");
! 663: Ask_Yes_or_No(ans);
! 664: outlog := (ans = 'y');
! 665: Driver_for_Input_Planes(file,m,p,k,target_planes);
! 666: Matrix_Indeterminates.Initialize_Symbols(m+p,p);
! 667: Add_t_Symbol;
! 668: Set_Parameters(file,report);
! 669: tstart(timer);
! 670: for i in index_poset(target_level)'range loop
! 671: declare
! 672: root : Node := index_poset(target_level)(i).all;
! 673: begin
! 674: Solve(file,m+p,k,deform_poset,root,planes,report,outlog,
! 675: npaths,timings);
! 676: end;
! 677: end loop;
! 678: tstop(timer);
! 679: new_line(file);
! 680: print_times(file,timer,"Solving along the deformation poset");
! 681: Write_Poset_Times(file,timer,
! 682: npaths(1..target_level),timings(1..target_level));
! 683: if deform_poset(target_level)(1) /= null
! 684: then new_line(file);
! 685: Solve_General_Target_System
! 686: (file,m,p,k,planes,target_planes,index_poset,deform_poset,
! 687: target_level,report);
! 688: end if;
! 689: end Solve_General_Deformation_Poset;
! 690:
! 691: -- HYPERSURFACE SOLVERS :
! 692:
! 693: procedure Create_Top_Hypersurface_Poset
! 694: ( file : in file_type; m,p : in natural ) is
! 695:
! 696: -- DESCRIPTION :
! 697: -- Create the poset by incrementing only top pivots.
! 698:
! 699: timer : Timing_Widget;
! 700: root : Node(p) := Trivial_Root(m,p);
! 701: lnkroot : Link_to_Node := new Node'(root);
! 702: level_poset : Array_of_Nodes(0..m*p);
! 703: index_poset : Array_of_Array_of_Nodes(0..m*p);
! 704:
! 705: begin
! 706: tstart(timer);
! 707: Top_Create(lnkroot,m+p);
! 708: put_line("The poset created from the top : ");
! 709: put_line(file,"The poset created from the top : ");
! 710: level_poset := Create_Leveled_Poset(lnkroot);
! 711: Count_Roots(level_poset);
! 712: index_poset := Create_Indexed_Poset(level_poset);
! 713: put(index_poset);
! 714: put(file,index_poset);
! 715: Solve_Hypersurface_Deformation_Poset(file,m,p,level_poset,index_poset);
! 716: tstop(timer);
! 717: new_line(file);
! 718: print_times(file,timer,
! 719: "Total time for Hypersurface Pieri Homotopy Algorithm");
! 720: end Create_Top_Hypersurface_Poset;
! 721:
! 722: procedure Create_Bottom_Hypersurface_Poset
! 723: ( file : in file_type; m,p : in natural ) is
! 724:
! 725: -- DESCRIPTION :
! 726: -- Create the poset by decrementing only bottom pivots.
! 727:
! 728: timer : Timing_Widget;
! 729: root : Node(p) := Trivial_Root(m,p);
! 730: lnkroot : Link_to_Node := new Node'(root);
! 731: level_poset : Array_of_Nodes(0..m*p);
! 732: index_poset : Array_of_Array_of_Nodes(0..m*p);
! 733:
! 734: begin
! 735: tstart(timer);
! 736: Bottom_Create(lnkroot);
! 737: put_line("The poset created from the bottom : ");
! 738: put_line(file,"The poset created from the bottom : ");
! 739: level_poset := Create_Leveled_Poset(lnkroot);
! 740: Count_Roots(level_poset);
! 741: index_poset := Create_Indexed_Poset(level_poset);
! 742: put(index_poset);
! 743: put(file,index_poset);
! 744: Solve_Hypersurface_Deformation_Poset(file,m,p,level_poset,index_poset);
! 745: tstop(timer);
! 746: new_line(file);
! 747: print_times(file,timer,
! 748: "Total time for Hypersurface Pieri Homotopy Algorithm");
! 749: end Create_Bottom_Hypersurface_Poset;
! 750:
! 751: procedure Create_Mixed_Hypersurface_Poset
! 752: ( file : in file_type; m,p : in natural ) is
! 753:
! 754: -- DESCRIPTION :
! 755: -- Create the poset by incrementing top and decrementing bottom pivots.
! 756:
! 757: timer : Timing_Widget;
! 758: root : Node(p) := Trivial_Root(m,p);
! 759: lnkroot : Link_to_Node := new Node'(root);
! 760: level_poset : Array_of_Nodes(0..m*p);
! 761: index_poset : Array_of_Array_of_Nodes(0..m*p);
! 762:
! 763: begin
! 764: tstart(timer);
! 765: Top_Bottom_Create(lnkroot,m+p);
! 766: put_line("The poset created in a mixed fashion : ");
! 767: put_line(file,"The poset created in a mixed fashion : ");
! 768: level_poset := Create_Leveled_Poset(lnkroot);
! 769: Count_Roots(level_poset);
! 770: index_poset := Create_Indexed_Poset(level_poset);
! 771: put(index_poset);
! 772: put(file,index_poset);
! 773: Solve_Hypersurface_Deformation_Poset(file,m,p,level_poset,index_poset);
! 774: tstop(timer);
! 775: new_line(file);
! 776: print_times(file,timer,
! 777: "Total time for Hypersurface Pieri Homotopy Algorithm");
! 778: end Create_Mixed_Hypersurface_Poset;
! 779:
! 780: -- GENERAL SOLVERS :
! 781:
! 782: procedure Create_Top_General_Poset
! 783: ( file : in file_type; m,p : in natural ) is
! 784:
! 785: -- DESCRIPTION :
! 786: -- Creates a poset for counting general subspace intersections,
! 787: -- by consistently incrementing the top pivots.
! 788:
! 789: timer : Timing_Widget;
! 790: root : Node(p) := Trivial_Root(m,p);
! 791: lnkroot : Link_to_Node := new Node'(root);
! 792: codim : constant Bracket := Read_Codimensions(m,p,0);
! 793: level_poset : Array_of_Nodes(0..m*p);
! 794: index_poset : Array_of_Array_of_Nodes(0..m*p);
! 795:
! 796: begin
! 797: put(file," k = "); put(file,codim); new_line(file);
! 798: tstart(timer);
! 799: Top_Create(lnkroot,codim,m+p);
! 800: put_line("The poset created from the top : ");
! 801: put_line(file,"The poset created from the top : ");
! 802: level_poset := Create_Leveled_Poset(lnkroot);
! 803: Count_Roots(level_poset);
! 804: index_poset := Create_Indexed_Poset(level_poset);
! 805: put(index_poset);
! 806: put(file,index_poset);
! 807: Solve_General_Deformation_Poset(file,m,p,codim,index_poset);
! 808: tstop(timer);
! 809: new_line(file);
! 810: print_times(file,timer,
! 811: "Total time for General Pieri Homotopy Algorithm");
! 812: end Create_Top_General_Poset;
! 813:
! 814: procedure Create_Bottom_General_Poset
! 815: ( file : in file_type; m,p : in natural ) is
! 816:
! 817: -- DESCRIPTION :
! 818: -- Creates a poset for counting general subspace intersections,
! 819: -- by consistently incrementing the top pivots.
! 820:
! 821: timer : Timing_Widget;
! 822: root : Node(p) := Trivial_Root(m,p);
! 823: lnkroot : Link_to_Node := new Node'(root);
! 824: codim : constant Bracket := Read_Codimensions(m,p,0);
! 825: level_poset : Array_of_Nodes(0..m*p);
! 826: index_poset : Array_of_Array_of_Nodes(0..m*p);
! 827:
! 828: begin
! 829: put(file," k = "); put(file,codim); new_line(file);
! 830: tstart(timer);
! 831: Bottom_Create(lnkroot,codim);
! 832: put_line("The poset created from the bottom : ");
! 833: put_line(file,"The poset created from the bottom : ");
! 834: level_poset := Create_Leveled_Poset(lnkroot);
! 835: Count_Roots(level_poset);
! 836: index_poset := Create_Indexed_Poset(level_poset);
! 837: put(index_poset);
! 838: put(file,index_poset);
! 839: Solve_General_Deformation_Poset(file,m,p,codim,index_poset);
! 840: tstop(timer);
! 841: new_line(file);
! 842: print_times(file,timer,
! 843: "Total time for General Pieri Homotopy Algorithm");
! 844: end Create_Bottom_General_Poset;
! 845:
! 846: procedure Create_Mixed_General_Poset
! 847: ( file : in file_type; m,p : in natural ) is
! 848:
! 849: -- DESCRIPTION :
! 850: -- Creates a poset for counting general subspace intersections,
! 851: -- by incrementing the top and decrementing the bottom pivots.
! 852:
! 853: timer : Timing_Widget;
! 854: root : Node(p) := Trivial_Root(m,p);
! 855: lnkroot : Link_to_Node := new Node'(root);
! 856: codim : constant Bracket := Read_Codimensions(m,p,0);
! 857: level_poset : Array_of_Nodes(0..m*p);
! 858: index_poset : Array_of_Array_of_Nodes(0..m*p);
! 859:
! 860: begin
! 861: put(file," k = "); put(file,codim); new_line(file);
! 862: tstart(timer);
! 863: Top_Bottom_Create(lnkroot,codim,m+p);
! 864: put_line("The poset created in a mixed fashion : ");
! 865: put_line(file,"The poset created in a mixed fashion : ");
! 866: level_poset := Create_Leveled_Poset(lnkroot);
! 867: Count_Roots(level_poset);
! 868: index_poset := Create_Indexed_Poset(level_poset);
! 869: put(index_poset);
! 870: put(file,index_poset);
! 871: Solve_General_Deformation_Poset(file,m,p,codim,index_poset);
! 872: tstop(timer);
! 873: new_line(file);
! 874: print_times(file,timer,
! 875: "Total time for General Pieri Homotopy Algorithm");
! 876: end Create_Mixed_General_Poset;
! 877:
! 878: procedure Main is
! 879:
! 880: p : constant natural := d;
! 881: m : constant natural := n-d;
! 882: ans : character;
! 883:
! 884: begin
! 885: new_line;
! 886: put_line("MENU for deforming p-planes in (m+p)-space, co-dimensions k : ");
! 887: put_line(" 1. k_i = 1 consistently incrementing the top pivots.");
! 888: put_line(" 2. consistently decrementing the bottom pivots.");
! 889: put_line(" 3. mixed top-bottom sequence for poset creation.");
! 890: put_line(" 4. k_i >= 1 consistently incrementing the top pivots.");
! 891: put_line(" 5. consistently decrementing the bottom pivots.");
! 892: put_line(" 6. mixed top-bottom sequence for poset creation.");
! 893: put("Type 1, 2, 3, 4, 5, or 6 to choose : ");
! 894: Ask_Alternative(ans,"123456");
! 895: new_line;
! 896: case ans is
! 897: when '1' => new_line(file); Create_Top_Hypersurface_Poset(file,m,p);
! 898: when '2' => new_line(file); Create_Bottom_Hypersurface_Poset(file,m,p);
! 899: when '3' => new_line(file); Create_Mixed_Hypersurface_Poset(file,m,p);
! 900: when '4' => Create_Top_General_Poset(file,m,p);
! 901: when '5' => Create_Bottom_General_Poset(file,m,p);
! 902: when '6' => Create_Mixed_General_Poset(file,m,p);
! 903: when others => put_line("Option not recognized. Please try again.");
! 904: end case;
! 905: end Main;
! 906:
! 907: begin
! 908: Main;
! 909: end Driver_for_Pieri_Homotopies;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to driver_for_pieri_homotopies.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