Return to mainpoco.adb CVS log

Up to [local] / OpenXM_contrib / PHC / Ada / Continuation

Annotation of OpenXM_contrib/PHC/Ada/Continuation/mainpoco.adb, Revision 1.1.1.1

1.1       maekawa     1: with text_io,integer_io;                 use text_io,integer_io;
                      2: with Communications_with_User;           use Communications_with_User;
                      3: with Standard_Floating_Numbers;          use Standard_Floating_Numbers;
                      4: with Standard_Complex_Numbers;           use Standard_Complex_Numbers;
                      5: with Numbers_io;                         use Numbers_io;
                      6: with Standard_Complex_Polynomials;       use Standard_Complex_Polynomials;
                      7: with Standard_Complex_Poly_Systems;      use Standard_Complex_Poly_Systems;
                      8: with Standard_Complex_Poly_Systems_io;   use Standard_Complex_Poly_Systems_io;
                      9: with Standard_Complex_Poly_SysFun;       use Standard_Complex_Poly_SysFun;
                     10: with Homotopy;
                     11: with Standard_Complex_Solutions;         use Standard_Complex_Solutions;
                     12: with Standard_Complex_Solutions_io;      use Standard_Complex_Solutions_io;
                     13: with Projective_Transformations;         use Projective_Transformations;
                     14: with Standard_Root_Refiners;             use Standard_Root_Refiners;
                     15: with Drivers_for_Poly_Continuation;      use Drivers_for_Poly_Continuation;
                     16: --with Bye_Bye_Message;
                     17:
                     18: procedure mainpoco ( infilename,outfilename : in string ) is
                     19:
                     20:   solsft,outft : file_type;
                     21:   lp : Link_to_Poly_Sys;
                     22:   sols,refsols : Solution_List;
                     23:   artificial,solsfile : boolean;
                     24:   k,len : natural;
                     25:   ans : character;
                     26:   tarre,tarim : double_float;
                     27:   target : Complex_Number;
                     28:
                     29:   procedure Read_System ( filename : in string ) is
                     30:
                     31:     file : file_type;
                     32:
                     33:   begin
                     34:     if filename /= ""
                     35:      then Open(file,in_file,filename);
                     36:           get(file,lp);
                     37:           Close(file);
                     38:     end if;
                     39:   exception
                     40:     when others =>
                     41:       new_line;
                     42:       put("Could not open file with name "); put_line(filename);
                     43:       lp := null; return;
                     44:   end Read_System;
                     45:
                     46: begin
                     47:   Read_System(infilename);
                     48:   if lp = null
                     49:    then new_line; get(lp);
                     50:   end if;
                     51:   Create_Output_File(outft,outfilename);
                     52:   put(outft,lp.all); new_line(outft);
                     53:   new_line;
                     54:   put("Do you want the solutions on separate file ? (y/n) ");
                     55:   Ask_Yes_or_No(ans);
                     56:   if ans = 'y'
                     57:    then
                     58:      put_line("Reading the name of the file to write the solutions on.");
                     59:      Read_Name_and_Create_File(solsft);
                     60:      solsfile := true;
                     61:    else
                     62:      solsfile := false;
                     63:   end if;
                     64:   artificial := (Number_of_Unknowns(lp(lp'first)) = lp'last);
                     65:   if artificial
                     66:    then Driver_for_Polynomial_Continuation(outft,lp.all,sols,target);
                     67:    else new_line;
                     68:         put("Give the index of the parameter : "); Read_Natural(k);
                     69:         new_line;
                     70:         put_line("Reading the target value of the continuation parameter.");
                     71:         put("Give the real part of the target : "); Read_Double_Float(tarre);
                     72:         put("Give the imaginary part of the target : ");
                     73:         Read_Double_Float(tarim);
                     74:         target := Create(tarre,tarim);
                     75:         Driver_for_Polynomial_Continuation(outft,lp.all,k,target,sols);
                     76:   end if;
                     77:   if Length_Of(sols) > 0
                     78:    then declare
                     79:           epsxa,epsfa,tolsing : constant double_float := 10.0**(-8);
                     80:           nb : natural := 0;
                     81:         begin
                     82:           if artificial
                     83:            then
                     84:              if not Is_Null(sols) and then Head_Of(sols).n > lp'last
                     85:               then Affine_Transformation(sols);
                     86:              end if;
                     87:              if target = Create(1.0)
                     88:               then
                     89:                 if solsfile
                     90:                  then Reporting_Root_Refiner
                     91:                         (outft,lp.all,sols,refsols,epsxa,epsfa,tolsing,
                     92:                          nb,5,false);
                     93:                  else Reporting_Root_Refiner
                     94:                         (outft,lp.all,sols,epsxa,epsfa,tolsing,nb,5,false);
                     95:                 end if;
                     96:               else
                     97:                 declare
                     98:                   pt : Poly_Sys(lp'range);
                     99:                 begin
                    100:                   pt := Homotopy.Eval(target);
                    101:                   if solsfile
                    102:                    then Reporting_Root_Refiner
                    103:                           (outft,pt,sols,refsols,epsxa,epsfa,tolsing,
                    104:                            nb,5,false);
                    105:                    else Reporting_Root_Refiner
                    106:                           (outft,pt,sols,epsxa,epsfa,tolsing,nb,5,false);
                    107:                   end if;
                    108:                   Clear(pt);
                    109:                 end;
                    110:              end if;
                    111:            else
                    112:              declare
                    113:                pt : Poly_Sys(lp'range);
                    114:              begin
                    115:                pt := Eval(lp.all,target,k);
                    116:                if solsfile
                    117:                 then Reporting_Root_Refiner
                    118:                        (outft,pt,sols,refsols,epsxa,epsfa,tolsing,nb,5,false);
                    119:                 else Reporting_Root_Refiner
                    120:                        (outft,pt,sols,epsxa,epsfa,tolsing,nb,5,false);
                    121:                end if;
                    122:                Clear(pt);
                    123:              end;
                    124:           end if;
                    125:         end;
                    126:   end if;
                    127:  -- put(outft,Bye_Bye_Message);
                    128:   Close(outft);
                    129:   if solsfile
                    130:    then len := Length_Of(refsols);
                    131:         if len > 0
                    132:          then put(solsft,len,Head_Of(refsols).n,refsols);
                    133:         end if;
                    134:         Close(solsft);
                    135:   end if;
                    136: end mainpoco;