00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015 #ifndef _VIENNACL_TOOLS_HPP_
00016 #define _VIENNACL_TOOLS_HPP_
00017
00022 #include <string>
00023 #include <fstream>
00024 #include <sstream>
00025 #include "viennacl/forwards.h"
00026 #include "viennacl/tools/adapter.hpp"
00027
00028
00029 #ifdef VIENNACL_HAVE_UBLAS
00030 #include <boost/numeric/ublas/matrix_sparse.hpp>
00031 #include <boost/numeric/ublas/matrix.hpp>
00032 #endif
00033
00034 #ifdef VIENNACL_HAVE_EIGEN
00035 #include <Eigen/Core>
00036 #include <Eigen/Sparse>
00037 #endif
00038
00039 #ifdef VIENNACL_HAVE_MTL4
00040 #include <boost/numeric/mtl/mtl.hpp>
00041 #endif
00042
00043 #include <vector>
00044 #include <map>
00045
00046 namespace viennacl
00047 {
00048 namespace tools
00049 {
00050 namespace traits
00051 {
00052
00053
00054 template <typename MatrixType>
00055 void resize(MatrixType & matrix, size_t rows, size_t cols)
00056 {
00057 matrix.resize(rows, cols);
00058 }
00059
00060 template <typename VectorType>
00061 void resize(VectorType & vec, size_t new_size)
00062 {
00063 vec.resize(new_size);
00064 }
00065
00066 #ifdef VIENNACL_HAVE_UBLAS
00067
00068 template <typename ScalarType>
00069 void resize(boost::numeric::ublas::compressed_matrix<ScalarType> & matrix,
00070 size_t rows,
00071 size_t cols)
00072 {
00073 matrix.resize(rows, cols, false);
00074 }
00075 #endif
00076
00077
00078 #ifdef VIENNACL_HAVE_MTL4
00079 template <typename ScalarType>
00080 void resize(mtl::compressed2D<ScalarType> & matrix,
00081 size_t rows,
00082 size_t cols)
00083 {
00084 matrix.change_dim(rows, cols);
00085 }
00086
00087 template <typename ScalarType>
00088 void resize(mtl::dense_vector<ScalarType> & vec,
00089 size_t new_size)
00090 {
00091 vec.change_dim(new_size);
00092 }
00093 #endif
00094
00095
00096 template <typename VectorType>
00097 void clear(VectorType & vec)
00098 {
00099 vec.clear();
00100 }
00101
00102 #ifdef VIENNACL_HAVE_EIGEN
00103 void clear(Eigen::VectorXf & vec) { for (int i=0; i<vec.size(); ++i) vec[i] = 0; }
00104 void clear(Eigen::VectorXd & vec) { for (int i=0; i<vec.size(); ++i) vec[i] = 0; }
00105 #endif
00106
00107 #ifdef VIENNACL_HAVE_MTL4
00108 template <typename ScalarType>
00109 void clear(mtl::dense_vector<ScalarType> & vec)
00110 {
00111 for (typename mtl::dense_vector<ScalarType>::size_type i=0;
00112 i<vec.used_memory();
00113 ++i)
00114 vec[i] = 0;
00115 }
00116 #endif
00117
00118
00119
00120 template <typename VectorType>
00121 unsigned int size(VectorType & vec)
00122 {
00123 return vec.size();
00124 }
00125
00126 #ifdef VIENNACL_HAVE_MTL4
00127 template <typename ScalarType>
00128 unsigned int size(mtl::dense_vector<ScalarType> const & vec) { return vec.used_memory(); }
00129 #endif
00130
00131 }
00132
00133
00134
00135 namespace result_of
00136 {
00137
00138 template <typename T>
00139 struct value_type
00140 {
00141 typedef typename T::value_type type;
00142 };
00143
00144 #ifdef VIENNACL_HAVE_EIGEN
00145 template <>
00146 struct value_type<Eigen::MatrixXf>
00147 {
00148 typedef Eigen::MatrixXf::RealScalar type;
00149 };
00150
00151 template <>
00152 struct value_type<Eigen::MatrixXd>
00153 {
00154 typedef Eigen::MatrixXd::RealScalar type;
00155 };
00156
00157 template <typename ScalarType, int option>
00158 struct value_type<Eigen::SparseMatrix<ScalarType, option> >
00159 {
00160 typedef ScalarType type;
00161 };
00162
00163 template <>
00164 struct value_type<Eigen::VectorXf>
00165 {
00166 typedef Eigen::VectorXf::RealScalar type;
00167 };
00168
00169 template <>
00170 struct value_type<Eigen::VectorXd>
00171 {
00172 typedef Eigen::VectorXd::RealScalar type;
00173 };
00174
00175 #endif
00176
00177 }
00178
00179
00180
00181
00183 template <class SCALARTYPE, typename F, unsigned int ALIGNMENT>
00184 struct MATRIX_ITERATOR_INCREMENTER<viennacl::row_iteration, viennacl::matrix<SCALARTYPE, F, ALIGNMENT> >
00185 {
00186 static void apply(const viennacl::matrix<SCALARTYPE, F, ALIGNMENT> & mat, unsigned int & row, unsigned int & col) { ++row; }
00187 };
00188
00189 template <class SCALARTYPE, typename F, unsigned int ALIGNMENT>
00190 struct MATRIX_ITERATOR_INCREMENTER<viennacl::col_iteration, viennacl::matrix<SCALARTYPE, F, ALIGNMENT> >
00191 {
00192 static void apply(const viennacl::matrix<SCALARTYPE, F, ALIGNMENT> & mat, unsigned int & row, unsigned int & col) { ++col; }
00193 };
00194
00195
00197 template <bool b, class T = void>
00198 struct enable_if
00199 {
00200 typedef T type;
00201 };
00202
00203 template <class T>
00204 struct enable_if<false, T> {};
00205
00206
00208 template <typename T>
00209 struct CHECK_SCALAR_TEMPLATE_ARGUMENT
00210 {
00211 typedef typename T::ERROR_SCALAR_MUST_HAVE_TEMPLATE_ARGUMENT_FLOAT_OR_DOUBLE ResultType;
00212 };
00213
00214 template <>
00215 struct CHECK_SCALAR_TEMPLATE_ARGUMENT<float>
00216 {
00217 typedef float ResultType;
00218 };
00219
00220 template <>
00221 struct CHECK_SCALAR_TEMPLATE_ARGUMENT<double>
00222 {
00223 typedef double ResultType;
00224 };
00225
00226
00227
00233 inline std::string readTextFromFile(const std::string & filename)
00234 {
00235 std::ifstream f(filename.c_str());
00236 if (!f) return std::string();
00237
00238 std::stringstream result;
00239 std::string tmp;
00240 while (std::getline(f, tmp))
00241 result << tmp << std::endl;
00242
00243 return result.str();
00244 }
00245
00253 inline std::string strReplace(const std::string & text, std::string to_search, std::string to_replace)
00254 {
00255 std::string::size_type pos = 0;
00256 std::string result;
00257 std::string::size_type found;
00258 while( (found = text.find(to_search, pos)) != std::string::npos )
00259 {
00260 result.append(text.substr(pos,found-pos));
00261 result.append(to_replace);
00262 pos = found + to_search.length();
00263 }
00264 if (pos < text.length())
00265 result.append(text.substr(pos));
00266 return result;
00267 }
00268
00276 template <class INT_TYPE>
00277 INT_TYPE roundUpToNextMultiple(INT_TYPE to_reach, INT_TYPE base)
00278 {
00279 if (to_reach % base == 0) return to_reach;
00280 return ((to_reach / base) + 1) * base;
00281 }
00282
00283
00290 inline std::string make_double_kernel(std::string const & source, std::string platform_info)
00291
00292 {
00293 std::stringstream ss;
00294 if (platform_info.compare(0, 8, "Advanced") == 0)
00295 ss << "#pragma OPENCL EXTENSION cl_amd_fp64 : enable\n\n";
00296 else
00297 ss << "#pragma OPENCL EXTENSION cl_khr_fp64 : enable\n\n";
00298
00299 std::string result = ss.str();
00300 result.append(strReplace(source, "float", "double"));
00301 return result;
00302 }
00303
00304
00306 template <typename T>
00307 struct CONST_REMOVER
00308 {
00309 typedef T ResultType;
00310 };
00311
00312 template <typename T>
00313 struct CONST_REMOVER<const T>
00314 {
00315 typedef T ResultType;
00316 };
00317
00318
00324 template <typename LHS, typename RHS>
00325 struct VECTOR_EXTRACTOR_IMPL
00326 {
00327 typedef typename LHS::ERROR_COULD_NOT_EXTRACT_VECTOR_INFORMATION_FROM_VECTOR_EXPRESSION ResultType;
00328 };
00329
00330 template <typename LHS, typename ScalarType, unsigned int A>
00331 struct VECTOR_EXTRACTOR_IMPL<LHS, viennacl::vector<ScalarType, A> >
00332 {
00333 typedef viennacl::vector<ScalarType, A> ResultType;
00334 };
00335
00336 template <typename RHS, typename ScalarType, unsigned int A>
00337 struct VECTOR_EXTRACTOR_IMPL<viennacl::vector<ScalarType, A>, RHS>
00338 {
00339 typedef viennacl::vector<ScalarType, A> ResultType;
00340 };
00341
00342
00343 template <typename ScalarType, unsigned int A>
00344 struct VECTOR_EXTRACTOR_IMPL<viennacl::vector<ScalarType, A>, viennacl::vector<ScalarType, A> >
00345 {
00346 typedef viennacl::vector<ScalarType, A> ResultType;
00347 };
00348
00349 template <typename LHS, typename RHS>
00350 struct VECTOR_EXTRACTOR
00351 {
00352 typedef typename VECTOR_EXTRACTOR_IMPL<typename CONST_REMOVER<LHS>::ResultType,
00353 typename CONST_REMOVER<RHS>::ResultType>::ResultType ResultType;
00354 };
00355
00362 template <typename LHS, typename RHS, typename OP>
00363 struct VECTOR_SIZE_DEDUCER
00364 {
00365
00366 };
00367
00368
00369 template <typename ScalarType, unsigned int A, typename RHS>
00370 struct VECTOR_SIZE_DEDUCER<const viennacl::vector<ScalarType, A>, RHS, viennacl::op_prod>
00371 {
00372 static size_t size(const viennacl::vector<ScalarType, A> & lhs,
00373 const RHS & rhs) { return lhs.size(); }
00374 };
00375
00376 template <typename ScalarType, unsigned int A, typename RHS>
00377 struct VECTOR_SIZE_DEDUCER<const viennacl::vector<ScalarType, A>, RHS, viennacl::op_div>
00378 {
00379 static size_t size(const viennacl::vector<ScalarType, A> & lhs,
00380 const RHS & rhs) { return lhs.size(); }
00381 };
00382
00383
00384 template <typename ScalarType, typename F, unsigned int Amat, unsigned int A>
00385 struct VECTOR_SIZE_DEDUCER<const viennacl::matrix<ScalarType, F, Amat>, const viennacl::vector<ScalarType, A>, viennacl::op_prod>
00386 {
00387 static size_t size(const viennacl::matrix<ScalarType, F, Amat> & lhs,
00388 const viennacl::vector<ScalarType, A> & rhs) { return lhs.size1(); }
00389 };
00390
00391 template <typename ScalarType, unsigned int Amat, unsigned int A>
00392 struct VECTOR_SIZE_DEDUCER<const viennacl::compressed_matrix<ScalarType, Amat>, const viennacl::vector<ScalarType, A>, viennacl::op_prod>
00393 {
00394 static size_t size(const viennacl::compressed_matrix<ScalarType, Amat> & lhs,
00395 const viennacl::vector<ScalarType, A> & rhs) { return lhs.size1(); }
00396 };
00397
00398 template <typename ScalarType, unsigned int Amat, unsigned int A>
00399 struct VECTOR_SIZE_DEDUCER<const viennacl::coordinate_matrix<ScalarType, Amat>, const viennacl::vector<ScalarType, A>, viennacl::op_prod>
00400 {
00401 static size_t size(const viennacl::coordinate_matrix<ScalarType, Amat> & lhs,
00402 const viennacl::vector<ScalarType, A> & rhs) { return lhs.size1(); }
00403 };
00404
00405
00406 template <typename ScalarType, typename F, unsigned int Amat, unsigned int A>
00407 struct VECTOR_SIZE_DEDUCER<const viennacl::matrix_expression< const viennacl::matrix<ScalarType, F, Amat>,
00408 const viennacl::matrix<ScalarType, F, Amat>,
00409 op_trans>,
00410 const viennacl::vector<ScalarType, A>,
00411 viennacl::op_prod>
00412 {
00413 static size_t size(const viennacl::matrix_expression< const viennacl::matrix<ScalarType, F, Amat>,
00414 const viennacl::matrix<ScalarType, F, Amat>,
00415 op_trans> & lhs,
00416 const viennacl::vector<ScalarType, A> & rhs) { return lhs.lhs().size2(); }
00417 };
00418
00419
00420
00421
00422
00427 template <typename T>
00428 struct CPU_SCALAR_TYPE_DEDUCER
00429 {
00430
00431
00432 };
00433
00434 template <>
00435 struct CPU_SCALAR_TYPE_DEDUCER< float >
00436 {
00437 typedef float ResultType;
00438 };
00439
00440 template <>
00441 struct CPU_SCALAR_TYPE_DEDUCER< double >
00442 {
00443 typedef double ResultType;
00444 };
00445
00446 template <typename T>
00447 struct CPU_SCALAR_TYPE_DEDUCER< viennacl::scalar<T> >
00448 {
00449 typedef T ResultType;
00450 };
00451
00452 template <typename T, unsigned int A>
00453 struct CPU_SCALAR_TYPE_DEDUCER< viennacl::vector<T, A> >
00454 {
00455 typedef T ResultType;
00456 };
00457
00458 template <typename T, typename F, unsigned int A>
00459 struct CPU_SCALAR_TYPE_DEDUCER< viennacl::matrix<T, F, A> >
00460 {
00461 typedef T ResultType;
00462 };
00463
00464
00465 template <typename T, typename F, unsigned int A>
00466 struct CPU_SCALAR_TYPE_DEDUCER< viennacl::matrix_expression<const matrix<T, F, A>, const matrix<T, F, A>, op_trans> >
00467 {
00468 typedef T ResultType;
00469 };
00470
00471
00472 }
00473 }
00474
00475
00476 #endif