TABLE OF CONTENTS
m_abi_linalg/abi_chpev [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_chpev
FUNCTION
INPUTS
SOURCE
120 subroutine abi_chpev(jobz,uplo,n,a,w,z,ldz) 121 122 !Arguments ------------------------------------ 123 character(len=1), intent(in) :: jobz 124 character(len=1), intent(in) :: uplo 125 integer, intent(in) :: n,ldz 126 complex(spc), intent(inout) :: a(:,:) 127 complex(spc), intent(out) :: z(:,:) 128 real(sp), intent(out) :: w(:) 129 130 !Local variables------------------------------- 131 integer :: info 132 real(sp),pointer :: rwork(:) 133 complex(spc),pointer :: work(:) 134 135 ! ********************************************************************* 136 137 ABI_CHECK(lapack_packed_storage,"BUG(1) in abi_chpev (storage)!") 138 ABI_CHECK(lapack_single_precision,"BUG(2) in abi_chpev (precision)!") 139 ABI_CHECK(n<=eigen_c_maxsize,"BUG(3) in abi_chpev (maxsize)!") 140 141 work => eigen_c_work ; rwork => eigen_c_rwork 142 143 !===== LAPACK 144 if (eigen_c_lwork==0) then 145 ABI_MALLOC(work,(2*n-1)) 146 end if 147 if (eigen_c_lrwork==0) then 148 ABI_MALLOC(rwork,(3*n-2)) 149 end if 150 call chpev(jobz,uplo,n,a,w,z,ldz,work,rwork,info) 151 if (eigen_c_lwork==0) then 152 ABI_FREE(work) 153 end if 154 if (eigen_c_lrwork==0) then 155 ABI_FREE(rwork) 156 end if 157 158 ABI_CHECK(info==0,"abi_chpev returned info!=0!") 159 160 end subroutine abi_chpev
m_abi_linalg/abi_dhpev [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_dhpev
FUNCTION
INPUTS
SOURCE
31 subroutine abi_dhpev(jobz,uplo,n,a,w,z,ldz,istwf_k,use_slk,use_gpu_elpa) 32 33 !Arguments ------------------------------------ 34 character(len=1), intent(in) :: jobz 35 character(len=1), intent(in) :: uplo 36 integer, intent(in) :: n,ldz 37 real(dp), intent(inout) :: a(:) 38 real(dp), intent(out) :: z(:,:) 39 real(dp), intent(out) :: w(:) 40 integer, optional, intent(in) :: istwf_k 41 integer, optional, intent(in) :: use_slk,use_gpu_elpa 42 43 !Local variables------------------------------- 44 integer :: info,use_slk_,use_gpu_elpa_,istwf_k_ 45 #ifdef HAVE_LINALG_SCALAPACK 46 type(matrix_scalapack) :: sca_a,sca_ev 47 real(dp),allocatable :: tmp_evec(:,:) 48 integer :: dim_evec1,ierr 49 #endif 50 51 ! ********************************************************************* 52 53 ABI_CHECK(lapack_packed_storage,"BUG(1) in abi_dhpev (storage)!") 54 ABI_CHECK(lapack_double_precision,"BUG(2) in abi_dhpev (precision)!") 55 ABI_CHECK(n<=eigen_d_maxsize,"BUG(3) in abi_dhpev (maxsize)!") 56 57 info = 0 ! to avoid unwanted warnings but scalapack doesn't check return 58 59 use_slk_ = 0; if (present(use_slk)) use_slk_ = use_slk 60 istwf_k_ = 1; if (present(istwf_k)) istwf_k_ = istwf_k 61 use_gpu_elpa_=0 62 #ifdef HAVE_LINALG_ELPA 63 if (present(use_gpu_elpa)) use_gpu_elpa_=use_gpu_elpa 64 #endif 65 66 !===== SCALAPACK 67 if (ABI_LINALG_SCALAPACK_ISON.and.use_slk_==1.and.n>slk_minsize) then 68 #if defined HAVE_LINALG_SCALAPACK 69 ! if istwfk=1, then dim_evec1=2*n and if istwfk=2, dim_evec1=n 70 dim_evec1= 2*n/istwf_k_ 71 ABI_MALLOC(tmp_evec,(dim_evec1,n)) 72 tmp_evec = zero 73 call sca_a%init(n,n,slk_processor,istwf_k_) 74 call sca_ev%init(n,n,slk_processor,istwf_k_) 75 #ifdef HAVE_LINALG_ELPA 76 call matrix_from_global_sym(sca_a,a,istwf_k_) 77 #else 78 call matrix_from_global(sca_a,a,istwf_k_) 79 #endif 80 call compute_eigen_problem(slk_processor,sca_a,sca_ev,w,slk_communicator,istwf_k_,& 81 & use_gpu_elpa=use_gpu_elpa_) 82 call matrix_to_global(sca_a,a,istwf_k_) 83 call matrix_to_reference(sca_ev,tmp_evec,istwf_k_) 84 call xmpi_sum(tmp_evec,z,dim_evec1*n,slk_communicator,ierr) 85 call sca_a%free() 86 call sca_ev%free() 87 ABI_FREE(tmp_evec) 88 #endif 89 90 !===== LAPACK 91 else 92 if (istwf_k_/=2) then 93 call zhpev(jobz,uplo,n,a,w,z,ldz,eigen_z_work,eigen_z_rwork,info) 94 else 95 call dspev(jobz,uplo,n,a,w,z,ldz,eigen_d_work,info) 96 end if 97 end if 98 99 ABI_CHECK(info==0,"dhpev returned info!=0") 100 101 #ifndef HAVE_LINALG_ELPA 102 ABI_UNUSED(use_gpu_elpa) 103 #endif 104 105 end subroutine abi_dhpev
m_abi_linalg/abi_xhpev [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_xhpev
FUNCTION
abi_xhpev is the generic function that compute all eigenvalues and, optionally, eigenvectors of a symmetric or hermitian matrix A in packed storage
COPYRIGHT
Copyright (C) 2001-2024 ABINIT group (LNguyen,FDahm,MT) This file is distributed under the terms of the GNU General Public License, see ~ABINIT/Infos/copyright or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
m_abi_linalg/abi_zhpev [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_zhpev
FUNCTION
INPUTS
SOURCE
175 subroutine abi_zhpev(jobz,uplo,n,a,w,z,ldz) 176 177 !Arguments ------------------------------------ 178 character(len=1), intent(in) :: jobz 179 character(len=1), intent(in) :: uplo 180 integer, intent(in) :: n,ldz 181 complex(dpc), intent(inout) :: a(:,:) 182 complex(dpc), intent(out) :: z(:,:) 183 real(dp), intent(out) :: w(:) 184 185 !Local variables------------------------------- 186 integer :: info 187 real(dp),pointer :: rwork(:) 188 complex(dpc),pointer :: work(:) 189 190 ! ********************************************************************* 191 192 ABI_CHECK(lapack_packed_storage,"BUG(1) in abi_zhpev (storage)!") 193 ABI_CHECK(lapack_double_precision,"BUG(2) in abi_zhpev (precision)!") 194 ABI_CHECK(n<=eigen_z_maxsize,"BUG(3) in abi_zhpev (maxsize)!") 195 196 work => eigen_z_work ; rwork => eigen_z_rwork 197 198 !===== LAPACK 199 if (eigen_z_lwork==0) then 200 ABI_MALLOC(work,(2*n-1)) 201 end if 202 if (eigen_z_lrwork==0) then 203 ABI_MALLOC(rwork,(3*n-2)) 204 end if 205 call zhpev(jobz,uplo,n,a,w,z,ldz,work,rwork,info) 206 if (eigen_z_lwork==0) then 207 ABI_FREE(work) 208 end if 209 if (eigen_z_lrwork==0) then 210 ABI_FREE(rwork) 211 end if 212 213 ABI_CHECK(info==0,"abi_zhpev returned info!=0!") 214 215 end subroutine abi_zhpev