uBLAS is, more or less, BLAS (see http://www.netlib.org/blas/). It provides the "building blocks" for storing vector and matricies and vector and matrix operations (scaling, addition, multiplication etc). uBLAS provides the C++ infrastructure on which safe and efficient linear algebra algorithms can be built.

uBLAS also has a set of functions for LU factorisation. These implement common LU factorisation algorithms directly in uBLAS. They are defined in the header file <boost/numeric/ublas/lu.hpp>. LU Matrix Inversion demonstrates how these functions can be used to invert matrices.

For example, for matrix inversion using an LU factorisation 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."

The `bindings' library is a very useful interface between Boost and uBLAS and traditional linear algebra libraries. It is being actively maintianed with the hope of eventual inclusion in Boost.

Development of the binding interface is located at Boost Sandbox CVS - http://www.boost.org/more/mailing_lists.htm#sandbox

Here you can also find newest, sometimes experimental,versions of the binding. You can access the header files directly at:

http://boost-sandbox.cvs.sourceforge.net/boost-sandbox/boost-sandbox/boost/numeric/bindings/

The starting point for documentation and the many examples is:

boost-sandbox/boost-sandbox/libs/numeric/bindings

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

You may also download ready-made snapshot releases of the bindings from http://news.tiker.net/software/boost-bindings.

From the 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/)

On the other hand, from the homepage of LAPACK++ v2.1.2: http://lapackpp.sourceforge.net/

However, they abandoned LAPACK in the year 2000 and stated: "Lapack++ is no longer actively supported. The successor to this project is that Template Numerical Toolkit (TNT), see http://math.nist.gov/tnt for details." Unfortunately, the project TNT never really took off. Therefore this fork from the original LAPACK++ has been started. There are a whole number of changes now in here. Most notably, this local copy has complex matrices enabled again by adding a custom copy of stdc++'s complex type (see include/lacomplex.h and include/lacomplex).

So my question, again: What about LAPACK++ ?

Blaisorblade answers:

(I'm not a Boost developer) LAPACK++ seems fairly complete, but it seems not written by C++ experts. Beyond that, it contains its own reimplementation of basic vector operation (indeed, its own reimplementation of BLAS for C++), but it is a really inferior one (it does not definitely use expression templates).

Say, in Lapack++ 2.5.2, LaGenMatDouble? class (a double dense matrix) has two identical methods, having these signatures:

LaGenMatDouble? operator()(const LaIndex?& I, const LaIndex?& J) ; LaGenMatDouble? operator()(const LaIndex?& I, const LaIndex?& J) const;the non-const overload returns a copy of the selected portion of the matrix, not a view. I wonder if they understand that the second overload is not needed.

In the overload of operator() to access a single element (operator()(int row, int column)), a reference to the element is returned; I wonder whether they'd need to return a proxy reference, since they can share the backing array.

However, I downloaded LAPACK++ yesterday, and gave just a look. I'll compare it with the bindings project, to compare their respective merits.

Back to Effective uBLAS

Disclaimer: This site not officially maintained by Boost Developers