Boost logo

Ublas :

From: Karl Meerbergen (karl.meerbergen_at_[hidden])
Date: 2007-09-07 02:33:27


Dear Feng Gao,

It surprises me that many people are trying to use ublas routines for
solving linear systems, while excellent external software is available:
LAPACK for dense matrices, UMFPACK, MUMPS for sparse.

In your example, you are trying to solve a 5x5 system of a sparse
matrix. I guess this is just a small test example to make sure the code
works. Suppose you want to solve a 20,000x20,0000 sparse system. Which
performance do you expect to obtain with lu_factorize? For a sparse
matrix, this is, in general, going to be a disaster.

I suggest to spend some time to develop a link with external software
rather than (sorry for being so rude) wasting your time using
lu_factorize. There are C++ bindings to UMFPACK in boost.sandbox. I am
currently working on MUMPS bindings, which I hope to make available
soon. You will have to spend some time installing packages and compiling
and linking, maybe a day or half a day. It may be a frustrating and hard
job, but once you have done this, your code will have a good performance.

Have a nice day,

Karl

pinacle wrote:
> Dear all,
>
> I've been using UBLAS for solving problems like Ax=b, where A is a complex
> dense matrix, and x and b are complex dense vectors. Today, I tried to
> change A to a sparse matrix, but I got compiling errors as:
> Error 10 error C2784: 'matrix_binary_scalar2_traits<E1,const
> T2,boost::numeric::ublas::scalar_divides<E1::value_type,T2>>::result_type
> boost::numeric::ublas::operator /(const
> boost::numeric::ublas::matrix_expression<E> &,const T2 &)' : could not
> deduce template argument for 'const
> boost::numeric::ublas::matrix_expression<E> &' from 'std::complex<double>'
> D:\boost\boost\numeric\ublas\lu.hpp 147
>
> But if I use real sparse A, no error appears. So the conclusion now is that
> lu_factorize works for real matrices no matter dense or sparse, and works
> for complex dense matrix, but doesn't work for complex sparse matrix.
>
> I do need a lu_factorize (and lu_substitute of course) working for complex
> sparse matrices. How could I solve the problem? BTW, I am using boost
> 1.34.0.
>
> Here is a sample code for testing:
>
> #include <complex>
> #include <iostream>
>
> #include <boost/numeric/ublas/storage_sparse.hpp>
> #include <boost/numeric/ublas/matrix.hpp>
> #include <boost/numeric/ublas/vector.hpp>
> #include <boost/numeric/ublas/matrix_sparse.hpp>
> #include <boost/numeric/ublas/io.hpp>
> #include <boost/numeric/ublas/triangular.hpp>
> #include <boost/numeric/ublas/lu.hpp>
>
> using namespace boost::numeric::ublas;
> using namespace std;
>
> int main(int argc, char *argv[])
> {
> const size_t size = 5;
>
> // error:
> //compressed_matrix<complex<double> > A(size, size);
> //compressed_matrix<complex<double> > ALU(size, size);
> // good:
> compressed_matrix<double> A(size, size);
> compressed_matrix<double> ALU(size, size);
>
> permutation_matrix<> pm(size);
>
> for (size_t i=0; i<A.size1(); ++i) {
> if (i>=1) {
> A.push_back(i,i-1,-2.0);
> }
> A.push_back(i,i,4);
> if (i+1<A.size2()) {
> A.push_back(i,i+1,-1.0);
> }
> }
> ALU = A;
> lu_factorize(ALU,pm); // LU factorized
>
> cout << "A = " << A << endl;
> cout << "ALU = " << ALU << endl;
>
> return EXIT_SUCCESS;
> };
>
> Thanks,
> Feng Gao
>
>
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://lists.boost.org/mailman/listinfo.cgi/ublas
>