31 #ifndef SHARK_LINALG_BLAS_KERNELS_DEFAULT_tpmv_HPP 32 #define SHARK_LINALG_BLAS_KERNELS_DEFAULT_tpmv_HPP 34 #include "../../matrix_proxy.hpp" 35 #include "../../vector_expression.hpp" 36 #include <boost/mpl/bool.hpp> 38 namespace shark{
namespace blas{
namespace bindings{
45 template<
class TriangularA,
class V>
47 matrix_expression<TriangularA>
const& A,
48 vector_expression<V>& b,
52 typedef typename V::value_type value_type;
53 typedef typename V::size_type size_type;
54 typedef typename TriangularA::const_row_iterator row_iterator;
55 size_type
size = A().size1();
57 for(size_type i = 0; i !=
size; ++i){
59 row_iterator end = A().row_end(i);
60 for(row_iterator pos = A().row_begin(i); pos != end; ++pos){
61 sum += *pos*b()(pos.index());
72 template<
class TriangularA,
class V>
74 matrix_expression<TriangularA>
const& A,
75 vector_expression<V>& b,
79 typedef typename V::value_type value_type;
80 typedef typename V::size_type size_type;
81 typedef typename TriangularA::const_row_iterator row_iterator;
82 size_type size = A().size1();
84 for(size_type irev = size; irev != 0; --irev){
87 row_iterator end = A().row_end(i);
88 for(row_iterator pos = A().row_begin(i); pos != end; ++pos){
89 sum += *pos*b()(pos.index());
98 template<
class TriangularA,
class V>
100 matrix_expression<TriangularA>
const& A,
101 vector_expression<V>& b,
105 typedef typename TriangularA::const_column_iterator column_iterator;
106 typedef typename V::size_type size_type;
107 typedef typename V::value_type value_type;
108 size_type size = A().size1();
109 for(size_type i = 0; i !=
size; ++i){
110 value_type bi = b()(i);
111 b()(i)= value_type();
112 column_iterator end = A().column_end(i);
113 for(column_iterator pos = A().column_begin(i); pos != end; ++pos){
114 b()(pos.index()) += *pos*bi;
123 template<
class TriangularA,
class V>
125 matrix_expression<TriangularA>
const& A,
126 vector_expression<V>& b,
130 typedef typename TriangularA::const_column_iterator column_iterator;
131 typedef typename V::size_type size_type;
132 typedef typename V::value_type value_type;
133 size_type size = A().size1();
135 for(size_type irev = size; irev != 0; --irev){
137 value_type bi = b()(i);
138 b()(i)= value_type();
139 column_iterator end = A().column_end(i);
140 for(column_iterator pos = A().column_begin(i); pos != end; ++pos){
141 b()(pos.index()) += *pos*bi;
147 template <
typename TriangularA,
typename V>
149 matrix_expression<TriangularA>
const& A,
150 vector_expression<V>& b,
157 typename TriangularA::orientation::orientation(),
158 typename TriangularA::orientation::triangular_type()