Return to ts_diclp.adb CVS log

Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Supports

Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Supports/ts_diclp.adb, Revision 1.1.1.1

1.1       maekawa     1: with text_io,integer_io;                 use text_io,integer_io;
                      2: with Standard_Floating_Numbers;          use Standard_Floating_Numbers;
                      3: with Standard_Floating_Numbers_io;       use Standard_Floating_Numbers_io;
                      4: with Standard_Integer_Vectors;
                      5: with Standard_Floating_Vectors;          use Standard_Floating_Vectors;
                      6: with Standard_Floating_Matrices;         use Standard_Floating_Matrices;
                      7: with Dictionaries,Linear_Programming;    use Dictionaries,Linear_Programming;
                      8:
                      9: procedure ts_diclp is
                     10:
                     11: -- DESCRIPTION :
                     12: --   This procedure can be used for solving the primal problem:
                     13: --
                     14: --      max <c,x>
                     15: --          <a,x> <= 0
                     16: --
                     17: --   where x = (1,x1,x2,..,xn)
                     18: --
                     19: --   and for solving the dual problem :
                     20: --
                     21: --      min <c,x>
                     22: --          <a,x> >= 0
                     23: --
                     24: --   where x = (1,x1,x2,..,xn).
                     25:
                     26:
                     27:   n,m : natural;
                     28:   ans : character;
                     29:   eps : constant double_float := 10.0**(-12);
                     30:
                     31: -- OUTPUT ROUTINES :
                     32:
                     33:   procedure Write ( file : in file_type; dic : in Matrix;
                     34:                     in_bas,out_bas : in Standard_Integer_Vectors.Vector;
                     35:                     fo,af,ex : in natural ) is
                     36:
                     37:   -- DESCRIPTION :
                     38:   --   This procedure writes the dictionary on standard output or on file.
                     39:
                     40:   -- ON ENTRY :
                     41:   --   dic          matrix for the dictionary;
                     42:   --   in_bas       unknowns in the basis;
                     43:   --   out_bas      unknowns out the basis;
                     44:   --   fo           the number of decimal literals before the dot;
                     45:   --   af           the number of decimal literals after the dot;
                     46:   --   ex           the number of decimal literals in the exponent.
                     47:
                     48:   begin
                     49:     new_line(file);
                     50:     put_line(file,"*****  THE DICTIONARY  *****");
                     51:     new_line(file);
                     52:     put_line(file,"The elements in the basis : ");
                     53:     for i in in_bas'range loop
                     54:       put(file," "); put(file,in_bas(i),1);
                     55:     end loop;
                     56:     new_line(file);
                     57:     put_line(file,"The dictionary : ");
                     58:     for i in dic'range(1) loop
                     59:       for j in dic'range(2) loop
                     60:         put(file,dic(i,j), fore => fo, aft => af, exp => ex );
                     61:       end loop;
                     62:       new_line(file);
                     63:     end loop;
                     64:   end Write;
                     65:
                     66:   procedure Report ( dic : in Matrix;
                     67:                      in_bas,out_bas : in Standard_Integer_Vectors.Vector ) is
                     68:   begin
                     69:     Write(Standard_Output,dic,in_bas,out_bas,3,3,3);
                     70:   end Report;
                     71:
                     72:   procedure plp is new Generic_Primal_Simplex(Report);
                     73:   procedure dlp is new Generic_Dual_Simplex(Report);
                     74:
                     75: begin
                     76:
                     77: -- INPUT OF THE DATA :
                     78:
                     79:   put("Give the number of unknowns : "); get(n);
                     80:   put("Give the number of constraints : "); get(m);
                     81:
                     82:   declare
                     83:     c : Standard_Floating_Vectors.vector(0..n);
                     84:     a : matrix(1..m,0..n);
                     85:     dic : matrix(0..m,0..n);
                     86:     inbas : Standard_Integer_Vectors.Vector(1..m) := (1..m => 0);
                     87:     outbas : Standard_Integer_Vectors.Vector(1..n) := (1..n => 0);
                     88:     primsol : Standard_Floating_Vectors.Vector(1..n);
                     89:     dualsol : Standard_Floating_Vectors.Vector(1..m);
                     90:     nit : natural := 0;
                     91:     ok : boolean := false;
                     92:     feasibound : boolean;
                     93:   begin
                     94:
                     95:     put_line("Give the coefficients of the target :");
                     96:     for i in c'range loop
                     97:       get(c(i));
                     98:     end loop;
                     99:     put_line("Give the coefficients of the constraints :");
                    100:     for i in a'range(1) loop
                    101:       for j in a'range(2) loop
                    102:        get(a(i,j));
                    103:       end loop;
                    104:     end loop;
                    105:     loop
                    106:       put("Primal of Dual simplex problem ? (p/d) "); get(ans);
                    107:       if ans = 'p'
                    108:        then Dictionaries.Primal_Init(c,a,dic,inbas,outbas); ok := true;
                    109:        elsif ans = 'd'
                    110:           then Dictionaries.Dual_Init(c,a,dic,inbas,outbas); ok := true;
                    111:           else put_line("try again ...");  ok := false;
                    112:       end if;
                    113:       exit when ok;
                    114:     end loop;
                    115:
                    116:     if ans = 'p'
                    117:      then plp(dic,eps,inbas,outbas,nit,feasibound);
                    118:      else dlp(dic,eps,inbas,outbas,nit,feasibound);
                    119:     end if;
                    120:
                    121:     primsol := Dictionaries.Primal_Solution(dic,inbas,outbas);
                    122:     new_line;
                    123:     put_line("THE PRIMAL SOLUTION :");
                    124:     for i in primsol'range loop
                    125:       put("  x"); put(i,1); put(" : ");
                    126:       put(primsol(i),3,3,3); new_line;
                    127:     end loop;
                    128:
                    129:     if ans = 'd'
                    130:      then dualsol := -Dictionaries.Dual_Solution(dic,inbas,outbas);
                    131:      else dualsol :=  Dictionaries.Dual_Solution(dic,inbas,outbas);
                    132:     end if;
                    133:     new_line;
                    134:     put_line("THE DUAL SOLUTION :");
                    135:     for i in dualsol'range loop
                    136:       put("  x"); put(i,1); put(" : ");
                    137:       put(dualsol(i),3,3,3); new_line;
                    138:     end loop;
                    139:
                    140:     new_line;
                    141:     put("OPTIMUM : "); put(Optimum(dic)); new_line;
                    142:     put("NUMBER OF ITERATIONS : "); put(nit,1); new_line;
                    143:     if ans = 'd'
                    144:      then if feasibound
                    145:            then put_line("THE PROBLEM IS FEASIBLE");
                    146:            else put_line("THE PROBLEM IS NOT FEASIBLE");
                    147:           end if;
                    148:      else if feasibound
                    149:            then put_line("THE PROBLEM IS UNBOUNDED");
                    150:            else put_line("THE PROBLEM IS NOT UNBOUNDED");
                    151:           end if;
                    152:     end if;
                    153:
                    154:   end;
                    155:
                    156: end ts_diclp;