All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
circulant_matrix.hpp
Go to the documentation of this file.
#ifndef VIENNACL_CIRCULANT_MATRIX_HPP
#define VIENNACL_CIRCULANT_MATRIX_HPP
/* =========================================================================
Copyright (c) 2010-2016, Institute for Microelectronics,
Institute for Analysis and Scientific Computing,
TU Wien.
Portions of this software are copyright by UChicago Argonne, LLC.
-----------------
ViennaCL - The Vienna Computing Library
-----------------
Project Head: Karl Rupp rupp@iue.tuwien.ac.at
(A list of authors and contributors can be found in the manual)
License: MIT (X11), see file LICENSE in the base directory
============================================================================= */
#include "viennacl/forwards.h"
#include "viennacl/vector.hpp"
#include "viennacl/ocl/backend.hpp"
#include "viennacl/linalg/circulant_matrix_operations.hpp"
#include "viennacl/fft.hpp"
namespace viennacl
{
template<class NumericT, unsigned int AlignmentV>
class circulant_matrix
{
public:
typedef viennacl::backend::mem_handle handle_type;
typedef scalar<typename viennacl::tools::CHECK_SCALAR_TEMPLATE_ARGUMENT<NumericT>::ResultType> value_type;
explicit circulant_matrix() {}
explicit circulant_matrix(vcl_size_t rows, vcl_size_t cols) : elements_(rows)
{
assert(rows == cols && bool("Circulant matrix must be square!"));
(void)cols; // avoid 'unused parameter' warning in optimized builds
}
void resize(vcl_size_t sz, bool preserve = true)
{
elements_.resize(sz, preserve);
}
handle_type const & handle() const { return elements_.handle(); }
viennacl::vector<NumericT, AlignmentV> & elements() { return elements_; }
viennacl::vector<NumericT, AlignmentV> const & elements() const { return elements_; }
vcl_size_t size1() const { return elements_.size(); }
vcl_size_t size2() const { return elements_.size(); }
vcl_size_t internal_size() const { return elements_.internal_size(); }
entry_proxy<NumericT> operator()(vcl_size_t row_index, vcl_size_t col_index)
{
long index = static_cast<long>(row_index) - static_cast<long>(col_index);
assert(row_index < size1() && col_index < size2() && bool("Invalid access"));
while (index < 0)
index += static_cast<long>(size1());
return elements_[static_cast<vcl_size_t>(index)];
}
circulant_matrix<NumericT, AlignmentV>& operator +=(circulant_matrix<NumericT, AlignmentV>& that)
{
elements_ += that.elements();
return *this;
}
private:
circulant_matrix(circulant_matrix const &) {}
circulant_matrix & operator=(circulant_matrix const & t);
viennacl::vector<NumericT, AlignmentV> elements_;
};
template<typename NumericT, unsigned int AlignmentV>
void copy(std::vector<NumericT>& cpu_vec, circulant_matrix<NumericT, AlignmentV>& gpu_mat)
{
assert( (gpu_mat.size1() == 0 || cpu_vec.size() == gpu_mat.size1()) && bool("Size mismatch"));
copy(cpu_vec, gpu_mat.elements());
}
template<typename NumericT, unsigned int AlignmentV>
void copy(circulant_matrix<NumericT, AlignmentV>& gpu_mat, std::vector<NumericT>& cpu_vec)
{
assert(cpu_vec.size() == gpu_mat.size1() && bool("Size mismatch"));
copy(gpu_mat.elements(), cpu_vec);
}
template<typename NumericT, unsigned int AlignmentV, typename MatrixT>
void copy(circulant_matrix<NumericT, AlignmentV>& circ_src, MatrixT& com_dst)
{
vcl_size_t size = circ_src.size1();
assert(size == viennacl::traits::size1(com_dst) && bool("Size mismatch"));
assert(size == viennacl::traits::size2(com_dst) && bool("Size mismatch"));
std::vector<NumericT> tmp(size);
copy(circ_src, tmp);
for (vcl_size_t i = 0; i < size; i++)
{
for (vcl_size_t j = 0; j < size; j++)
{
long index = static_cast<long>(i) - static_cast<long>(j);
if (index < 0)
index += static_cast<long>(size);
com_dst(i, j) = tmp[static_cast<vcl_size_t>(index)];
}
}
}
template<typename NumericT, unsigned int AlignmentV, typename MatrixT>
void copy(MatrixT& com_src, circulant_matrix<NumericT, AlignmentV>& circ_dst)
{
assert( (circ_dst.size1() == 0 || circ_dst.size1() == viennacl::traits::size1(com_src)) && bool("Size mismatch"));
assert( (circ_dst.size2() == 0 || circ_dst.size2() == viennacl::traits::size2(com_src)) && bool("Size mismatch"));
vcl_size_t size = viennacl::traits::size1(com_src);
std::vector<NumericT> tmp(size);
for (vcl_size_t i = 0; i < size; i++) tmp[i] = com_src(i, 0);
copy(tmp, circ_dst);
}
/*namespace linalg
{
template<typename NumericT, unsigned int AlignmentV, unsigned int VECTOR_AlignmentV>
void prod_impl(circulant_matrix<NumericT, AlignmentV> const & mat,
vector<NumericT, VECTOR_AlignmentV> const & vec,
vector<NumericT, VECTOR_AlignmentV>& result) {
viennacl::vector<NumericT, VECTOR_AlignmentV> circ(mat.elements().size() * 2);
fft::real_to_complex(mat.elements(), circ, mat.elements().size());
viennacl::vector<NumericT, VECTOR_AlignmentV> tmp(vec.size() * 2);
viennacl::vector<NumericT, VECTOR_AlignmentV> tmp2(vec.size() * 2);
fft::real_to_complex(vec, tmp, vec.size());
fft::convolve(circ, tmp, tmp2);
fft::complex_to_real(tmp2, result, vec.size());
}
}*/
template<class NumericT, unsigned int AlignmentV>
std::ostream & operator<<(std::ostream& s, circulant_matrix<NumericT, AlignmentV>& gpu_matrix)
{
vcl_size_t size = gpu_matrix.size1();
std::vector<NumericT> tmp(size);
copy(gpu_matrix, tmp);
s << "[" << size << "," << size << "](";
for (vcl_size_t i = 0; i < size; i++)
{
s << "(";
for (vcl_size_t j = 0; j < size; j++)
{
long index = static_cast<long>(i) - static_cast<long>(j);
if (index < 0) index = static_cast<long>(size) + index;
s << tmp[vcl_size_t(index)];
//s << index;
if (j < (size - 1)) s << ",";
}
s << ")";
}
s << ")";
return s;
}
//
// Specify available operations:
//
namespace linalg
{
namespace detail
{
// x = A * y
template<typename T, unsigned int A>
struct op_executor<vector_base<T>, op_assign, vector_expression<const circulant_matrix<T, A>, const vector_base<T>, op_prod> >
{
static void apply(vector_base<T> & lhs, vector_expression<const circulant_matrix<T, A>, const vector_base<T>, op_prod> const & rhs)
{
// check for the special case x = A * x
if (viennacl::traits::handle(lhs) == viennacl::traits::handle(rhs.rhs()))
{
viennacl::vector<T> temp(lhs);
viennacl::linalg::prod_impl(rhs.lhs(), rhs.rhs(), temp);
lhs = temp;
}
else
viennacl::linalg::prod_impl(rhs.lhs(), rhs.rhs(), lhs);
}
};
template<typename T, unsigned int A>
struct op_executor<vector_base<T>, op_inplace_add, vector_expression<const circulant_matrix<T, A>, const vector_base<T>, op_prod> >
{
static void apply(vector_base<T> & lhs, vector_expression<const circulant_matrix<T, A>, const vector_base<T>, op_prod> const & rhs)
{
viennacl::vector<T> temp(lhs);
viennacl::linalg::prod_impl(rhs.lhs(), rhs.rhs(), temp);
lhs += temp;
}
};
template<typename T, unsigned int A>
struct op_executor<vector_base<T>, op_inplace_sub, vector_expression<const circulant_matrix<T, A>, const vector_base<T>, op_prod> >
{
static void apply(vector_base<T> & lhs, vector_expression<const circulant_matrix<T, A>, const vector_base<T>, op_prod> const & rhs)
{
viennacl::vector<T> temp(lhs);
viennacl::linalg::prod_impl(rhs.lhs(), rhs.rhs(), temp);
lhs -= temp;
}
};
// x = A * vec_op
template<typename T, unsigned int A, typename LHS, typename RHS, typename OP>
struct op_executor<vector_base<T>, op_assign, vector_expression<const circulant_matrix<T, A>, const vector_expression<const LHS, const RHS, OP>, op_prod> >
{
static void apply(vector_base<T> & lhs, vector_expression<const circulant_matrix<T, A>, const vector_expression<const LHS, const RHS, OP>, op_prod> const & rhs)
{
viennacl::vector<T> temp(rhs.rhs());
viennacl::linalg::prod_impl(rhs.lhs(), temp, lhs);
}
};
// x = A * vec_op
template<typename T, unsigned int A, typename LHS, typename RHS, typename OP>
struct op_executor<vector_base<T>, op_inplace_add, vector_expression<const circulant_matrix<T, A>, vector_expression<const LHS, const RHS, OP>, op_prod> >
{
static void apply(vector_base<T> & lhs, vector_expression<const circulant_matrix<T, A>, vector_expression<const LHS, const RHS, OP>, op_prod> const & rhs)
{
viennacl::vector<T> temp(rhs.rhs());
viennacl::vector<T> temp_result(lhs);
viennacl::linalg::prod_impl(rhs.lhs(), temp, temp_result);
lhs += temp_result;
}
};
// x = A * vec_op
template<typename T, unsigned int A, typename LHS, typename RHS, typename OP>
struct op_executor<vector_base<T>, op_inplace_sub, vector_expression<const circulant_matrix<T, A>, const vector_expression<const LHS, const RHS, OP>, op_prod> >
{
static void apply(vector_base<T> & lhs, vector_expression<const circulant_matrix<T, A>, const vector_expression<const LHS, const RHS, OP>, op_prod> const & rhs)
{
viennacl::vector<T> temp(rhs.rhs());
viennacl::vector<T> temp_result(lhs);
viennacl::linalg::prod_impl(rhs.lhs(), temp, temp_result);
lhs -= temp_result;
}
};
} // namespace detail
} // namespace linalg
}
#endif // VIENNACL_CIRCULANT_MATRIX_HPP
viennacl::circulant_matrix::elements
viennacl::vector< NumericT, AlignmentV > const & elements() const
Definition: circulant_matrix.hpp:87
viennacl::scalar
This class represents a single scalar value on the GPU and behaves mostly like a built-in scalar type...
Definition: forwards.h:227
viennacl::circulant_matrix::elements
viennacl::vector< NumericT, AlignmentV > & elements()
Returns an internal viennacl::vector, which represents a circulant matrix elements.
Definition: circulant_matrix.hpp:86
circulant_matrix_operations.hpp
Implementations of operations using circulant_matrix. Experimental.
viennacl::circulant_matrix::handle_type
viennacl::backend::mem_handle handle_type
Definition: circulant_matrix.hpp:44
viennacl::traits::size1
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
viennacl::circulant_matrix::value_type
scalar< typename viennacl::tools::CHECK_SCALAR_TEMPLATE_ARGUMENT< NumericT >::ResultType > value_type
Definition: circulant_matrix.hpp:45
forwards.h
This file provides the forward declarations for the main types used within ViennaCL.
viennacl::circulant_matrix::resize
void resize(vcl_size_t sz, bool preserve=true)
Resizes the matrix. Existing entries can be preserved.
Definition: circulant_matrix.hpp:71
viennacl::traits::size2
result_of::size_type< MatrixType >::type size2(MatrixType const &mat)
Generic routine for obtaining the number of columns of a matrix (ViennaCL, uBLAS, etc...
Definition: size.hpp:201
viennacl::traits::size
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.)
Definition: size.hpp:239
viennacl::circulant_matrix::internal_size
vcl_size_t internal_size() const
Returns the internal size of matrix representtion. Usually required for launching OpenCL kernels only...
Definition: circulant_matrix.hpp:104
viennacl::circulant_matrix::size1
vcl_size_t size1() const
Returns the number of rows of the matrix.
Definition: circulant_matrix.hpp:92
viennacl::vector::resize
void resize(size_type new_size, bool preserve=true)
Resizes the allocated memory for the vector. Pads the memory to be a multiple of 'AlignmentV'.
Definition: vector.hpp:1046
viennacl::vcl_size_t
std::size_t vcl_size_t
Definition: forwards.h:75
viennacl::vector< NumericT, AlignmentV >
viennacl::circulant_matrix::operator()
entry_proxy< NumericT > operator()(vcl_size_t row_index, vcl_size_t col_index)
Read-write access to a single element of the matrix.
Definition: circulant_matrix.hpp:113
viennacl::circulant_matrix::circulant_matrix
circulant_matrix()
The default constructor. Does not allocate any memory.
Definition: circulant_matrix.hpp:51
backend.hpp
Implementations of the OpenCL backend, where all contexts are stored in.
viennacl::circulant_matrix::handle
handle_type const & handle() const
Returns the OpenCL handle.
Definition: circulant_matrix.hpp:80
vector.hpp
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
viennacl::copy
void copy(std::vector< NumericT > &cpu_vec, circulant_matrix< NumericT, AlignmentV > &gpu_mat)
Copies a circulant matrix from the std::vector to the OpenCL device (either GPU or multi-core CPU) ...
Definition: circulant_matrix.hpp:150
viennacl::vector_base::size
size_type size() const
Returns the length of the vector (cf. std::vector)
Definition: vector_def.hpp:118
fft.hpp
All routines related to the Fast Fourier Transform. Experimental.
viennacl::backend::mem_handle
Main abstraction class for multiple memory domains. Represents a buffer in either main RAM...
Definition: mem_handle.hpp:89
viennacl::circulant_matrix
A Circulant matrix class.
Definition: circulant_matrix.hpp:41
viennacl::circulant_matrix::circulant_matrix
circulant_matrix(vcl_size_t rows, vcl_size_t cols)
Creates the matrix with the given size.
Definition: circulant_matrix.hpp:59
viennacl::circulant_matrix::size2
vcl_size_t size2() const
Returns the number of columns of the matrix.
Definition: circulant_matrix.hpp:97
viennacl::linalg::prod_impl
void prod_impl(const matrix_base< NumericT > &mat, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication.
Definition: matrix_operations.hpp:438
viennacl::vector_base::internal_size
size_type internal_size() const
Returns the internal length of the vector, which is given by size() plus the extra memory due to padd...
Definition: vector_def.hpp:120
viennacl::traits::handle
viennacl::backend::mem_handle & handle(T &obj)
Returns the generic memory handle of an object. Non-const version.
Definition: handle.hpp:41
viennacl::entry_proxy
A proxy class for a single element of a vector or matrix. This proxy should not be noticed by end-use...
Definition: forwards.h:233
op_assign
Definition: self_assign.cpp:132
viennacl::vector_base::handle
const handle_type & handle() const
Returns the memory handle.
Definition: vector_def.hpp:128
viennacl::circulant_matrix::operator+=
circulant_matrix< NumericT, AlignmentV > & operator+=(circulant_matrix< NumericT, AlignmentV > &that)
+= operation for circulant matrices
Definition: circulant_matrix.hpp:130