From: Piotr Stanczyk (piotr-s_at_[hidden])
Date: 2006-05-19 04:54:09
So, the solve function is only for LU systems? I was hoping that this
would be done in the solve function and an error would be generated on
encountering singular matrices.
So I have another question, where in the code should I look for taking a
SVD approach?
Thanks,
Piotr
Piotr Stanczyk wrote:
> Hi all,
>
> this a newbie question so please bear with me.
>
> I am simply trying to solve a bog standard linear system: Ax=b
>
> I have found this in the doc pages:
>
> Examples
>
> #include <boost/numeric/ublas/triangular.hpp>
> #include <boost/numeric/ublas/io.hpp>
>
> int main () {
> using namespace boost::numeric::ublas;
> matrix<double> m (3, 3);
> vector<double> v (3);
> for (unsigned i = 0; i < std::min (m.size1 (), v.size ()); ++ i) {
> for (unsigned j = 0; j <= i; ++ j)
> m (i, j) = 3 * i + j + 1;
> v (i) = i;
> }
>
> std::cout << solve (m, v, lower_tag ()) << std::endl;
> std::cout << solve (v, m, lower_tag ()) << std::endl;
> }
>
> but I have no idea what the args to solve are, looking through the hpp
> is a little daunting to say the least, and I am not getting the correct
> values back on most of my test cases ....
>
> Here is a trivial example, 1st solve is fine
>
> #include <boost/numeric/ublas/triangular.hpp>
> #include <boost/numeric/ublas/io.hpp>
>
> int main () {
> using namespace boost::numeric::ublas;
>
> matrix<double> m (3, 3);
> vector<double> v (3);
> vector<double> v_1 (3);
> vector<double> v_2 (3);
> vector<double> v_3 (3);
>
>
> m(0,0) = 0.158564; m(0,1) = 0.0780387; m(0,2) = 0.0534636;
> m(1,0) = 0.498959; m(1,1) = 0.28208; m(1,2) = 0.206461;
> m(2,0) = 0.112965; m(2,1) = 0.18; m(2,2) = 0.311839;
> std::cout << "m =\n" << m << std::endl;
>
> v(0) = m(0,0); v(1) = m(1,0); v(2) = m(2,0);
> std::cout << "v =\n" << v << std::endl;
>
> v_1=solve (m, v, lower_tag ());
> std::cout << "v1:\n" << v_1 << std::endl;
> std::cout << "Check Solve:" << std::endl;
>
> v_3 = prod (m, v_1);
> std::cout <<v - v_3 << std::endl;
>
>
> v(0) = m(0,1); v(1) = m(1,1); v(2) = m(2,1);
> std::cout << "v =\n" << v << std::endl;
>
> v_1=solve (m, v, lower_tag ());
> std::cout << "v1:\n" << v_1 << std::endl;
> std::cout << "Check Solve:" << std::endl;
>
> v_3 = prod (m, v_1);
> std::cout <<v - v_3 << std::endl;
>
>
> }
>
> Can anyone help me?
>
> Thanks
>
>
> Piotr
>
>
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://listarchives.boost.org/mailman/listinfo.cgi/ublas
>
-- Dr Piotr Stanczyk 2D R&D - Digital Film The Moving Picture Company London, UK