File: [local] / OpenXM_contrib / PHC / Ada / Root_Counts / Stalift / floating_polyhedral_continuation.adb (download)
Revision 1.1.1.1 (vendor branch), Sun Oct 29 17:45:29 2000 UTC (23 years, 10 months ago) by maekawa
Branch: PHC, MAIN
CVS Tags: v2, maekawa-ipv6, RELEASE_1_2_3, RELEASE_1_2_2_KNOPPIX_b, RELEASE_1_2_2_KNOPPIX, RELEASE_1_2_2, RELEASE_1_2_1, HEAD Changes since 1.1: +0 -0
lines
Import the second public release of PHCpack.
OKed by Jan Verschelde.
|
with integer_io; use integer_io;
with Standard_Floating_Numbers; use Standard_Floating_Numbers;
with Standard_Floating_Numbers_io; use Standard_Floating_Numbers_io;
with Standard_Mathematical_Functions; use Standard_Mathematical_Functions;
with Standard_Complex_Numbers; use Standard_Complex_Numbers;
with Standard_Floating_Vectors_io; use Standard_Floating_Vectors_io;
with Standard_Complex_Vectors;
with Standard_Complex_Vectors_io; use Standard_Complex_Vectors_io;
with Standard_Complex_Norms_Equals; use Standard_Complex_Norms_Equals;
with Standard_Integer_VecVecs;
with Standard_Floating_VecVecs;
with Standard_Complex_Matrices; use Standard_Complex_Matrices;
with Standard_Complex_Polynomials;
with Standard_Complex_Laur_Polys;
with Standard_Complex_Laur_Functions; use Standard_Complex_Laur_Functions;
with Standard_Complex_Poly_Systems; use Standard_Complex_Poly_Systems;
with Standard_Laur_Poly_Convertors; use Standard_Laur_Poly_Convertors;
with Standard_Complex_Substitutors; use Standard_Complex_Substitutors;
with Power_Lists; use Power_Lists;
with Lists_of_Floating_Vectors; use Lists_of_Floating_Vectors;
with Floating_Lifting_Utilities; use Floating_Lifting_Utilities;
with Floating_Integer_Convertors; use Floating_Integer_Convertors;
with Transforming_Integer_Vector_Lists; use Transforming_Integer_Vector_Lists;
with Transforming_Laurent_Systems; use Transforming_Laurent_Systems;
with Continuation_Parameters;
with Increment_and_Fix_Continuation; use Increment_and_Fix_Continuation;
with Fewnomial_System_Solvers; use Fewnomial_System_Solvers;
with BKK_Bound_Computations; use BKK_Bound_Computations;
package body Floating_Polyhedral_Continuation is
-- UTILITIES FOR POLYHEDRAL COEFFICIENT HOMOTOPY :
procedure put ( file : in file_type;
v : in Standard_Floating_Vectors.Vector;
fore,aft,exp : in natural ) is
begin
for i in v'range loop
put(file,v(i),fore,aft,exp);
end loop;
end put;
procedure put ( file : in file_type;
v : in Standard_Complex_Vectors.Vector;
fore,aft,exp : in natural ) is
begin
for i in v'range loop
put(file,REAL_PART(v(i)),fore,aft,exp);
put(file,IMAG_PART(v(i)),fore,aft,exp);
end loop;
end put;
function Power ( x,m : double_float ) return double_float is
-- DESCRIPTION :
-- Computes x**m for high powers of m to avoid RANGE_ERROR.
-- intm : integer := integer(m);
-- fltm : double_float;
begin
-- if m < 10.0
-- then
return (x**m); -- GNAT is better at this
-- else if double_float(intm) > m
-- then intm := intm-1;
-- end if;
-- fltm := m - double_float(intm);
-- return ((x**intm)*(x**fltm));
-- end if;
end Power;
function Minimum ( v : Standard_Floating_Vectors.Vector )
return double_float is
-- DESCRIPTION :
-- Returns the minimal element (/= 0) in the vector v.
tol : constant double_float := 10.0**(-7);
min : double_float := 0.0;
tmp : double_float;
begin
for i in v'range loop
tmp := abs(v(i));
if tmp > tol
then if (min = 0.0) or else (tmp < min)
then min := tmp;
end if;
end if;
end loop;
return min;
end Minimum;
function Minimum ( v : Standard_Floating_VecVecs.VecVec )
return double_float is
-- DESCRIPTION :
-- Returns the minimal nonzero element in the vectors of v.
tol : constant double_float := 10.0**(-7);
min : double_float := 0.0;
tmp : double_float;
begin
for i in v'range loop
tmp := Minimum(v(i).all);
if abs(tmp) > tol
then if (min = 0.0) or else (tmp < min)
then min := tmp;
end if;
end if;
end loop;
return min;
end Minimum;
function Scale ( v : Standard_Floating_Vectors.Vector; s : double_float )
return Standard_Floating_Vectors.Vector is
-- DESCRIPTION :
-- Returns the scaled vector, by dividing every element by s.
res : Standard_Floating_Vectors.Vector(v'range);
begin
for i in res'range loop
res(i) := v(i)/s;
end loop;
return res;
end Scale;
function Scale ( v : Standard_Floating_VecVecs.VecVec )
return Standard_Floating_VecVecs.VecVec is
-- DESCRIPTION :
-- Returns an array of scaled vectors.
res : Standard_Floating_VecVecs.VecVec(v'range);
min : constant double_float := Minimum(v);
tol : constant double_float := 10.0**(-8);
begin
if (abs(min) > tol) and (min /= 1.0)
then
for i in v'range loop
res(i) := new Standard_Floating_Vectors.Vector'(Scale(v(i).all,min));
end loop;
else
for i in v'range loop
res(i) := new Standard_Floating_Vectors.Vector'(v(i).all);
end loop;
end if;
return res;
end Scale;
procedure Shift ( v : in out Standard_Floating_Vectors.Vector ) is
-- DESCRIPTION :
-- Shifts the elements in v, such that the minimal element equals zero.
min : double_float := v(v'first);
begin
for i in v'first+1..v'last loop
if v(i) < min
then min := v(i);
end if;
end loop;
if min /= 0.0
then for i in v'range loop
v(i) := v(i) - min;
end loop;
end if;
end Shift;
function Create ( e : Standard_Integer_VecVecs.VecVec;
l : List; normal : Standard_Floating_Vectors.Vector )
return Standard_Floating_Vectors.Vector is
-- DESCRIPTION :
-- Returns a vector with all inner products of the normal with
-- the exponents in the list, such that minimal value equals zero.
res : Standard_Floating_Vectors.Vector(e'range);
use Standard_Floating_Vectors;
found : boolean;
lif : double_float;
begin
for i in e'range loop
declare
fei : constant Standard_Floating_Vectors.Vector := Convert(e(i).all);
begin
Search_Lifting(l,fei,found,lif);
if not found
then res(i) := 1000.0;
else res(i) := fei*normal(fei'range) + lif*normal(normal'last);
end if;
end;
end loop;
Shift(res);
return res;
end Create;
function Create ( e : Exponent_Vectors_Array;
l : Array_of_Lists; mix : Standard_Integer_Vectors.Vector;
normal : Standard_Floating_Vectors.Vector )
return Standard_Floating_VecVecs.VecVec is
res : Standard_Floating_VecVecs.VecVec(normal'first..normal'last-1);
cnt : natural := res'first;
begin
for i in mix'range loop
declare
rescnt : constant Standard_Floating_Vectors.Vector
:= Create(e(cnt).all,l(i),normal);
begin
res(cnt) := new Standard_Floating_Vectors.Vector(rescnt'range);
for j in rescnt'range loop
res(cnt)(j) := rescnt(j);
end loop;
end;
Shift(res(cnt).all);
for k in 1..(mix(i)-1) loop
res(cnt+k) := new Standard_Floating_Vectors.Vector(res(cnt)'range);
for j in res(cnt)'range loop
res(cnt+k)(j) := res(cnt)(j);
end loop;
end loop;
cnt := cnt + mix(i);
end loop;
return res;
end Create;
procedure Eval ( c : in Standard_Complex_Vectors.Vector;
t : in double_float; m : in Standard_Floating_Vectors.Vector;
ctm : in out Standard_Complex_Vectors.Vector ) is
-- DESCRIPTION : ctm = c*t**m.
begin
for i in ctm'range loop
-- ctm(i) := c(i)*Create((t**m(i)));
ctm(i) := c(i)*Create(Power(t,m(i)));
end loop;
end Eval;
procedure Eval ( c : in Standard_Complex_VecVecs.VecVec;
t : in double_float; m : in Standard_Floating_VecVecs.VecVec;
ctm : in out Standard_Complex_VecVecs.VecVec ) is
-- DESCRIPTION : ctm = c*t**m.
begin
for i in ctm'range loop
Eval(c(i).all,t,m(i).all,ctm(i).all);
end loop;
end Eval;
-- USEFUL PRIMITIVES :
function Is_In ( l : List;
d : Standard_Complex_Laur_Polys.Degrees )
return boolean is
-- DESCRIPTION :
-- Returns true if the degrees belong to the list l.
tmp : List := l;
pt : Standard_Floating_Vectors.Link_to_Vector;
fld : Standard_Floating_Vectors.Vector(d'range);
equ : boolean;
begin
for i in fld'range loop
fld(i) := double_float(d(i));
end loop;
while not Is_Null(tmp) loop
pt := Head_Of(tmp);
equ := true;
for i in fld'range loop
if pt(i) /= fld(i)
then equ := false;
end if;
exit when not equ;
end loop;
if equ
then return true;
else tmp := Tail_Of(tmp);
end if;
end loop;
return false;
end Is_In;
function Select_Terms ( p : Standard_Complex_Laur_Polys.Poly; l : List )
return Standard_Complex_Laur_Polys.Poly is
use Standard_Complex_Laur_Polys;
res : Poly := Null_Poly;
procedure Visit_Term ( t : in Term; cont : out boolean ) is
begin
if Is_In(l,t.dg)
then Add(res,t);
end if;
cont := true;
end Visit_Term;
procedure Visit_Terms is new Visiting_Iterator(Visit_Term);
begin
Visit_Terms(p);
return res;
end Select_Terms;
function Select_Subsystem
( p : Laur_Sys; mix : Standard_Integer_Vectors.Vector;
mic : Mixed_Cell ) return Laur_Sys is
-- DESCRIPTION :
-- Given a Laurent polynomial system and a mixed cell,
-- the corresponding subsystem will be returned.
-- ON ENTRY :
-- p a Laurent polynomial system;
-- mix type of mixture: occurencies of the supports;
-- mic a mixed cell.
-- REQUIRED :
-- The polynomials in p must be ordered according to the type of mixture.
res : Laur_Sys(p'range);
cnt : natural := 0;
begin
for k in mix'range loop
for l in 1..mix(k) loop
cnt := cnt + 1;
res(cnt) := Select_Terms(p(cnt),mic.pts(k));
end loop;
end loop;
return res;
end Select_Subsystem;
procedure Extract_Regular ( sols : in out Solution_List ) is
function To_Be_Removed ( flag : in natural ) return boolean is
begin
return ( flag /= 1 );
end To_Be_Removed;
procedure Extract_Regular_Solutions is new
Standard_Complex_Solutions.Delete(To_Be_Removed);
begin
Extract_Regular_Solutions(sols);
end Extract_Regular;
procedure Refine ( file : in file_type; p : in Laur_Sys;
sols : in out Solution_List ) is
-- DESCRIPTION :
-- Given a polyhedral homotopy p and a list of solution for t=1,
-- this list of solutions will be refined.
pp : Poly_Sys(p'range) := Laurent_to_Polynomial_System(p);
n : constant natural := p'length;
eps : constant double_float := 10.0**(-12);
tolsing : constant double_float := 10.0**(-8);
max : constant natural := 3;
numb : natural := 0;
begin
pp := Laurent_to_Polynomial_System(p);
Substitute(n+1,Create(1.0),pp);
-- Reporting_Root_Refiner(file,pp,sols,eps,eps,tolsing,numb,max,false);
Clear(pp); Extract_Regular(sols);
end Refine;
-- FIRST LAYER OF CONTINUATION ROUTINES :
procedure Mixed_Continuation
( mix : in Standard_Integer_Vectors.Vector;
lifted : in Array_of_Lists; h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
normal : in Standard_Floating_Vectors.Vector;
sols : in out Solution_List ) is
pow : Standard_Floating_VecVecs.VecVec(c'range)
:= Create(e,lifted,mix,normal);
scapow : Standard_Floating_VecVecs.VecVec(c'range) := Scale(pow);
ctm : Standard_Complex_VecVecs.VecVec(c'range);
use Standard_Floating_Vectors;
use Standard_Complex_Laur_Polys;
function Eval ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
return Standard_Complex_Vectors.Vector is
begin
Eval(c,REAL_PART(t),scapow,ctm);
return Eval(h,ctm,x);
end Eval;
function dHt ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
return Standard_Complex_Vectors.Vector is
res : Standard_Complex_Vectors.Vector(h'range);
xtl : constant integer := x'last+1;
begin
Eval(c,REAL_PART(t),scapow,ctm);
for i in res'range loop
res(i) := Eval(j(i,xtl),m(i,xtl).all,ctm(i).all,x);
end loop;
return res;
end dHt;
function dHx ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
return matrix is
mt : Matrix(x'range,x'range);
begin
Eval(c,REAL_PART(t),scapow,ctm);
for k in mt'range(1) loop
for l in mt'range(2) loop
mt(k,l) := Eval(j(k,l),m(k,l).all,ctm(k).all,x);
end loop;
end loop;
return mt;
end dHx;
procedure Laur_Cont is new Silent_Continue(Max_Norm,Eval,dHt,dHx);
begin
for i in c'range loop
ctm(i)
:= new Standard_Complex_Vectors.Vector'(c(i).all'range => Create(0.0));
end loop;
Laur_Cont(sols,false);
Standard_Floating_VecVecs.Clear(pow);
Standard_Floating_VecVecs.Clear(scapow);
Standard_Complex_VecVecs.Clear(ctm);
Extract_Regular(sols);
end Mixed_Continuation;
procedure Mixed_Continuation
( file : in file_type;
mix : in Standard_Integer_Vectors.Vector;
lifted : in Array_of_Lists; h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
normal : in Standard_Floating_Vectors.Vector;
sols : in out Solution_List ) is
pow : Standard_Floating_VecVecs.VecVec(c'range)
:= Create(e,lifted,mix,normal);
scapow : Standard_Floating_VecVecs.VecVec(c'range) := Scale(pow);
ctm : Standard_Complex_VecVecs.VecVec(c'range);
use Standard_Floating_Vectors;
use Standard_Complex_Laur_Polys;
function Eval ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
return Standard_Complex_Vectors.Vector is
begin
Eval(c,REAL_PART(t),scapow,ctm);
return Eval(h,ctm,x);
end Eval;
function dHt ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
return Standard_Complex_Vectors.Vector is
res : Standard_Complex_Vectors.Vector(h'range);
xtl : constant integer := x'last+1;
begin
Eval(c,REAL_PART(t),scapow,ctm);
for i in res'range loop
res(i) := Eval(j(i,xtl),m(i,xtl).all,ctm(i).all,x);
end loop;
return res;
end dHt;
function dHx ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
return Matrix is
mt : Matrix(x'range,x'range);
begin
Eval(c,REAL_PART(t),scapow,ctm);
for k in m'range(1) loop
for l in m'range(1) loop
mt(k,l) := Eval(j(k,l),m(k,l).all,ctm(k).all,x);
end loop;
end loop;
return mt;
end dHx;
procedure Laur_Cont is new Reporting_Continue(Max_Norm,Eval,dHt,dHx);
begin
-- put_line(file,"The coefficient vectors :" );
-- for i in c'range loop
-- put(file,c(i).all,3,3,3); new_line(file);
-- end loop;
for i in c'range loop
ctm(i)
:= new Standard_Complex_Vectors.Vector'(c(i).all'range => Create(0.0));
end loop;
-- put(file,"The normal : "); put(file,normal,3,3,3); new_line(file);
-- put_line(file,"The exponent vector : ");
-- for i in pow'range loop
-- put(file,pow(i).all,3,3,3); new_line(file);
-- end loop;
-- put_line(file,"The scaled exponent vector : ");
-- for i in pow'range loop
-- put(file,scapow(i).all,3,3,3); new_line(file);
-- end loop;
Laur_Cont(file,sols,false);
Standard_Floating_VecVecs.Clear(pow);
Standard_Floating_VecVecs.Clear(scapow);
Standard_Complex_VecVecs.Clear(ctm);
Extract_Regular(sols);
end Mixed_Continuation;
-- UTILITIES FOR SECOND LAYER :
function Remove_Lifting ( l : List ) return List is
-- DESCRIPTION :
-- Removes the lifting value from the vectors.
tmp,res,res_last : List;
begin
tmp := l;
while not Is_Null(tmp) loop
declare
d1 : constant Standard_Floating_Vectors.Vector := Head_Of(tmp).all;
d2 : constant Standard_Floating_Vectors.Vector
:= d1(d1'first..d1'last-1);
begin
Append(res,res_last,d2);
end;
tmp := Tail_Of(tmp);
end loop;
return res;
end Remove_Lifting;
function Sub_Lifting ( q : Laur_Sys; mix : Standard_Integer_Vectors.Vector;
mic : Mixed_Cell ) return Array_of_Lists is
-- DESCRIPTION :
-- Returns the lifting used to subdivide the cell.
res : Array_of_Lists(mix'range);
sup : Array_of_Lists(q'range);
n : constant natural := q'last;
cnt : natural := sup'first;
begin
for i in mic.pts'range loop
sup(cnt) := Remove_Lifting(mic.pts(i));
for j in 1..(mix(i)-1) loop
Copy(sup(cnt),sup(cnt+j));
end loop;
cnt := cnt + mix(i);
end loop;
res := Induced_Lifting(n,mix,sup,mic.sub.all);
Deep_Clear(sup);
return res;
end Sub_Lifting;
function Sub_Polyhedral_Homotopy
( l : List; e : Standard_Integer_VecVecs.VecVec;
c : Standard_Complex_Vectors.Vector )
return Standard_Complex_Vectors.Vector is
-- DESCRIPTION :
-- For every vector in e that does not belong to l, the corresponding
-- index in c will be set to zero, otherwise it is copied to the result.
res : Standard_Complex_Vectors.Vector(c'range);
found : boolean;
lif : double_float;
begin
for i in e'range loop
declare
fei : constant Standard_Floating_Vectors.Vector := Convert(e(i).all);
begin
Search_Lifting(l,fei,found,lif);
if not found
then res(i) := Create(0.0);
else res(i) := c(i);
end if;
end;
end loop;
return res;
end Sub_Polyhedral_Homotopy;
function Sub_Polyhedral_Homotopy
( mix : Standard_Integer_Vectors.Vector; mic : Mixed_Cell;
e : Exponent_Vectors_Array;
c : Standard_Complex_VecVecs.VecVec )
return Standard_Complex_VecVecs.VecVec is
-- DESCRIPTION :
-- Given a subsystem q of p and the coefficient vector of p, the
-- vector on return will have only nonzero entries for coefficients
-- that belong to q.
res : Standard_Complex_VecVecs.VecVec(c'range);
begin
for i in mix'range loop
declare
cri : constant Standard_Complex_Vectors.Vector
:= Sub_Polyhedral_Homotopy(mic.pts(i),e(i).all,c(i).all);
begin
res(i) := new Standard_Complex_Vectors.Vector'(cri);
for j in 1..(mix(i)-1) loop
declare
crj : constant Standard_Complex_Vectors.Vector
:= Sub_Polyhedral_Homotopy(mic.pts(i),e(i+j).all,c(i+j).all);
begin
res(i+j) := new Standard_Complex_Vectors.Vector'(crj);
end;
end loop;
end;
end loop;
return res;
end Sub_Polyhedral_Homotopy;
procedure Refined_Mixed_Solve
( q : in Laur_Sys; mix : in Standard_Integer_Vectors.Vector;
mic : in Mixed_Cell; h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
qsols : in out Solution_List ) is
-- DESCRIPTION :
-- Polyhedral coeffient-homotopy for subsystem q.
-- REQUIRED : mic.sub /= null.
lif : Array_of_Lists(mix'range) := Sub_Lifting(q,mix,mic);
cq : Standard_Complex_VecVecs.VecVec(c'range)
:= Sub_Polyhedral_Homotopy(mix,mic,e,c);
begin
Mixed_Solve(q,lif,h,cq,e,j,m,mix,mic.sub.all,qsols);
Standard_Complex_VecVecs.Clear(cq); Deep_Clear(lif);
end Refined_Mixed_Solve;
procedure Refined_Mixed_Solve
( file : in file_type; q : in Laur_Sys;
mix : in Standard_Integer_Vectors.Vector;
mic : in Mixed_Cell; h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
qsols : in out Solution_List ) is
-- DESCRIPTION :
-- Polyhedral coeffient-homotopy for subsystem q.
-- REQUIRED : mic.sub /= null.
lif : Array_of_Lists(mix'range) := Sub_Lifting(q,mix,mic);
cq : Standard_Complex_VecVecs.VecVec(c'range)
:= Sub_Polyhedral_Homotopy(mix,mic,e,c);
begin
Mixed_Solve(file,q,lif,h,cq,e,j,m,mix,mic.sub.all,qsols);
Standard_Complex_VecVecs.Clear(cq); Deep_Clear(lif);
end Refined_Mixed_Solve;
-- SECOND LAYER :
procedure Mixed_Solve
( p : in Laur_Sys; lifted : in Array_of_Lists;
h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
mix : in Standard_Integer_Vectors.Vector;
mic : in Mixed_Cell;
sols,sols_last : in out Solution_List ) is
q : Laur_Sys(p'range) := Select_Subsystem(p,mix,mic);
sq : Laur_Sys(q'range) := Shift(q);
pq : Poly_Sys(q'range);
qsols : Solution_List;
len : natural := 0;
fail : boolean;
begin
Solve(sq,qsols,fail);
if fail
then if mic.sub = null
then pq := Laurent_to_Polynomial_System(sq);
qsols := Solve_by_Static_Lifting(pq);
Clear(pq);
else Refined_Mixed_Solve(q,mix,mic,h,c,e,j,m,qsols);
end if;
Set_Continuation_Parameter(qsols,Create(0.0));
end if;
len := Length_Of(qsols);
if len > 0
then Mixed_Continuation(mix,lifted,h,c,e,j,m,mic.nor.all,qsols);
Concat(sols,sols_last,qsols);
end if;
Clear(sq); Clear(q); Clear(qsols);
end Mixed_Solve;
procedure Mixed_Solve
( file : in file_type;
p : in Laur_Sys; lifted : in Array_of_Lists;
h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
mix : in Standard_Integer_Vectors.Vector;
mic : in Mixed_Cell;
sols,sols_last : in out Solution_List ) is
q : Laur_Sys(p'range) := Select_Subsystem(p,mix,mic);
sq : Laur_Sys(q'range) := Shift(q);
pq : Poly_Sys(q'range);
qsols : Solution_List;
len : natural := 0;
fail : boolean;
begin
Solve(sq,qsols,fail);
if not fail
then put_line(file,"It is a fewnomial system.");
else put_line(file,"No fewnomial system.");
if mic.sub = null
then put_line(file,"Calling the black box solver.");
pq := Laurent_to_Polynomial_System(sq);
qsols := Solve_by_Static_Lifting(file,pq);
Clear(pq);
else put_line(file,"Using the refinement of the cell.");
Refined_Mixed_Solve(file,q,mix,mic,h,c,e,j,m,qsols);
end if;
Set_Continuation_Parameter(qsols,Create(0.0));
end if;
len := Length_Of(qsols);
put(file,len,1); put_line(file," solutions found.");
if len > 0
then
Mixed_Continuation(file,mix,lifted,h,c,e,j,m,mic.nor.all,qsols);
Concat(sols,sols_last,qsols);
end if;
Clear(sq); Clear(q); Clear(qsols);
end Mixed_Solve;
-- THIRD LAYER :
procedure Mixed_Solve
( p : in Laur_Sys; lifted : in Array_of_Lists;
h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
mix : in Standard_Integer_Vectors.Vector;
mixsub : in Mixed_Subdivision;
sols : in out Solution_List ) is
tmp : Mixed_Subdivision := mixsub;
mic : Mixed_Cell;
sols_last : Solution_List;
begin
while not Is_Null(tmp) loop
mic := Head_Of(tmp);
Mixed_Solve(p,lifted,h,c,e,j,m,mix,mic,sols,sols_last);
tmp := Tail_Of(tmp);
end loop;
end Mixed_Solve;
procedure Mixed_Solve
( file : in file_type;
p : in Laur_Sys; lifted : in Array_of_Lists;
h : in Eval_Coeff_Laur_Sys;
c : in Standard_Complex_VecVecs.VecVec;
e : in Exponent_Vectors_Array;
j : in Eval_Coeff_Jaco_Mat; m : in Mult_Factors;
mix : in Standard_Integer_Vectors.Vector;
mixsub : in Mixed_Subdivision;
sols : in out Solution_List ) is
tmp : Mixed_Subdivision := mixsub;
mic : Mixed_Cell;
sols_last : Solution_List;
cnt : natural := 0;
begin
while not Is_Null(tmp) loop
mic := Head_Of(tmp);
cnt := cnt + 1;
new_line(file);
put(file,"*** PROCESSING SUBSYSTEM "); put(file,cnt,1);
put_line(file," ***");
new_line(file);
Mixed_Solve(file,p,lifted,h,c,e,j,m,mix,mic,sols,sols_last);
tmp := Tail_Of(tmp);
end loop;
end Mixed_Solve;
end Floating_Polyhedral_Continuation;