diff --git a/Makefile.am b/Makefile.am index 0cacb141573495a60b933c9ebc129daa63c0e028..41bfe0c336a3dcb445645889a1a446edba45a7b2 100644 --- a/Makefile.am +++ b/Makefile.am @@ -382,7 +382,8 @@ dist_files_DATA = \ test/Fortran/test_invert_trm_real.F90 \ test/Fortran/test_cholesky_complex.F90 \ test/Fortran/test_invert_trm_complex.F90 \ - test/Fortran/test_new_interface.F90 \ + test/Fortran/test_new_interface_real.F90 \ + test/Fortran/test_new_interface_complex.F90 \ test/Fortran/elpa_tests.F90 \ src/elpa2/elpa2_print_kernels.F90 @@ -426,7 +427,8 @@ noinst_PROGRAMS = \ elpa2_test_real_c_version@SUFFIX@ \ elpa2_test_complex_c_version@SUFFIX@ \ elpa_driver_real_c_version@SUFFIX@ \ - elpa_test_new_interface@SUFFIX@ \ + elpa_test_new_interface_real@SUFFIX@ \ + elpa_test_new_interface_complex@SUFFIX@ \ elpa_driver_complex_c_version@SUFFIX@ if WANT_SINGLE_PRECISION_COMPLEX @@ -439,7 +441,8 @@ noinst_PROGRAMS += \ elpa1_complex_cholesky_single_precision@SUFFIX@ \ elpa1_complex_invert_trm_single_precision@SUFFIX@ \ elpa2_test_complex_api_single_precision@SUFFIX@ \ - elpa_driver_complex_c_version_single_precision@SUFFIX@ + elpa_driver_complex_c_version_single_precision@SUFFIX@ \ + elpa_test_new_interface_complex_single@SUFFIX@ endif if WANT_SINGLE_PRECISION_REAL @@ -454,7 +457,8 @@ noinst_PROGRAMS += \ elpa1_real_cholesky_single_precision@SUFFIX@ \ elpa1_real_invert_trm_single_precision@SUFFIX@ \ elpa_driver_real_c_version_single_precision@SUFFIX@ \ - elpa1_real_toeplitz_single_precision@SUFFIX@ + elpa1_real_toeplitz_single_precision@SUFFIX@ \ + elpa_test_new_interface_real_single@SUFFIX@ endif if WITH_GPU_VERSION @@ -501,11 +505,15 @@ libelpatest@SUFFIX@_la_SOURCES += \ test/shared/redir.c \ test/shared/redirect.F90 endif -elpa_test_new_interface@SUFFIX@_SOURCES = test/Fortran/test_new_interface.F90 -elpa_test_new_interface@SUFFIX@_LDADD = $(build_lib) $(FCLIBS) -elpa_test_new_interface@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules -EXTRA_elpa_test_new_interface@SUFFIX@_DEPENDENCIES = test/Fortran/elpa_print_headers.X90 +elpa_test_new_interface_real@SUFFIX@_SOURCES = test/Fortran/test_new_interface_real.F90 +elpa_test_new_interface_real@SUFFIX@_LDADD = $(build_lib) $(FCLIBS) +elpa_test_new_interface_real@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules +EXTRA_elpa_test_new_interface_real@SUFFIX@_DEPENDENCIES = test/Fortran/elpa_print_headers.X90 +elpa_test_new_interface_complex@SUFFIX@_SOURCES = test/Fortran/test_new_interface_complex.F90 +elpa_test_new_interface_complex@SUFFIX@_LDADD = $(build_lib) $(FCLIBS) +elpa_test_new_interface_complex@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules +EXTRA_elpa_test_new_interface_complex@SUFFIX@_DEPENDENCIES = test/Fortran/elpa_print_headers.X90 elpa1_test_real_c_version@SUFFIX@_SOURCES = test/C/elpa1_test_real_c_version.c elpa1_test_real_c_version@SUFFIX@_LDADD = $(build_lib) $(FCLIBS) @@ -653,6 +661,11 @@ elpa_tests@SUFFIX@_LDADD = $(build_lib) elpa_tests@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules if WANT_SINGLE_PRECISION_REAL +elpa_test_new_interface_real_single@SUFFIX@_SOURCES = test/Fortran/test_new_interface_real_single.F90 +elpa_test_new_interface_real_single@SUFFIX@_LDADD = $(build_lib) $(FCLIBS) +elpa_test_new_interface_real_single@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules +EXTRA_elpa_test_new_interface_real_single@SUFFIX@_DEPENDENCIES = test/Fortran/elpa_print_headers.X90 + elpa1_test_real_single_precision@SUFFIX@_SOURCES = test/Fortran/test_real_single.F90 elpa1_test_real_single_precision@SUFFIX@_LDADD = $(build_lib) elpa1_test_real_single_precision@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules @@ -713,6 +726,11 @@ EXTRA_elpa2_test_real_api_single_precision@SUFFIX@_DEPENDENCIES = test/Fortran/e endif if WANT_SINGLE_PRECISION_COMPLEX +elpa_test_new_interface_complex_single@SUFFIX@_SOURCES = test/Fortran/test_new_interface_complex_single.F90 +elpa_test_new_interface_complex_single@SUFFIX@_LDADD = $(build_lib) $(FCLIBS) +elpa_test_new_interface_complex_single@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules +EXTRA_elpa_test_new_interface_complex_single@SUFFIX@_DEPENDENCIES = test/Fortran/elpa_print_headers.X90 + elpa1_test_complex_single_precision@SUFFIX@_SOURCES = test/Fortran/test_complex_single.F90 elpa1_test_complex_single_precision@SUFFIX@_LDADD = $(build_lib) elpa1_test_complex_single_precision@SUFFIX@_FCFLAGS = $(AM_FCFLAGS) @FC_MODOUT@private_modules @FC_MODINC@private_modules @@ -836,7 +854,8 @@ check_SCRIPTS = \ elpa2_test_real_c_version@SUFFIX@.sh \ elpa2_test_complex_c_version@SUFFIX@.sh \ elpa_driver_real_c_version@SUFFIX@.sh \ - elpa_test_new_interface@SUFFIX@.sh \ + elpa_test_new_interface_real@SUFFIX@.sh \ + elpa_test_new_interface_complex@SUFFIX@.sh \ elpa_driver_complex_c_version@SUFFIX@.sh if WANT_SINGLE_PRECISION_REAL @@ -850,7 +869,8 @@ check_SCRIPTS += \ elpa1_real_invert_trm_single_precision@SUFFIX@.sh \ elpa1_real_transpose_multiply_single_precision@SUFFIX@.sh \ elpa2_test_real_api_single_precision@SUFFIX@.sh \ - elpa_driver_real_c_version_single_precision@SUFFIX@.sh + elpa_driver_real_c_version_single_precision@SUFFIX@.sh \ + elpa_test_new_interface_real_single@SUFFIX@.sh endif if WANT_SINGLE_PRECISION_COMPLEX @@ -863,7 +883,8 @@ check_SCRIPTS += \ elpa1_complex_invert_trm_single_precision@SUFFIX@.sh \ elpa1_complex_transpose_multiply_single_precision@SUFFIX@.sh \ elpa2_test_complex_api_single_precision@SUFFIX@.sh \ - elpa_driver_complex_c_version_single_precision@SUFFIX@.sh + elpa_driver_complex_c_version_single_precision@SUFFIX@.sh \ + elpa_test_new_interface_complex_single@SUFFIX@.sh endif if WITH_GPU_VERSION diff --git a/src/elpa_t.F90 b/src/elpa_t.F90 index fa9db9743fe9af37200d8f7ee708b4b8761794dc..ee9c0a8e6a7e36c24f36ebc2bab5856827d7e4fe 100644 --- a/src/elpa_t.F90 +++ b/src/elpa_t.F90 @@ -64,7 +64,10 @@ module elpa_type procedure, public :: get_communicators => get_communicators generic, public :: solve => elpa_solve_real_double, & - elpa_solve_complex_double + elpa_solve_real_single, & + elpa_solve_complex_double, & + elpa_solve_complex_single + procedure, public :: destroy => elpa_destroy @@ -72,7 +75,9 @@ module elpa_type procedure, private :: elpa_set_integer procedure, private :: elpa_get_integer procedure, private :: elpa_solve_real_double + procedure, private :: elpa_solve_real_single procedure, private :: elpa_solve_complex_double + procedure, private :: elpa_solve_complex_single end type elpa_t @@ -274,6 +279,40 @@ module elpa_type end subroutine + subroutine elpa_solve_real_single(self, a, ev, q, success) + use elpa + use elpa_utilities, only : error_unit + + use iso_c_binding + implicit none + class(elpa_t) :: self + + real(kind=c_float) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols), & + ev(self%na) + integer, optional :: success + logical :: success_l + +#ifdef WANT_SINGLE_PRECISION_REAL + success_l = elpa_solve_evp_real_single(self%na, self%nev, a, self%local_nrows, ev, q, & + self%local_nrows, self%nblk, self%local_ncols, & + self%mpi_comm_rows, self%mpi_comm_cols, & + self%mpi_comm_parent) + + if (present(success)) then + if (success_l) then + success = ELPA_OK + else + success = ELPA_ERROR + endif + else if (.not. success_l) then + write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!" + endif +#else + success = ELPA_ERROR +#endif + + end subroutine + subroutine elpa_solve_complex_double(self, a, ev, q, success) use elpa @@ -307,6 +346,43 @@ module elpa_type end subroutine + subroutine elpa_solve_complex_single(self, a, ev, q, success) + use elpa + use elpa_utilities, only : error_unit + + use iso_c_binding + implicit none + class(elpa_t) :: self + + complex(kind=c_float_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols) + real(kind=c_float) :: ev(self%na) + + integer, optional :: success + logical :: success_l + +#ifdef WANT_SINGLE_PRECISION_COMPLEX + success_l = elpa_solve_evp_complex_single(self%na, self%nev, a, self%local_nrows, ev, q, & + self%local_nrows, self%nblk, self%local_ncols, & + self%mpi_comm_rows, self%mpi_comm_cols, & + self%mpi_comm_parent) + + if (present(success)) then + if (success_l) then + success = ELPA_OK + else + success = ELPA_ERROR + endif + else if (.not. success_l) then + write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!" + endif +#else + success = ELPA_ERROR +#endif + + + end subroutine + + subroutine elpa_destroy(self) class(elpa_t) :: self call elpa_free_options(self%options) diff --git a/test/Fortran/test_new_interface_complex.F90 b/test/Fortran/test_new_interface_complex.F90 new file mode 100644 index 0000000000000000000000000000000000000000..8d249f77e0cea579783ba41bf3be7da6f10cce31 --- /dev/null +++ b/test/Fortran/test_new_interface_complex.F90 @@ -0,0 +1,173 @@ +! 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" + +#define stringify_(x) "x" +#define stringify(x) stringify_(x) +#define assert(x) call x_assert(x, stringify(x), __FILE__, __LINE__) + +module assert + implicit none + contains + subroutine x_assert(condition, condition_string, file, line) + use elpa_utilities, only : error_unit + logical, intent(in) :: condition + character(len=*), intent(in) :: condition_string + character(len=*), intent(in) :: file + integer, intent(in) :: line + + if (.not. condition) then + write(error_unit,'(a,i0)') "Assertion failed:" // condition_string // " at " // file // ":", line + end if + end subroutine +end module + +program test_interface + use precision + use assert + use mod_setup_mpi + use elpa_mpi + use elpa_type + use mod_prepare_matrix + use mod_read_input_parameters + use mod_blacs_infrastructure + use mod_check_correctness + + implicit none + + ! matrix dimensions + integer :: na, nev, nblk + + ! mpi + integer :: myid, nprocs + integer :: na_cols, na_rows ! local matrix size + integer :: np_cols, np_rows ! number of MPI processes per column/row + integer :: my_prow, my_pcol ! local MPI task position (my_prow, my_pcol) in the grid (0..np_cols -1, 0..np_rows -1) + integer :: mpierr + + ! blacs + integer :: my_blacs_ctxt, sc_desc(9), info, nprow, npcol + + ! The Matrix + real(kind=C_DOUBLE_COMPLEX), allocatable :: a(:,:), as(:,:) + ! eigenvectors + real(kind=C_DOUBLE_COMPLEX), allocatable :: z(:,:) + ! eigenvalues + real(kind=C_DOUBLE), allocatable :: ev(:) + + integer :: success, status + + integer(kind=c_int) :: solver + integer(kind=c_int) :: qr + + type(output_t) :: write_to_file + type(elpa_t) :: e + + call read_input_parameters(na, nev, nblk, write_to_file) + call setup_mpi(myid, nprocs) + + do np_cols = NINT(SQRT(REAL(nprocs))),2,-1 + if(mod(nprocs,np_cols) == 0 ) exit + enddo + + np_rows = nprocs/np_cols + + my_prow = mod(myid, np_cols) + my_pcol = myid / np_cols + + call set_up_blacsgrid(mpi_comm_world, my_blacs_ctxt, np_rows, np_cols, & + nprow, npcol, my_prow, my_pcol) + + call set_up_blacs_descriptor(na, nblk, my_prow, my_pcol, np_rows, np_cols, & + na_rows, na_cols, sc_desc, my_blacs_ctxt, info) + + allocate(a (na_rows,na_cols), as(na_rows,na_cols)) + allocate(z (na_rows,na_cols)) + allocate(ev(na)) + + a(:,:) = 0.0 + z(:,:) = 0.0 + ev(:) = 0.0 + + call prepare_matrix_double(na, myid, sc_desc, a, z, as) + + if (elpa_init(20170403) /= ELPA_OK) then + error stop "ELPA API version not supported" + endif + + e = elpa_create(na, nev, na_rows, na_cols, nblk, mpi_comm_world, my_prow, my_pcol, success) + assert(success == ELPA_OK) + + qr = e%get("qr", success) + print *, "qr =", qr + assert(success == ELPA_OK) + + solver = e%get("solver", success) + print *, "solver =", solver + assert(success == ELPA_OK) + + call e%set("solver", ELPA_SOLVER_2STAGE, success) + assert(success == ELPA_OK) + + call e%set("complex_kernel", ELPA_2STAGE_COMPLEX_GENERIC, success) + assert(success == ELPA_OK) + + call e%solve(a, ev, z, success) + assert(success == ELPA_OK) + + call e%destroy() + + call elpa_uninit() + + status = check_correctness_double(na, nev, as, z, ev, sc_desc, myid) + + deallocate(a) + deallocate(as) + deallocate(z) + deallocate(ev) + +#ifdef WITH_MPI + call mpi_finalize(mpierr) +#endif + +end program diff --git a/test/Fortran/test_new_interface_complex_single.F90 b/test/Fortran/test_new_interface_complex_single.F90 new file mode 100644 index 0000000000000000000000000000000000000000..86a39fc4a2632cf455ea4dc6680713fe2845a373 --- /dev/null +++ b/test/Fortran/test_new_interface_complex_single.F90 @@ -0,0 +1,173 @@ +! 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" + +#define stringify_(x) "x" +#define stringify(x) stringify_(x) +#define assert(x) call x_assert(x, stringify(x), __FILE__, __LINE__) + +module assert + implicit none + contains + subroutine x_assert(condition, condition_string, file, line) + use elpa_utilities, only : error_unit + logical, intent(in) :: condition + character(len=*), intent(in) :: condition_string + character(len=*), intent(in) :: file + integer, intent(in) :: line + + if (.not. condition) then + write(error_unit,'(a,i0)') "Assertion failed:" // condition_string // " at " // file // ":", line + end if + end subroutine +end module + +program test_interface + use precision + use assert + use mod_setup_mpi + use elpa_mpi + use elpa_type + use mod_prepare_matrix + use mod_read_input_parameters + use mod_blacs_infrastructure + use mod_check_correctness + + implicit none + + ! matrix dimensions + integer :: na, nev, nblk + + ! mpi + integer :: myid, nprocs + integer :: na_cols, na_rows ! local matrix size + integer :: np_cols, np_rows ! number of MPI processes per column/row + integer :: my_prow, my_pcol ! local MPI task position (my_prow, my_pcol) in the grid (0..np_cols -1, 0..np_rows -1) + integer :: mpierr + + ! blacs + integer :: my_blacs_ctxt, sc_desc(9), info, nprow, npcol + + ! The Matrix + real(kind=C_FLOAT_COMPLEX), allocatable :: a(:,:), as(:,:) + ! eigenvectors + real(kind=C_FLOAT_COMPLEX), allocatable :: z(:,:) + ! eigenvalues + real(kind=C_FLOAT), allocatable :: ev(:) + + integer :: success, status + + integer(kind=c_int) :: solver + integer(kind=c_int) :: qr + + type(output_t) :: write_to_file + type(elpa_t) :: e + + call read_input_parameters(na, nev, nblk, write_to_file) + call setup_mpi(myid, nprocs) + + do np_cols = NINT(SQRT(REAL(nprocs))),2,-1 + if(mod(nprocs,np_cols) == 0 ) exit + enddo + + np_rows = nprocs/np_cols + + my_prow = mod(myid, np_cols) + my_pcol = myid / np_cols + + call set_up_blacsgrid(mpi_comm_world, my_blacs_ctxt, np_rows, np_cols, & + nprow, npcol, my_prow, my_pcol) + + call set_up_blacs_descriptor(na, nblk, my_prow, my_pcol, np_rows, np_cols, & + na_rows, na_cols, sc_desc, my_blacs_ctxt, info) + + allocate(a (na_rows,na_cols), as(na_rows,na_cols)) + allocate(z (na_rows,na_cols)) + allocate(ev(na)) + + a(:,:) = 0.0 + z(:,:) = 0.0 + ev(:) = 0.0 + + call prepare_matrix_single(na, myid, sc_desc, a, z, as) + + if (elpa_init(20170403) /= ELPA_OK) then + error stop "ELPA API version not supported" + endif + + e = elpa_create(na, nev, na_rows, na_cols, nblk, mpi_comm_world, my_prow, my_pcol, success) + assert(success == ELPA_OK) + + qr = e%get("qr", success) + print *, "qr =", qr + assert(success == ELPA_OK) + + solver = e%get("solver", success) + print *, "solver =", solver + assert(success == ELPA_OK) + + call e%set("solver", ELPA_SOLVER_2STAGE, success) + assert(success == ELPA_OK) + + call e%set("complex_kernel", ELPA_2STAGE_COMPLEX_GENERIC, success) + assert(success == ELPA_OK) + + call e%solve(a, ev, z, success) + assert(success == ELPA_OK) + + call e%destroy() + + call elpa_uninit() + + status = check_correctness_single(na, nev, as, z, ev, sc_desc, myid) + + deallocate(a) + deallocate(as) + deallocate(z) + deallocate(ev) + +#ifdef WITH_MPI + call mpi_finalize(mpierr) +#endif + +end program diff --git a/test/Fortran/test_new_interface.F90 b/test/Fortran/test_new_interface_real.F90 similarity index 100% rename from test/Fortran/test_new_interface.F90 rename to test/Fortran/test_new_interface_real.F90 diff --git a/test/Fortran/test_new_interface_real_single.F90 b/test/Fortran/test_new_interface_real_single.F90 new file mode 100644 index 0000000000000000000000000000000000000000..8afece0a9bf849f408681f204d5d8c462d381a8b --- /dev/null +++ b/test/Fortran/test_new_interface_real_single.F90 @@ -0,0 +1,176 @@ +! 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" + +#define stringify_(x) "x" +#define stringify(x) stringify_(x) +#define assert(x) call x_assert(x, stringify(x), __FILE__, __LINE__) + +module assert + implicit none + contains + subroutine x_assert(condition, condition_string, file, line) + use elpa_utilities, only : error_unit + logical, intent(in) :: condition + character(len=*), intent(in) :: condition_string + character(len=*), intent(in) :: file + integer, intent(in) :: line + + if (.not. condition) then + write(error_unit,'(a,i0)') "Assertion failed:" // condition_string // " at " // file // ":", line + end if + end subroutine +end module + +program test_interface + use precision + use assert + use mod_setup_mpi + use elpa_mpi + use elpa_type + use mod_prepare_matrix + use mod_read_input_parameters + use mod_blacs_infrastructure + use mod_check_correctness + + implicit none + + ! matrix dimensions + integer :: na, nev, nblk + + ! mpi + integer :: myid, nprocs + integer :: na_cols, na_rows ! local matrix size + integer :: np_cols, np_rows ! number of MPI processes per column/row + integer :: my_prow, my_pcol ! local MPI task position (my_prow, my_pcol) in the grid (0..np_cols -1, 0..np_rows -1) + integer :: mpierr + + ! blacs + integer :: my_blacs_ctxt, sc_desc(9), info, nprow, npcol + + ! The Matrix + real(kind=C_FLOAT), allocatable :: a(:,:), as(:,:) + ! eigenvectors + real(kind=C_FLOAT), allocatable :: z(:,:) + ! eigenvalues + real(kind=C_FLOAT), allocatable :: ev(:) + + integer :: success, status + + integer(kind=c_int) :: solver + integer(kind=c_int) :: qr + + type(output_t) :: write_to_file + type(elpa_t) :: e + + call read_input_parameters(na, nev, nblk, write_to_file) + call setup_mpi(myid, nprocs) + + do np_cols = NINT(SQRT(REAL(nprocs))),2,-1 + if(mod(nprocs,np_cols) == 0 ) exit + enddo + + np_rows = nprocs/np_cols + + my_prow = mod(myid, np_cols) + my_pcol = myid / np_cols + + call set_up_blacsgrid(mpi_comm_world, my_blacs_ctxt, np_rows, np_cols, & + nprow, npcol, my_prow, my_pcol) + + call set_up_blacs_descriptor(na, nblk, my_prow, my_pcol, np_rows, np_cols, & + na_rows, na_cols, sc_desc, my_blacs_ctxt, info) + + allocate(a (na_rows,na_cols), as(na_rows,na_cols)) + allocate(z (na_rows,na_cols)) + allocate(ev(na)) + + a(:,:) = 0.0 + z(:,:) = 0.0 + ev(:) = 0.0 + + call prepare_matrix_single(na, myid, sc_desc, a, z, as) + + if (elpa_init(20170403) /= ELPA_OK) then + error stop "ELPA API version not supported" + endif + + e = elpa_create(na, nev, na_rows, na_cols, nblk, mpi_comm_world, my_prow, my_pcol, success) + assert(success == ELPA_OK) + + qr = e%get("qr", success) + print *, "qr =", qr + assert(success == ELPA_OK) + + solver = e%get("solver", success) + print *, "solver =", solver + assert(success == ELPA_OK) + + call e%set("solver", ELPA_SOLVER_2STAGE, success) + assert(success == ELPA_OK) + + call e%set("real_kernel", ELPA_2STAGE_REAL_GENERIC, success) + assert(success == ELPA_OK) + + call e%set("complex_kernel", ELPA_2STAGE_COMPLEX_GENERIC, success) + assert(success == ELPA_OK) + + call e%solve(a, ev, z, success) + assert(success == ELPA_OK) + + call e%destroy() + + call elpa_uninit() + + status = check_correctness_single(na, nev, as, z, ev, sc_desc, myid) + + deallocate(a) + deallocate(as) + deallocate(z) + deallocate(ev) + +#ifdef WITH_MPI + call mpi_finalize(mpierr) +#endif + +end program