Return to minor_computations.adb CVS log

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

Annotation of OpenXM_contrib/PHC/Ada/Schubert/minor_computations.adb, Revision 1.1

1.1     ! maekawa     1: with Standard_Floating_Numbers_io;       use Standard_Floating_Numbers_io;
        !             2: with Standard_Natural_Vectors;
        !             3: with Standard_Natural_Vectors_io;        use Standard_Natural_Vectors_io;
        !             4: with Standard_Floating_Linear_Solvers;   use Standard_Floating_Linear_Solvers;
        !             5:
        !             6: package body Minor_Computations is
        !             7:
        !             8:   function Determinant
        !             9:               ( mat : Matrix; rows,cols : Standard_Natural_Vectors.Vector )
        !            10:               return double_float is
        !            11:
        !            12:   -- DESCRIPTION :
        !            13:   --   Computes the determinant of the matrix selecting d rows and columns,
        !            14:   --   where d = rows'length = cols'length.
        !            15:
        !            16:     res : double_float := 1.0;
        !            17:     sqm : Matrix(rows'range,rows'range);
        !            18:     piv : Standard_Natural_Vectors.Vector(rows'range);
        !            19:     inf : natural;
        !            20:
        !            21:   begin
        !            22:     for i in rows'range loop
        !            23:       piv(i) := i;
        !            24:       for j in cols'range loop
        !            25:         sqm(i,j) := mat(rows(i),cols(j));
        !            26:       end loop;
        !            27:     end loop;
        !            28:     lufac(sqm,rows'last,piv,inf);
        !            29:     for i in rows'range loop
        !            30:       res := res*sqm(i,i);
        !            31:     end loop;
        !            32:     for i in piv'range loop
        !            33:       if piv(i) > i
        !            34:        then res := -res;
        !            35:       end if;
        !            36:     end loop;
        !            37:     return res;
        !            38:   end Determinant;
        !            39:
        !            40:   procedure Minors ( file : in file_type;
        !            41:                      n,m : in natural; mat : in Matrix ) is
        !            42:
        !            43:     procedure Minors ( d : in natural ) is
        !            44:
        !            45:       rows,cols : Standard_Natural_Vectors.Vector(1..d);
        !            46:
        !            47:       procedure Select_Rows ( k,start : in natural ) is
        !            48:
        !            49:         det : double_float;
        !            50:
        !            51:       begin
        !            52:         if k > d
        !            53:          then det := Determinant(mat,rows,cols);
        !            54:               put(file,"Minor ");
        !            55:               put(file,rows); put(file," X"); put(file,cols);
        !            56:               put(file," equals ");
        !            57:               put(file,det); new_line(file);
        !            58:          else for j in start..n loop
        !            59:                 rows(k) := j;
        !            60:                 Select_Rows(k+1,j+1);
        !            61:               end loop;
        !            62:         end if;
        !            63:       end Select_Rows;
        !            64:
        !            65:       procedure Select_Columns ( k,start : in natural ) is
        !            66:       begin
        !            67:         if k > d
        !            68:          then Select_Rows(1,1);
        !            69:          else for j in start..m loop
        !            70:                 cols(k) := j;
        !            71:                 Select_Columns(k+1,j+1);
        !            72:               end loop;
        !            73:         end if;
        !            74:       end Select_Columns;
        !            75:
        !            76:     begin
        !            77:       Select_Columns(1,1);
        !            78:     end Minors;
        !            79:
        !            80:   begin
        !            81:     for d in 2..m loop
        !            82:       Minors(d);
        !            83:     end loop;
        !            84:   end Minors;
        !            85:
        !            86:   function Number_of_Minors ( n,m : natural ) return natural is
        !            87:
        !            88:     res : natural := 0;
        !            89:
        !            90:     procedure Count_Minors ( d : in natural ) is
        !            91:
        !            92:       procedure Select_Rows ( k,start : in natural ) is
        !            93:       begin
        !            94:         if k > d
        !            95:          then res := res + 1;
        !            96:          else for j in start..n loop
        !            97:                 Select_Rows(k+1,j+1);
        !            98:               end loop;
        !            99:         end if;
        !           100:       end Select_Rows;
        !           101:
        !           102:       procedure Select_Columns ( k,start : in natural ) is
        !           103:       begin
        !           104:         if k > d
        !           105:          then Select_Rows(1,1);
        !           106:          else for j in start..m loop
        !           107:                 Select_Columns(k+1,j+1);
        !           108:               end loop;
        !           109:         end if;
        !           110:       end Select_Columns;
        !           111:
        !           112:     begin
        !           113:       Select_Columns(1,1);
        !           114:     end Count_Minors;
        !           115:
        !           116:   begin
        !           117:     for d in 2..m loop
        !           118:       Count_Minors(d);
        !           119:     end loop;
        !           120:     return res;
        !           121:   end Number_of_Minors;
        !           122:
        !           123:   function Sign_of_Minors ( n,m : natural; mat : Matrix ) return Vector is
        !           124:
        !           125:     res : Vector(1..Number_of_Minors(n,m));
        !           126:     tol : constant double_float := 10.0**(-10);
        !           127:     cnt : natural := 0;
        !           128:
        !           129:     procedure Compute_Minors ( d : in natural ) is
        !           130:
        !           131:       rows,cols : Standard_Natural_Vectors.Vector(1..d);
        !           132:
        !           133:       procedure Select_Rows ( k,start : in natural ) is
        !           134:
        !           135:         det : double_float;
        !           136:
        !           137:       begin
        !           138:         if k > d
        !           139:          then det := Determinant(mat,rows,cols);
        !           140:               cnt := cnt + 1;
        !           141:               if abs(det) < tol
        !           142:                then res(cnt) := 0;
        !           143:                elsif det > 0.0
        !           144:                    then res(cnt) := +1;
        !           145:                    else res(cnt) := -1;
        !           146:               end if;
        !           147:          else for j in start..n loop
        !           148:                 rows(k) := j;
        !           149:                 Select_Rows(k+1,j+1);
        !           150:               end loop;
        !           151:         end if;
        !           152:       end Select_Rows;
        !           153:
        !           154:       procedure Select_Columns ( k,start : in natural ) is
        !           155:       begin
        !           156:         if k > d
        !           157:          then Select_Rows(1,1);
        !           158:          else for j in start..m loop
        !           159:                 cols(k) := j;
        !           160:                 Select_Columns(k+1,j+1);
        !           161:               end loop;
        !           162:         end if;
        !           163:       end Select_Columns;
        !           164:
        !           165:     begin
        !           166:       Select_Columns(1,1);
        !           167:     end Compute_Minors;
        !           168:
        !           169:   begin
        !           170:     for d in 2..m loop
        !           171:       Compute_Minors(d);
        !           172:     end loop;
        !           173:     return res;
        !           174:   end Sign_of_Minors;
        !           175:
        !           176: end Minor_Computations;