Etienne Grossmann > Meschach > PDL::Meschach

Download:
Meschach-0.03.tar.gz

Dependencies

Annotate this POD

View/Report Bugs
Module Version: 0.03   Source  

NAME ^

PDL::Meschach - Link PDL to meschach 1.2 matrix library

Version 0.03 alpha

PDL::Meschach Etienne Grossmann, 11/11/1996 [etienne@isr.ist.utl.pt] meschach 1.2 David E. Stewart [david.stewart@anu.edu.au] Zbigniew Leyk [zbigniew.leyk@anu.edu.au]

DESCRIPTION ^

 PDL::Meschach links PDL to a few matrix functions from meschach 1.2 : 

    Diagonal,upper,lower triangle extraction ...

    Matrix exponentiation ...

    LU , Cholesky , QR Factorisation and associated 

    Linear equation solvers.

    Symmetric matrix eigenvector/eigenvalue extraction.

    Singular value decomposition.

INTRODUCTION ^

 The functionalities are available either through "friendly"
 functions, e.g. :

   # Solve A.x == b using LU decomposition

   $x = lusolve($A,$b);

 or through "raw" functions :

   $LU = $A + 0 ;             # Copy $A
   lufac_($LU,$Perm)          # LU decomposition overwrites $LU
   lusolve_($x,$b,$LU,$Perm) 

 That may be more efficient in terms of memory allocation.

 All "raw" function names have a trailing underscore.

CAVEAT ^

 * Dimensions are expressed 

        ALWAYS as COLUMN, ROW 

 instead of the usual matrix row, column.



 * pdl members "incs" and "offs" are not yet handled correctly, and
   neither can they be detected reliably. The pdls passed to Meschach
   routines are assumed to hold data consecutively arranged in memory.

USAGE ^

To load Meschach in a script : either :

  use PDL::Meschach;           # Only "Friendly" functions.

 or   

  use PDL::Meschach qw( :Raw ) # "Raw" functions too.

To use Meschach in perldl, : insert

        eval "use Meschach qw( :All )";
        if($@ ne ""){
          print "Meschach NOT AVAILABLE : \n$@\n";
        } else { print "Meschach found\n";}


  somewhere in it, e.g. before the line :

        eval "use Term::ReadLine"; 

AVAILABLE FUNCTIONS ^

ut, lt
    $T = ut( $A );

 Puts the upper triangle of $A in $T.

    $T = ut( $Col, $Row );

 Sets to 1 the upper-triangle of T, to 0 its strictly lower triangle.

 Raw function :
    
    ut_($T,$A);
    ut_($T,$Col,$Row);


  The output is the same type as th input.

  For lower triangle, use  lt_, lt. 
diag
    $Vec = diag( $Mat );

    diag_( $Vec, $Mat);

 Both put into $Vec the diagonal of $Mat. output is PDL_D (!).


    $Mat = diag ( $Vec [,$cols, $rows] );

    diag_($Vec,$Mat);

 Make $Mat a diagonal matrix with $Vec as diagonal. $Mat may be
 of arbitrary size if $cols, $rows are used :

    $Mat = diag ( $Vec,$c );     # $c x $c 
    $Mat = diag ( $Vec,$c, $r ); # $c x $m

 output is PDL_F (!)
ident
  $Mat = ident($Cols [,$Rows = $Cols ] )

  Returns a PDL_D Identity Matrix.
mpow
    $Out = mpow( $In, $Pow );

    mpow_( $Out, $In, $Pow [,$Coerce] );

 Integer (negative or positive) powers of a matrix,  If $Coerce is
 true (default), the result is Real typed. Otherwise, $$Out{Datatype}
 is unchanged. 
inv
                $Out = inv( $In );

    inv_( $Out, $In [,$Coerce] );

 Inverse of a matrix.   
lusolve
 Solve $A x $x == $b , by LU Factorization, with pivoting.

    ($LU,$Perm,$x) = lusolve( $b, $A );

 $LU,$Perm describe the factorization. They may be re-used (which
 spares some computation). $LU is a matrix the same size as $A. $Perm
 is an integer vector describing the pivoting used in the
 factorization.

    $x = lusolve($b, $LU, $Perm );

 The third argument is what decides lusolve not to factor. In all
 cases ($LU,$Perm,$x) is returned. 

 An estimate of the conditioning (the ratio of the greatest and
 smallest eigenvalues) is returned by :

                lucond($LU,$Perm);


 Raw method :

    $LU = $A + 0 ;                   # Copy $A. 
    lufac_( $LU, $Perm );            # Factorize. 
    lusolve_( $x, $b, $LU, $Perm );  # Solve.
qrsolve
 Solve $A x $x == $b , by QR Factorization 

  (A = Q.R where Q is orthogonal, and R upper triangular).
  The "QR" functions are almost equivalent to the "LU" ones.

    ($QR,$V,$x) = qrsolve( $b, $A );

    $x = qrsolve($b, $QR, $V );

 $QR,$V describe the factorization. The "R" matrix is represented by
 the strictly higher triangle in $QR, and its diagonal is $V. $Q is
 represented by the lower triangle.

 The third argument is what decides qrsolve not to factor. In all
 cases ($QR,$V,$x) is returned. 

 Raw method :

    $QR = $A + 0 ;                   # Copy $A. 
    qrfac_( $QR, $V );               # Factorize. 
    qrsolve_( $x, $b, $QR, $V );     # Solve.


 An estimate of the conditioning (the ratio of the greatest and
 smallest eigenvalues) is returned by :

                qrcond($QR);                                                                             # only QR is passed.
chsolve
 If A is positive definite, solving  $A x $x == $b 
 by Cholesky Factorization may be more efficient :

    ($CH,$x) = chsolve( $b, $A );

 $CH may be re-used :

    $x = chsolve ( $b , $CH , 1 );

 If the third argument is true, chsolve considers that it is already
 in factored form. 
    

 Raw method :
    
    $CH = $A + 0 ;            # Copy
    chfac_($CH);              # Factor
    chsolve_($x,$b,$CH);      # Solve
symmeig
 Symmetric Matrix Eigenvalues/vectors 

    ( $Mat, $Vec ) = symmeig( $A );

 $Mat contains the eigenvectors of $A,
 $Vec contains the eigenvalues  of $A,

 Raw method :

    symmeig_( $Mat, $Vec, $A );

 symmeig_ returns true upon success.   
svd
 Singular Value Decomposition of a  ncol,nrow matrix :

    ($U,$V,$l) = svd( $A ) ;

 $U Contains the left "singular vectors" (nrow,nrow matrix).
 $V Contains the right "singular vectors" (ncol,ncol matrix).
 $l Contains the "singular values"     (min(ncol,nrow) vector).

 Raw method :
   
    svd_( $U, $V, $l, $A );

 or

    svd_( $l, $A );        # Only the "singular values".  

 svd_ returns true upon success.   

TYPE CONVERSION ^

The equivalence between types is done in the BOOT: part of Meschach.xs

The results of most functions are "Real"-typed pdls.

Sometimes (e.g. the $Perm argument in lufac_ ) it is an integer-typed pdls.

"Raw" functions may conserve the type of their return-argument when they accept a $coerce argument.

syntax highlighting: