linalg.texi

来自「math library from gnu」· TEXI 代码 · 共 1,201 行 · 第 1/4 页

@end deftypefun

@node Tridiagonal Systems
@section Tridiagonal Systems
@cindex tridiagonal systems

The functions described in this section efficiently solve symmetric,
non-symmetric and cyclic tridiagonal systems with minimal storage.
Note that the current implementations of these functions use a variant
of Cholesky decomposition, so the tridiagonal matrix must be positive
definite.  For non-positive definite matrices, the functions return
the error code @code{GSL_ESING}.

@deftypefun int gsl_linalg_solve_tridiag (const gsl_vector * @var{diag}, const gsl_vector * @var{e}, const gsl_vector * @var{f}, const gsl_vector * @var{b}, gsl_vector * @var{x})
This function solves the general @math{N}-by-@math{N} system @math{A x =
b} where @var{A} is tridiagonal (@c{$N\geq 2$}
@math{N >= 2}). The super-diagonal and
sub-diagonal vectors @var{e} and @var{f} must be one element shorter
than the diagonal vector @var{diag}.  The form of @var{A} for the 4-by-4
case is shown below,
@tex
\beforedisplay
$$
A = \pmatrix{d_0&e_0&  0& 0\cr
             f_0&d_1&e_1& 0\cr
             0  &f_1&d_2&e_2\cr 
             0  &0  &f_2&d_3\cr}
$$
\afterdisplay
@end tex
@ifinfo

@example
A = ( d_0 e_0  0   0  )
    ( f_0 d_1 e_1  0  )
    (  0  f_1 d_2 e_2 )
    (  0   0  f_2 d_3 )
@end example
@end ifinfo
@noindent
@end deftypefun

@deftypefun int gsl_linalg_solve_symm_tridiag (const gsl_vector * @var{diag}, const gsl_vector * @var{e}, const gsl_vector * @var{b}, gsl_vector * @var{x})
This function solves the general @math{N}-by-@math{N} system @math{A x =
b} where @var{A} is symmetric tridiagonal (@c{$N\geq 2$}
@math{N >= 2}).  The off-diagonal vector
@var{e} must be one element shorter than the diagonal vector @var{diag}.
The form of @var{A} for the 4-by-4 case is shown below,
@tex
\beforedisplay
$$
A = \pmatrix{d_0&e_0&  0& 0\cr
             e_0&d_1&e_1& 0\cr
             0  &e_1&d_2&e_2\cr 
             0  &0  &e_2&d_3\cr}
$$
\afterdisplay
@end tex
@ifinfo

@example
A = ( d_0 e_0  0   0  )
    ( e_0 d_1 e_1  0  )
    (  0  e_1 d_2 e_2 )
    (  0   0  e_2 d_3 )
@end example
@end ifinfo
@end deftypefun

@deftypefun int gsl_linalg_solve_cyc_tridiag (const gsl_vector * @var{diag}, const gsl_vector * @var{e}, const gsl_vector * @var{f}, const gsl_vector * @var{b}, gsl_vector * @var{x})
This function solves the general @math{N}-by-@math{N} system @math{A x =
b} where @var{A} is cyclic tridiagonal (@c{$N\geq 3$}
@math{N >= 3}).  The cyclic super-diagonal and
sub-diagonal vectors @var{e} and @var{f} must have the same number of
elements as the diagonal vector @var{diag}.  The form of @var{A} for the
4-by-4 case is shown below,
@tex
\beforedisplay
$$
A = \pmatrix{d_0&e_0& 0 &f_3\cr
             f_0&d_1&e_1& 0 \cr
              0 &f_1&d_2&e_2\cr 
             e_3& 0 &f_2&d_3\cr}
$$
\afterdisplay
@end tex
@ifinfo

@example
A = ( d_0 e_0  0  f_3 )
    ( f_0 d_1 e_1  0  )
    (  0  f_1 d_2 e_2 )
    ( e_3  0  f_2 d_3 )
@end example
@end ifinfo
@end deftypefun


@deftypefun int gsl_linalg_solve_symm_cyc_tridiag (const gsl_vector * @var{diag}, const gsl_vector * @var{e}, const gsl_vector * @var{b}, gsl_vector * @var{x})
This function solves the general @math{N}-by-@math{N} system @math{A x =
b} where @var{A} is symmetric cyclic tridiagonal (@c{$N\geq 3$}
@math{N >= 3}).  The cyclic
off-diagonal vector @var{e} must have the same number of elements as the
diagonal vector @var{diag}.  The form of @var{A} for the 4-by-4 case is
shown below,
@tex
\beforedisplay
$$
A = \pmatrix{d_0&e_0& 0 &e_3\cr
             e_0&d_1&e_1& 0 \cr
              0 &e_1&d_2&e_2\cr 
             e_3& 0 &e_2&d_3\cr}
$$
\afterdisplay
@end tex
@ifinfo

