Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Dynlift/dynamic_polyhedral_continuation.adb, Revision 1.1.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>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to dynamic_polyhedral_continuation.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Root_Counts / Dynlift