diff options
author | 2012-07-29 18:53:26 +0200 | |
---|---|---|
committer | 2012-07-29 18:53:26 +0200 | |
commit | 637d5057499669550edd1a3e14bf49664808198f (patch) | |
tree | 591638b56138b6274aec69005b26c8cdf5774968 | |
parent | Inserted patched BTL into numbench (diff) | |
download | auto-numerical-bench-637d5057499669550edd1a3e14bf49664808198f.tar.gz auto-numerical-bench-637d5057499669550edd1a3e14bf49664808198f.tar.bz2 auto-numerical-bench-637d5057499669550edd1a3e14bf49664808198f.zip |
New BTL module lapacke.
-rw-r--r-- | btl/libs/LAPACK/lapack_.hh | 2 | ||||
-rw-r--r-- | btl/libs/LAPACK/lapack_interface.hh | 33 | ||||
-rw-r--r-- | btl/libs/LAPACK/lapack_interface_impl.hh | 4 | ||||
-rw-r--r-- | btl/libs/LAPACK/lapacke_interface_impl.hh | 101 | ||||
-rw-r--r-- | btl/libs/LAPACK/main.cpp | 86 |
5 files changed, 160 insertions, 66 deletions
diff --git a/btl/libs/LAPACK/lapack_.hh b/btl/libs/LAPACK/lapack_.hh index e7c8c3f..bc56dda 100644 --- a/btl/libs/LAPACK/lapack_.hh +++ b/btl/libs/LAPACK/lapack_.hh @@ -32,6 +32,8 @@ void LAPACKFUNC(gesvd)(const char*, const char*, const int*, const int*, SCALAR* void LAPACKFUNC(syev)(const char*, const char*, const int*, SCALAR*, const int*, SCALAR*, SCALAR*, const int*, int*); void LAPACKFUNC(stev)(const char*, const int*, SCALAR*, SCALAR*, SCALAR*, const int*, SCALAR*, int*); +void LAPACKFUNC(gels)(char*, int*, int*, int*, SCALAR*, int*, SCALAR*, int*, SCALAR*, int*, int*); + #ifdef __cplusplus } #endif diff --git a/btl/libs/LAPACK/lapack_interface.hh b/btl/libs/LAPACK/lapack_interface.hh index a0ec52b..87b3cbb 100644 --- a/btl/libs/LAPACK/lapack_interface.hh +++ b/btl/libs/LAPACK/lapack_interface.hh @@ -20,26 +20,7 @@ #include <../BLAS/c_interface_base.h> #include <complex> -#include "lapack.hh" - -extern "C" { -#include "../BLAS/blas.h" - -//void sgesv_(int*, int*, float *, int*, int*, float *, int*, int*); -//void dgesv_(int*, int*, double*, int*, int*, double*, int*, int*); - -void sgels_(char*, int*, int*, int*, float *, int*, float *, int*, float *, int*, int*); -void dgels_(char*, int*, int*, int*, double*, int*, double*, int*, double*, int*, int*); - -//void sgetrf_(int*, int*, float *, int*, int*, int*); -//void dgetrf_(int*, int*, double*, int*, int*, int*); - -//void spotrf_(char*, int*, float *, int*, int*); -//void dpotrf_(char*, int*, double*, int*, int*); - -//void ssyev_(char*, char*, int*, float *, int*, float *, float *, int*, int*); -//void dsyev_(char*, char*, int*, double*, int*, double*, double*, int*, int*); -} +//#include "lapack.hh" #define MAKE_STRING2(S) #S @@ -65,13 +46,21 @@ static int zeroint = 0; #define SCALAR float #define SCALAR_PREFIX s -#include "lapack_interface_impl.hh" +#ifdef LAPACKE_INTERFACE +# include "lapacke_interface_impl.hh" +#else +# include "lapack_interface_impl.hh" +#endif #undef SCALAR #undef SCALAR_PREFIX #define SCALAR double #define SCALAR_PREFIX d -#include "lapack_interface_impl.hh" +#ifdef LAPACKE_INTERFACE +# include "lapacke_interface_impl.hh" +#else +# include "lapack_interface_impl.hh" +#endif #undef SCALAR #undef SCALAR_PREFIX diff --git a/btl/libs/LAPACK/lapack_interface_impl.hh b/btl/libs/LAPACK/lapack_interface_impl.hh index a95bc5b..523e53a 100644 --- a/btl/libs/LAPACK/lapack_interface_impl.hh +++ b/btl/libs/LAPACK/lapack_interface_impl.hh @@ -17,6 +17,9 @@ // #define LAPACKFUNC(NAME) CAT(CAT(SCALAR_PREFIX,NAME),_) +#include "lapack_.hh" +#include "blas.h" + template<> class lapack_interface<SCALAR> : public c_interface_base<SCALAR> { public: @@ -114,7 +117,6 @@ public: static inline void stev(const gene_vector& D, const gene_vector& E, gene_vector& W, gene_matrix& V, const int& N) { int N0 = N; - int N1 = N-1; LAPACKFUNC(copy)(&N0, D, &intone, W, &intone); stl_vector E_(E, E+N-1), work(max(1, 2*N-2)); diff --git a/btl/libs/LAPACK/lapacke_interface_impl.hh b/btl/libs/LAPACK/lapacke_interface_impl.hh new file mode 100644 index 0000000..eef3638 --- /dev/null +++ b/btl/libs/LAPACK/lapacke_interface_impl.hh @@ -0,0 +1,101 @@ +//===================================================== +// Copyright (C) 2012 Andrea Arteaga <andyspiros@gmail.com> +//===================================================== +// +// This program is free software; you can redistribute it and/or +// modify it under the terms of the GNU General Public License +// as published by the Free Software Foundation; either version 2 +// of the License, or (at your option) any later version. +// +// This program is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// You should have received a copy of the GNU General Public License +// along with this program; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. +// +#define LAPACKEFUNC(NAME) CAT(LAPACKE_,CAT(SCALAR_PREFIX,NAME)) +#define BLASFUNC(NAME) CAT(CAT(SCALAR_PREFIX,NAME),_) + +#include "lapacke.h" + + +// Define BLAS functions used by LAPACKE interface +extern "C" { + int BLASFUNC(copy) (int *, SCALAR*, int *, SCALAR*, int *); +} + +template<> class lapack_interface<SCALAR> : public c_interface_base<SCALAR> +{ +public: + + static inline std::string name() + { + return MAKE_STRING(LAPACKNAME); + } + + static inline void general_solve(gene_matrix& A, gene_vector& b, gene_vector& x, int N) + { + std::vector<int> ipiv(N); + BLASFUNC(copy)(&N, b, &intone, x, &intone); + LAPACKEFUNC(gesv)(LAPACK_COL_MAJOR, N, 1, A, N, &ipiv[0], x, N); + } + + static inline void least_squares(gene_matrix& A, gene_vector& b, gene_vector& x, int N) + { + BLASFUNC(copy)(&N, b, &intone, x, &intone); + LAPACKEFUNC(gels)(LAPACK_COL_MAJOR, 'N', N, N, 1, A, N, x, N); + } + + static inline void lu_decomp(const gene_matrix& X, gene_matrix& C, int N) + { + int N2 = N*N; + std::vector<int> ipiv(N); + BLASFUNC(copy)(&N2, X, &intone, C, &intone); + LAPACKEFUNC(getrf)(LAPACK_COL_MAJOR, N, N, C, N, &ipiv[0]); + } + + static inline void cholesky(const gene_matrix& X, gene_matrix& C, int N) + { + int N2 = N*N; + BLASFUNC(copy)(&N2, X, &intone, C, &intone); + LAPACKEFUNC(potrf)(LAPACK_COL_MAJOR, 'L', N, C, N); + } + + static inline void qr_decomp(const gene_matrix& X, gene_matrix& QR, gene_vector& tau, const int& N) + { + int N2 = N*N; + BLASFUNC(copy)(&N2, X, &intone, QR, &intone); + LAPACKEFUNC(geqrf)(LAPACK_COL_MAJOR, N, N, QR, N, tau); + } + + static inline void svd_decomp(const gene_matrix& X, gene_matrix& U, gene_vector& S, gene_matrix& VT, const int& N) + { + int N2 = N*N; + stl_vector Xcopy(N2), superb(N-1); + BLASFUNC(copy)(&N2, X, &intone, &Xcopy[0], &intone); + LAPACKEFUNC(gesvd)(LAPACK_COL_MAJOR, 'A', 'A', N, N, &Xcopy[0], N, S, U, N, VT, N, &superb[0]); + } + + static inline void syev(const gene_matrix& X, gene_matrix& V, gene_vector& W, const int& N) + { + int N2 = N*N; + BLASFUNC(copy)(&N2, X, &intone, V, &intone); + LAPACKEFUNC(syev)(LAPACK_COL_MAJOR, 'V', 'U', N, V, N, W); + } + + /* Size of D, W: N; size of E: N-1, size of V: NxN */ + static inline void stev(const gene_vector& D, const gene_vector& E, gene_vector& W, gene_matrix& V, int N) + { + stl_vector E_(E, E+N-1); + BLASFUNC(copy)(&N, D, &intone, W, &intone); + LAPACKEFUNC(stev)(LAPACK_COL_MAJOR, 'V', N, W, &E_[0], V, N); + } + + static inline void symm_ev(const gene_matrix& X, gene_vector& W, int N) + { + LAPACKEFUNC(syev)(LAPACK_COL_MAJOR, 'N', 'L', N, X, N, W); + } + +}; diff --git a/btl/libs/LAPACK/main.cpp b/btl/libs/LAPACK/main.cpp index e70c161..78f314a 100644 --- a/btl/libs/LAPACK/main.cpp +++ b/btl/libs/LAPACK/main.cpp @@ -35,61 +35,61 @@ BTL_MAIN; int main(int argc, char **argv) { - bool - general_solve=false, least_squares=false, lu_decomp=false, cholesky=false, qr_decomp=false, svd_decomp=false, - syev=false, stev=false, symm_ev=false - ; - int N = 100; + bool + general_solve=false, least_squares=false, lu_decomp=false, cholesky=false, + qr_decomp=false, svd_decomp=false, syev=false, stev=false, + symm_ev=false; + int N = 100; - for (int i = 1; i < argc; ++i) { - std::string arg = argv[i]; - if (arg == "general_solve") general_solve = true; - else if (arg == "least_squares") least_squares = true; - else if (arg == "lu_decomp") lu_decomp = true; - else if (arg == "cholesky") cholesky = true; - else if (arg == "qr_decomp") qr_decomp = true; - else if (arg == "svd_decomp") svd_decomp = true; - else if (arg == "syev") syev = true; - else if (arg == "stev") stev = true; - else if (arg == "symm_ev") symm_ev = true; + for (int i = 1; i < argc; ++i) { + std::string arg = argv[i]; + if (arg == "general_solve") general_solve = true; + else if (arg == "least_squares") least_squares = true; + else if (arg == "lu_decomp") lu_decomp = true; + else if (arg == "cholesky") cholesky = true; + else if (arg == "qr_decomp") qr_decomp = true; + else if (arg == "svd_decomp") svd_decomp = true; + else if (arg == "syev") syev = true; + else if (arg == "stev") stev = true; + else if (arg == "symm_ev") symm_ev = true; - // Check switch -N - else if (arg[0] == '-' && arg[1] == 'N') { - if (arg[2] != '\0') - N = atoi(arg.c_str()+2); - else - N = atoi(argv[++i]); - } - } + // Check switch -N + else if (arg[0] == '-' && arg[1] == 'N') { + if (arg[2] != '\0') + N = atoi(arg.c_str()+2); + else + N = atoi(argv[++i]); + } + } - if (general_solve) - bench<Action_general_solve<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); + if (general_solve) + bench<Action_general_solve<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); - if (least_squares) - bench<Action_least_squares<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); + if (least_squares) + bench<Action_least_squares<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); - if (lu_decomp) - bench<Action_lu_decomp<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); + if (lu_decomp) + bench<Action_lu_decomp<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); - if (cholesky) - bench<Action_cholesky<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); + if (cholesky) + bench<Action_cholesky<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); - if (qr_decomp) - bench<Action_qr_decomp<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); + if (qr_decomp) + bench<Action_qr_decomp<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); - if (svd_decomp) - bench<Action_svd_decomp<lapack_interface<REAL_TYPE> > >(MIN_MM,750,N); + if (svd_decomp) + bench<Action_svd_decomp<lapack_interface<REAL_TYPE> > >(MIN_MM,750,N); - if (syev) - bench<Action_syev<lapack_interface<REAL_TYPE> > >(MIN_MM,750,N); + if (syev) + bench<Action_syev<lapack_interface<REAL_TYPE> > >(MIN_MM,750,N); - if (stev) - bench<Action_stev<lapack_interface<REAL_TYPE> > >(MIN_MM,1000,N); + if (stev) + bench<Action_stev<lapack_interface<REAL_TYPE> > >(MIN_MM,1000,N); - if (symm_ev) - bench<Action_symm_ev<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); + if (symm_ev) + bench<Action_symm_ev<lapack_interface<REAL_TYPE> > >(MIN_MM,MAX_MM,N); - return 0; + return 0; } |