File: [local] / OpenXM_contrib / PHC / Ada / Homotopy / drivers_for_reduction.adb (download)
Revision 1.1.1.1 (vendor branch), Sun Oct 29 17:45:23 2000 UTC (23 years, 11 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 Communications_with_User; use Communications_with_User;
with Timing_Package; use Timing_Package;
with Numbers_io; use Numbers_io;
with Standard_Natural_Vectors;
with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
with Standard_Complex_Poly_Systems_io; use Standard_Complex_Poly_Systems_io;
with Reduction_of_Polynomial_Systems; use Reduction_of_Polynomial_Systems;
with Reduction_of_Nonsquare_Systems; use Reduction_of_Nonsquare_Systems;
package body Drivers_for_Reduction is
procedure Display_Info is
i : array(1..12)of string(1..65);
begin
i( 1):="The goal of reduction is to rewrite the system into an equivalent";
i( 2):="one (i.e.: the same finite solutions) that has a lower total";
i( 3):="degree, so that fewer solution paths need to be followed.";
i( 4):="Sometimes reduction can already detect whether a system has no";
i( 5):="solutions or an infinite number of solutions. ";
i( 6):=" We distinguish between linear and nonlinear reduction. The";
i( 7):="first type performs row-reduction on the coefficient matrix of";
i( 8):="the system. By nonlinear reduction, highest-degree monomials are";
i( 9):="eliminated by replacing a polynomial in the system by a";
i(10):="Subtraction-polynomial. This second type is more powerful, but";
i(11):="also more expensive. Bounds have to be set to limit the";
i(12):="combinatorial enumeration. ";
for k in i'range loop
put_line(i(k));
end loop;
end Display_Info;
procedure Display_Menu ( exit_opt : in boolean; ans : in out character ) is
m : array(0..3) of string(1..65);
begin
m(0):=" 0 : No Reduction : leave the menu ";
m(1):=" 1 : Linear Reduction : triangulate coefficient matrix ";
m(2):=" 2 : Sparse Linear Reduction : diagonalize coefficient matrix ";
m(3):=" 3 : Nonlinear Reduction : S-polynomial combinations ";
loop
new_line;
put_line("MENU for Reducing Polynomial Systems :");
if exit_opt
then for i in m'range loop
put_line(m(i));
end loop;
put("Type 0, 1, 2, or 3 to select reduction, or i for info : ");
Ask_Alternative(ans,"0123i");
else for i in 1..m'last loop
put_line(m(i));
end loop;
put("Type 1 , 2, or 3 to select reduction, or i for info : ");
Ask_Alternative(ans,"123i");
end if;
if ans = 'i'
then new_line; Display_Info;
end if;
exit when ans /= 'i';
end loop;
end Display_Menu;
procedure Rationalize ( p : in out Poly_Sys ) is
-- DESCRIPTION :
-- Shortens the exponent vectors when an unknown disappears.
n : natural := p'length;
to_rationalize : boolean;
procedure rationalize ( p : in out Poly; k : in natural ) is
procedure rat_term ( t : in out Term; cont : out boolean ) is
tmp : Degrees := new Standard_Natural_Vectors.Vector'(1..n-1 => 0);
begin
for i in t.dg'range loop
if i < k
then tmp(i) := t.dg(i);
elsif i > k
then tmp(i-1) := t.dg(i);
end if;
end loop;
Standard_Natural_Vectors.Clear
(Standard_Natural_Vectors.Link_to_Vector(t.dg));
t.dg := tmp;
cont := true;
end rat_term;
procedure rat_terms is new Changing_Iterator(rat_term);
begin
rat_terms(p);
end rationalize;
begin
for i in reverse 1..n loop
to_rationalize := true; -- suppose the i-th unknown may disappear
for j in 1..n loop
if Degree(p(j),i) > 0
then to_rationalize := false;
end if;
exit when not to_rationalize;
end loop;
if to_rationalize
then for j in 1..n loop
rationalize(p(j),i);
end loop;
end if;
end loop;
end Rationalize;
procedure Write_Diagnostics
( file : in file_type; p : in Poly_Sys;
diagonal,inconsistent,infinite : in boolean;
d : out natural ) is
-- DESCRIPTION :
-- Writes diagnostics after reduction.
b : natural;
begin
if not diagonal
then if inconsistent
then put_line(" Inconsistent system: no solutions.");
put_line(file," Inconsistent system: no solutions.");
elsif infinite
then put(" Probably an infinite number");
put_line(" of solutions.");
put(file," Probably an infinite number");
put_line(file," of solutions.");
else put(" The total degree is ");
put(file," The total degree is ");
b := Total_Degree(p);
put(b,1); put(file,b,1);
put_line("."); put_line(file,".");
d := b;
end if;
else put_line(" No initial terms could be eliminated.");
put_line(file," No initial terms could be eliminated.");
end if;
end Write_Diagnostics;
procedure Write_Results ( file : in file_type; p : in Poly_Sys;
timer : in Timing_Widget; banner : in string ) is
begin
new_line(file);
put_line(file,"THE REDUCED SYSTEM :");
put(file,p);
new_line(file);
print_times(file,timer,banner);
end Write_Results;
procedure Driver_for_Linear_Reduction
( file : in file_type; p : in out Poly_Sys; d : out natural ) is
timer : Timing_Widget;
diagonal,inconsistent,infinite : boolean := false;
begin
new_line(file);
put_line(file,"LINEAR REDUCTION : ");
tstart(timer);
Reduce(p,diagonal,inconsistent,infinite);
tstop(timer);
Write_Diagnostics(file,p,diagonal,inconsistent,infinite,d);
Write_Results(file,p,timer,"Linear Reduction");
end Driver_for_Linear_Reduction;
procedure Driver_for_Sparse_Linear_Reduction
( file : in file_type; p : in out Poly_Sys; d : out natural ) is
timer : Timing_Widget;
diagonal,inconsistent,infinite : boolean := false;
begin
new_line(file);
put_line(file,"SPARSE LINEAR REDUCTION : ");
tstart(timer);
Sparse_Reduce(p,inconsistent,infinite);
tstop(timer);
Write_Diagnostics(file,p,diagonal,inconsistent,infinite,d);
Write_Results(file,p,timer,"Sparse Reduction");
end Driver_for_Sparse_Linear_Reduction;
procedure Driver_for_Nonlinear_Reduction
( file : in file_type;
p : in out Poly_Sys; d : out natural ) is
res : Poly_Sys(p'range);
sparse : boolean;
cnt_eq,cnt_sp,cnt_rp,max_eq,max_sp,max_rp : natural;
timer : Timing_Widget;
b : natural;
begin
Rationalize(p);
Clear(res); Copy(p,res);
new_line(file);
put_line(file,"NONLINEAR REDUCTION :");
put(" Give the limit on #equal degree replacements : ");
Read_Natural(max_eq); cnt_eq := 0;
put(" Give the limit on #computed S-polynomials : ");
Read_Natural(max_sp); cnt_sp := 0;
put(" Give the limit on #computed R-polynomials : ");
Read_Natural(max_rp); cnt_rp := 0;
put(file," The limit on #equal degree replacements : ");
put(file,max_eq,1); new_line(file);
put(file," The limit on #computed S-polynomials : ");
put(file,max_sp,1); new_line(file);
put(file," The limit on #computed R-polynomials : ");
put(file,max_rp,1); new_line(file);
sparse := false;
tstart(timer);
if sparse
then Sparse_Reduce(p,res,cnt_eq,max_eq);
else Reduce(p,res,cnt_eq,max_eq,cnt_sp,max_sp,cnt_rp,max_rp);
end if;
tstop(timer);
Clear(p); Copy(res,p); Clear(res);
b := Total_Degree(p);
put("The total degree is "); put(b,1);
put(file,"The total degree is "); put(file,b,1);
put_line("."); put_line(file,".");
d := b;
new_line(file);
put_line(file,"Amount of arithmetic work");
put(file," #equal replacements : ");
put(file,cnt_eq,4); new_line(file);
put(file," #computed S-polynomials : ");
put(file,cnt_sp,4); new_line(file);
put(file," #computed R-polynomials : ");
put(file,cnt_rp,4); new_line(file);
new_line(file);
put_line(file,"The reduced system :");
put(file,p'length,p);
new_line(file);
print_times(file,timer,"Nonlinear Reduction");
new_line(file);
end Driver_for_Nonlinear_Reduction;
procedure Driver_for_Overconstrained_Reduction
( file : in file_type; p : in out Poly_Sys ) is
ans : character;
n : constant natural := p'length;
m : constant natural := Number_of_Unknowns(p(p'first));
res : Poly_Sys(1..m);
begin
new_line;
put_line("MENU for Overconstrained Reduction :");
put_line(" 0. No reduction : leave the menu.");
put_line(" 1. Add a random combination of the remainder to the first.");
put_line(" 2. Reduce the first equations with the remainder.");
put("Type 0, 1, or 2 to select reduction : ");
Ask_Alternative(ans,"012");
case ans is
when '1' => res := Random_Square(p);
when '2' => res := Reduced_Square(p);
when others => null;
end case;
if ans /= '0'
then for i in 1..m loop
Copy(res(i),p(i)); Clear(res(i));
end loop;
for i in m+1..n loop
Clear(p(i));
end loop;
end if;
end Driver_for_Overconstrained_Reduction;
procedure Driver_for_Reduction
( file : in file_type; p : in out Poly_Sys; d : out natural;
exit_option : in boolean ) is
n : constant natural := p'length;
ans : character := '0';
begin
Display_Menu(exit_option,ans);
case ans is
when '1' => Driver_for_Linear_Reduction(file,p,d);
when '2' => Driver_for_Sparse_Linear_Reduction(file,p,d);
when '3' => Driver_for_Sparse_Linear_Reduction(file,p,d);
Driver_for_Nonlinear_Reduction(file,p,d);
when others => null;
end case;
if Number_of_Unknowns(p(p'first)) < n
then Driver_for_Overconstrained_Reduction(file,p);
end if;
end Driver_for_Reduction;
end Drivers_for_Reduction;