From c9f2e983686c181e368b4dcb7b4974cc18a13dcb Mon Sep 17 00:00:00 2001 From: Daniel Wortmann Date: Wed, 21 Dec 2016 11:51:26 +0100 Subject: [PATCH] Fixed ELPA for version 2016.05.003 --- cmake/tests/test_ELPA.cmake | 4 ++-- diagonalization/eigen_diag.F90 | 4 ++-- diagonalization/elpa.F90 | 17 ++++++++++------- diagonalization/elpa_cpp.F90 | 35 ++++++++++++++++++++++++++++------ 4 files changed, 43 insertions(+), 17 deletions(-) diff --git a/cmake/tests/test_ELPA.cmake b/cmake/tests/test_ELPA.cmake index e14a25d7..b5ee0d1f 100644 --- a/cmake/tests/test_ELPA.cmake +++ b/cmake/tests/test_ELPA.cmake @@ -6,5 +6,5 @@ LINK_LIBRARIES ${FLEUR_LIBRARIES} message("ELPA Library found:${FLEUR_USE_ELPA}") if (FLEUR_USE_ELPA) - set(FLEUR_MPI_DEFINITIONS ${FLEUR_MPI_DEFINITIONS} "CPP_ELPA") -endif() \ No newline at end of file + set(FLEUR_MPI_DEFINITIONS ${FLEUR_MPI_DEFINITIONS} "CPP_ELPA" "CPP_ELPA_201605003") +endif() diff --git a/diagonalization/eigen_diag.F90 b/diagonalization/eigen_diag.F90 index b387d15a..6e4f8826 100644 --- a/diagonalization/eigen_diag.F90 +++ b/diagonalization/eigen_diag.F90 @@ -163,9 +163,9 @@ CONTAINS #ifdef CPP_ELPA CASE (diag_elpa) IF (hamovlp%l_real) THEN - CALL elpa(lapw%nmat,n,SUB_COMM,hamOvlp%a_r,hamOvlp%b_r,zMat%z_r,eig,ne_found) + CALL elpa_diag(lapw%nmat,n,SUB_COMM,hamOvlp%a_r,hamOvlp%b_r,zMat%z_r,eig,ne_found) ELSE - CALL elpa(lapw%nmat,n,SUB_COMM,hamOvlp%a_c,hamOvlp%b_c,zMat%z_c,eig,ne_found) + CALL elpa_diag(lapw%nmat,n,SUB_COMM,hamOvlp%a_c,hamOvlp%b_c,zMat%z_c,eig,ne_found) ENDIF #endif #ifdef CPP_ELEMENTAL diff --git a/diagonalization/elpa.F90 b/diagonalization/elpa.F90 index 507e791a..de7a393a 100644 --- a/diagonalization/elpa.F90 +++ b/diagonalization/elpa.F90 @@ -5,20 +5,21 @@ !-------------------------------------------------------------------------------- MODULE m_elpa PRIVATE - - interface elpa +#ifdef CPP_ELPA + interface elpa_diag module procedure elpa_r,elpa_c - end interface elpa + end interface elpa_diag !Module to call elpa library for parallel diagonalization !uses ssubredist1/2 for redistribution - PUBLIC elpa + PUBLIC elpa_diag CONTAINS ! First the real version of the code #define CPP_CHOLESKY cholesky_real #define CPP_invert_trm invert_trm_real +#define CPP_solve_evp_1stage solve_evp_real_1stage #define CPP_solve_evp solve_evp_real #define CPP_solve_evp_2stage solve_evp_real_2stage #define CPP_REALCOMPLEX real @@ -55,6 +56,7 @@ CONTAINS ! Now the complex version #define CPP_CHOLESKY cholesky_complex #define CPP_invert_trm invert_trm_complex +#define CPP_solve_evp_1stage solve_evp_complex_1stage #define CPP_solve_evp solve_evp_complex #define CPP_solve_evp_2stage solve_evp_complex_2stage #define CPP_REALCOMPLEX complex @@ -190,7 +192,7 @@ CONTAINS !Create communicators for ELPA !print *,"creating ELPA comms" -#ifdef CPP_ELPA_NEW +#if (defined CPP_ELPA_201605003)||defined(CPP_ELPA_NEW) ierr=get_elpa_row_col_comms(mpi_subcom, myrowblacs, mycolblacs,mpi_comm_rows, mpi_comm_cols) #else CALL get_elpa_row_col_comms(mpi_subcom, myrowblacs, mycolblacs,mpi_comm_rows, mpi_comm_cols) @@ -202,7 +204,8 @@ CONTAINS IF (ierr /=0 ) CALL juDFT_error('descinit1 failed',calledby='elpa') END SUBROUTINE priv_create_blacsgrid - - +#else + LOGICAL :: elpa_used=.false. +#endif end MODULE m_elpa diff --git a/diagonalization/elpa_cpp.F90 b/diagonalization/elpa_cpp.F90 index 860db613..a956f290 100644 --- a/diagonalization/elpa_cpp.F90 +++ b/diagonalization/elpa_cpp.F90 @@ -4,9 +4,9 @@ ! of the MIT license as expressed in the LICENSE file in more detail. !-------------------------------------------------------------------------------- USE m_juDFT - USE elpa1 USE m_subredist1 USE m_subredist2 + USE elpa1 #ifdef CPP_ELPA2 USE elpa2 #endif @@ -84,7 +84,10 @@ ! Please note: cholesky_complex/invert_trm_complex are not trimmed for speed. ! The only reason having them is that the Scalapack counterpart ! PDPOTRF very often fails on higher processor numbers for unknown reasons! -#ifdef CPP_ELPA_NEW +#ifdef CPP_ELPA_201605003 + ok=CPP_CHOLESKY (m,bsca,SIZE(bsca,1),nb,mycolssca,mpi_comm_rows,mpi_comm_cols,.false.) + ok=CPP_invert_trm(m,bsca,SIZE(bsca,1),nb,mycolssca,mpi_comm_rows,mpi_comm_cols,.false.) +#elif defined CPP_ELPA_NEW CALL CPP_CHOLESKY (m,bsca,SIZE(bsca,1),nb,mycolssca,mpi_comm_rows,mpi_comm_cols,.false.,ok) CALL CPP_invert_trm(m,bsca,SIZE(bsca,1),nb,mycolssca,mpi_comm_rows,mpi_comm_cols,.false.,ok) #else @@ -106,15 +109,23 @@ n_row = numroc (n_col, nb, myrow, 0, nprow) asca(n_row+1:myrowssca,i) = eigvec(n_row+1:myrowssca,i) ENDDO - +#ifdef CPP_ELPA_201605003 + ok=CPP_mult ('U', 'L',m, m,bsca,myrowssca,asca,SIZE(asca,1),nb, mpi_comm_rows, mpi_comm_cols,eigvec,myrowssca) +#else CALL CPP_mult ('U', 'L',m, m,bsca,myrowssca,asca,SIZE(asca,1),nb, mpi_comm_rows, mpi_comm_cols,eigvec,myrowssca) - +#endif ! 2b. tmp2 = eigvec**T CALL CPP_transpose(m,m,CPP_ONE,eigvec,1,1,sc_desc,CPP_ZERO,tmp2,1,1,sc_desc) ! 2c. A = U**-T * tmp2 ( = U**-T * Aorig * U**-1 ) +#ifdef CPP_ELPA_201605003 + ok=CPP_mult ('U', 'U', m, m, bsca, SIZE(bsca,1), tmp2,& + SIZE(tmp2,1),nb, mpi_comm_rows, mpi_comm_cols, asca, SIZE(asca,1)) +#else CALL CPP_mult ('U', 'U', m, m, bsca, SIZE(bsca,1), tmp2,& SIZE(tmp2,1),nb, mpi_comm_rows, mpi_comm_cols, asca, SIZE(asca,1)) +#endif + ! A is only set in the upper half, solve_evp_real needs a full matrix ! Set lower half from upper half @@ -133,7 +144,15 @@ ! 3. Calculate eigenvalues/eigenvectors of U**-T * A * U**-1 ! Eigenvectors go to eigvec num2=num -#ifdef CPP_ELPA_NEW +#ifdef CPP_ELPA_201605003 +#ifdef CPP_ELPA2 + ok=CPP_solve_evp_2stage(m,num2,asca,SIZE(asca,1),& + eig2,eigvec,SIZE(asca,1), nb,mycolssca, mpi_comm_rows, mpi_comm_cols,sub_comm) +#else + ok=CPP_solve_evp_1stage(m, num2,asca,SIZE(asca,1),& + eig2,eigvec, SIZE(asca,1), nb,mycolssca,mpi_comm_rows, mpi_comm_cols) +#endif +#elif defined CPP_ELPA_NEW #ifdef CPP_ELPA2 err=CPP_solve_evp_2stage(m,num2,asca,SIZE(asca,1),& eig2,eigvec,SIZE(asca,1), nb,mycolssca, mpi_comm_rows, mpi_comm_cols,sub_comm) @@ -156,9 +175,13 @@ ! mult_ah_b_complex needs the transpose of U**-1, thus tmp2 = (U**-1)**T CALL CPP_transpose(m,m,CPP_ONE,bsca,1,1,sc_desc,CPP_ZERO,tmp2,1,1,sc_desc) +#ifdef CPP_ELPA_201605003 + ok= CPP_mult ('L', 'N',m, num2, tmp2, SIZE(asca,1),& + eigvec, SIZE(asca,1),nb,mpi_comm_rows, mpi_comm_cols, asca, SIZE(asca,1)) +#else CALL CPP_mult ('L', 'N',m, num2, tmp2, SIZE(asca,1),& eigvec, SIZE(asca,1),nb,mpi_comm_rows, mpi_comm_cols, asca, SIZE(asca,1)) - +#endif ! END of ELPA stuff CALL BLACS_GRIDEXIT(ictextblacs,ierr) CALL MPI_COMM_FREE(mpi_comm_rows,ierr) -- GitLab