@example
A = ( d_0 e_0  0  e_3 )
    ( e_0 d_1 e_1  0  )
    (  0  e_1 d_2 e_2 )
    ( e_3  0  e_2 d_3 )
@end example
@end ifinfo
@end deftypefun

@node Balancing
@section Balancing
@cindex balancing matrices

The process of balancing a matrix applies similarity transformations
to make the rows and columns have comparable norms. This is
useful, for example, to reduce roundoff errors in the solution
of eigenvalue problems. Balancing a matrix @math{A} consists
of replacing @math{A} with a similar matrix
@tex
\beforedisplay
$$
A' = D^{-1} A D
$$
\afterdisplay
@end tex
@ifinfo

@example
A' = D^(-1) A D
@end example

@end ifinfo
where @math{D} is a diagonal matrix whose entries are powers
of the floating point radix.

@deftypefun int gsl_linalg_balance_matrix (gsl_matrix * @var{A}, gsl_vector * @var{D})
This function replaces the matrix @var{A} with its balanced counterpart
and stores the diagonal elements of the similarity transformation
into the vector @var{D}.
@end deftypefun

@node Linear Algebra Examples
@section Examples

The following program solves the linear system @math{A x = b}. The
system to be solved is,
@tex
\beforedisplay
$$
\left(
\matrix{0.18& 0.60& 0.57& 0.96\cr
0.41& 0.24& 0.99& 0.58\cr
0.14& 0.30& 0.97& 0.66\cr
0.51& 0.13& 0.19& 0.85}
\right)
\left(
\matrix{x_0\cr
x_1\cr
x_2\cr
x_3}
\right)
=
\left(
\matrix{1.0\cr
2.0\cr
3.0\cr
4.0}
\right)
$$
\afterdisplay
@end tex
@ifinfo

@example
[ 0.18 0.60 0.57 0.96 ] [x0]   [1.0]
[ 0.41 0.24 0.99 0.58 ] [x1] = [2.0]
[ 0.14 0.30 0.97 0.66 ] [x2]   [3.0]
[ 0.51 0.13 0.19 0.85 ] [x3]   [4.0]
@end example

@end ifinfo
@noindent
and the solution is found using LU decomposition of the matrix @math{A}.

@example
@verbatiminclude examples/linalglu.c
@end example

@noindent
Here is the output from the program,

@example
@verbatiminclude examples/linalglu.out
@end example

@noindent
This can be verified by multiplying the solution @math{x} by the
original matrix @math{A} using @sc{gnu octave},

@example
octave> A = [ 0.18, 0.60, 0.57, 0.96;
              0.41, 0.24, 0.99, 0.58; 
              0.14, 0.30, 0.97, 0.66; 
              0.51, 0.13, 0.19, 0.85 ];

octave> x = [ -4.05205; -12.6056; 1.66091; 8.69377];

octave> A * x
ans =
  1.0000
  2.0000
  3.0000
  4.0000
@end example

@noindent
This reproduces the original right-hand side vector, @math{b}, in
accordance with the equation @math{A x = b}.

@node Linear Algebra References and Further Reading
@section References and Further Reading

Further information on the algorithms described in this section can be
found in the following book,

@itemize @asis
@item
G. H. Golub, C. F. Van Loan, @cite{Matrix Computations} (3rd Ed, 1996),
Johns Hopkins University Press, ISBN 0-8018-5414-8.
@end itemize

@noindent
The @sc{lapack} library is described in the following manual,

@itemize @asis
@item
@cite{LAPACK Users' Guide} (Third Edition, 1999), Published by SIAM,
ISBN 0-89871-447-8.

@uref{http://www.netlib.org/lapack} 
@end itemize

@noindent
The @sc{lapack} source code can be found at the website above, along
with an online copy of the users guide.

@noindent
The Modified Golub-Reinsch algorithm is described in the following paper,

@itemize @asis
@item
T.F. Chan, ``An Improved Algorithm for Computing the Singular Value
Decomposition'', @cite{ACM Transactions on Mathematical Software}, 8
(1982), pp 72--83.
@end itemize

@noindent
The Jacobi algorithm for singular value decomposition is described in
the following papers,

@itemize @asis
@item
J.C. Nash, ``A one-sided transformation method for the singular value
decomposition and algebraic eigenproblem'', @cite{Computer Journal},
Volume 18, Number 1 (1975), p 74--76

@item
J.C. Nash and S. Shlien ``Simple algorithms for the partial singular
value decomposition'', @cite{Computer Journal}, Volume 30 (1987), p
268--275.

@item
James Demmel, Kresimir Veselic, ``Jacobi's Method is more accurate than
QR'', @cite{Lapack Working Note 15} (LAWN-15), October 1989. Available
from netlib, @uref{http://www.netlib.org/lapack/} in the @code{lawns} or
@code{lawnspdf} directories.
@end itemize