|
|
Ublas : |
From: Max Weinberg (max.weinberg_at_[hidden])
Date: 2006-08-06 15:44:36
Hi,
I tried to solve a system A.x = b where A is full rank but over-determined.
This triggered a UBLAS_CHECK assertion in lu_factorize. I commented that
assertion out and the result of lu_factorize looks correct. But
lu_substitute doesn't like the input. Looks like it can't handle a
non-square matrix. The code below demonstrates the problem.
Best, Max
#include <boost/numeric/ublas/vector.hpp>
#include <boost/numeric/ublas/vector_proxy.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/triangular.hpp>
#include <boost/numeric/ublas/lu.hpp>
#include <boost/numeric/ublas/io.hpp>
namespace ublas = boost::numeric::ublas;
int main( int /*argc*/, char* /*argv*/[] ) {
typedef ublas::permutation_matrix<std::size_t> pmatrix;
ublas::matrix<double> mat_A(3,2);
ublas::vector<double> vec_b(3);
//right hand side
vec_b(0) = 6; vec_b(1) = 12; vec_b(2) = 12;
//matrix, rank 2, one equation duplicated
mat_A(0,0) = 1; mat_A(0,1) = 2;
mat_A(1,0) = 1; mat_A(1,1) = 8;
mat_A(2,0) = 1; mat_A(2,1) = 8;
//create a permutation matrix for the LU-factorization
pmatrix pm(mat_A.size1());
lu_factorize(mat_A,pm);
ublas::vector<double> solution(vec_b);
lu_substitute(mat_A, pm, solution);
//correct result would be solution = [4; 1]
}
___________________________________________________________________________
mymail - der unschlagbare und kostenlose E-Mail-Dienst der Schweiz!
http://mymail.ch/?redirect=9999
Unser Produkte-Berater bringt Sie rasch zum Ziel. Nur bei www.arp.com!
http://clk.tradedoubler.com/click?p=18634&a=1036965&g=346289