Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Dynlift/dynamic_polyhedral_continuation.adb, Revision 1.1
1.1 ! maekawa 1: with integer_io; use integer_io;
! 2: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
! 3: with Standard_Complex_Laur_Polys; use Standard_Complex_Laur_Polys;
! 4: with Standard_Complex_Laur_Systems; use Standard_Complex_Laur_Systems;
! 5: with Standard_Poly_Laur_Convertors; use Standard_Poly_Laur_Convertors;
! 6: with Fewnomial_System_Solvers; use Fewnomial_System_Solvers;
! 7: with Transforming_Laurent_Systems;
! 8: with Simplices; use Simplices;
! 9: with Dynamic_Triangulations; use Dynamic_Triangulations;
! 10: with Cayley_Trick; use Cayley_Trick;
! 11: with Integer_Lifting_Utilities; use Integer_Lifting_Utilities;
! 12: with Unfolding_Subdivisions; use Unfolding_Subdivisions;
! 13: with Integer_Mixed_Subdivisions_io; use Integer_Mixed_Subdivisions_io;
! 14: with Integer_Polyhedral_Continuation; use Integer_Polyhedral_Continuation;
! 15:
! 16: package body Dynamic_Polyhedral_Continuation is
! 17:
! 18: -- CAUTION : PATCH FOR FEWNOMIALS => NO SHIFT !!!
! 19:
! 20: -- AUXILIAIRIES :
! 21:
! 22: procedure Flatten ( t : in out Term ) is
! 23:
! 24: -- DESCRIPTION :
! 25: -- Flattens the Laurent term, i.e., the last exponent of t becomes zero.
! 26:
! 27: begin
! 28: t.dg(t.dg'last) := 0;
! 29: end Flatten;
! 30:
! 31: procedure Flatten ( p : in out Poly ) is
! 32:
! 33: -- DESCRIPTION :
! 34: -- Flattens the Laurent polynomial,
! 35: -- i.e., the last exponent of every monomial becomes zero.
! 36:
! 37: procedure Flatten_Term ( t : in out Term; cont : out boolean ) is
! 38: begin
! 39: Flatten(t); cont := true;
! 40: end Flatten_Term;
! 41: procedure Flatten_Terms is new Changing_Iterator(Flatten_Term);
! 42:
! 43: begin
! 44: Flatten_Terms(p);
! 45: end Flatten;
! 46:
! 47: procedure Flatten ( p : in out Laur_Sys ) is
! 48:
! 49: -- DESCRIPTION :
! 50: -- Flattens the Laurent polynomial system,
! 51: -- i.e., the last exponent of every monomial becomes zero.
! 52:
! 53: begin
! 54: for i in p'range loop
! 55: Flatten(p(i));
! 56: end loop;
! 57: end Flatten;
! 58:
! 59: function Non_Flattened_Points ( l : List ) return List is
! 60:
! 61: -- DESCRIPTION :
! 62: -- Returns the list of points with last coordinate /= 0.
! 63:
! 64: tmp,res,res_last : List;
! 65: pt : Link_to_Vector;
! 66:
! 67: begin
! 68: tmp := l;
! 69: while not Is_Null(tmp) loop
! 70: pt := Head_Of(tmp);
! 71: if pt(pt'last) /= 0
! 72: then Append(res,res_last,pt);
! 73: end if;
! 74: tmp := Tail_Of(tmp);
! 75: end loop;
! 76: return res;
! 77: end Non_Flattened_Points;
! 78:
! 79: function Non_Flattened_Points ( l : Array_of_Lists )
! 80: return Array_of_Lists is
! 81:
! 82: -- DESCRIPTION :
! 83: -- Returns the points with last coordinate /= 0.
! 84:
! 85: res : Array_of_Lists(l'range);
! 86:
! 87: begin
! 88: for i in l'range loop
! 89: res(i) := Non_Flattened_Points(l(i));
! 90: end loop;
! 91: return res;
! 92: end Non_Flattened_Points;
! 93:
! 94: function Non_Flattened_Points ( fs : Face_Structures )
! 95: return Array_of_Lists is
! 96:
! 97: -- DESCRIPTION :
! 98: -- Returns the points with last coordinate /= 0.
! 99:
! 100: res : Array_of_Lists(fs'range);
! 101:
! 102: begin
! 103: for i in fs'range loop
! 104: res(i) := Non_Flattened_Points(fs(i).l);
! 105: end loop;
! 106: return res;
! 107: end Non_Flattened_Points;
! 108:
! 109: function Non_Flattened_Points
! 110: ( mix : Vector; mixsub : Mixed_Subdivision )
! 111: return Array_of_Lists is
! 112:
! 113: -- DESCRIPTION :
! 114: -- Returns the points with last coordinate /= 0.
! 115:
! 116: res,res_last : Array_of_Lists(mix'range);
! 117: tmp : Mixed_Subdivision := mixsub;
! 118:
! 119: begin
! 120: while not Is_Null(tmp) loop
! 121: declare
! 122: mic : Mixed_Cell := Head_Of(tmp);
! 123: begin
! 124: for i in mic.pts'range loop
! 125: Deep_Concat_Diff(res(i),res_last(i),Non_Flattened_Points(mic.pts(i)));
! 126: end loop;
! 127: end;
! 128: tmp := Tail_Of(tmp);
! 129: end loop;
! 130: return res;
! 131: end Non_Flattened_Points;
! 132:
! 133: function Convert ( s : Simplex ) return Mixed_Cell is
! 134:
! 135: res : Mixed_Cell;
! 136:
! 137: begin
! 138: res.nor := new vector'(Normal(s));
! 139: res.pts := new Array_of_Lists(1..1);
! 140: res.pts(1) := Shallow_Create(Vertices(s));
! 141: return res;
! 142: end Convert;
! 143:
! 144: function Convert ( normal : in vector; s : Simplex ) return Mixed_Cell is
! 145:
! 146: res : Mixed_Cell;
! 147:
! 148: begin
! 149: res.nor := new vector'(normal);
! 150: res.pts := new Array_of_Lists(1..1);
! 151: res.pts(1) := Shallow_Create(Vertices(s));
! 152: return res;
! 153: end Convert;
! 154:
! 155: function Non_Flattened_Cells
! 156: ( flatnor : Link_to_Vector; mixsub : Mixed_Subdivision )
! 157: return Mixed_Subdivision is
! 158:
! 159: -- DESCRIPTION :
! 160: -- Returns the cells which are not flattened.
! 161:
! 162: tmp,res,res_last : Mixed_Subdivision;
! 163: mic : Mixed_Cell;
! 164:
! 165: begin
! 166: tmp := mixsub;
! 167: while not Is_Null(tmp) loop
! 168: mic := Head_Of(tmp);
! 169: if mic.nor.all /= flatnor.all
! 170: then Append(res,res_last,mic);
! 171: end if;
! 172: tmp := Tail_Of(tmp);
! 173: end loop;
! 174: return res;
! 175: end Non_Flattened_Cells;
! 176:
! 177: function Convert ( t : Triangulation ) return Mixed_Subdivision is
! 178:
! 179: res,res_last : Mixed_Subdivision;
! 180: tmp : Triangulation := t;
! 181:
! 182: begin
! 183: while not Is_Null(tmp) loop
! 184: declare
! 185: mic : Mixed_Cell := Convert(Head_Of(tmp));
! 186: begin
! 187: Append(res,res_last,mic);
! 188: end;
! 189: tmp := Tail_Of(tmp);
! 190: end loop;
! 191: return res;
! 192: end Convert;
! 193:
! 194: function Non_Flattened_Cells
! 195: ( flatnor : Link_to_Vector; t : Triangulation )
! 196: return Mixed_Subdivision is
! 197:
! 198: tmp : Triangulation;
! 199: res,res_last : Mixed_Subdivision;
! 200:
! 201: begin
! 202: tmp := t;
! 203: while not Is_Null(tmp) loop
! 204: declare
! 205: s : Simplex := Head_Of(tmp);
! 206: nor : constant vector := Normal(s);
! 207: begin
! 208: if nor /= flatnor.all
! 209: then declare
! 210: mic : Mixed_Cell := Convert(nor,s);
! 211: begin
! 212: Append(res,res_last,mic);
! 213: end;
! 214: end if;
! 215: end;
! 216: tmp := Tail_Of(tmp);
! 217: end loop;
! 218: return res;
! 219: end Non_Flattened_Cells;
! 220:
! 221: -- UTILITIES :
! 222:
! 223: function Pointer_to_Last ( l : Solution_List ) return Solution_List is
! 224:
! 225: -- DESCRIPTION :
! 226: -- Returns a pointer to the last element of l.
! 227:
! 228: begin
! 229: if Is_Null(l)
! 230: then return l;
! 231: elsif Is_Null(Tail_Of(l))
! 232: then return l;
! 233: else return Pointer_to_Last(Tail_Of(l));
! 234: end if;
! 235: end Pointer_to_Last;
! 236:
! 237: function Initialize_Polyhedral_Homotopy
! 238: ( n : natural; mix : Vector; fs : Face_Structures;
! 239: p : Laur_Sys ) return Laur_Sys is
! 240:
! 241: -- DESCRIPTION :
! 242: -- Initializes the polyhedral homotopy w.r.t. to the already lifted
! 243: -- points in the face structure.
! 244:
! 245: res : Laur_Sys(p'range);
! 246:
! 247: begin
! 248: if not Is_Null(fs(fs'first).l)
! 249: then declare
! 250: lifted : Array_of_Lists(fs'range);
! 251: begin
! 252: for i in lifted'range loop
! 253: lifted(i) := fs(i).l;
! 254: end loop;
! 255: res := Perform_Lifting(n,mix,lifted,p);
! 256: end;
! 257: end if;
! 258: return res;
! 259: end Initialize_Polyhedral_Homotopy;
! 260:
! 261: function Initialize_Polyhedral_Homotopy
! 262: ( n : natural; l : List; p : Laur_Sys ) return Laur_Sys is
! 263:
! 264: -- DESCRIPTION :
! 265: -- Initializes the polyhedral homotopy w.r.t. to the already lifted
! 266: -- points in the list.
! 267:
! 268: res : Laur_Sys(p'range);
! 269:
! 270: begin
! 271: if not Is_Null(l)
! 272: then declare
! 273: lifted : Array_of_Lists(p'range);
! 274: mix : Vector(1..1) := (1..1 => n);
! 275: begin
! 276: for i in lifted'range loop
! 277: lifted(i) := l;
! 278: end loop;
! 279: res := Perform_Lifting(n,mix,lifted,p);
! 280: end;
! 281: end if;
! 282: return res;
! 283: end Initialize_Polyhedral_Homotopy;
! 284:
! 285: function Select_Terms ( p : Poly; l : List ) return Poly is
! 286:
! 287: -- DESCRIPTION :
! 288: -- Given a tuple of lifted points, selected terms of p, whose exponents
! 289: -- occur in l, will be returned.
! 290:
! 291: res : Poly := Null_Poly;
! 292: tmp : List := l;
! 293: point : Link_to_Vector;
! 294:
! 295: begin
! 296: while not Is_Null(tmp) loop
! 297: point := Head_Of(tmp);
! 298: declare
! 299: d : degrees := new Vector(point'first..point'last-1);
! 300: t : Term;
! 301: begin
! 302: d.all := point(point'first..point'last-1);
! 303: t.cf := Coeff(p,d);
! 304: t.dg := d;
! 305: Add(res,t);
! 306: Standard_Integer_Vectors.Clear
! 307: (Standard_Integer_Vectors.Link_to_Vector(d));
! 308: end;
! 309: tmp := Tail_Of(tmp);
! 310: end loop;
! 311: return res;
! 312: end Select_Terms;
! 313:
! 314: function Select_Terms ( p : Laur_Sys; mix : Vector; l : Array_of_Lists )
! 315: return Laur_Sys is
! 316:
! 317: -- DESCRIPTION :
! 318: -- Given a tuple of lifted points, selected terms of p, whose exponents
! 319: -- occur in l, will be returned.
! 320:
! 321: res : Laur_Sys(p'range);
! 322: cnt : natural := p'first;
! 323:
! 324: begin
! 325: for i in mix'range loop
! 326: for j in 1..mix(i) loop
! 327: res(cnt) := Select_Terms(p(cnt),l(i));
! 328: cnt := cnt + 1;
! 329: end loop;
! 330: end loop;
! 331: return res;
! 332: end Select_Terms;
! 333:
! 334: procedure Update_Polyhedral_Homotopy
! 335: ( p : in Laur_Sys; point : in Vector; i : in integer;
! 336: hom : in out Laur_Sys ) is
! 337:
! 338: -- DESCRIPTION :
! 339: -- Given the lifted point of the ith support list,
! 340: -- the ith polynomial of hom will be updated with a new term.
! 341:
! 342: d : degrees
! 343: := new Standard_Integer_Vectors.Vector(point'first..point'last-1);
! 344: t : Term;
! 345:
! 346: begin
! 347: d.all := point(point'first..point'last-1);
! 348: t.cf := Coeff(p(i),d);
! 349: t.dg := new Standard_Integer_Vectors.Vector'(point);
! 350: Add(hom(i),t);
! 351: Standard_Integer_Vectors.Clear(Standard_Integer_Vectors.Link_to_Vector(d));
! 352: Clear(t);
! 353: end Update_Polyhedral_homotopy;
! 354:
! 355: procedure Update_Polyhedral_Homotopy
! 356: ( p : in Laur_Sys; l : in List; i : in integer;
! 357: hom : in out Laur_Sys ) is
! 358:
! 359: -- DESCRIPTION :
! 360: -- Given a list of lifted points of the ith support list,
! 361: -- the ith polynomial of hom will be updated with new terms.
! 362:
! 363: tmp : List := l;
! 364:
! 365: begin
! 366: while not Is_Null(tmp) loop
! 367: Update_Polyhedral_Homotopy(p,Head_Of(tmp).all,i,hom);
! 368: tmp := Tail_Of(tmp);
! 369: end loop;
! 370: end Update_Polyhedral_Homotopy;
! 371:
! 372: procedure Update_Polyhedral_Homotopy
! 373: ( p : in Laur_Sys; lifted : in Array_of_Lists;
! 374: mix : in Vector; hom : in out Laur_Sys ) is
! 375:
! 376: -- DESCRIPTION :
! 377: -- Given a lists of lifted points, the polyhedral homotopy hom
! 378: -- will be updated with new terms.
! 379:
! 380: cnt : natural := hom'first;
! 381:
! 382: begin
! 383: for i in mix'range loop
! 384: for j in 1..mix(i) loop
! 385: Update_Polyhedral_Homotopy(p,lifted(i),cnt,hom);
! 386: cnt := cnt + 1;
! 387: end loop;
! 388: end loop;
! 389: end Update_Polyhedral_Homotopy;
! 390:
! 391: procedure Update_Polyhedral_Homotopy
! 392: ( p : in Laur_Sys; lifted : in List; hom : in out Laur_Sys ) is
! 393:
! 394: -- DESCRIPTION :
! 395: -- Given a list of lifted points, the unmixed polyhedral homotopy hom
! 396: -- will be updated with new terms.
! 397:
! 398: begin
! 399: for i in hom'range loop
! 400: Update_Polyhedral_Homotopy(p,lifted,i,hom);
! 401: end loop;
! 402: end Update_Polyhedral_Homotopy;
! 403:
! 404: procedure Solve_by_Unfolding
! 405: ( file : in file_type; p,hom : in Laur_Sys; n : in natural;
! 406: mix : in Vector; nor : in Link_to_Vector;
! 407: mixsub : in Mixed_Subdivision;
! 408: sols : in out Solution_List ) is
! 409:
! 410: -- DESCRIPTION :
! 411: -- All cells in the given mixed subdivision have the same normal.
! 412: -- The solutions which correspond to these cells will be computed
! 413: -- by means of the given polyhedral homotopy.
! 414:
! 415: -- ON ENTRY :
! 416: -- file to write intermediate results on;
! 417: -- p randomized system, without lifting;
! 418: -- hom the polyhedral homotopy;
! 419: -- n dimension before lifting;
! 420: -- mix type of mixture;
! 421: -- nor the same normal to all cells in the subdivision;
! 422: -- mixsub the mixed subdivision all with the same normal nor.
! 423:
! 424: -- ON RETURN :
! 425: -- sols the solutions of p, w.r.t. the cells in mixsub.
! 426:
! 427: work : Mixed_Subdivision;
! 428: mv : natural;
! 429: homsub : Laur_Sys(p'range);
! 430: first : boolean := true;
! 431: flatnor : Link_to_Vector;
! 432: sols_last : Solution_List;
! 433:
! 434: procedure Process ( mic : in Mixed_Cell; newpts : in Array_of_Lists ) is
! 435:
! 436: subsys : Laur_Sys(p'range) := Select_Terms(p,mix,mic.pts.all);
! 437: subsols : Solution_List;
! 438: fail : boolean;
! 439:
! 440: begin
! 441: Transforming_Laurent_Systems.Shift(subsys); -- patch for Fewnomials !!
! 442: Solve(subsys,subsols,fail);
! 443: put(file,"Number of solutions of subsystem : ");
! 444: put(file,Length_Of(subsols),1); new_line(file);
! 445: if first
! 446: then homsub := Perform_Lifting(n,mix,mic.pts.all,p);
! 447: sols := subsols;
! 448: first := false;
! 449: else Update_Polyhedral_Homotopy(p,newpts,mix,homsub);
! 450: Mixed_Continuation(file,homsub,mic.nor.all,subsols);
! 451: Set_Continuation_Parameter(sols,Create(0.0));
! 452: Mixed_Continuation(file,homsub,flatnor.all,sols);
! 453: Flatten(homsub);
! 454: sols_last := Pointer_to_Last(sols);
! 455: Concat(sols,sols_last,subsols);
! 456: Shallow_Clear(subsols);
! 457: end if;
! 458: Clear(subsys);
! 459: end Process;
! 460: procedure Solve_Subsystem is new Unfolding(Process);
! 461:
! 462: begin
! 463: new_line(file);
! 464: put_line(file,"**** SOLVING BY UNFOLDING ****");
! 465: new_line(file);
! 466: flatnor := new vector(1..n+1);
! 467: flatnor(1..n) := (1..n => 0);
! 468: flatnor(n+1) := 1;
! 469: Copy(mixsub,work);
! 470: put(file,n,mix,work,mv);
! 471: Solve_Subsystem(work);
! 472: --Deep_Clear(work);
! 473: Clear(flatnor); Clear(homsub);
! 474: Set_Continuation_Parameter(sols,Create(0.0));
! 475: Mixed_Continuation(file,hom,nor.all,sols);
! 476: end Solve_by_Unfolding;
! 477:
! 478: procedure Solve_with_Unfolding
! 479: ( file : in file_type; p,hom : in Laur_Sys; n : in natural;
! 480: mix : in Vector; mixsub : in Mixed_Subdivision;
! 481: sols : in out Solution_List ) is
! 482:
! 483: -- DESCRIPTION :
! 484: -- Computes the new solutions corresponding to the cells in the
! 485: -- mixed subdivision.
! 486:
! 487: sols_last : Solution_List;
! 488: newnor : List := Different_Normals(mixsub);
! 489: tmp : List := newnor;
! 490: normal : Link_to_Vector;
! 491:
! 492: begin
! 493: while not Is_Null(tmp) loop
! 494: normal := Head_Of(tmp);
! 495: declare
! 496: newcells : Mixed_Subdivision := Extract(normal.all,mixsub);
! 497: bigcell : Mixed_Subdivision;
! 498: newsols : Solution_List;
! 499: begin
! 500: -- put_line(file,"THE LIST OF NEW CELLS"); put(file,n,mix,newcells);
! 501: if Length_Of(newcells) = 1
! 502: then Mixed_Solve(file,hom,mix,newcells,newsols);
! 503: -- else Solve_by_Unfolding(file,p,hom,n,mix,normal,newcells,newsols);
! 504: else bigcell := Merge_Same_Normal(newcells);
! 505: -- put_line(file,"THE BIG CELL"); put(file,n,mix,bigcell);
! 506: Mixed_Solve(file,hom,mix,bigcell,newsols);
! 507: -- Deep_Clear(bigcell);
! 508: Shallow_Clear(bigcell);
! 509: end if;
! 510: sols_last := Pointer_to_Last(sols);
! 511: Concat(sols,sols_last,newsols);
! 512: Shallow_Clear(newsols);
! 513: Shallow_Clear(newcells);
! 514: end;
! 515: tmp := Tail_Of(tmp);
! 516: end loop;
! 517: Deep_Clear(newnor);
! 518: end Solve_with_Unfolding;
! 519:
! 520: -- DYNAMIC LIFTING FOR UNMIXED SYSTEMS :
! 521:
! 522: procedure Dynamic_Unmixed_Solve
! 523: ( file : in file_type; n : in natural;
! 524: l : in List; order,inter : in boolean; maxli : in natural;
! 525: lifted,lifted_last : in out List; t : in out Triangulation;
! 526: q : in Poly_Sys; qsols : in out Solution_List ) is
! 527:
! 528: p : Laur_Sys(q'range) := Polynomial_to_Laurent_System(q);
! 529: qt : Laur_Sys(q'range); -- the polyhedral homotopy
! 530: firstflat : boolean := true; -- first time flattening
! 531: flatnor : Link_to_Vector; -- the flat normal
! 532: qsols_last : Solution_List;
! 533: mix : Vector(1..1) := (1..1 => n);
! 534:
! 535: procedure Solve_Before_Flattening
! 536: ( t : in Triangulation; lifted : in List ) is
! 537:
! 538: -- DESCRIPTION :
! 539: -- Computes the new solutions, right before the subdivision
! 540: -- is flattened.
! 541:
! 542: newpts : List;
! 543: newcells : Mixed_Subdivision;
! 544:
! 545: begin
! 546: if firstflat
! 547: then newpts := lifted;
! 548: newcells := Convert(t);
! 549: else newcells := Non_Flattened_Cells(flatnor,t);
! 550: newpts := Non_Flattened_Points(lifted);
! 551: end if;
! 552: Update_Polyhedral_Homotopy(p,newpts,qt);
! 553: if not firstflat
! 554: then new_line(file);
! 555: put_line(file,"*** EXTENDING THE SOLUTIONS ***");
! 556: new_line(file);
! 557: Set_Continuation_Parameter(qsols,Create(0.0));
! 558: Mixed_Continuation(file,qt,flatnor.all,qsols);
! 559: end if;
! 560: Solve_with_Unfolding(file,p,qt,n,mix,newcells,qsols);
! 561: if not firstflat
! 562: then Shallow_Clear(newpts); Shallow_Clear(newcells);
! 563: else firstflat := false;
! 564: end if;
! 565: Flatten(qt);
! 566: end Solve_Before_Flattening;
! 567:
! 568: procedure S_Dynamic_Lifting is
! 569: new Dynamic_Triangulations.Dynamic_Lifting_with_Flat
! 570: ( Before_Flattening => Solve_Before_Flattening );
! 571:
! 572: begin
! 573: qt := Initialize_Polyhedral_Homotopy(n,lifted,p);
! 574: flatnor := new vector(1..n+1);
! 575: flatnor(1..n) := (1..n => 0);
! 576: flatnor(n+1) := 1;
! 577: S_Dynamic_Lifting(l,order,inter,maxli,lifted,lifted_last,t);
! 578: Solve_Before_Flattening(t,lifted);
! 579: Clear(flatnor); Clear(mix);
! 580: Clear(qt); Clear(p);
! 581: end Dynamic_Unmixed_Solve;
! 582:
! 583: -- DYNAMIC LIFTING FOR THE CAYLEY TRICK :
! 584:
! 585: procedure Dynamic_Cayley_Solve
! 586: ( file : in file_type; n : in natural; mix : in Vector;
! 587: supports : in Array_of_Lists; order,inter : in boolean;
! 588: maxli : in natural; lifted : in out Array_of_Lists;
! 589: mixsub : in out Mixed_Subdivision; numtri : out natural;
! 590: q : in Poly_Sys; qsols : in out Solution_List ) is
! 591:
! 592: p : Laur_Sys(q'range) := Polynomial_to_Laurent_System(q);
! 593: qt : Laur_Sys(q'range); -- the polyhedral homotopy
! 594: firstflat : boolean := true; -- first time flattening
! 595: flatnor : Link_to_Vector; -- the flat normal
! 596: qsols_last : Solution_List;
! 597:
! 598: procedure Solve_Before_Flattening
! 599: ( mixsub : in out Mixed_Subdivision;
! 600: lifted : in Array_of_Lists ) is
! 601:
! 602: -- DESCRIPTION :
! 603: -- Computes the new solutions, right before the subdivision
! 604: -- is flattened.
! 605:
! 606: newpts : Array_of_Lists(lifted'range);
! 607: newcells : Mixed_Subdivision;
! 608:
! 609: begin
! 610: if firstflat
! 611: then newpts := lifted;
! 612: newcells := mixsub;
! 613: else newcells := Non_Flattened_Cells(flatnor,mixsub);
! 614: newpts := Non_Flattened_Points(lifted);
! 615: end if;
! 616: Update_Polyhedral_Homotopy(p,newpts,mix,qt);
! 617: if not Is_Null(qsols)
! 618: then new_line(file);
! 619: put_line(file,"*** EXTENDING THE SOLUTIONS ***");
! 620: new_line(file);
! 621: Set_Continuation_Parameter(qsols,Create(0.0));
! 622: Mixed_Continuation(file,qt,flatnor.all,qsols);
! 623: end if;
! 624: if not Is_Null(newcells)
! 625: then
! 626: Solve_with_Unfolding(file,p,qt,n,mix,newcells,qsols);
! 627: if not firstflat
! 628: then Shallow_Clear(newpts); Shallow_Clear(newcells);
! 629: else firstflat := false;
! 630: end if;
! 631: end if;
! 632: Flatten(qt);
! 633: end Solve_Before_Flattening;
! 634:
! 635: procedure S_Dynamic_Lifting is
! 636: new Cayley_Trick.Dynamic_Cayley_with_Flat
! 637: ( Before_Flattening => Solve_Before_Flattening );
! 638:
! 639: begin
! 640: flatnor := new vector(1..n+1);
! 641: flatnor(1..n) := (1..n => 0);
! 642: flatnor(n+1) := 1;
! 643: S_Dynamic_Lifting(n,mix,supports,order,inter,maxli,lifted,mixsub,numtri);
! 644: Solve_Before_Flattening(mixsub,lifted);
! 645: Clear(flatnor);
! 646: end Dynamic_Cayley_Solve;
! 647:
! 648: -- DYNAMIC LIFTING FOR SEMI-MIXED SYSTEMS :
! 649:
! 650: procedure Dynamic_Mixed_Solve
! 651: ( file : in file_type; n : in natural;
! 652: mix : in Vector; supports : in Array_of_Lists;
! 653: order,inter,conmv : in boolean; maxli : in natural;
! 654: mixsub : in out Mixed_Subdivision;
! 655: fs : in out Face_Structures;
! 656: nbsucc,nbfail : in out Standard_Floating_Vectors.Vector;
! 657: q : in Poly_Sys; qsols : in out Solution_List ) is
! 658:
! 659: p : Laur_Sys(q'range) := Polynomial_to_Laurent_System(q);
! 660: qt : Laur_Sys(q'range); -- the polyhedral homotopy
! 661: firstflat : boolean := true; -- first time flattening
! 662: flatnor : Link_to_Vector; -- the flat normal
! 663: qsols_last : Solution_List;
! 664:
! 665: procedure Solve_Before_Flattening
! 666: ( mixsub : in out Mixed_Subdivision;
! 667: fs : in Face_Structures ) is
! 668:
! 669: -- DESCRIPTION :
! 670: -- Computes the new solutions, right before the subdivision
! 671: -- is flattened.
! 672:
! 673: newpts : Array_of_Lists(fs'range);
! 674: newcells : Mixed_Subdivision;
! 675: newsols : Solution_List;
! 676:
! 677: begin
! 678: if firstflat
! 679: then for i in fs'range loop
! 680: newpts(i) := fs(i).l;
! 681: end loop;
! 682: newcells := mixsub;
! 683: else newcells := Non_Flattened_Cells(flatnor,mixsub);
! 684: newpts := Non_Flattened_Points(fs);
! 685: end if;
! 686: Update_Polyhedral_Homotopy(p,newpts,mix,qt);
! 687: if not firstflat
! 688: then new_line(file);
! 689: put_line(file,"*** EXTENDING THE SOLUTIONS ***");
! 690: new_line(file);
! 691: Set_Continuation_Parameter(qsols,Create(0.0));
! 692: Mixed_Continuation(file,qt,flatnor.all,qsols);
! 693: end if;
! 694: Mixed_Solve(file,qt,mix,newcells,newsols);
! 695: qsols_last := Pointer_to_Last(qsols);
! 696: Concat(qsols,qsols_last,newsols);
! 697: Shallow_Clear(newsols);
! 698: if not firstflat
! 699: then Shallow_Clear(newpts); Shallow_Clear(newcells);
! 700: else firstflat := false;
! 701: end if;
! 702: Flatten(qt);
! 703: end Solve_Before_Flattening;
! 704:
! 705: procedure S_Dynamic_Lifting is
! 706: new Dynamic_Mixed_Subdivisions.Dynamic_Lifting_with_Flat
! 707: ( Before_Flattening => Solve_Before_Flattening );
! 708:
! 709: begin
! 710: qt := Initialize_Polyhedral_Homotopy(n,mix,fs,p);
! 711: flatnor := new vector(1..n+1);
! 712: flatnor(1..n) := (1..n => 0);
! 713: flatnor(n+1) := 1;
! 714: S_Dynamic_Lifting(n,mix,supports,order,inter,conmv,maxli,mixsub,fs,
! 715: nbsucc,nbfail);
! 716: Solve_Before_Flattening(mixsub,fs);
! 717: Clear(flatnor);
! 718: Clear(qt); Clear(p);
! 719: end Dynamic_Mixed_Solve;
! 720:
! 721: end Dynamic_Polyhedral_Continuation;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>