ViennaCL - The Vienna Computing Library  1.7.0
Free open-source GPU-accelerated linear algebra and solver library.
row_scaling.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_LINALG_ROW_SCALING_HPP_
2 #define VIENNACL_LINALG_ROW_SCALING_HPP_
3 
4 /* =========================================================================
5  Copyright (c) 2010-2015, Institute for Microelectronics,
6  Institute for Analysis and Scientific Computing,
7  TU Wien.
8  Portions of this software are copyright by UChicago Argonne, LLC.
9 
10  -----------------
11  ViennaCL - The Vienna Computing Library
12  -----------------
13 
14  Project Head: Karl Rupp rupp@iue.tuwien.ac.at
15 
16  (A list of authors and contributors can be found in the manual)
17 
18  License: MIT (X11), see file LICENSE in the base directory
19 ============================================================================= */
20 
25 #include <vector>
26 #include <cmath>
27 #include "viennacl/forwards.h"
28 #include "viennacl/vector.hpp"
30 #include "viennacl/tools/tools.hpp"
31 
32 #include <map>
33 
34 namespace viennacl
35 {
36  namespace linalg
37  {
38 
41  {
42  public:
47  row_scaling_tag(unsigned int p = 2) : norm_(p) {}
48 
50  unsigned int norm() const { return norm_; }
51 
52  private:
53  unsigned int norm_;
54  };
55 
56 
58  namespace detail
59  {
60  template<typename T>
61  struct row_scaling_for_viennacl
62  {
63  enum { value = false };
64  };
65 
66  template<typename ScalarType, unsigned int ALIGNMENT>
67  struct row_scaling_for_viennacl< viennacl::compressed_matrix<ScalarType, ALIGNMENT> >
68  {
69  enum { value = true };
70  };
71 
72  template<typename ScalarType, unsigned int ALIGNMENT>
73  struct row_scaling_for_viennacl< viennacl::coordinate_matrix<ScalarType, ALIGNMENT> >
74  {
75  enum { value = true };
76  };
77  }
85  template<typename MatrixType,
86  bool is_viennacl = detail::row_scaling_for_viennacl<MatrixType>::value >
88  {
89  typedef typename MatrixType::value_type ScalarType;
90 
91  public:
97  row_scaling(MatrixType const & mat, row_scaling_tag const & tag) : diag_M(viennacl::traits::size1(mat))
98  {
99  assert(mat.size1() == mat.size2() && bool("Size mismatch"));
100  init(mat, tag);
101  }
102 
103  void init(MatrixType const & mat, row_scaling_tag const & tag)
104  {
105  diag_M.resize(mat.size1()); //resize without preserving values
106 
107  for (typename MatrixType::const_iterator1 row_it = mat.begin1();
108  row_it != mat.end1();
109  ++row_it)
110  {
111  for (typename MatrixType::const_iterator2 col_it = row_it.begin();
112  col_it != row_it.end();
113  ++col_it)
114  {
115  if (tag.norm() == 0)
116  diag_M[col_it.index1()] = std::max<ScalarType>(diag_M[col_it.index1()], std::fabs(*col_it));
117  else if (tag.norm() == 1)
118  diag_M[col_it.index1()] += std::fabs(*col_it);
119  else if (tag.norm() == 2)
120  diag_M[col_it.index1()] += (*col_it) * (*col_it);
121  }
122  if (!diag_M[row_it.index1()])
123  throw zero_on_diagonal_exception("ViennaCL: Zero row encountered while setting up row scaling preconditioner!");
124 
125  if (tag.norm() == 2)
126  diag_M[row_it.index1()] = std::sqrt(diag_M[row_it.index1()]);
127  }
128  }
129 
130 
132  template<typename VectorType>
133  void apply(VectorType & vec) const
134  {
135  assert(vec.size() == diag_M.size() && bool("Size mismatch"));
136  for (vcl_size_t i=0; i<vec.size(); ++i)
137  vec[i] /= diag_M[i];
138  }
139 
140  private:
141  std::vector<ScalarType> diag_M;
142  };
143 
144 
149  template<typename MatrixType>
150  class row_scaling< MatrixType, true>
151  {
153 
154 
155  public:
161  row_scaling(MatrixType const & mat, row_scaling_tag const & tag) : diag_M(mat.size1(), viennacl::traits::context(mat))
162  {
163  init(mat, tag);
164  }
165 
166  void init(MatrixType const & mat, row_scaling_tag const & tag)
167  {
168  switch (tag.norm())
169  {
170  case 0:
172  break;
173  case 1:
175  break;
176  case 2:
178  break;
179  default:
180  throw unknown_norm_exception("Unknown norm when initializing row_scaling preconditioner!");
181  }
182  }
183 
184  template<unsigned int ALIGNMENT>
186  {
187  assert(viennacl::traits::size(diag_M) == viennacl::traits::size(vec) && bool("Size mismatch"));
188  vec = element_div(vec, diag_M);
189  }
190 
191  private:
193  };
194 
195  }
196 }
197 
198 
199 
200 
201 #endif
202 
203 
204 
viennacl::vector_expression< const vector_base< T >, const vector_base< T >, op_element_binary< op_div > > element_div(vector_base< T > const &v1, vector_base< T > const &v2)
row_scaling(MatrixType const &mat, row_scaling_tag const &tag)
Constructor for the preconditioner.
Definition: row_scaling.hpp:97
Various little tools used here and there in ViennaCL.
vcl_size_t size1(MatrixType const &mat)
Generic routine for obtaining the number of rows of a matrix (ViennaCL, uBLAS, etc.)
Definition: size.hpp:163
This file provides the forward declarations for the main types used within ViennaCL.
unsigned int norm() const
Returns the index p of the l^p-norm (0 ... ||x||_sup, 1... sum(abs(x)), 2... sqrt(sum(x_i^2))). Currently only p=0, p=1, and p=2 supported.
Definition: row_scaling.hpp:50
Jacobi-type preconditioner class, can be supplied to solve()-routines. This is a diagonal preconditio...
Definition: row_scaling.hpp:87
Represents a generic 'context' similar to an OpenCL context, but is backend-agnostic and thus also su...
Definition: context.hpp:39
Main namespace in ViennaCL. Holds all the basic types such as vector, matrix, etc. and defines operations upon them.
Definition: cpu_ram.hpp:34
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.)
Definition: size.hpp:235
Definition: blas3.hpp:36
void init(MatrixType const &mat, row_scaling_tag const &tag)
Implementation of the compressed_matrix class.
A tag for a row scaling preconditioner which merely normalizes the equation system such that each row...
Definition: row_scaling.hpp:40
std::size_t vcl_size_t
Definition: forwards.h:75
T::ERROR_CANNOT_DEDUCE_CPU_SCALAR_TYPE_FOR_T type
Definition: result_of.hpp:271
row_scaling(MatrixType const &mat, row_scaling_tag const &tag)
Constructor for the preconditioner.
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
void apply(VectorType &vec) const
Apply to res = b - Ax, i.e. row applied vec (right hand side),.
viennacl::enable_if< viennacl::is_any_sparse_matrix< SparseMatrixType >::value >::type row_info(SparseMatrixType const &mat, vector< SCALARTYPE, VEC_ALIGNMENT > &vec, row_info_types info_selector)
row_scaling_tag(unsigned int p=2)
Constructor.
Definition: row_scaling.hpp:47
void apply(viennacl::vector< ScalarType, ALIGNMENT > &vec) const
void init(MatrixType const &mat, row_scaling_tag const &tag)