Annotation of OpenXM_contrib/PHC/Ada/Homotopy/drivers_for_reduction.adb, Revision 1.1.1.1
1.1 maekawa 1: with integer_io; use integer_io;
2: with Communications_with_User; use Communications_with_User;
3: with Timing_Package; use Timing_Package;
4: with Numbers_io; use Numbers_io;
5: with Standard_Natural_Vectors;
6: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
7: with Standard_Complex_Poly_Systems_io; use Standard_Complex_Poly_Systems_io;
8: with Reduction_of_Polynomial_Systems; use Reduction_of_Polynomial_Systems;
9: with Reduction_of_Nonsquare_Systems; use Reduction_of_Nonsquare_Systems;
10:
11: package body Drivers_for_Reduction is
12:
13: procedure Display_Info is
14:
15: i : array(1..12)of string(1..65);
16:
17: begin
18: i( 1):="The goal of reduction is to rewrite the system into an equivalent";
19: i( 2):="one (i.e.: the same finite solutions) that has a lower total";
20: i( 3):="degree, so that fewer solution paths need to be followed.";
21: i( 4):="Sometimes reduction can already detect whether a system has no";
22: i( 5):="solutions or an infinite number of solutions. ";
23: i( 6):=" We distinguish between linear and nonlinear reduction. The";
24: i( 7):="first type performs row-reduction on the coefficient matrix of";
25: i( 8):="the system. By nonlinear reduction, highest-degree monomials are";
26: i( 9):="eliminated by replacing a polynomial in the system by a";
27: i(10):="Subtraction-polynomial. This second type is more powerful, but";
28: i(11):="also more expensive. Bounds have to be set to limit the";
29: i(12):="combinatorial enumeration. ";
30: for k in i'range loop
31: put_line(i(k));
32: end loop;
33: end Display_Info;
34:
35: procedure Display_Menu ( exit_opt : in boolean; ans : in out character ) is
36:
37: m : array(0..3) of string(1..65);
38:
39: begin
40: m(0):=" 0 : No Reduction : leave the menu ";
41: m(1):=" 1 : Linear Reduction : triangulate coefficient matrix ";
42: m(2):=" 2 : Sparse Linear Reduction : diagonalize coefficient matrix ";
43: m(3):=" 3 : Nonlinear Reduction : S-polynomial combinations ";
44: loop
45: new_line;
46: put_line("MENU for Reducing Polynomial Systems :");
47: if exit_opt
48: then for i in m'range loop
49: put_line(m(i));
50: end loop;
51: put("Type 0, 1, 2, or 3 to select reduction, or i for info : ");
52: Ask_Alternative(ans,"0123i");
53: else for i in 1..m'last loop
54: put_line(m(i));
55: end loop;
56: put("Type 1 , 2, or 3 to select reduction, or i for info : ");
57: Ask_Alternative(ans,"123i");
58: end if;
59: if ans = 'i'
60: then new_line; Display_Info;
61: end if;
62: exit when ans /= 'i';
63: end loop;
64: end Display_Menu;
65:
66: procedure Rationalize ( p : in out Poly_Sys ) is
67:
68: -- DESCRIPTION :
69: -- Shortens the exponent vectors when an unknown disappears.
70:
71: n : natural := p'length;
72: to_rationalize : boolean;
73:
74: procedure rationalize ( p : in out Poly; k : in natural ) is
75:
76: procedure rat_term ( t : in out Term; cont : out boolean ) is
77:
78: tmp : Degrees := new Standard_Natural_Vectors.Vector'(1..n-1 => 0);
79:
80: begin
81: for i in t.dg'range loop
82: if i < k
83: then tmp(i) := t.dg(i);
84: elsif i > k
85: then tmp(i-1) := t.dg(i);
86: end if;
87: end loop;
88: Standard_Natural_Vectors.Clear
89: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
90: t.dg := tmp;
91: cont := true;
92: end rat_term;
93: procedure rat_terms is new Changing_Iterator(rat_term);
94:
95: begin
96: rat_terms(p);
97: end rationalize;
98:
99: begin
100: for i in reverse 1..n loop
101: to_rationalize := true; -- suppose the i-th unknown may disappear
102: for j in 1..n loop
103: if Degree(p(j),i) > 0
104: then to_rationalize := false;
105: end if;
106: exit when not to_rationalize;
107: end loop;
108: if to_rationalize
109: then for j in 1..n loop
110: rationalize(p(j),i);
111: end loop;
112: end if;
113: end loop;
114: end Rationalize;
115:
116: procedure Write_Diagnostics
117: ( file : in file_type; p : in Poly_Sys;
118: diagonal,inconsistent,infinite : in boolean;
119: d : out natural ) is
120:
121: -- DESCRIPTION :
122: -- Writes diagnostics after reduction.
123:
124: b : natural;
125:
126: begin
127: if not diagonal
128: then if inconsistent
129: then put_line(" Inconsistent system: no solutions.");
130: put_line(file," Inconsistent system: no solutions.");
131: elsif infinite
132: then put(" Probably an infinite number");
133: put_line(" of solutions.");
134: put(file," Probably an infinite number");
135: put_line(file," of solutions.");
136: else put(" The total degree is ");
137: put(file," The total degree is ");
138: b := Total_Degree(p);
139: put(b,1); put(file,b,1);
140: put_line("."); put_line(file,".");
141: d := b;
142: end if;
143: else put_line(" No initial terms could be eliminated.");
144: put_line(file," No initial terms could be eliminated.");
145: end if;
146: end Write_Diagnostics;
147:
148: procedure Write_Results ( file : in file_type; p : in Poly_Sys;
149: timer : in Timing_Widget; banner : in string ) is
150: begin
151: new_line(file);
152: put_line(file,"THE REDUCED SYSTEM :");
153: put(file,p);
154: new_line(file);
155: print_times(file,timer,banner);
156: end Write_Results;
157:
158: procedure Driver_for_Linear_Reduction
159: ( file : in file_type; p : in out Poly_Sys; d : out natural ) is
160:
161: timer : Timing_Widget;
162: diagonal,inconsistent,infinite : boolean := false;
163:
164: begin
165: new_line(file);
166: put_line(file,"LINEAR REDUCTION : ");
167: tstart(timer);
168: Reduce(p,diagonal,inconsistent,infinite);
169: tstop(timer);
170: Write_Diagnostics(file,p,diagonal,inconsistent,infinite,d);
171: Write_Results(file,p,timer,"Linear Reduction");
172: end Driver_for_Linear_Reduction;
173:
174: procedure Driver_for_Sparse_Linear_Reduction
175: ( file : in file_type; p : in out Poly_Sys; d : out natural ) is
176:
177: timer : Timing_Widget;
178: diagonal,inconsistent,infinite : boolean := false;
179:
180: begin
181: new_line(file);
182: put_line(file,"SPARSE LINEAR REDUCTION : ");
183: tstart(timer);
184: Sparse_Reduce(p,inconsistent,infinite);
185: tstop(timer);
186: Write_Diagnostics(file,p,diagonal,inconsistent,infinite,d);
187: Write_Results(file,p,timer,"Sparse Reduction");
188: end Driver_for_Sparse_Linear_Reduction;
189:
190: procedure Driver_for_Nonlinear_Reduction
191: ( file : in file_type;
192: p : in out Poly_Sys; d : out natural ) is
193:
194: res : Poly_Sys(p'range);
195: sparse : boolean;
196: cnt_eq,cnt_sp,cnt_rp,max_eq,max_sp,max_rp : natural;
197:
198: timer : Timing_Widget;
199: b : natural;
200:
201: begin
202: Rationalize(p);
203: Clear(res); Copy(p,res);
204:
205: new_line(file);
206: put_line(file,"NONLINEAR REDUCTION :");
207: put(" Give the limit on #equal degree replacements : ");
208: Read_Natural(max_eq); cnt_eq := 0;
209: put(" Give the limit on #computed S-polynomials : ");
210: Read_Natural(max_sp); cnt_sp := 0;
211: put(" Give the limit on #computed R-polynomials : ");
212: Read_Natural(max_rp); cnt_rp := 0;
213: put(file," The limit on #equal degree replacements : ");
214: put(file,max_eq,1); new_line(file);
215: put(file," The limit on #computed S-polynomials : ");
216: put(file,max_sp,1); new_line(file);
217: put(file," The limit on #computed R-polynomials : ");
218: put(file,max_rp,1); new_line(file);
219:
220: sparse := false;
221: tstart(timer);
222: if sparse
223: then Sparse_Reduce(p,res,cnt_eq,max_eq);
224: else Reduce(p,res,cnt_eq,max_eq,cnt_sp,max_sp,cnt_rp,max_rp);
225: end if;
226: tstop(timer);
227: Clear(p); Copy(res,p); Clear(res);
228: b := Total_Degree(p);
229: put("The total degree is "); put(b,1);
230: put(file,"The total degree is "); put(file,b,1);
231: put_line("."); put_line(file,".");
232: d := b;
233:
234: new_line(file);
235: put_line(file,"Amount of arithmetic work");
236: put(file," #equal replacements : ");
237: put(file,cnt_eq,4); new_line(file);
238: put(file," #computed S-polynomials : ");
239: put(file,cnt_sp,4); new_line(file);
240: put(file," #computed R-polynomials : ");
241: put(file,cnt_rp,4); new_line(file);
242: new_line(file);
243:
244: put_line(file,"The reduced system :");
245: put(file,p'length,p);
246: new_line(file);
247: print_times(file,timer,"Nonlinear Reduction");
248: new_line(file);
249: end Driver_for_Nonlinear_Reduction;
250:
251: procedure Driver_for_Overconstrained_Reduction
252: ( file : in file_type; p : in out Poly_Sys ) is
253:
254: ans : character;
255: n : constant natural := p'length;
256: m : constant natural := Number_of_Unknowns(p(p'first));
257: res : Poly_Sys(1..m);
258:
259: begin
260: new_line;
261: put_line("MENU for Overconstrained Reduction :");
262: put_line(" 0. No reduction : leave the menu.");
263: put_line(" 1. Add a random combination of the remainder to the first.");
264: put_line(" 2. Reduce the first equations with the remainder.");
265: put("Type 0, 1, or 2 to select reduction : ");
266: Ask_Alternative(ans,"012");
267: case ans is
268: when '1' => res := Random_Square(p);
269: when '2' => res := Reduced_Square(p);
270: when others => null;
271: end case;
272: if ans /= '0'
273: then for i in 1..m loop
274: Copy(res(i),p(i)); Clear(res(i));
275: end loop;
276: for i in m+1..n loop
277: Clear(p(i));
278: end loop;
279: end if;
280: end Driver_for_Overconstrained_Reduction;
281:
282: procedure Driver_for_Reduction
283: ( file : in file_type; p : in out Poly_Sys; d : out natural;
284: exit_option : in boolean ) is
285:
286: n : constant natural := p'length;
287: ans : character := '0';
288:
289: begin
290: Display_Menu(exit_option,ans);
291: case ans is
292: when '1' => Driver_for_Linear_Reduction(file,p,d);
293: when '2' => Driver_for_Sparse_Linear_Reduction(file,p,d);
294: when '3' => Driver_for_Sparse_Linear_Reduction(file,p,d);
295: Driver_for_Nonlinear_Reduction(file,p,d);
296: when others => null;
297: end case;
298: if Number_of_Unknowns(p(p'first)) < n
299: then Driver_for_Overconstrained_Reduction(file,p);
300: end if;
301: end Driver_for_Reduction;
302:
303: end Drivers_for_Reduction;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to drivers_for_reduction.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Homotopy