Commit cc0bf8f5 authored by Pavel Kus's avatar Pavel Kus Committed by Andreas Marek

templating test_scalapack

adding variants for single precision and complex math datatype
for the scalapack test
parent 347e6a87
......@@ -670,6 +670,11 @@ EXTRA_DIST += \
src/elpa1/legacy_interface/elpa_solve_tridi.F90
endif
if WITH_SCALAPACK_TESTS
EXTRA_DIST += \
test/shared/test_scalapack_template.F90
endif
LIBTOOL_DEPS = @LIBTOOL_DEPS@
libtool: $(LIBTOOL_DEPS)
$(SHELL) ./config.status libtool
......
......@@ -50,7 +50,7 @@ for m, g, t, p, d, s, l in product(
if(m == "analytic" and (g == 1 or t != "eigenvectors")):
continue
if(s == "scalapack_all" and (g == 1 or t != "eigenvectors" or p == "single" or d == "complex" or m != "analytic")):
if(s == "scalapack_all" and (g == 1 or t != "eigenvectors" or m != "analytic")):
continue
if (t == "solve_tridiagonal" and (s == "2stage" or d == "complex")):
......
......@@ -3,6 +3,8 @@
#undef SINGLE_PRECSION_REAL
#undef MATH_DATATYPE
#undef BLAS_TRANS_OR_CONJ
#undef BLAS_CHAR
#undef BLAS_CHAR_AND_SY_OR_HE
#undef PRECISION
#undef SPECIAL_COMPLEX_DATATYPE
#undef PRECISION_STR
......@@ -65,6 +67,8 @@
#define PRECISION_STR 'double'
#define PRECISION_SUFFIX "_double"
#define REAL_DATATYPE rk8
#define BLAS_CHAR D
#define BLAS_CHAR_AND_SY_OR_HE DSY
#define SPECIAL_COMPLEX_DATATYPE ck8
#define PRECISION_TRTRI DTRTRI
......@@ -121,6 +125,8 @@
#define PRECISION_STR 'single'
#define PRECISION_SUFFIX "_single"
#define REAL_DATATYPE rk4
#define BLAS_CHAR S
#define BLAS_CHAR_AND_SY_OR_HE SSY
#define SPECIAL_COMPLEX_DATATYPE ck4
#define PRECISION_TRTRI STRTRI
......@@ -177,6 +183,8 @@
#undef SINGLE_PRECISION_COMPLEX
#undef MATH_DATATYPE
#undef BLAS_TRANS_OR_CONJ
#undef BLAS_CHAR
#undef BLAS_CHAR_AND_SY_OR_HE
#undef PRECISION
#undef COMPLEX_DATATYPE
/* in the complex case also sometime real valued variables are needed */
......@@ -241,7 +249,6 @@
/* General definitions needed in single and double case */
#define MATH_DATATYPE complex
#define BLAS_TRANS_OR_CONJ 'C'
#ifdef DOUBLE_PRECISION
#define DOUBLE_PRECISION_COMPLEX
......@@ -249,6 +256,8 @@
#define PRECISION_STR 'double'
#define PRECISION_SUFFIX "_double"
#define COMPLEX_DATATYPE CK8
#define BLAS_CHAR Z
#define BLAS_CHAR_AND_SY_OR_HE ZHE
#define REAL_DATATYPE RK8
#define PRECISION_TRTRI ZTRTRI
......@@ -311,6 +320,8 @@
#define PRECISION_STR 'single'
#define PRECISION_SUFFIX "_single"
#define COMPLEX_DATATYPE CK4
#define BLAS_CHAR C
#define BLAS_CHAR_AND_SY_OR_HE CHE
#define REAL_DATATYPE RK4
#define PRECISION_TRTRI CTRTRI
......
......@@ -48,37 +48,52 @@ module test_scalapack
interface solve_scalapack_all
module procedure solve_pdsyevd
end interface
module procedure solve_pzheevd
#ifdef WANT_SINGLE_PRECISION_REAL
module procedure solve_pssyevd
#endif
#ifdef WANT_SINGLE_PRECISION_COMPLEX
module procedure solve_pcheevd
#endif
end interface
contains
subroutine solve_pdsyevd(na, a, sc_desc, ev, z)
implicit none
integer(kind=ik), intent(in) :: na
real(kind=rk8), intent(in) :: a(:,:)
real(kind=rk8), intent(inout) :: z(:,:), ev(:)
integer(kind=ik), intent(in) :: sc_desc(:)
integer(kind=ik) :: info, lwork, liwork
real(kind=rk8), allocatable :: work(:)
integer, allocatable :: iwork(:)
#define COMPLEXCASE 1
#define DOUBLE_PRECISION 1
#include "../../src/general/precision_macros.h"
#include "test_scalapack_template.F90"
#undef DOUBLE_PRECISION
#undef COMPLEXCASE
allocate(work(1), iwork(1))
#ifdef WANT_SINGLE_PRECISION_COMPLEX
! query for required workspace
call pdsyevd('V', 'U', na, a, 1, 1, sc_desc, ev, z, 1, 1, sc_desc, work, -1, iwork, -1, info)
! write(*,*) "computed sizes", lwork, liwork, "required sizes ", work(1), iwork(1)
lwork = work(1)
liwork = iwork(1)
#define COMPLEXCASE 1
#define SINGLE_PRECISION 1
#include "../../src/general/precision_macros.h"
#include "test_scalapack_template.F90"
#undef SINGLE_PRECISION
#undef COMPLEXCASE
deallocate(work, iwork)
allocate(work(lwork), stat = info)
allocate(iwork(liwork), stat = info)
#endif /* WANT_SINGLE_PRECISION_COMPLEX */
! the actuall call to the method
call pdsyevd('V', 'U', na, a, 1, 1, sc_desc, ev, z, 1, 1, sc_desc, work, lwork, iwork, liwork, info)
#define REALCASE 1
#define DOUBLE_PRECISION 1
#include "../../src/general/precision_macros.h"
#include "test_scalapack_template.F90"
#undef DOUBLE_PRECISION
#undef REALCASE
#ifdef WANT_SINGLE_PRECISION_REAL
#define REALCASE 1
#define SINGLE_PRECISION 1
#include "../../src/general/precision_macros.h"
#include "test_scalapack_template.F90"
#undef SINGLE_PRECISION
#undef REALCASE
#endif /* WANT_SINGLE_PRECISION_REAL */
deallocate(iwork)
deallocate(work)
end subroutine
end module
! (c) Copyright Pavel Kus, 2017, MPCDF
!
! This file is part of ELPA.
!
! The ELPA library was originally created by the ELPA consortium,
! consisting of the following organizations:
!
! - Max Planck Computing and Data Facility (MPCDF), formerly known as
! Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
! - Bergische Universität Wuppertal, Lehrstuhl für angewandte
! Informatik,
! - Technische Universität München, Lehrstuhl für Informatik mit
! Schwerpunkt Wissenschaftliches Rechnen ,
! - Fritz-Haber-Institut, Berlin, Abt. Theorie,
! - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
! Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
! and
! - IBM Deutschland GmbH
!
!
! More information can be found here:
! http://elpa.mpcdf.mpg.de/
!
! ELPA is free software: you can redistribute it and/or modify
! it under the terms of the version 3 of the license of the
! GNU Lesser General Public License as published by the Free
! Software Foundation.
!
! ELPA 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 Lesser General Public License for more details.
!
! You should have received a copy of the GNU Lesser General Public License
! along with ELPA. If not, see <http://www.gnu.org/licenses/>
!
! ELPA reflects a substantial effort on the part of the original
! ELPA consortium, and we ask you to respect the spirit of the
! license that we chose: i.e., please contribute any changes you
! may have back to the original ELPA library distribution, and keep
! any derivatives of ELPA under the same license that we chose for
! the original distribution, the GNU Lesser General Public License.
subroutine solve_p&
&BLAS_CHAR_AND_SY_OR_HE&
&evd(na, a, sc_desc, ev, z)
implicit none
#include "../../src/general/precision_kinds.F90"
integer(kind=ik), intent(in) :: na
MATH_DATATYPE(kind=rck), intent(in) :: a(:,:)
MATH_DATATYPE(kind=rck), intent(inout) :: z(:,:)
real(kind=rk), intent(inout) :: ev(:)
integer(kind=ik), intent(in) :: sc_desc(:)
integer(kind=ik) :: info, lwork, liwork, lrwork
MATH_DATATYPE(kind=rck), allocatable :: work(:)
real(kind=rk), allocatable :: rwork(:)
integer, allocatable :: iwork(:)
allocate(work(1), iwork(1), rwork(1))
! query for required workspace
#ifdef REALCASE
call p&
&BLAS_CHAR&
&syevd('V', 'U', na, a, 1, 1, sc_desc, ev, z, 1, 1, sc_desc, work, -1, iwork, -1, info)
#endif
#ifdef COMPLEXCASE
call p&
&BLAS_CHAR&
&heevd('V', 'U', na, a, 1, 1, sc_desc, ev, z, 1, 1, sc_desc, work, -1, rwork, -1, iwork, -1, info)
#endif
! write(*,*) "computed sizes", lwork, liwork, "required sizes ", work(1), iwork(1)
lwork = work(1)
liwork = iwork(1)
deallocate(work, iwork)
allocate(work(lwork), stat = info)
allocate(iwork(liwork), stat = info)
#ifdef COMPLEXCASE
lrwork = rwork(1)
deallocate(rwork)
allocate(rwork(lrwork), stat = info)
#endif
! the actuall call to the method
#ifdef REALCASE
call p&
&BLAS_CHAR&
&syevd('V', 'U', na, a, 1, 1, sc_desc, ev, z, 1, 1, sc_desc, work, lwork, iwork, liwork, info)
#endif
#ifdef COMPLEXCASE
call p&
&BLAS_CHAR&
&heevd('V', 'U', na, a, 1, 1, sc_desc, ev, z, 1, 1, sc_desc, work, lwork, rwork, lrwork, iwork, liwork, info)
#endif
deallocate(iwork, work, rwork)
end subroutine
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment