Commit e865804e authored by Victor Yu's avatar Victor Yu
Browse files

Check every error code

And better message in case of errors.
parent 377f7bcc
......@@ -7,7 +7,7 @@ SET(elsi_URL "http://elsi-interchange.org")
SET(elsi_EMAIL "elsi-team@duke.edu")
SET(elsi_LICENSE "BSD 3")
SET(elsi_DESCRIPTION "Electronic Structure Infrastructure")
SET(elsi_DATESTAMP "20200516")
SET(elsi_DATESTAMP "20200524")
### CMake modules ###
LIST(APPEND CMAKE_MODULE_PATH ${PROJECT_SOURCE_DIR}/cmake)
......
......@@ -11,9 +11,9 @@ module ELSI_BSEPACK
use ELSI_DATATYPE, only: elsi_param_t,elsi_basic_t
use ELSI_MALLOC, only: elsi_allocate,elsi_deallocate
use ELSI_MPI, only: elsi_stop
use ELSI_OUTPUT, only: elsi_say,elsi_get_time
use ELSI_PRECISION, only: r8,i4
use ELSI_UTIL, only: elsi_check_err
implicit none
......@@ -79,14 +79,11 @@ subroutine elsi_solve_bsepack_real(ph,bh,mat_a,mat_b,eval,evec)
call pdbseig(0,ph%n_basis,mat_a,1,1,bh%desc,mat_b,1,1,bh%desc,eval,evec,1,1,&
ph%bse_desc,work,lwork,iwork,liwork,ierr)
call elsi_check_err(bh,"BSEPACK eigensolver",ierr,caller)
call elsi_deallocate(bh,work,"work")
call elsi_deallocate(bh,iwork,"iwork")
if(ierr /= 0) then
write(msg,"(A)") "BSEPACK eigensolver failed"
call elsi_stop(bh,msg,caller)
end if
call elsi_get_time(t1)
write(msg,"(A)") "Finished solving BSE eigenproblem"
......@@ -154,15 +151,12 @@ subroutine elsi_solve_bsepack_cmplx(ph,bh,mat_a,mat_b,eval,evec)
call pzbseig(0,ph%n_basis,mat_a,1,1,bh%desc,mat_b,1,1,bh%desc,eval,evec,1,1,&
ph%bse_desc,work,lwork,rwork,lrwork,iwork,liwork,ierr)
call elsi_check_err(bh,"BSEPACK eigensolver",ierr,caller)
call elsi_deallocate(bh,work,"work")
call elsi_deallocate(bh,rwork,"rwork")
call elsi_deallocate(bh,iwork,"iwork")
if(ierr /= 0) then
write(msg,"(A)") "BSEPACK eigensolver failed"
call elsi_stop(bh,msg,caller)
end if
call elsi_get_time(t1)
write(msg,"(A)") "Finished solving BSE eigenproblem"
......
......@@ -13,9 +13,10 @@ module ELSI_DECISION
use ELSI_CONSTANT, only: AUTO_SOLVER,ELPA_SOLVER,PEXSI_SOLVER,NTPOLY_SOLVER,&
UNSET
use ELSI_DATATYPE, only: elsi_param_t,elsi_basic_t
use ELSI_MPI, only: elsi_check_mpi,MPI_SUM,MPI_INTEGER4,MPI_REAL8
use ELSI_MPI, only: MPI_SUM,MPI_INTEGER4,MPI_REAL8
use ELSI_OUTPUT, only: elsi_say
use ELSI_PRECISION, only: r8,i4
use ELSI_UTIL, only: elsi_check_err
implicit none
......@@ -79,14 +80,14 @@ subroutine elsi_decide_dm_real(ph,bh,mat)
call MPI_Allreduce(nnz_l,nnz_g,1,MPI_INTEGER4,MPI_SUM,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
sparsity = 1.0_r8-(1.0_r8*nnz_g/ph%n_basis/ph%n_basis)
end if
call MPI_Bcast(sparsity,1,MPI_REAL8,0,bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
call elsi_decide_dm_smart(ph,bh,sparsity)
end if
......@@ -117,14 +118,14 @@ subroutine elsi_decide_dm_cmplx(ph,bh,mat)
call MPI_Allreduce(nnz_l,nnz_g,1,MPI_INTEGER4,MPI_SUM,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
sparsity = 1.0_r8-(1.0_r8*nnz_g/ph%n_basis/ph%n_basis)
end if
call MPI_Bcast(sparsity,1,MPI_REAL8,0,bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
call elsi_decide_dm_smart(ph,bh,sparsity)
end if
......@@ -151,7 +152,7 @@ subroutine elsi_decide_dm_sparse(ph,bh)
call MPI_Bcast(sparsity,1,MPI_REAL8,0,bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
call elsi_decide_dm_smart(ph,bh,sparsity)
end if
......
......@@ -14,10 +14,11 @@ module ELSI_EIGENEXA
use ELSI_ELPA, only: elsi_factor_ovlp_elpa,elsi_reduce_evp_elpa,&
elsi_back_ev_elpa
use ELSI_MALLOC, only: elsi_allocate,elsi_deallocate
use ELSI_MPI, only: elsi_check_mpi,MPI_SUM,MPI_INTEGER4
use ELSI_MPI, only: MPI_SUM,MPI_INTEGER4
use ELSI_OUTPUT, only: elsi_say,elsi_get_time
use ELSI_PRECISION, only: r8,i4
use ELSI_REDIST, only: elsi_blacs_to_eigenexa_h,elsi_eigenexa_to_blacs_ev
use ELSI_UTIL, only: elsi_check_err
use EIGEN_LIBS_MOD, only: eigen_init,eigen_get_procs,eigen_get_id,&
eigen_get_matdims,eigen_s,eigen_sx,eigen_free
......@@ -94,7 +95,7 @@ subroutine elsi_solve_eigenexa_real(ph,bh,ham,ovlp,eval,evec)
call MPI_Allreduce(bh%nnz_l,bh%nnz_g,1,MPI_INTEGER4,MPI_SUM,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
end if
! Transform to standard form
......
This diff is collapsed.
......@@ -147,7 +147,7 @@ subroutine elsi_solve_magma_real(ph,bh,ham,ovlp,eval,evec)
end select
end if
if(ierr /= 0) then
if(ierr /= 0 .or. n_solved(1) < ph%n_states) then
write(msg,"(A)") "MAGMA eigensolver failed"
call elsi_stop(bh,msg,caller)
end if
......
......@@ -5,7 +5,7 @@
! which may be found in the LICENSE file in the ELSI root directory.
!>
!! Provide a collection of MPI utilities.
!! MPI interface.
!!
module ELSI_MPI
......@@ -16,7 +16,6 @@ module ELSI_MPI
implicit none
public :: elsi_stop
public :: elsi_check_mpi
contains
......@@ -55,26 +54,4 @@ subroutine elsi_stop(bh,info,caller)
end subroutine
!>
!! Checks if an MPI call is successful.
!!
subroutine elsi_check_mpi(bh,routine,ierr,caller)
implicit none
type(elsi_basic_t), intent(in) :: bh
character(len=*), intent(in) :: routine
integer(kind=i4), intent(in) :: ierr
character(len=*), intent(in) :: caller
character(len=200) :: msg
if(ierr /= MPI_SUCCESS) then
write(msg,"(2A)") trim(routine)," failed"
call elsi_stop(bh,msg,caller)
end if
end subroutine
end module ELSI_MPI
......@@ -5,7 +5,7 @@
! which may be found in the LICENSE file in the ELSI root directory.
!>
!! Provide a collection of MPI utilities.
!! MPI interface.
!!
module ELSI_MPI
......@@ -17,7 +17,6 @@ module ELSI_MPI
include "mpif.h"
public :: elsi_stop
public :: elsi_check_mpi
contains
......@@ -56,26 +55,4 @@ subroutine elsi_stop(bh,info,caller)
end subroutine
!>
!! Checks if an MPI call is successful.
!!
subroutine elsi_check_mpi(bh,routine,ierr,caller)
implicit none
type(elsi_basic_t), intent(in) :: bh
character(len=*), intent(in) :: routine
integer(kind=i4), intent(in) :: ierr
character(len=*), intent(in) :: caller
character(len=200) :: msg
if(ierr /= MPI_SUCCESS) then
write(msg,"(2A)") trim(routine)," failed"
call elsi_stop(bh,msg,caller)
end if
end subroutine
end module ELSI_MPI
......@@ -12,9 +12,10 @@ module ELSI_NTPOLY
use ELSI_CONSTANT, only: NTPOLY_SOLVER,NTPOLY_PM,NTPOLY_TRS2,NTPOLY_TRS4,&
NTPOLY_HPCP,EXTRA_FACTOR,EXTRA_TRS2
use ELSI_DATATYPE, only: elsi_param_t,elsi_basic_t
use ELSI_MPI, only: elsi_check_mpi,MPI_LOGICAL
use ELSI_MPI, only: MPI_LOGICAL
use ELSI_OUTPUT, only: elsi_say,elsi_get_time
use ELSI_PRECISION, only: r8,i4
use ELSI_UTIL, only: elsi_check_err
use NTPOLY, only: PM,TRS2,TRS4,HPCP,EnergyDensityMatrix,LowdinExtrapolate,&
PurificationExtrapolate,Matrix_ps,ConstructEmptyMatrix,DestructMatrix,&
CopyMatrix,ScaleMatrix,FillMatrixFromTripletList,GetMatrixTripletList,&
......@@ -92,7 +93,7 @@ subroutine elsi_init_ntpoly(ph,bh)
call MPI_Bcast(ph%nt_output,1,MPI_LOGICAL,0,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
ph%nt_filter = max(ph%nt_filter,bh%def0)
......
......@@ -13,10 +13,11 @@ module ELSI_OCC
SQRT_PI,INVERT_SQRT_PI
use ELSI_DATATYPE, only: elsi_handle,elsi_param_t,elsi_basic_t
use ELSI_MALLOC, only: elsi_allocate,elsi_deallocate
use ELSI_MPI, only: elsi_stop,elsi_check_mpi,MPI_SUM,MPI_REAL8
use ELSI_MPI, only: elsi_stop,MPI_SUM,MPI_REAL8
use ELSI_OUTPUT, only: elsi_say
use ELSI_PRECISION, only: r8,i4
use ELSI_SORT, only: elsi_heapsort,elsi_permute,elsi_unpermute
use ELSI_UTIL, only: elsi_check_err
implicit none
......@@ -620,7 +621,7 @@ subroutine elsi_get_occ_for_dm(ph,bh,eval,occ)
call MPI_Allreduce(tmp1,k_wt,ph%n_kpts,MPI_REAL8,MPI_SUM,bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
call elsi_deallocate(bh,tmp1,"tmp")
else
......@@ -637,7 +638,7 @@ subroutine elsi_get_occ_for_dm(ph,bh,eval,occ)
call MPI_Allreduce(tmp2,eval_all,ph%n_states*ph%n_spins*ph%n_kpts,&
MPI_REAL8,MPI_SUM,bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
call elsi_deallocate(bh,tmp2,"tmp")
else
......
......@@ -11,10 +11,10 @@ module ELSI_OMM
use ELSI_CONSTANT, only: UNSET
use ELSI_DATATYPE, only: elsi_param_t,elsi_basic_t
use ELSI_ELPA, only: elsi_elpa_cholesky,elsi_elpa_invert
use ELSI_MPI, only: elsi_check_mpi,MPI_SUM,MPI_INTEGER4
use ELSI_MPI, only: MPI_SUM,MPI_INTEGER4
use ELSI_OUTPUT, only: elsi_say,elsi_get_time
use ELSI_PRECISION, only: r8,i4
use ELSI_UTIL, only: elsi_check_err
use MATRIXSWITCH, only: matrix,m_register_pdbc,ms_scalapack_setup,&
m_deallocate
......@@ -115,7 +115,7 @@ subroutine elsi_solve_omm_real(ph,bh,ham,ovlp,coeff,dm)
call MPI_Allreduce(bh%nnz_l,bh%nnz_g,1,MPI_INTEGER4,MPI_SUM,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
end if
write(msg,"(A)") "Starting libOMM density matrix solver"
......@@ -131,7 +131,9 @@ subroutine elsi_solve_omm_real(ph,bh,ham,ovlp,coeff,dm)
if(ph%omm_first .and. ph%elpa_first) then
call elsi_get_time(t0)
call elsi_elpa_cholesky(ph,bh,ovlp)
call ph%elpa_aux%cholesky(ovlp,ierr)
call elsi_check_err(bh,"ELPA Cholesky factorization",ierr,caller)
call elsi_get_time(t1)
......@@ -142,7 +144,9 @@ subroutine elsi_solve_omm_real(ph,bh,ham,ovlp,coeff,dm)
end if
if(ph%omm_first) then
call elsi_elpa_invert(ph,bh,ovlp)
call ph%elpa_aux%invert_triangular(ovlp,ierr)
call elsi_check_err(bh,"ELPA matrix inversion",ierr,caller)
end if
end if
end if
......@@ -273,7 +277,7 @@ subroutine elsi_solve_omm_cmplx(ph,bh,ham,ovlp,coeff,dm)
call MPI_Allreduce(bh%nnz_l,bh%nnz_g,1,MPI_INTEGER4,MPI_SUM,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
end if
write(msg,"(A)") "Starting libOMM density matrix solver"
......@@ -289,7 +293,9 @@ subroutine elsi_solve_omm_cmplx(ph,bh,ham,ovlp,coeff,dm)
if(ph%omm_first .and. ph%elpa_first) then
call elsi_get_time(t0)
call elsi_elpa_cholesky(ph,bh,ovlp)
call ph%elpa_aux%cholesky(ovlp,ierr)
call elsi_check_err(bh,"ELPA Cholesky factorization",ierr,caller)
call elsi_get_time(t1)
......@@ -300,7 +306,9 @@ subroutine elsi_solve_omm_cmplx(ph,bh,ham,ovlp,coeff,dm)
end if
if(ph%omm_first) then
call elsi_elpa_invert(ph,bh,ovlp)
call ph%elpa_aux%invert_triangular(ovlp,ierr)
call elsi_check_err(bh,"ELPA matrix inversion",ierr,caller)
end if
end if
end if
......
......@@ -12,9 +12,10 @@ module ELSI_PEXSI
use ELSI_CONSTANT, only: UNSET,PEXSI_CSC,PEXSI_DM,PEXSI_EDM,PEXSI_FDM
use ELSI_DATATYPE, only: elsi_param_t,elsi_basic_t
use ELSI_MALLOC, only: elsi_allocate,elsi_deallocate
use ELSI_MPI, only: elsi_stop,elsi_check_mpi,MPI_SUM,MPI_REAL8,MPI_COMPLEX16
use ELSI_MPI, only: MPI_SUM,MPI_REAL8,MPI_COMPLEX16
use ELSI_OUTPUT, only: elsi_say,elsi_get_time
use ELSI_PRECISION, only: r8,i4
use ELSI_UTIL, only: elsi_check_err
use F_PPEXSI_INTERFACE, only: f_ppexsi_plan_initialize,&
f_ppexsi_plan_finalize,f_ppexsi_set_default_options,&
f_ppexsi_load_real_hs_matrix,f_ppexsi_load_complex_hs_matrix,&
......@@ -70,7 +71,6 @@ subroutine elsi_init_pexsi(ph,bh)
integer(kind=i4) :: j
integer(kind=i4) :: log_id
integer(kind=i4) :: ierr
character(len=200) :: msg
character(len=*), parameter :: caller = "elsi_init_pexsi"
......@@ -113,17 +113,17 @@ subroutine elsi_init_pexsi(ph,bh)
call MPI_Comm_split(bh%comm,ph%pexsi_my_prow,ph%pexsi_my_pcol,&
ph%pexsi_comm_intra_pole,ierr)
call elsi_check_mpi(bh,"MPI_Comm_split",ierr,caller)
call elsi_check_err(bh,"MPI_Comm_split",ierr,caller)
call MPI_Comm_split(bh%comm,ph%pexsi_my_pcol,ph%pexsi_my_prow,&
ph%pexsi_comm_inter_pole,ierr)
call elsi_check_mpi(bh,"MPI_Comm_split",ierr,caller)
call elsi_check_err(bh,"MPI_Comm_split",ierr,caller)
call MPI_Comm_split(bh%comm,ph%pexsi_myid_point,ph%pexsi_my_point,&
ph%pexsi_comm_inter_point,ierr)
call elsi_check_mpi(bh,"MPI_Comm_split",ierr,caller)
call elsi_check_err(bh,"MPI_Comm_split",ierr,caller)
if(.not. bh%pexsi_csc_ready) then
! Set up 1D block distribution
......@@ -145,10 +145,7 @@ subroutine elsi_init_pexsi(ph,bh)
ph%pexsi_plan = f_ppexsi_plan_initialize(bh%comm,ph%pexsi_n_prow,&
ph%pexsi_n_pcol,log_id,ierr)
if(ierr /= 0) then
write(msg,"(A)") "Initialization failed"
call elsi_stop(bh,msg,caller)
end if
call elsi_check_err(bh,"PEXSI initialization",ierr,caller)
if(bh%n_lcol_sp == UNSET) then
bh%n_lcol_sp = bh%n_lcol_sp1
......@@ -230,10 +227,7 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
0,ovlp,ierr)
end if
if(ierr /= 0) then
write(msg,"(A)") "Failed to load matrices"
call elsi_stop(bh,msg,caller)
end if
call elsi_check_err(bh,"PEXSI load matrices",ierr,caller)
call elsi_get_time(t1)
......@@ -249,9 +243,13 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call f_ppexsi_symbolic_factorize_real_symmetric_matrix(ph%pexsi_plan,&
ph%pexsi_options,ierr)
call elsi_check_err(bh,"PEXSI symbolic factorization",ierr,caller)
call f_ppexsi_symbolic_factorize_complex_symmetric_matrix(ph%pexsi_plan,&
ph%pexsi_options,ierr)
call elsi_check_err(bh,"PEXSI symbolic factorization",ierr,caller)
call elsi_get_time(t1)
write(msg,"(A)") "Finished symbolic factorization"
......@@ -260,11 +258,6 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call elsi_say(bh,msg)
end if
if(ierr /= 0) then
write(msg,"(A)") "Symbolic factorization failed"
call elsi_stop(bh,msg,caller)
end if
! Inertia counting
mu_range = ph%pexsi_options%muMax0-ph%pexsi_options%muMin0
......@@ -292,6 +285,8 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call f_ppexsi_inertia_count_real_matrix(ph%pexsi_plan,&
ph%pexsi_options,n_shift,shifts,inertias,ierr)
call elsi_check_err(bh,"PEXSI inertia counting",ierr,caller)
inertias(:) = inertias*ph%spin_degen*ph%i_wt
! Get global inertias
......@@ -307,7 +302,7 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call MPI_Allreduce(send_buf,inertias,n_shift,MPI_REAL8,MPI_SUM,&
bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
call elsi_deallocate(bh,send_buf,"send_buf")
end if
......@@ -355,11 +350,6 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call elsi_say(bh,msg)
write(msg,"(A,F10.3,A)") "| Time :",t1-t0," s"
call elsi_say(bh,msg)
if(ierr /= 0) then
write(msg,"(A)") "Inertia counting failed"
call elsi_stop(bh,msg,caller)
end if
end if
! Save chemical potential bounds
......@@ -377,6 +367,8 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
if(ph%pexsi_my_point == i-1) then
call f_ppexsi_calculate_fermi_operator_real3(ph%pexsi_plan,&
ph%pexsi_options,ph%n_electrons,ph%mu,ph%pexsi_ne,ierr)
call elsi_check_err(bh,"PEXSI Fermi operator expansion",ierr,caller)
end if
end do
......@@ -387,7 +379,7 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call MPI_Allreduce(send_buf,ne_vec,ph%pexsi_options%nPoints,MPI_REAL8,&
MPI_SUM,ph%pexsi_comm_inter_point,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
! Get global number of electrons
if(ph%n_spins*ph%n_kpts > 1) then
......@@ -400,7 +392,7 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call MPI_Allreduce(send_buf,ne_vec,ph%pexsi_options%nPoints,MPI_REAL8,&
MPI_SUM,bh%comm_all,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
end if
call elsi_deallocate(bh,send_buf,"send_buf")
......@@ -412,11 +404,6 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
write(msg,"(A,F10.3,A)") "| Time :",t1-t0," s"
call elsi_say(bh,msg)
if(ierr /= 0) then
write(msg,"(A)") "Fermi operator calculation failed"
call elsi_stop(bh,msg,caller)
end if
if(ph%pexsi_options%method /= 2) then
! Get free energy density matrix
call elsi_get_time(t0)
......@@ -430,12 +417,12 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call MPI_Reduce(local_energy,free_energy,1,MPI_REAL8,MPI_SUM,0,&
ph%pexsi_comm_intra_pole,ierr)
call elsi_check_mpi(bh,"MPI_Reduce",ierr,caller)
call elsi_check_err(bh,"MPI_Reduce",ierr,caller)
end if
call MPI_Bcast(free_energy,1,MPI_REAL8,0,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
end if
! Get density matrix
......@@ -450,12 +437,12 @@ subroutine elsi_solve_pexsi_real(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call MPI_Reduce(local_energy,ph%ebs,1,MPI_REAL8,MPI_SUM,0,&
ph%pexsi_comm_intra_pole,ierr)
call elsi_check_mpi(bh,"MPI_Reduce",ierr,caller)
call elsi_check_err(bh,"MPI_Reduce",ierr,caller)
end if
call MPI_Bcast(ph%ebs,1,MPI_REAL8,0,bh%comm,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
! Compute entropy
if(ph%pexsi_options%method /= 2) then
......@@ -529,7 +516,6 @@ subroutine elsi_retrieve_dm_pexsi_real(ph,bh,which,ne_vec,dm)
integer(kind=i4) :: i
integer(kind=i4) :: ierr
logical :: converged
character(len=200) :: msg
real(kind=r8), allocatable :: tmp(:)
real(kind=r8), allocatable :: send_buf(:)
......@@ -543,29 +529,16 @@ subroutine elsi_retrieve_dm_pexsi_real(ph,bh,which,ne_vec,dm)
select case(which)
case(PEXSI_DM)
call f_ppexsi_retrieve_real_dm(ph%pexsi_plan,tmp,local_energy,ierr)
if(ierr /= 0) then
write(msg,"(A)") "Failed to get density matirx"
call elsi_stop(bh,msg,caller)
end if
case(PEXSI_EDM)
call f_ppexsi_retrieve_real_edm(ph%pexsi_plan,ph%pexsi_options,tmp,&
local_energy,ierr)
if(ierr /= 0) then
write(msg,"(A)") "Failed to get energy density matirx"
call elsi_stop(bh,msg,caller)
end if
case(PEXSI_FDM)
call f_ppexsi_retrieve_real_fdm(ph%pexsi_plan,ph%pexsi_options,tmp,&
local_energy,ierr)
if(ierr /= 0) then
write(msg,"(A)") "Failed to get free energy density matirx"
call elsi_stop(bh,msg,caller)
end if
end select
call elsi_check_err(bh,"PEXSI get density matrix",ierr,caller)
! Check convergence
mu_range = ph%pexsi_mu_max-ph%pexsi_mu_min
shift_width = mu_range/(ph%pexsi_options%nPoints+1)
......@@ -627,7 +600,7 @@ subroutine elsi_retrieve_dm_pexsi_real(ph,bh,which,ne_vec,dm)
call MPI_Bcast(tmp,bh%nnz_l_sp1,MPI_REAL8,i-1,&
ph%pexsi_comm_inter_point,ierr)
call elsi_check_mpi(bh,"MPI_Bcast",ierr,caller)
call elsi_check_err(bh,"MPI_Bcast",ierr,caller)
exit
end if
......@@ -657,7 +630,7 @@ subroutine elsi_retrieve_dm_pexsi_real(ph,bh,which,ne_vec,dm)
call MPI_Allreduce(send_buf,tmp,bh%nnz_l_sp1,MPI_REAL8,MPI_SUM,&
ph%pexsi_comm_inter_point,ierr)
call elsi_check_mpi(bh,"MPI_Allreduce",ierr,caller)
call elsi_check_err(bh,"MPI_Allreduce",ierr,caller)
call elsi_deallocate(bh,send_buf,"send_buf")
end if
......@@ -730,10 +703,7 @@ subroutine elsi_solve_pexsi_cmplx(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
0,ovlp,ierr)
end if
if(ierr /= 0) then
write(msg,"(A)") "Failed to load matirces"
call elsi_stop(bh,msg,caller)
end if
call elsi_check_err(bh,"PEXSI load matrices",ierr,caller)
call elsi_get_time(t1)
......@@ -749,9 +719,13 @@ subroutine elsi_solve_pexsi_cmplx(ph,bh,row_ind,col_ptr,ne_vec,ham,ovlp,dm)
call f_ppexsi_symbolic_factorize_complex_symmetric_matrix(ph%pexsi_plan,&
ph%pexsi_options,ierr