diff --git a/Makefile.am b/Makefile.am index ac14263308380199b5be545ca52a8516554efe8f..2752f4f55ad95564ae6a111722c84405f418dad8 100644 --- a/Makefile.am +++ b/Makefile.am @@ -462,6 +462,7 @@ libelpatest@SUFFIX@_la_SOURCES = \ test/shared/test_blacs_infrastructure.F90 \ test/shared/test_prepare_matrix.F90 \ test/shared/test_analytic.F90 \ + test/shared/test_scalapack.F90 \ test/shared/test_output_type.F90 if HAVE_REDIRECT diff --git a/generate_automake_test_programs.py b/generate_automake_test_programs.py index 3731ef0664d168305da8b7fcda300b0ad40d5b4f..e7826d9e014ca04484c5307386c4cd189121bf34 100755 --- a/generate_automake_test_programs.py +++ b/generate_automake_test_programs.py @@ -13,6 +13,7 @@ prec_flag = { solver_flag = { "1stage" : "-DTEST_SOLVER_1STAGE", "2stage" : "-DTEST_SOLVER_2STAGE", + "scalapack_all" : "-DTEST_SCALAPACK_ALL", } gpu_flag = { 0 : "-DTEST_GPU=0", @@ -47,6 +48,9 @@ for m, g, t, p, d, s, l in product( if(m == "analytic" and (g == 1 or t != "eigenvectors" or p == "single" or d == "complex")): continue + if(s == "scalapack_all" and (g == 1 or t != "eigenvectors" or p == "single" or d == "complex" or m != "analytic")): + continue + if (t == "solve_tridiagonal" and (s == "2stage" or d == "complex")): continue diff --git a/test/Fortran/test.F90 b/test/Fortran/test.F90 index f3d7a6c85ebbeec47aec1756682f5a6fff307072..48dd7aa4f4eced1371539619ab643939bd3e430b 100644 --- a/test/Fortran/test.F90 +++ b/test/Fortran/test.F90 @@ -48,16 +48,16 @@ ! Define TEST_GPU \in [0, 1] ! Define either TEST_ALL_KERNELS or a TEST_KERNEL \in [any valid kernel] -#if !(defined(TEST_REAL) ^ defined(TEST_COMPLEX)) -error: define exactly one of TEST_REAL or TEST_COMPLEX +#if !(defined(TEST_REAL) ^ defined(TEST_COMPLEX) ) +error: define exactly one of TEST_REAL or TEST_COMPLEX #endif #if !(defined(TEST_SINGLE) ^ defined(TEST_DOUBLE)) error: define exactly one of TEST_SINGLE or TEST_DOUBLE #endif -#if !(defined(TEST_SOLVER_1STAGE) ^ defined(TEST_SOLVER_2STAGE)) -error: define exactly one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE +#if !(defined(TEST_SOLVER_1STAGE) ^ defined(TEST_SOLVER_2STAGE) ^ defined(TEST_SCALAPACK_ALL)) +error: define exactly one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE or TEST_SCALAPACK_ALL #endif #ifdef TEST_SOLVER_1STAGE @@ -111,6 +111,7 @@ program test use test_blacs_infrastructure use test_check_correctness use test_analytic + use test_scalapack #ifdef HAVE_REDIRECT use test_redirect @@ -331,7 +332,11 @@ program test ! The actual solve step #ifdef TEST_EIGENVECTORS call e%timer_start("e%eigenvectors()") +#ifdef TEST_SCALAPACK_ALL + call solve_scalapack_all(na, a, sc_desc, ev, z) +#else call e%eigenvectors(a, ev, z, error) +#endif call e%timer_stop("e%eigenvectors()") #endif @@ -348,6 +353,7 @@ program test ev(:) = d(:) #endif + assert_elpa_ok(error) #ifdef TEST_ALL_KERNELS diff --git a/test/shared/test_scalapack.F90 b/test/shared/test_scalapack.F90 new file mode 100644 index 0000000000000000000000000000000000000000..57110bb8763f9ebf518cefabb5b99a56d6d79297 --- /dev/null +++ b/test/shared/test_scalapack.F90 @@ -0,0 +1,84 @@ +! (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 +! +! 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. + +#include "config-f90.h" + +module test_scalapack + use test_util + + interface solve_scalapack_all + module procedure solve_pdsyevd + 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(:) + + allocate(work(1), iwork(1)) + + ! 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) + + deallocate(work, iwork) + allocate(work(lwork), stat = info) + allocate(iwork(liwork), stat = info) + + ! 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) + + deallocate(iwork) + deallocate(work) + end subroutine + +end module