• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 //          Copyright Gunter Winkler 2004 - 2009.
2 // Distributed under the Boost Software License, Version 1.0.
3 //    (See accompanying file LICENSE_1_0.txt or copy at
4 //          http://www.boost.org/LICENSE_1_0.txt)
5 
6 
7 #include <iostream>
8 
9 #include <boost/numeric/ublas/matrix.hpp>
10 #include <boost/numeric/ublas/triangular.hpp>
11 
12 #include <boost/numeric/ublas/io.hpp>
13 
14 using std::cout;
15 using std::endl;
16 
17 
18 
19 namespace ublas = boost::numeric::ublas;
20 
21 
main(int argc,char * argv[])22 int main(int argc, char * argv[] ) {
23 
24   ublas::matrix<double> M (3, 3);
25   for (std::size_t i=0; i < M.data().size(); ++i) { M.data()[i] = 1+i ; }
26 
27   std::cout << "full         M  = " << M << "\n" ;
28 
29   ublas::triangular_matrix<double, ublas::lower> L;
30   ublas::triangular_matrix<double, ublas::unit_lower> UL;
31   ublas::triangular_matrix<double, ublas::strict_lower> SL;
32 
33   L  = ublas::triangular_adaptor<ublas::matrix<double>, ublas::lower> (M);
34   SL = ublas::triangular_adaptor<ublas::matrix<double>, ublas::strict_lower> (M);
35   UL = ublas::triangular_adaptor<ublas::matrix<double>, ublas::unit_lower> (M);
36 
37   std::cout << "lower        L  = " << L << "\n"
38             << "strict lower SL = " << SL << "\n"
39             << "unit lower   UL = " << UL << "\n" ;
40 
41   ublas::triangular_matrix<double, ublas::upper> U;
42   ublas::triangular_matrix<double, ublas::unit_upper> UU;
43   ublas::triangular_matrix<double, ublas::strict_upper> SU;
44 
45   U =  ublas::triangular_adaptor<ublas::matrix<double>, ublas::upper> (M);
46   SU = ublas::triangular_adaptor<ublas::matrix<double>, ublas::strict_upper> (M);
47   UU = ublas::triangular_adaptor<ublas::matrix<double>, ublas::unit_upper> (M);
48 
49   std::cout << "upper        U  = " << U << "\n"
50             << "strict upper SU = " << SU << "\n"
51             << "unit upper   UU = " << UU << "\n" ;
52 
53   std::cout << "M = L + SU ? " << ((norm_inf( M - (L + SU) ) == 0.0)?"ok":"failed") << "\n";
54   std::cout << "M = U + SL ? " << ((norm_inf( M - (U + SL) ) == 0.0)?"ok":"failed") << "\n";
55 
56 }
57 
58 
59