Annotation of OpenXM_contrib/PHC/Ada/Continuation/drivers_for_poly_continuation.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 File_Scanning; use File_Scanning;
5: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
6: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
7: with Standard_Complex_Numbers_io; use Standard_Complex_Numbers_io;
8: with Numbers_io; use Numbers_io;
9: with Standard_Floating_Vectors; use Standard_Floating_Vectors;
10: with Standard_Floating_VecVecs; use Standard_Floating_VecVecs;
11: with Standard_Complex_Norms_Equals; use Standard_Complex_Norms_Equals;
12: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
13: with Symbol_Table,Symbol_Table_io; use Symbol_Table;
14: with Standard_Complex_Poly_Systems_io; use Standard_Complex_Poly_Systems_io;
15: with Standard_Complex_Solutions_io; use Standard_Complex_Solutions_io;
16: with Homotopy;
17: with Projective_Transformations; use Projective_Transformations;
18: with Drivers_for_Homotopy_Creation; use Drivers_for_Homotopy_Creation;
19: with Continuation_Parameters;
20: with Continuation_Parameters_io;
21: with Increment_and_Fix_Continuation; use Increment_and_Fix_Continuation;
22: with Process_io; use Process_io;
23: with Drivers_for_Path_Directions; use Drivers_for_Path_Directions;
24:
25: package body Drivers_for_Poly_Continuation is
26:
27: -- AUXILIARIES :
28:
29: procedure Continue ( file : in file_type; sols : in out Solution_List;
30: proj,report : in boolean;
31: target : in Complex_Number ) is
32:
33: -- DESCRIPTION :
34: -- Instantiates the path-trackers.
35:
36: timer : timing_widget;
37:
38: procedure Sil_Cont is
39: new Silent_Continue(Max_Norm,
40: Homotopy.Eval,Homotopy.Diff,Homotopy.Diff);
41: procedure Rep_Cont is
42: new Reporting_Continue(Max_Norm,
43: Homotopy.Eval,Homotopy.Diff,Homotopy.Diff);
44:
45: begin
46: tstart(timer);
47: if report
48: then Rep_Cont(file,sols,proj,target);
49: else Sil_Cont(sols,proj,target);
50: end if;
51: tstop(timer);
52: new_line(file); print_times(file,timer,"continuation");
53: end Continue;
54:
55: procedure Ask_Symbol is
56:
57: -- DESCRIPTION :
58: -- This procedure asks for the symbol to display the additional unknown.
59:
60: sb : Symbol;
61:
62: begin
63: put("Give symbol to display additional unknown : ");
64: sb := (sb'range => ' ');
65: Symbol_Table.Enlarge(1);
66: Symbol_Table_io.Get(sb);
67: Symbol_Table.Add(sb);
68: end Ask_Symbol;
69:
70: -- TARGET ROUTINES :
71:
72: procedure Driver_for_Process_io ( file : in file_type; oc : out natural ) is
73:
74: ans : character;
75: m : array(0..8) of string(1..65);
76:
77: begin
78: put_line("MENU for Output Information during Continuation : ");
79: m(0):=" 0 : no intermediate output information during continuation ";
80: m(1):=" 1 : only the final solutions at the end of the paths ";
81: m(2):=" 2 : intermediate solutions at each step along the paths ";
82: m(3):=" 3 : information of the predictor: t and step length ";
83: m(4):=" 4 : information of the corrector: corrections and residuals ";
84: m(5):=" 5 : intermediate solutions and information of the predictor ";
85: m(6):=" 6 : intermediate solutions and information of the corrector ";
86: m(7):=" 7 : information of predictor and corrector ";
87: m(8):=" 8 : intermediate solutions, info of predictor and corrector ";
88: for i in m'range loop
89: put_line(m(i));
90: end loop;
91: put("Type a number between 0 and 8 to select output information : ");
92: Ask_Alternative(ans,"012345678");
93: new_line(file);
94: put_line(file,"OUTPUT INFORMATION DURING CONTINUATION :");
95: case ans is
96: when '0' => Set_output_code(nil); oc := 0; put_line(file,m(0));
97: when '1' => Set_output_code(nil); oc := 1; put_line(file,m(1));
98: when '2' => Set_output_code(s); oc := 2; put_line(file,m(2));
99: when '3' => Set_output_code(p); oc := 3; put_line(file,m(3));
100: when '4' => Set_output_code(c); oc := 4; put_line(file,m(4));
101: when '5' => Set_output_code(sp); oc := 5; put_line(file,m(5));
102: when '6' => Set_output_code(sc); oc := 6; put_line(file,m(6));
103: when '7' => Set_output_code(pc); oc := 7; put_line(file,m(7));
104: when '8' => Set_output_code(spc); oc := 8; put_line(file,m(8));
105: when others => null;
106: end case;
107: end Driver_for_Process_io;
108:
109: procedure Driver_for_Continuation_Parameters ( file : in file_type ) is
110:
111: nb : natural := 0;
112:
113: procedure Begin_Banner ( file : in file_type ) is
114: begin
115: put_line(file,"****************** CURRENT CONTINUATION PARAMETERS "
116: & "*****************");
117: end Begin_Banner;
118:
119: procedure End_Banner ( file : in file_type ) is
120: begin
121: put_line(file,"***************************************************"
122: & "*****************");
123: end End_Banner;
124:
125: begin
126: loop
127: Begin_Banner(Standard_Output);
128: Continuation_Parameters_io.put;
129: End_Banner(Standard_Output);
130: Continuation_Parameters_io.get(nb);
131: exit when (nb = 0);
132: end loop;
133: new_line(file);
134: Begin_Banner(file);
135: Continuation_Parameters_io.put(file);
136: End_Banner(file);
137: end Driver_for_Continuation_Parameters;
138:
139: procedure Check_Continuation_Parameter
140: ( sols : in out Solution_List ) is
141:
142: ans : character;
143: tre,tim : double_float;
144:
145: begin
146: if not Is_Null(sols)
147: then
148: if Head_Of(sols).t = Create(1.0)
149: then put_line("The first solution has continuation parameter t = 1.0.");
150: put("Do you want to change t ? (y/n) "); Ask_Yes_or_No(ans);
151: if ans = 'y'
152: then put("Give real part of t : "); Read_Double_Float(tre);
153: put("Give imaginary part of t : "); Read_Double_Float(tim);
154: Set_Continuation_Parameter(sols,Create(tre,tim));
155: end if;
156: end if;
157: end if;
158: end Check_Continuation_Parameter;
159:
160: procedure Driver_for_Polynomial_Continuation
161: ( file : in file_type; p : in Poly_Sys;
162: sols : out Solution_List; target : out Complex_Number ) is
163:
164: infile : file_type;
165: pp,q : Poly_Sys(p'range);
166: t : Complex_Number;
167: qsols : Solution_List;
168: found,proj : boolean;
169:
170: procedure Read_Start_System is
171:
172: firstfail : boolean := true;
173: n : natural;
174:
175: begin
176: put_line("Reading the name of the file for start system.");
177: Read_Name_and_Open_File(infile);
178: get(infile,n,q);
179: exception
180: when others => put("The system on the file is not correct.");
181: if firstfail
182: then put_line(" Try again..."); Close(infile);
183: firstfail := false;
184: Read_Start_System;
185: end if;
186: end Read_Start_System;
187:
188: begin
189: new_line;
190: Read_Start_System;
191: put_line(file,"THE START SYSTEM : "); put(file,q); new_line(file);
192: Scan_and_Skip(infile,"SOLUTIONS",found);
193: if found
194: then get(infile,qsols);
195: else new_line; Read(qsols);
196: end if;
197: Close(infile);
198: Check_Continuation_Parameter(qsols);
199: put_line(file,"THE START SOLUTIONS : ");
200: put(file,Length_Of(qsols),Head_Of(qsols).n,qsols); new_line(file);
201: Copy(p,pp);
202: Driver_for_Homotopy_Construction(file,pp,q,qsols,t);
203: proj := (Number_of_Unknowns(q(q'first)) > q'last);
204: if proj
205: then Ask_Symbol;
206: end if;
207: new_line;
208: Driver_for_Polynomial_Continuation(file,qsols,proj,t);
209: -- Homotopy.Clear; --> clearing here creates difficulties for root refiner
210: sols := qsols;
211: target := t;
212: end Driver_for_Polynomial_Continuation;
213:
214: procedure Driver_for_Polynomial_Continuation
215: ( file : in file_type; p : in Poly_Sys; k : in natural;
216: target : in Complex_Number; sols : out Solution_list ) is
217:
218: qsols : Solution_List;
219:
220: begin
221: new_line; Read(qsols);
222: put_line(file,"THE START SOLUTIONS :");
223: put(file,Length_Of(qsols),Head_Of(qsols).n,qsols); new_line(file);
224: Homotopy.Create(p,k);
225: put_line(file,"HOMOTOPY PARAMETERS :");
226: put(file," k : "); put(file,k,2); new_line(file);
227: put(file," a : "); put(file,target); new_line(file);
228: Driver_for_Polynomial_Continuation(file,qsols,false,target);
229: -- Homotopy.Clear; --> clearing here creates difficulties for root refiner
230: sols := qsols;
231: end Driver_for_Polynomial_Continuation;
232:
233: procedure Driver_for_Polynomial_Continuation
234: ( file : in file_type; sols : in out Solution_List;
235: proj : in boolean;
236: target : Complex_Number := Create(1.0) ) is
237:
238: oc : natural;
239: timer : timing_widget;
240: report : boolean;
241: n : constant natural := Head_Of(sols).n;
242: nv : constant natural := Length_Of(sols);
243: v : Link_to_VecVec;
244: errv : Link_to_Vector;
245:
246: begin
247: new_line;
248: Driver_for_Continuation_Parameters(file);
249: if Continuation_Parameters.endext_order > 0
250: then Init_Path_Directions(n,nv,v,errv);
251: end if;
252: new_line;
253: Driver_for_Process_io(file,oc);
254: report := (oc /= 0);
255: new_line;
256: put_line("No more input expected. See output file for results.");
257: new_line;
258: if Continuation_Parameters.endext_order > 0
259: then Toric_Continue(file,sols,proj,report,v.all,errv.all,target);
260: Write_Directions(file,v.all,errv.all);
261: else Continue(file,sols,proj,report,target);
262: end if;
263: end Driver_for_Polynomial_Continuation;
264:
265: end Drivers_for_Poly_Continuation;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>