Return to drivers_for_path_directions.adb CVS log

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

Annotation of OpenXM_contrib/PHC/Ada/Continuation/drivers_for_path_directions.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 Standard_Floating_Numbers;          use Standard_Floating_Numbers;
                      5: with Standard_Floating_Numbers_io;       use Standard_Floating_Numbers_io;
                      6: with Numbers_io;                         use Numbers_io;
                      7: with Standard_Complex_Vectors;
                      8: with Standard_Complex_Norms_Equals;      use Standard_Complex_Norms_Equals;
                      9: with Standard_Complex_Matrices;          use Standard_Complex_Matrices;
                     10: with Standard_Complex_Poly_Systems;      use Standard_Complex_Poly_Systems;
                     11: with Standard_Complex_Poly_Systems_io;   use Standard_Complex_Poly_Systems_io;
                     12: with Homotopy;
                     13: with Projective_Transformations;         use Projective_Transformations;
                     14: with Continuation_Data;                  use Continuation_Data;
                     15: with Continuation_Parameters;            use Continuation_Parameters;
                     16: with Increment_and_Fix_Continuation;     use Increment_and_Fix_Continuation;
                     17:
                     18: package body Drivers_for_Path_Directions is
                     19:
                     20: -- AUXILIARIES TO INSTANTIATE :
                     21:
                     22:   function Scale ( v : Vector ) return Vector is
                     23:
                     24:   -- DESCRIPTION :
                     25:   --   Divides the vector by its largest component.
                     26:
                     27:     res : Vector(v'range);
                     28:     tol : constant double_float := 10.0**(-12);
                     29:     ind : integer := v'first;
                     30:     max : double_float := abs(v(ind));
                     31:
                     32:   begin
                     33:     for i in v'range loop
                     34:       if abs(v(i)) > max
                     35:        then max := abs(v(i));
                     36:             ind := i;
                     37:       end if;
                     38:     end loop;
                     39:     if max > tol
                     40:      then for i in v'range loop
                     41:             res(i) := v(i)/max;
                     42:           end loop;
                     43:     end if;
                     44:     return res;
                     45:   end Scale;
                     46:
                     47:   function Toric_Evaluator ( x : Standard_Complex_Vectors.Vector;
                     48:                              t : Complex_Number )
                     49:                            return Standard_Complex_Vectors.Vector is
                     50:   begin
                     51:     return Homotopy.Eval(x,t);
                     52:   end Toric_Evaluator;
                     53:
                     54:   function Toric_Differentiator
                     55:                ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
                     56:                return Standard_Complex_Vectors.Vector is
                     57:   begin
                     58:     return Homotopy.Eval(x,t);
                     59:   end Toric_Differentiator;
                     60:
                     61:   function Toric_Differentiator
                     62:                ( x : Standard_Complex_Vectors.Vector; t : Complex_Number )
                     63:                return Standard_Complex_Matrices.Matrix is
                     64:   begin
                     65:     return Homotopy.Diff(x,t);
                     66:   end Toric_Differentiator;
                     67:
                     68: -- TARGET ROUTINES :
                     69:
                     70:   procedure Init_Path_Directions
                     71:                ( n,nv : in natural; v : in out Link_to_VecVec;
                     72:                  errv : in out Link_to_Vector ) is
                     73:
                     74:   begin
                     75:     v := new VecVec(1..nv);
                     76:     for i in v'range loop
                     77:       v(i) := new Vector'(1..n => 0.0);
                     78:     end loop;
                     79:     errv := new Vector'(1..nv => 1.0);
                     80:   end Init_Path_Directions;
                     81:
                     82:   procedure Toric_Continue
                     83:                ( file : in file_type; sols : in out Solution_List;
                     84:                  proj,report : in boolean; v : in out VecVec;
                     85:                  errv : in out Vector; target : in Complex_Number ) is
                     86:
                     87:   -- DESCRIPTION :
                     88:   --   Performs the continuation with online toric compactifications.
                     89:
                     90:     timer : timing_widget;
                     91:     h : constant Poly_Sys := Homotopy.Eval(target); --Homotopy.Homotopy_System;
                     92:     n : constant natural := h'length;
                     93:     hh : Poly_Sys(h'range);
                     94:
                     95:     procedure Sil_Cont is
                     96:       new Silent_Toric_Continue(Max_Norm,Toric_Evaluator,
                     97:                                 Toric_Differentiator,Toric_Differentiator);
                     98:     procedure Rep_Cont is
                     99:       new Reporting_Toric_Continue(Max_Norm,Toric_Evaluator,
                    100:                                    Toric_Differentiator,Toric_Differentiator);
                    101:   begin
                    102:     tstart(timer);
                    103:     if report
                    104:      then Rep_Cont(file,sols,proj,v,errv,target);
                    105:      else Sil_Cont(sols,proj,v,errv,target);
                    106:     end if;
                    107:     tstop(timer);
                    108:     new_line(file);
                    109:     print_times(file,timer,"toric continuation"); new_line(file);
                    110:   end Toric_Continue;
                    111:
                    112:   procedure Write_Directions
                    113:                ( file : in file_type; v : in VecVec; errv : in Vector ) is
                    114:
                    115:     procedure Write ( file : in file_type; v : in Vector ) is
                    116:     begin
                    117:       for i in v'range loop
                    118:         put(file,v(i)); new_line(file);
                    119:       end loop;
                    120:     end Write;
                    121:
                    122:     procedure Write_Direction
                    123:                  ( file : in file_type;
                    124:                    v : in Vector; error : in double_float; i : natural ) is
                    125:
                    126:     begin
                    127:       put(file,"Computed direction of path ");
                    128:       put(file,i,1); put_line(file," :"); Write(file,v);
                    129:       put(file,"with error : "); put(file,error); new_line(file);
                    130:     end Write_Direction;
                    131:
                    132:   begin
                    133:     for i in v'range loop
                    134:       Write_Direction(file,v(i).all,errv(i),i);
                    135:     end loop;
                    136:   end Write_Directions;
                    137:
                    138: end Drivers_for_Path_Directions;