Re: [Boost-Users] uBLAS - Matrix Inversion

From: Kresimir Fresl (fresl_at_[hidden])
Date: 2003-02-25 11:41:59

Thomas Reisinger wrote:

> I'm new in the field of numeric computing and I would like to learn more 
> about it. To solve my problem I need to invert a matix and use the 
> intermediate result as the input for following calculations.
> 
> I browsed the uBLAS library and played a little bit with it but I didn't 
> find a routine for matrix inversion. 

As Jon Agiato already replied, uBLAS doesn't have matrix
inversion.

uBLAS is, more or less, BLAS, that is, it provides elementary
vector and matrix operation (scaling, addition, multiplication etc.).
One day we will hopefully have uLAPACK -- solution of linear
systems, matrix inversion, eigenvalues & eigenvectors...

(Roland aka speedsnaii mentioned uBLAS's solve() functions
in his reply.  But those functions are `triangular' solvers, i.e.
matrix A must be triangular -- corresponding BLAS trsm() functions
are used in LAPACK for `backward' substitutions after LU or
Cholesky factorizations.)

> Is there any possibility to preform 
> it in a single step? Is this package the wrong tool for doing this and I 
> should use LAPACK++ or similar stuff?

 From  LAPACK++ home page:
``[NOTE: This package is being superseded by the Template Numerical
Toolkit (TNT), which utilizes new features of the ANSI C++ specification.
TNT is a newer design, and will integrate the functionlaity of Lapack++,
IML++, SparseLib++, and MV++.]''
(TNT home page is: http://math.nist.gov/tnt/)

But I'd recommend the `bindings' library (as one of the co-authors
I am, of course, biased ;o) which *can be used with uBLAS
vectors and matrices*.  Idea is similar to the LAPACK++, that is,
library provides interfaces to other, well-known numerical
libraries (written in Fortran or C).

For example, for matrix inversion you can use ATLAS
bindings:
=========================================
#include <boost/numeric/bindings/atlas/clapack.hpp>
#include <boost/numeric/bindings/traits/ublas_matrix.hpp>
#include <boost/numeric/bindings/traits/std_vector.hpp>

namespace ublas = boost::numeric::ublas;
namespace atlas = boost::numeric::bindings::atlas;

int main() {
     std::size_t n = 5;
     ublas::matrix<double> A (n, n);
     // fill matrix A

     std::vector<int> ipiv (n);   // pivot vector
     atlas::lu_factor (a, ipiv);  // alias for getrf()
     atlas::lu_invert (a, ipiv);  // alias for getri()
}
=========================================
(Of course, you must install ATLAS first:

        http://math-atlas.sourceforge.net/

``The ATLAS (Automatically Tuned Linear Algebra Software)
project is an ongoing research effort focusing on applying
empirical techniques in order to provide portable performance.
At present, it provides C and Fortran77 interfaces to a portably
efficient BLAS implementation, as well as a few routines from
LAPACK.'')

`bindings' library is still in its infancy -- it hasn't its home
page. Older version can be found in boost sandbox, and
you can download the latest version from the files section
of ublas-dev mailing list:

       http://groups.yahoo.com/group/ublas-dev/files/

File is `bindings_2003_01_30.zip'.

There you can also find newest, sometimes experimental,
versions of uBLAS. Also, if you are interested, you can join
the ublas-dev mailing list:

     ublas-dev-subscribe_at_[hidden]

Its message archive contains some discussions about
`bindings' , which can serve as `poor man's documentation'. ;o)
Namely, there's almost no documentation about the library.
For ATLAS bindings there are some examples (with comments)
in the sub-directory

      /libs/numeric/bindings/atlas

Examples for matrix inversion are in `ublas_getri.cc' and,
if your matrix is SPD, in `ublas_potri.cc'.

If you have any problems and/or questions, don't hesitate
to ask.

Regards,

fres