[BACK]Return to multprec_complex_linear_solvers.ads CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Matrices

Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Matrices/multprec_complex_linear_solvers.ads, Revision 1.1.1.1

1.1       maekawa     1: with Standard_Natural_Vectors;
                      2: with Multprec_Floating_Numbers;          use Multprec_Floating_Numbers;
                      3: with Multprec_Complex_Vectors;           use Multprec_Complex_Vectors;
                      4: with Multprec_Complex_Matrices;          use Multprec_Complex_Matrices;
                      5:
                      6: package Multprec_Complex_Linear_Solvers is
                      7:
                      8: -- DESCRIPTION :
                      9: --   This package offers a few routines to solve linear systems of equations.
                     10: --   The code for lufac, lufco and lusolve is a literal translation from the
                     11: --   f77-linpack code.
                     12:
                     13:   procedure Scale ( a : in out Matrix; b : in out Vector );
                     14:
                     15:   -- DESCRIPTION :
                     16:   --   Divides the ith equation in the system a*x = b by the largest
                     17:   --   element on the ith row of a, for i in a'range(1).
                     18:
                     19:   -- REQUIRED : a'range(1) = b'range(1).
                     20:
                     21:   procedure lufac ( a : in out Matrix; n : in integer;
                     22:                     ipvt : out Standard_Natural_Vectors.Vector;
                     23:                     info : out natural );
                     24:
                     25:   -- DESCRIPTION :
                     26:   --   lufac factors a complex matrix by gaussian elimination
                     27:
                     28:   --   lufac is usually called by lufco, but it can be called
                     29:   --   directly with a saving of time if rcond is not needed.
                     30:   --   (time for lufco) = (1 + 9/n)*(time for lufac).
                     31:
                     32:   -- ON ENTRY :
                     33:   --   a       complex matrix(1..n,1..n) to be factored
                     34:   --   n       the dimension of the matrix a
                     35:
                     36:   -- ON RETURN :
                     37:   --   a       an upper triangular matrix and the multipliers
                     38:   --           which were used to obtain it.
                     39:   --           The factorization can be written a = l*u where
                     40:   --           l is a product of permutation and unit lower
                     41:   --           triangular matrices and u is upper triangular.
                     42:   --   ipvt    an integer vector of pivot indices
                     43:   --   info    = 0  normal value
                     44:   --           = k  if u(k,k) = 0.0.
                     45:   --                This is not an error for this routine,
                     46:   --                but it does indicate that lusolve will
                     47:   --                divide by zero if called.  Use rcond in
                     48:   --                lufco for a reliable indication of singularity.
                     49:
                     50:   procedure lufco ( a : in out Matrix; n : in integer;
                     51:                     ipvt : out Standard_Natural_Vectors.Vector;
                     52:                     rcond : out Floating_Number );
                     53:
                     54:   -- DESCRIPTION :
                     55:   --   lufco factors a complex matrix by gaussian elimination
                     56:   --   and estimates the condition of the matrix.
                     57:
                     58:   -- If rcond is not needed, lufac is slightly faster.
                     59:   -- To solve a*x = b, follow lufco by lusolve.
                     60:
                     61:   -- ON ENTRY :
                     62:   --   a       complex matrix(1..n,1..n) to be factored
                     63:   --   n       the dimension of the matrix a
                     64:
                     65:   -- ON RETURN :
                     66:   --   a       an upper triangular matrix and the multipliers
                     67:   --           which are used to obtain it.
                     68:   --           The factorization can be written a = l*u, where
                     69:   --           l is a product of permutation and unit lower triangular
                     70:   --           matrices and u is upper triangular.
                     71:   --   ipvt    an integer vector of pivot indices
                     72:   --   rcond   an estimate of the reciprocal condition of a.
                     73:   --           For the system a*x = b, relative perturbations
                     74:   --           in a and b of size epsilon may cause relative
                     75:   --           perturbations in x of size epsilon/rcond.
                     76:   --           If rcond is so small that the logical expression
                     77:   --                  1.0 + rcond = 1.0
                     78:   --           is true, than a may be singular to working precision.
                     79:   --           In particular, rcond is zero if exact singularity is
                     80:   --           detected or the estimate underflows.
                     81:
                     82:   procedure lusolve ( a : in Matrix; n : in integer;
                     83:                       ipvt : in Standard_Natural_Vectors.Vector;
                     84:                       b : in out Vector );
                     85:
                     86:   -- DESCRIPTION :
                     87:   --   lusolve solves the complex system a*x = b using the factors
                     88:   --   computed by lufac or lufco
                     89:
                     90:   -- ON ENTRY :
                     91:   --   a       a complex matrix(1..n,1..n), the output from
                     92:   --           lufac or lufco
                     93:   --   n       the dimension of the matrix a
                     94:   --   ipvt    the pivot vector from lufac or lufco
                     95:   --   b       the right hand side vector
                     96:
                     97:   -- ON RETURN :
                     98:   --   b       the solution vector x
                     99:
                    100:   procedure Triangulate ( a : in out Matrix; n,m : in integer );
                    101:
                    102:   -- DESCRIPTION :
                    103:   --   triangulate makes the n*m complex matrix a triangular using
                    104:   --   Gaussian elimination.
                    105:
                    106:   -- ON ENTRY :
                    107:   --   a       a complex matrix(1..n,1..m)
                    108:   --   n       the number of rows of a
                    109:   --   m       the number of columns of a
                    110:
                    111:   -- ON RETURN :
                    112:   --   a       the triangulated matrix
                    113:
                    114:   procedure Diagonalize ( a : in out Matrix; n,m : in integer );
                    115:
                    116:   -- DESCRIPTION :
                    117:   --   diagonalize makes the n*m complex matrix a diagonal using
                    118:   --   Gauss-Jordan.
                    119:
                    120:   -- ON ENTRY :
                    121:   --   a       a complex matrix(1..n,1..m)
                    122:   --   n       the number of rows of a
                    123:   --   m       the number of columns of a
                    124:
                    125:   -- ON RETURN :
                    126:   --   a       the diagonalized matrix
                    127:
                    128: end Multprec_Complex_Linear_Solvers;

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>