[BACK]Return to dynamic_polyhedral_continuation.adb CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / PHC / Ada / Root_Counts / Dynlift

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>