Return to standard_integer_linear_equalities.adb CVS log

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

Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Matrices/standard_integer_linear_equalities.adb, Revision 1.1

1.1     ! maekawa     1: with Standard_Common_Divisors;           use Standard_Common_Divisors;
        !             2: with Standard_Integer_Linear_Solvers;    use Standard_Integer_Linear_Solvers;
        !             3:
        !             4: package body Standard_Integer_Linear_Equalities is
        !             5:
        !             6:   procedure Triangulate ( l : in natural; m : in matrix;
        !             7:                           index : in natural; ineq : in out matrix ) is
        !             8:
        !             9:     a,b,lcmab,faca,facb : integer;
        !            10:
        !            11:   begin
        !            12:     if ineq(index,l) /= 0
        !            13:      then a := m(l,l); b := ineq(index,l);
        !            14:           lcmab := lcm(a,b);
        !            15:           if lcmab < 0 then lcmab := -lcmab; end if;
        !            16:           faca := lcmab/a;  facb := lcmab/b;
        !            17:           if facb < 0
        !            18:            then facb := -facb; faca := -faca;
        !            19:           end if;
        !            20:           for j in l..ineq'last(2) loop
        !            21:             ineq(index,j) := facb*ineq(index,j) - faca*m(l,j);
        !            22:           end loop;
        !            23:     end if;
        !            24:   end Triangulate;
        !            25:
        !            26:   procedure Triangulate ( l : in natural; m : in matrix;
        !            27:                           first,last : in natural; ineq : in out matrix ) is
        !            28:   begin
        !            29:     for i in first..last loop
        !            30:       if ineq(i,l) /= 0
        !            31:        then Triangulate(l,m,i,ineq);
        !            32:       end if;
        !            33:     end loop;
        !            34:   end Triangulate;
        !            35:
        !            36:   procedure Triangulate ( index,start : in natural; ineq : in out matrix ) is
        !            37:
        !            38:     column : natural := start;  -- current column counter
        !            39:     firstineq : natural := ineq'first(1);
        !            40:     found : boolean;
        !            41:     ind2 : natural;
        !            42:     a,b,lcmab,faca,facb : integer;
        !            43:
        !            44:   begin
        !            45:     --put_line("INEQUALTITY");
        !            46:     --put(" BEFORE UPDATE : ");
        !            47:     --for l in ineq'range(2) loop
        !            48:     --  put(ineq(index,l),1); put(' ');
        !            49:     --end loop;
        !            50:     --new_line;
        !            51:     loop
        !            52:      -- SEARCH FOR FIRST NONZERO ENTRY IN CURRENT INEQUALITY :
        !            53:       while column < ineq'last(2) and then ineq(index,column) = 0 loop
        !            54:         column := column + 1;
        !            55:       end loop;
        !            56:       exit when ((ineq(index,column) = 0) or else (column = ineq'last(2)));
        !            57:                                         -- nothing to eliminate
        !            58:      -- SEARCH FOR INEQUALITY,
        !            59:      --    WITH SAME FIRST NONZERO COLUMN, BUT WITH OPPOSITE SIGN :
        !            60:       found := false;
        !            61:       for k in firstineq..(index-1) loop
        !            62:         if ineq(index,column)*ineq(k,column) < 0 -- check for sign
        !            63:          then found := true;
        !            64:               for l in start..column-1 loop      -- check if same zero pattern
        !            65:                -- if ineq(index,l) = 0
        !            66:                --  then
        !            67:                 found := (ineq(k,l) = 0);
        !            68:                -- end if;
        !            69:                 exit when not found;
        !            70:               end loop;
        !            71:               if found then ind2 := k; end if;
        !            72:         end if;
        !            73:         exit when found;
        !            74:       end loop;
        !            75:       exit when not found;  -- no possibility for elimination
        !            76:      -- if found
        !            77:      --  then
        !            78:         -- ELIMINATE BY MAKING A POSITIVE LINEAR COMBINATION :
        !            79:          a := ineq(index,column);  b := ineq(ind2,column);
        !            80:          if a < 0
        !            81:           then lcmab := lcm(-a,b);
        !            82:                faca := lcmab/(-a); facb := lcmab/b;
        !            83:           else lcmab := lcm(a,-b);
        !            84:                facb := lcmab/(-b); faca := lcmab/a;
        !            85:          end if;
        !            86:          if ineq(index,ineq'last(2)) >= 0
        !            87:            or else  -- PRESERVE SIGN OF ineq(index,ineq'last(2)) !!!!
        !            88:                 (faca*ineq(index,ineq'last(2))
        !            89:                             + facb*ineq(ind2,ineq'last(2)) < 0)
        !            90:           then
        !            91:             for l in start..ineq'last(2) loop
        !            92:               ineq(index,l) := faca*ineq(index,l) + facb*ineq(ind2,l);
        !            93:             end loop;
        !            94:            -- PROCEED :
        !            95:             column := column + 1;
        !            96:             firstineq := ineq'first(1);
        !            97:           else
        !            98:            -- TRY TO ELIMINATE WITH OTHER INEQUALITIES :
        !            99:             firstineq := ind2 + 1;
        !           100:          end if;
        !           101:          if (firstineq >= index)
        !           102:            -- impossible to eliminate with sign preservation
        !           103:           then firstineq := ineq'first(2);
        !           104:                column := column + 1;
        !           105:          end if;
        !           106:      --  else
        !           107:      --    column := column + 1;
        !           108:      --    firstineq := ineq'first(2);
        !           109:      -- end if;
        !           110:       exit when (column >= ineq'last(2)-1);
        !           111:     end loop;
        !           112:     --put(" AFTER UPDATE : ");
        !           113:     --for l in ineq'range(2) loop
        !           114:     --  put(ineq(index,l),1); put(' ');
        !           115:     --end loop;
        !           116:     --new_line;
        !           117:   end Triangulate;
        !           118:
        !           119: end Standard_Integer_Linear_Equalities;