All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
qr-method-common.hpp
Go to the documentation of this file.
#ifndef VIENNACL_LINALG_QR_METHOD_COMMON_HPP
#define VIENNACL_LINALG_QR_METHOD_COMMON_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 <cmath>
#ifdef VIENNACL_WITH_OPENCL
#include "viennacl/ocl/device.hpp"
#include "viennacl/ocl/handle.hpp"
#include "viennacl/ocl/kernel.hpp"
#include "viennacl/linalg/opencl/kernels/svd.hpp"
#endif
#ifdef VIENNACL_WITH_CUDA
#include "viennacl/linalg/cuda/matrix_operations.hpp"
#endif
#include "viennacl/meta/result_of.hpp"
#include "viennacl/vector.hpp"
#include "viennacl/matrix.hpp"
//#include <boost/numeric/ublas/vector.hpp>
//#include <boost/numeric/ublas/io.hpp>
namespace viennacl
{
namespace linalg
{
const std::string SVD_HOUSEHOLDER_UPDATE_QR_KERNEL = "house_update_QR";
const std::string SVD_MATRIX_TRANSPOSE_KERNEL = "transpose_inplace";
const std::string SVD_INVERSE_SIGNS_KERNEL = "inverse_signs";
const std::string SVD_GIVENS_PREV_KERNEL = "givens_prev";
const std::string SVD_FINAL_ITER_UPDATE_KERNEL = "final_iter_update";
const std::string SVD_UPDATE_QR_COLUMN_KERNEL = "update_qr_column";
const std::string SVD_HOUSEHOLDER_UPDATE_A_LEFT_KERNEL = "house_update_A_left";
const std::string SVD_HOUSEHOLDER_UPDATE_A_RIGHT_KERNEL = "house_update_A_right";
const std::string SVD_HOUSEHOLDER_UPDATE_QL_KERNEL = "house_update_QL";
namespace detail
{
static const double EPS = 1e-10;
static const vcl_size_t ITER_MAX = 50;
template <typename SCALARTYPE>
SCALARTYPE pythag(SCALARTYPE a, SCALARTYPE b)
{
return std::sqrt(a*a + b*b);
}
template <typename SCALARTYPE>
SCALARTYPE sign(SCALARTYPE val)
{
return (val >= 0) ? SCALARTYPE(1) : SCALARTYPE(-1);
}
// DEPRECATED: Replace with viennacl::linalg::norm_2
template <typename VectorType>
typename VectorType::value_type norm_lcl(VectorType const & x, vcl_size_t size)
{
typename VectorType::value_type x_norm = 0.0;
for(vcl_size_t i = 0; i < size; i++)
x_norm += std::pow(x[i], 2);
return std::sqrt(x_norm);
}
template <typename VectorType>
void normalize(VectorType & x, vcl_size_t size)
{
typename VectorType::value_type x_norm = norm_lcl(x, size);
for(vcl_size_t i = 0; i < size; i++)
x[i] /= x_norm;
}
template <typename VectorType>
void householder_vector(VectorType & v, vcl_size_t start)
{
typedef typename VectorType::value_type ScalarType;
ScalarType x_norm = norm_lcl(v, v.size());
ScalarType alpha = -sign(v[start]) * x_norm;
v[start] += alpha;
normalize(v, v.size());
}
template <typename SCALARTYPE>
void transpose(matrix_base<SCALARTYPE> & A)
{
(void)A;
#ifdef VIENNACL_WITH_OPENCL
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(A).context());
if(A.row_major())
{
viennacl::linalg::opencl::kernels::svd<SCALARTYPE, row_major>::init(ctx);
viennacl::ocl::kernel & kernel = viennacl::ocl::get_kernel(viennacl::linalg::opencl::kernels::svd<SCALARTYPE, row_major>::program_name(), SVD_MATRIX_TRANSPOSE_KERNEL);
viennacl::ocl::enqueue(kernel(A,
static_cast<cl_uint>(A.internal_size1()),
static_cast<cl_uint>(A.internal_size2())
)
);
}
else
{
viennacl::linalg::opencl::kernels::svd<SCALARTYPE, row_major>::init(ctx);
viennacl::ocl::kernel & kernel = viennacl::ocl::get_kernel(viennacl::linalg::opencl::kernels::svd<SCALARTYPE, column_major>::program_name(), SVD_MATRIX_TRANSPOSE_KERNEL);
viennacl::ocl::enqueue(kernel(A,
static_cast<cl_uint>(A.internal_size1()),
static_cast<cl_uint>(A.internal_size2())
)
);
}
#endif
}
template <typename T>
void cdiv(T xr, T xi, T yr, T yi, T& cdivr, T& cdivi)
{
// Complex scalar division.
T r;
T d;
if (std::fabs(yr) > std::fabs(yi))
{
r = yi / yr;
d = yr + r * yi;
cdivr = (xr + r * xi) / d;
cdivi = (xi - r * xr) / d;
}
else
{
r = yr / yi;
d = yi + r * yr;
cdivr = (r * xr + xi) / d;
cdivi = (r * xi - xr) / d;
}
}
template<typename SCALARTYPE>
void prepare_householder_vector(
matrix_base<SCALARTYPE>& A,
vector_base<SCALARTYPE>& D,
vcl_size_t size,
vcl_size_t row_start,
vcl_size_t col_start,
vcl_size_t start,
bool is_column
)
{
//boost::numeric::ublas::vector<SCALARTYPE> tmp = boost::numeric::ublas::scalar_vector<SCALARTYPE>(size, 0);
std::vector<SCALARTYPE> tmp(size);
copy_vec(A, D, row_start, col_start, is_column);
fast_copy(D.begin(), D.begin() + vcl_ptrdiff_t(size - start), tmp.begin() + vcl_ptrdiff_t(start));
detail::householder_vector(tmp, start);
fast_copy(tmp, D);
}
} //detail
}
}
#endif
EPS
#define EPS
Definition: bisect.cpp:38
device.hpp
Represents an OpenCL device within ViennaCL.
viennacl::ocl::get_kernel
viennacl::ocl::kernel & get_kernel(std::string const &prog_name, std::string const &kernel_name)
Convenience function for getting the kernel for a particular program from the current active context...
Definition: backend.hpp:339
viennacl::linalg::SVD_INVERSE_SIGNS_KERNEL
const std::string SVD_INVERSE_SIGNS_KERNEL
Definition: qr-method-common.hpp:51
viennacl::linalg::copy_vec
void copy_vec(matrix_base< SCALARTYPE > &A, vector_base< SCALARTYPE > &V, vcl_size_t row_start, vcl_size_t col_start, bool copy_col)
This function copies a row or a column from a matrix to a vector.
Definition: matrix_operations.hpp:942
viennacl::ocl::kernel
Represents an OpenCL kernel within ViennaCL.
Definition: kernel.hpp:58
matrix.hpp
Implementation of the dense matrix class.
viennacl::ocl::context
Manages an OpenCL context and provides the respective convenience functions for creating buffers...
Definition: context.hpp:55
viennacl::linalg::SVD_UPDATE_QR_COLUMN_KERNEL
const std::string SVD_UPDATE_QR_COLUMN_KERNEL
Definition: qr-method-common.hpp:54
viennacl::matrix_base
Definition: matrix_def.hpp:103
viennacl::linalg::SVD_HOUSEHOLDER_UPDATE_A_RIGHT_KERNEL
const std::string SVD_HOUSEHOLDER_UPDATE_A_RIGHT_KERNEL
Definition: qr-method-common.hpp:56
viennacl::linalg::detail::cdiv
void cdiv(T xr, T xi, T yr, T yi, T &cdivr, T &cdivi)
Definition: qr-method-common.hpp:142
viennacl::linalg::SVD_HOUSEHOLDER_UPDATE_QR_KERNEL
const std::string SVD_HOUSEHOLDER_UPDATE_QR_KERNEL
Definition: qr-method-common.hpp:49
viennacl::linalg::detail::transpose
void transpose(matrix_base< SCALARTYPE > &A)
Definition: qr-method-common.hpp:107
viennacl::linalg::opencl::kernels::svd::init
static void init(viennacl::ocl::context &ctx)
Definition: svd.hpp:652
viennacl::linalg::detail::prepare_householder_vector
void prepare_householder_vector(matrix_base< SCALARTYPE > &A, vector_base< SCALARTYPE > &D, vcl_size_t size, vcl_size_t row_start, vcl_size_t col_start, vcl_size_t start, bool is_column)
Definition: qr-method-common.hpp:165
viennacl::linalg::detail::normalize
void normalize(VectorType &x, vcl_size_t size)
Definition: qr-method-common.hpp:87
svd.hpp
OpenCL kernel file for singular value decomposition.
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::linalg::opencl::kernels::svd
Main kernel class for generating OpenCL kernels for singular value decomposition of dense matrices...
Definition: svd.hpp:644
viennacl::vector_base::begin
iterator begin()
Returns an iterator pointing to the beginning of the vector (STL like)
Definition: vector.hpp:841
handle.hpp
Implementation of a smart-pointer-like class for handling OpenCL handles.
viennacl::traits::start
result_of::size_type< T >::type start(T const &obj)
Definition: start.hpp:44
viennacl::vector_base
Common base class for dense vectors, vector ranges, and vector slices.
Definition: vector_def.hpp:104
viennacl::linalg::detail::householder_vector
void householder_vector(VectorType &v, vcl_size_t start)
Definition: qr-method-common.hpp:97
viennacl::linalg::SVD_HOUSEHOLDER_UPDATE_A_LEFT_KERNEL
const std::string SVD_HOUSEHOLDER_UPDATE_A_LEFT_KERNEL
Definition: qr-method-common.hpp:55
viennacl::vcl_size_t
std::size_t vcl_size_t
Definition: forwards.h:75
viennacl::linalg::detail::pythag
SCALARTYPE pythag(SCALARTYPE a, SCALARTYPE b)
Definition: qr-method-common.hpp:65
matrix_operations.hpp
Implementations of dense matrix related operations, including matrix-vector products, using CUDA.
viennacl::ocl::enqueue
void enqueue(KernelType &k, viennacl::ocl::command_queue const &queue)
Enqueues a kernel in the provided queue.
Definition: enqueue.hpp:50
kernel.hpp
Representation of an OpenCL kernel in ViennaCL.
vector.hpp
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
viennacl::matrix_base::row_major
bool row_major() const
Definition: matrix_def.hpp:248
viennacl::linalg::SVD_HOUSEHOLDER_UPDATE_QL_KERNEL
const std::string SVD_HOUSEHOLDER_UPDATE_QL_KERNEL
Definition: qr-method-common.hpp:57
ScalarType
float ScalarType
Definition: fft_1d.cpp:42
viennacl::linalg::SVD_FINAL_ITER_UPDATE_KERNEL
const std::string SVD_FINAL_ITER_UPDATE_KERNEL
Definition: qr-method-common.hpp:53
viennacl::linalg::SVD_MATRIX_TRANSPOSE_KERNEL
const std::string SVD_MATRIX_TRANSPOSE_KERNEL
Definition: qr-method-common.hpp:50
viennacl::matrix_base::internal_size2
size_type internal_size2() const
Returns the internal number of columns. Usually required for launching OpenCL kernels only...
Definition: matrix_def.hpp:240
viennacl::matrix_base::internal_size1
size_type internal_size1() const
Returns the internal number of rows. Usually required for launching OpenCL kernels only...
Definition: matrix_def.hpp:238
viennacl::linalg::detail::norm_lcl
VectorType::value_type norm_lcl(VectorType const &x, vcl_size_t size)
Definition: qr-method-common.hpp:78
result_of.hpp
A collection of compile time type deductions.
viennacl::linalg::SVD_GIVENS_PREV_KERNEL
const std::string SVD_GIVENS_PREV_KERNEL
Definition: qr-method-common.hpp:52
viennacl::vcl_ptrdiff_t
std::ptrdiff_t vcl_ptrdiff_t
Definition: forwards.h:76
viennacl::linalg::detail::sign
SCALARTYPE sign(SCALARTYPE val)
Definition: qr-method-common.hpp:71
viennacl::fast_copy
void fast_copy(const const_vector_iterator< SCALARTYPE, ALIGNMENT > &gpu_begin, const const_vector_iterator< SCALARTYPE, ALIGNMENT > &gpu_end, CPU_ITERATOR cpu_begin)