All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
fft_operations.hpp
Go to the documentation of this file.
#ifndef VIENNACL_LINALG_OPENCL_FFT_OPERATIONS_HPP_
#define VIENNACL_LINALG_OPENCL_FFT_OPERATIONS_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/ocl/device.hpp"
#include "viennacl/ocl/kernel.hpp"
#include "viennacl/traits/handle.hpp"
#include "viennacl/traits/stride.hpp"
#include "viennacl/linalg/host_based/fft_operations.hpp"
#include "viennacl/linalg/opencl/kernels/fft.hpp"
#include "viennacl/linalg/opencl/kernels/matrix.hpp"
#include <viennacl/vector.hpp>
#include <viennacl/matrix.hpp>
#include <cmath>
#include <stdexcept>
namespace viennacl
{
namespace linalg
{
namespace detail
{
namespace fft
{
const vcl_size_t MAX_LOCAL_POINTS_NUM = 512;
inline vcl_size_t num_bits(vcl_size_t size)
{
vcl_size_t bits_datasize = 0;
vcl_size_t ds = 1;
while (ds < size)
{
ds = ds << 1;
bits_datasize++;
}
return bits_datasize;
}
inline vcl_size_t next_power_2(vcl_size_t n)
{
n = n - 1;
vcl_size_t power = 1;
while (power < sizeof(vcl_size_t) * 8)
{
n = n | (n >> power);
power *= 2;
}
return n + 1;
}
} //namespce fft
} //namespace detail
namespace opencl
{
template<typename NumericT>
void direct(viennacl::ocl::handle<cl_mem> const & in,
viennacl::ocl::handle<cl_mem> const & out,
vcl_size_t size, vcl_size_t stride, vcl_size_t batch_num, NumericT sign = NumericT(-1),
viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER data_order = viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(in.context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
std::string program_string = viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, row_major>::program_name();
if (data_order == viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::COL_MAJOR)
{
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, column_major>::init(ctx);
program_string =
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, column_major>::program_name();
} else
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, row_major>::init(ctx);
viennacl::ocl::kernel & k = ctx.get_kernel(program_string, "fft_direct");
viennacl::ocl::enqueue(k(in, out,
static_cast<cl_uint>(size),
static_cast<cl_uint>(stride),
static_cast<cl_uint>(batch_num),
sign)
);
}
/*
* This function performs reorder of input data. Indexes are sorted in bit-reversal order.
* Such reordering should be done before in-place FFT.
*/
template<typename NumericT>
void reorder(viennacl::ocl::handle<cl_mem> const & in,
vcl_size_t size, vcl_size_t stride,
vcl_size_t bits_datasize, vcl_size_t batch_num,
viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER data_order = viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(in.context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
std::string program_string = viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, row_major>::program_name();
if (data_order == viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::COL_MAJOR)
{
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, column_major>::init(ctx);
program_string = viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, column_major>::program_name();
} else
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, row_major>::init(ctx);
viennacl::ocl::kernel& k = ctx.get_kernel(program_string, "fft_reorder");
viennacl::ocl::enqueue(k(in,
static_cast<cl_uint>(bits_datasize), static_cast<cl_uint>(size),
static_cast<cl_uint>(stride), static_cast<cl_uint>(batch_num))
);
}
template<typename NumericT>
void radix2(viennacl::ocl::handle<cl_mem> const & in,
vcl_size_t size, vcl_size_t stride,
vcl_size_t batch_num, NumericT sign = NumericT(-1),
viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER data_order = viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(in.context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
assert(batch_num != 0 && bool("batch_num must be larger than 0"));
std::string program_string = viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, row_major>::program_name();
if (data_order == viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::COL_MAJOR)
{
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, column_major>::init(ctx);
program_string = viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, column_major>::program_name();
} else
viennacl::linalg::opencl::kernels::matrix_legacy<NumericT, row_major>::init(ctx);
vcl_size_t bits_datasize = viennacl::linalg::detail::fft::num_bits(size);
if (size <= viennacl::linalg::detail::fft::MAX_LOCAL_POINTS_NUM)
{
viennacl::ocl::kernel & k = ctx.get_kernel(program_string, "fft_radix2_local");
viennacl::ocl::enqueue(k(in,
viennacl::ocl::local_mem((size * 4) * sizeof(NumericT)),
static_cast<cl_uint>(bits_datasize), static_cast<cl_uint>(size),
static_cast<cl_uint>(stride), static_cast<cl_uint>(batch_num), sign));
}
else
{
viennacl::linalg::opencl::reorder<NumericT>(in, size, stride, bits_datasize, batch_num);
for (vcl_size_t step = 0; step < bits_datasize; step++)
{
viennacl::ocl::kernel & k = ctx.get_kernel(program_string, "fft_radix2");
viennacl::ocl::enqueue(k(in,
static_cast<cl_uint>(step), static_cast<cl_uint>(bits_datasize),
static_cast<cl_uint>(size), static_cast<cl_uint>(stride),
static_cast<cl_uint>(batch_num), sign));
}
}
}
template<typename NumericT, unsigned int AlignmentV>
void bluestein(viennacl::vector<NumericT, AlignmentV>& in,
viennacl::vector<NumericT, AlignmentV>& out, vcl_size_t /*batch_num*/)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(in).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
vcl_size_t size = in.size() >> 1;
vcl_size_t ext_size = viennacl::linalg::detail::fft::next_power_2(2 * size - 1);
viennacl::vector<NumericT, AlignmentV> A(ext_size << 1);
viennacl::vector<NumericT, AlignmentV> B(ext_size << 1);
viennacl::vector<NumericT, AlignmentV> Z(ext_size << 1);
{
viennacl::ocl::kernel & k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "zero2");
viennacl::ocl::enqueue(k(A, B, static_cast<cl_uint>(ext_size)));
}
{
viennacl::ocl::kernel & k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "bluestein_pre");
viennacl::ocl::enqueue(k(in, A, B, static_cast<cl_uint>(size), static_cast<cl_uint>(ext_size)));
}
viennacl::linalg::convolve_i(A, B, Z);
{
viennacl::ocl::kernel& k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "bluestein_post");
viennacl::ocl::enqueue(k(Z, out, static_cast<cl_uint>(size)));
}
}
template<typename NumericT, unsigned int AlignmentV>
void multiply_complex(viennacl::vector<NumericT, AlignmentV> const & input1,
viennacl::vector<NumericT, AlignmentV> const & input2,
viennacl::vector<NumericT, AlignmentV> & output)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(input1).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
vcl_size_t size = input1.size() >> 1;
viennacl::ocl::kernel & k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "fft_mult_vec");
viennacl::ocl::enqueue(k(input1, input2, output, static_cast<cl_uint>(size)));
}
template<typename NumericT, unsigned int AlignmentV>
void normalize(viennacl::vector<NumericT, AlignmentV> & input)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(input).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
viennacl::ocl::kernel & k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "fft_div_vec_scalar");
vcl_size_t size = input.size() >> 1;
NumericT norm_factor = static_cast<NumericT>(size);
viennacl::ocl::enqueue(k(input, static_cast<cl_uint>(size), norm_factor));
}
template<typename NumericT, unsigned int AlignmentV>
void transpose(viennacl::matrix<NumericT, viennacl::row_major, AlignmentV> & input)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(input).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
viennacl::ocl::kernel& k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "transpose_inplace");
viennacl::ocl::enqueue(k(input, static_cast<cl_uint>(input.internal_size1() >> 1),
static_cast<cl_uint>(input.internal_size2()) >> 1));
}
template<typename NumericT, unsigned int AlignmentV>
void transpose(viennacl::matrix<NumericT, viennacl::row_major, AlignmentV> const & input,
viennacl::matrix<NumericT, viennacl::row_major, AlignmentV> & output)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(input).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
viennacl::ocl::kernel& k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "transpose");
viennacl::ocl::enqueue(k(input, output, static_cast<cl_uint>(input.internal_size1() >> 1),
static_cast<cl_uint>(input.internal_size2() >> 1)));
}
template<typename NumericT>
void real_to_complex(viennacl::vector_base<NumericT> const & in,
viennacl::vector_base<NumericT> & out, vcl_size_t size)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(in).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
viennacl::ocl::kernel & k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "real_to_complex");
viennacl::ocl::enqueue(k(in, out, static_cast<cl_uint>(size)));
}
template<typename NumericT>
void complex_to_real(viennacl::vector_base<NumericT> const & in,
viennacl::vector_base<NumericT> & out, vcl_size_t size)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(in).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
viennacl::ocl::kernel& k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "complex_to_real");
viennacl::ocl::enqueue(k(in, out, static_cast<cl_uint>(size)));
}
template<typename NumericT>
void reverse(viennacl::vector_base<NumericT>& in)
{
viennacl::ocl::context & ctx = const_cast<viennacl::ocl::context &>(viennacl::traits::opencl_handle(in).context());
viennacl::linalg::opencl::kernels::fft<NumericT>::init(ctx);
vcl_size_t size = in.size();
viennacl::ocl::kernel& k = ctx.get_kernel(viennacl::linalg::opencl::kernels::fft<NumericT>::program_name(), "reverse_inplace");
viennacl::ocl::enqueue(k(in, static_cast<cl_uint>(size)));
}
} //namespace opencl
} //namespace linalg
} //namespace viennacl
#endif /* FFT_OPERATIONS_HPP_ */
viennacl::linalg::opencl::reorder
void reorder(viennacl::ocl::handle< cl_mem > const &in, vcl_size_t size, vcl_size_t stride, vcl_size_t bits_datasize, vcl_size_t batch_num, viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER data_order=viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR)
Definition: fft_operations.hpp:130
fft.hpp
OpenCL kernel file for FFT operations.
viennacl::linalg::detail::fft::next_power_2
vcl_size_t next_power_2(vcl_size_t n)
Find next power of two.
Definition: fft_operations.hpp:71
device.hpp
Represents an OpenCL device within ViennaCL.
viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER
DATA_ORDER
Definition: fft_operations.hpp:49
viennacl::linalg::opencl::complex_to_real
void complex_to_real(viennacl::vector_base< NumericT > const &in, viennacl::vector_base< NumericT > &out, vcl_size_t size)
Create real vector from complex vector (even elements(2*k) = real part, odd elements(2*k+1) = imagina...
Definition: fft_operations.hpp:320
viennacl::ocl::kernel
Represents an OpenCL kernel within ViennaCL.
Definition: kernel.hpp:58
matrix.hpp
Implementation of the dense matrix class.
viennacl::linalg::opencl::kernels::fft
Main kernel class for generating OpenCL kernels for the fast Fourier transform.
Definition: fft.hpp:260
viennacl::ocl::context
Manages an OpenCL context and provides the respective convenience functions for creating buffers...
Definition: context.hpp:55
viennacl::traits::stride
result_of::size_type< viennacl::vector_base< T > >::type stride(viennacl::vector_base< T > const &s)
Definition: stride.hpp:45
forwards.h
This file provides the forward declarations for the main types used within ViennaCL.
step
endcode *Final step
Definition: least-squares.cpp:145
viennacl::matrix
A dense matrix class.
Definition: forwards.h:375
stride.hpp
Determines row and column increments for matrices and matrix proxies.
viennacl::linalg::opencl::transpose
void transpose(viennacl::matrix< NumericT, viennacl::row_major, AlignmentV > &input)
Inplace_transpose matrix.
Definition: fft_operations.hpp:277
viennacl::ocl::handle::context
viennacl::ocl::context const & context() const
Definition: handle.hpp:193
NumericT
float NumericT
Definition: bisect.cpp:40
viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR
Definition: fft_operations.hpp:51
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::ocl::local_mem
A class representing local (shared) OpenCL memory. Typically used as kernel argument.
Definition: local_mem.hpp:33
viennacl::ocl::context::get_kernel
viennacl::ocl::kernel & get_kernel(std::string const &program_name, std::string const &kernel_name)
Convenience function for retrieving the kernel of a program directly from the context.
Definition: context.hpp:605
viennacl::linalg::opencl::kernels::matrix_legacy::init
static void init(viennacl::ocl::context &ctx)
Definition: matrix.hpp:1025
viennacl::vector_base< NumericT >
viennacl::linalg::opencl::radix2
void radix2(viennacl::ocl::handle< cl_mem > const &in, vcl_size_t size, vcl_size_t stride, vcl_size_t batch_num, NumericT sign=NumericT(-1), viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER data_order=viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR)
Radix-2 algorithm for computing Fourier transformation.
Definition: fft_operations.hpp:161
viennacl::linalg::opencl::reverse
void reverse(viennacl::vector_base< NumericT > &in)
Reverse vector to oposite order and save it in input vector.
Definition: fft_operations.hpp:334
viennacl::linalg::opencl::real_to_complex
void real_to_complex(viennacl::vector_base< NumericT > const &in, viennacl::vector_base< NumericT > &out, vcl_size_t size)
Create complex vector from real vector (even elements(2*k) = real part, odd elements(2*k+1) = imagina...
Definition: fft_operations.hpp:306
viennacl::vcl_size_t
std::size_t vcl_size_t
Definition: forwards.h:75
viennacl::linalg::opencl::kernels::matrix_legacy::program_name
static std::string program_name()
Definition: matrix.hpp:1020
viennacl::vector< NumericT, AlignmentV >
viennacl::linalg::convolve_i
void convolve_i(viennacl::vector< SCALARTYPE, ALIGNMENT > &input1, viennacl::vector< SCALARTYPE, ALIGNMENT > &input2, viennacl::vector< SCALARTYPE, ALIGNMENT > &output)
viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::COL_MAJOR
Definition: fft_operations.hpp:51
viennacl::linalg::detail::fft::MAX_LOCAL_POINTS_NUM
const vcl_size_t MAX_LOCAL_POINTS_NUM
Definition: fft_operations.hpp:49
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::linalg::opencl::multiply_complex
void multiply_complex(viennacl::vector< NumericT, AlignmentV > const &input1, viennacl::vector< NumericT, AlignmentV > const &input2, viennacl::vector< NumericT, AlignmentV > &output)
Mutiply two complex vectors and store result in output.
Definition: fft_operations.hpp:246
viennacl::vector_base::size
size_type size() const
Returns the length of the vector (cf. std::vector)
Definition: vector_def.hpp:118
viennacl::linalg::opencl::direct
void direct(viennacl::ocl::handle< cl_mem > const &in, viennacl::ocl::handle< cl_mem > const &out, vcl_size_t size, vcl_size_t stride, vcl_size_t batch_num, NumericT sign=NumericT(-1), viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::DATA_ORDER data_order=viennacl::linalg::host_based::detail::fft::FFT_DATA_ORDER::ROW_MAJOR)
Direct algorithm for computing Fourier transformation.
Definition: fft_operations.hpp:99
fft_operations.hpp
Implementations of Fast Furier Transformation using a plain single-threaded or OpenMP-enabled executi...
viennacl::linalg::detail::fft::num_bits
vcl_size_t num_bits(vcl_size_t size)
Get number of bits.
Definition: fft_operations.hpp:54
viennacl::matrix_base< NumericT >::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< NumericT >::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
handle.hpp
Extracts the underlying OpenCL handle from a vector, a matrix, an expression etc. ...
viennacl::linalg::opencl::kernels::fft::init
static void init(viennacl::ocl::context &ctx)
Definition: fft.hpp:267
viennacl::ocl::handle< cl_mem >
viennacl::linalg::opencl::normalize
void normalize(viennacl::vector< NumericT, AlignmentV > &input)
Normalize vector on with his own size.
Definition: fft_operations.hpp:261
viennacl::linalg::opencl::bluestein
void bluestein(viennacl::vector< NumericT, AlignmentV > &in, viennacl::vector< NumericT, AlignmentV > &out, vcl_size_t)
Bluestein's algorithm for computing Fourier transformation.
Definition: fft_operations.hpp:212
matrix.hpp
Runtime generation of OpenCL kernels for matrix operations.
fft
ScalarType fft(std::vector< ScalarType > &in, std::vector< ScalarType > &out, unsigned int, unsigned int, unsigned int batch_size)
Definition: fft_1d.cpp:719
viennacl::linalg::detail::sign
SCALARTYPE sign(SCALARTYPE val)
Definition: qr-method-common.hpp:71