alineso.F90 12.2 KB
Newer Older
 Markus Betzinger committed Apr 26, 2016 1 2 3 4 5 6 7 8 ``````MODULE m_alineso USE m_juDFT !---------------------------------------------------------------------- ! Set up SO-hamiltonian for 2nd variation (hsohelp and hsoham) and ! diagonalize by lapack routines. ! Eigenvalues and vectors (eig_so and zso) are returned !---------------------------------------------------------------------- CONTAINS `````` Daniel Wortmann committed Feb 26, 2018 9 10 `````` SUBROUTINE alineso(eig_id,lapw,& mpi,DIMENSION,atoms,sym,kpts,& `````` 11 12 `````` input,noco,cell,oneD, nk, usdus,rsoc,& nsize,nmat, eig_so,zso) `````` Markus Betzinger committed Apr 26, 2016 13 14 `````` #include"cpp_double.h" `````` Gustav Bihlmayer committed Sep 20, 2018 15 `````` USE m_types `````` Markus Betzinger committed Apr 26, 2016 16 17 18 19 20 `````` USE m_hsohelp USE m_hsoham USE m_eig66_io, ONLY : read_eig IMPLICIT NONE TYPE(t_mpi),INTENT(IN) :: mpi `````` Daniel Wortmann committed Feb 26, 2018 21 `````` TYPE(t_lapw),INTENT(IN) :: lapw `````` Markus Betzinger committed Apr 26, 2016 22 23 24 25 26 27 28 `````` TYPE(t_dimension),INTENT(IN) :: DIMENSION TYPE(t_oneD),INTENT(IN) :: oneD TYPE(t_input),INTENT(IN) :: input TYPE(t_noco),INTENT(IN) :: noco TYPE(t_sym),INTENT(IN) :: sym TYPE(t_cell),INTENT(IN) :: cell TYPE(t_atoms),INTENT(IN) :: atoms `````` Daniel Wortmann committed Feb 26, 2018 29 `````` TYPE(t_kpts),INTENT(IN) :: kpts `````` Markus Betzinger committed Apr 26, 2016 30 `````` TYPE(t_usdus),INTENT(IN) :: usdus `````` 31 `````` TYPE(t_rsoc),INTENT(IN) :: rsoc `````` Markus Betzinger committed Apr 26, 2016 32 33 34 35 `````` ! .. ! .. Scalar Arguments .. INTEGER, INTENT (IN) :: eig_id INTEGER, INTENT (IN) :: nk `````` Daniel Wortmann committed May 02, 2016 36 `````` INTEGER, INTENT (OUT):: nsize,nmat `````` Markus Betzinger committed Apr 26, 2016 37 38 39 `````` ! .. ! .. Array Arguments .. COMPLEX, INTENT (OUT) :: zso(:,:,:)!(dimension%nbasfcn,2*dimension%neigd,wannierspin) `````` 40 `````` REAL, INTENT (OUT) :: eig_so(2*DIMENSION%neigd) `````` Markus Betzinger committed Apr 26, 2016 41 42 43 44 45 46 `````` !-odim !+odim ! .. ! .. Local Scalars .. REAL r2 INTEGER i,i1 ,j,jsp,jsp1,k,ne,nn,nn1,nrec,info `````` Gustav Bihlmayer committed Sep 20, 2018 47 `````` INTEGER idim_c,idim_r,jsp2,nbas,j1,ierr `````` Markus Betzinger committed Apr 26, 2016 48 `````` CHARACTER vectors `````` Gregor Michalicek committed Nov 29, 2017 49 `````` LOGICAL l_socvec,l_qsgw,l_open,l_real `````` Markus Betzinger committed Apr 26, 2016 50 `````` INTEGER irec,irecl_qsgw `````` Gustav Bihlmayer committed Sep 20, 2018 51 `````` INTEGER nat_l, extra, nat_start, nat_stop `````` Markus Betzinger committed Apr 26, 2016 52 53 54 55 `````` COMPLEX cdum ! .. ! .. Local Arrays .. INTEGER :: nsz(2) `````` 56 `````` REAL :: eig(DIMENSION%neigd,DIMENSION%jspd),s(3) `````` Markus Betzinger committed Apr 26, 2016 57 58 `````` REAL, ALLOCATABLE :: rwork(:) COMPLEX,ALLOCATABLE :: cwork(:),chelp(:,:,:,:,:) `````` Gustav Bihlmayer committed Sep 20, 2018 59 `````` COMPLEX,ALLOCATABLE :: ahelp(:,:,:,:),bhelp(:,:,:,:) `````` Markus Betzinger committed Apr 26, 2016 60 61 62 `````` COMPLEX,ALLOCATABLE :: zhelp1(:,:),zhelp2(:,:) COMPLEX,ALLOCATABLE :: hso(:,:),hsomtx(:,:,:,:) COMPLEX,ALLOCATABLE :: sigma_xc_apw(:,:),sigma_xc(:,:) `````` Daniel Wortmann committed Feb 09, 2017 63 `````` `````` Gregor Michalicek committed Jun 04, 2018 64 `````` TYPE(t_mat)::zmat(dimension%jspd) `````` Markus Betzinger committed Apr 26, 2016 65 66 67 68 69 70 71 72 73 74 75 76 `````` ! .. ! .. External Subroutines .. EXTERNAL CPP_LAPACK_cheev ! .. ! .. External Functions .. COMPLEX CPP_BLAS_cdotu,CPP_BLAS_cdotc EXTERNAL CPP_BLAS_cdotu,CPP_BLAS_cdotc ! .. ! read from eigenvalue and -vector file ! `````` Daniel Wortmann committed Nov 29, 2017 77 `````` l_real=sym%invs.and..not.noco%l_noco.and..not.(noco%l_soc.and.atoms%n_u>0) `````` Daniel Wortmann committed Feb 09, 2017 78 `````` zmat%l_real=l_real `````` Gregor Michalicek committed Jun 04, 2018 79 80 `````` zMat(1:dimension%jspd)%matsize1=lapw%nv(1:dimension%jspd)+atoms%nlotot zmat%matsize2=dimension%neigd `````` Daniel Wortmann committed Feb 09, 2017 81 `````` `````` Markus Betzinger committed Apr 26, 2016 82 83 `````` INQUIRE (4649,opened=l_socvec) INQUIRE (file='fleur.qsgw',exist=l_qsgw) `````` Daniel Wortmann committed Oct 07, 2016 84 `````` if (l_real) THEN `````` Gregor Michalicek committed Jun 04, 2018 85 86 `````` ALLOCATE (zmat(1)%data_r(zmat(1)%matsize1,DIMENSION%neigd) ) zmat(1)%data_r(:,:)= 0. `````` Daniel Wortmann committed Feb 09, 2017 87 `````` if (size(zmat)==2)THEN `````` Gregor Michalicek committed Jun 04, 2018 88 89 `````` ALLOCATE(zmat(2)%data_r(zmat(2)%matsize1,DIMENSION%neigd) ) zmat(2)%data_r=0.0 `````` Daniel Wortmann committed Feb 09, 2017 90 `````` ENDIF `````` Daniel Wortmann committed Oct 07, 2016 91 `````` else `````` Gregor Michalicek committed Jun 04, 2018 92 93 `````` ALLOCATE (zmat(1)%data_c(zmat(1)%matsize1,DIMENSION%neigd) ) zmat(1)%data_c(:,:)= 0. `````` Daniel Wortmann committed Feb 09, 2017 94 `````` if (size(zmat)==2)THEN `````` Gregor Michalicek committed Jun 04, 2018 95 96 `````` ALLOCATE(zmat(2)%data_c(zmat(2)%matsize1,DIMENSION%neigd) ) zmat(2)%data_c=0.0 `````` Daniel Wortmann committed Feb 09, 2017 97 `````` ENDIF `````` Daniel Wortmann committed Oct 07, 2016 98 `````` endif `````` Markus Betzinger committed Apr 26, 2016 99 100 101 102 `````` zso(:,:,:)= CMPLX(0.,0.) DO jsp = 1,input%jspins CALL read_eig(& `````` 103 `````` eig_id,nk,jsp, neig=ne,eig=eig(:,jsp)) `````` Daniel Wortmann committed Feb 09, 2017 104 `````` CALL read_eig(& `````` Daniel Wortmann committed Oct 07, 2016 105 106 `````` eig_id,nk,jsp,& n_start=1,n_end=ne,& `````` Daniel Wortmann committed Feb 09, 2017 107 `````` zmat=zmat(jsp)) `````` Markus Betzinger committed Apr 26, 2016 108 109 110 `````` ! write(*,*) 'process',irank,' reads ',nk `````` Daniel Wortmann committed Feb 26, 2018 111 112 113 114 115 ``````!!\$ DO i = 1, lapw%nv(1) !!\$ s = lapw%bkpt +lapw%gvec(:,i,1) !!\$ r2 = DOT_PRODUCT(s,MATMUL(s,cell%bbmat)) !!\$ lapw%rk(i,1) = SQRT(r2) !!\$ ENDDO `````` Markus Betzinger committed Apr 26, 2016 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 `````` IF (ne.GT.DIMENSION%neigd) THEN WRITE (6,'(a13,i4,a8,i4)') 'alineso: ne=',ne,' > dimension%neigd=',DIMENSION%neigd CALL juDFT_error("alineso: ne > neigd",calledby="alineso") ENDIF nsz(jsp) = ne ENDDO ! ! set up size of e.v. problem in second variation: nsize ! nsize = 0 DO jsp = 1,input%jspins IF (input%jspins.EQ.1) THEN nsize = 2*nsz(jsp) nsz(2) = nsz(1) ELSE nsize = nsize + nsz(jsp) ENDIF ENDDO ! `````` Gustav Bihlmayer committed Sep 20, 2018 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 `````` ! distribution of (abc)cof over atoms ! ! ! in case of ev-parallelization, now distribute the atoms: ! IF (mpi%n_size > 1) THEN nat_l = FLOOR(real(atoms%nat)/mpi%n_size) extra = atoms%nat - nat_l*mpi%n_size nat_start = mpi%n_rank*nat_l + 1 + extra nat_stop = (mpi%n_rank+1)*nat_l + extra IF (mpi%n_rank < extra) THEN nat_start = nat_start - (extra - mpi%n_rank) nat_stop = nat_stop - (extra - mpi%n_rank - 1) ENDIF ELSE nat_start = 1 nat_stop = atoms%nat ENDIF nat_l = nat_stop - nat_start + 1 ! `````` Markus Betzinger committed Apr 26, 2016 156 157 `````` ! set up A and B coefficients ! `````` Gustav Bihlmayer committed Sep 20, 2018 158 159 160 `````` ALLOCATE ( ahelp(atoms%lmaxd*(atoms%lmaxd+2),nat_l,DIMENSION%neigd,DIMENSION%jspd) ) ALLOCATE ( bhelp(atoms%lmaxd*(atoms%lmaxd+2),nat_l,DIMENSION%neigd,DIMENSION%jspd) ) ALLOCATE ( chelp(-atoms%llod :atoms%llod, DIMENSION%neigd,atoms%nlod,nat_l,DIMENSION%jspd) ) `````` Markus Betzinger committed Apr 26, 2016 161 `````` CALL timestart("alineso SOC: -help") `````` Gustav Bihlmayer committed Sep 20, 2018 162 `````` write(*,*) nat_start,nat_stop,nat_l `````` Markus Betzinger committed Apr 26, 2016 163 164 165 `````` CALL hsohelp(& & DIMENSION,atoms,sym,& & input,lapw,nsz,& `````` Daniel Wortmann committed Feb 26, 2018 166 `````` & cell,& `````` Daniel Wortmann committed Feb 09, 2017 167 `````` & zmat,usdus,& `````` Markus Betzinger committed Apr 26, 2016 168 `````` & zso,noco,oneD,& `````` Gustav Bihlmayer committed Sep 20, 2018 169 `````` & nat_start,nat_stop,nat_l,& `````` Markus Betzinger committed Apr 26, 2016 170 `````` & ahelp,bhelp,chelp) `````` Gustav Bihlmayer committed Sep 20, 2018 171 `````` CALL timestop("alineso SOC: -help") `````` Markus Betzinger committed Apr 26, 2016 172 173 174 `````` ! ! set up hamilton matrix ! `````` Daniel Wortmann committed Apr 29, 2016 175 `````` CALL timestart("alineso SOC: -ham") `````` Gustav Bihlmayer committed Sep 20, 2018 176 177 178 179 180 181 182 ``````#ifdef CPP_MPI CALL MPI_BARRIER(mpi%MPI_COMM,ierr) #endif ALLOCATE ( hsomtx(DIMENSION%neigd,DIMENSION%neigd,2,2) ) CALL hsoham(atoms,noco,input,nsz,chelp,rsoc,ahelp,bhelp,& nat_start,nat_stop,mpi%n_rank,mpi%n_size,mpi%SUB_COMM,& hsomtx) `````` Markus Betzinger committed Apr 26, 2016 183 184 `````` DEALLOCATE ( ahelp,bhelp,chelp ) CALL timestop("alineso SOC: -ham") `````` Gustav Bihlmayer committed Sep 20, 2018 185 `````` IF (mpi%n_rank==0) THEN `````` Markus Betzinger committed Apr 26, 2016 186 187 188 189 190 191 192 `````` ! ! add e.v. on diagonal ! ! write(*,*) '!!!!!!!!!!! remove SOC !!!!!!!!!!!!!!' ! hsomtx = 0 !!!!!!!!!!!! DO jsp = 1,input%jspins DO i = 1,nsz(jsp) `````` Gustav Bihlmayer committed Sep 20, 2018 193 `````` hsomtx(i,i,jsp,jsp) = hsomtx(i,i,jsp,jsp) +& `````` Markus Betzinger committed Apr 26, 2016 194 195 `````` & CMPLX(eig(i,jsp),0.) IF (input%jspins.EQ.1) THEN `````` Gustav Bihlmayer committed Sep 20, 2018 196 `````` hsomtx(i,i,2,2) = hsomtx(i,i,2,2) +& `````` Markus Betzinger committed Apr 26, 2016 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 `````` & CMPLX(eig(i,jsp),0.) ENDIF ENDDO ENDDO ! ! resort H-matrix ! ALLOCATE ( hso(2*DIMENSION%neigd,2*DIMENSION%neigd) ) DO jsp = 1,2 DO jsp1 = 1,2 IF (jsp.EQ.1) nn = 0 IF (jsp1.EQ.1) nn1 = 0 IF (jsp.EQ.2) nn = nsz(1) IF (jsp1.EQ.2) nn1 = nsz(1) ! `````` Gustav Bihlmayer committed Sep 20, 2018 213 `````` !write(3333,'(2i3,4e15.8)') jsp,jsp1,hsomtx(jsp,jsp1,8,8),hsomtx(jsp,jsp1,32,109) `````` Markus Betzinger committed Apr 26, 2016 214 215 `````` DO i = 1,nsz(jsp) DO j = 1,nsz(jsp1) `````` Gustav Bihlmayer committed Sep 20, 2018 216 `````` hso(i+nn,j+nn1) = hsomtx(i,j,jsp,jsp1) `````` Markus Betzinger committed Apr 26, 2016 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 `````` ENDDO ENDDO ! ENDDO ENDDO DEALLOCATE ( hsomtx ) ! ! add Sigma-vxc (QSGW) ! IF( l_qsgw ) THEN nbas = lapw%nv(1) + atoms%nlotot WRITE(*,'(A,I3,A,I5,A)') 'Read fleur.qsgw (',nk,',',nbas,')' IF( DIMENSION%jspd .EQ. 2 ) STOP 'alineso: GW+noco not implemented.' ALLOCATE ( sigma_xc(2*nsz(1),2*nsz(1)) ) ALLOCATE ( sigma_xc_apw(nbas,nbas) ) INQUIRE(667,opened=l_open) IF( .NOT.l_open ) THEN IF( nk.NE.1 ) STOP 'unit 667 not opened but not at 1st k' OPEN(667,file='fleur.qsgw',form='unformatted') ELSE IF( nk.EQ.1) THEN REWIND(667) ENDIF DO jsp1 = 1,2 DO jsp2 = 1,jsp1 IF(jsp1.EQ.jsp2) THEN READ(667) ((sigma_xc_apw(i,j),j=1,i),i=1,nbas) DO i = 1,nbas DO j = 1,i-1 sigma_xc_apw(j,i) = CONJG(sigma_xc_apw(i,j)) ENDDO ENDDO ELSE READ(667) sigma_xc_apw ENDIF ! write(*,*) 'lo part set to zero!' ! sigma_xc_apw(nv+1:,nv+1:) = 0 i = nsz(1) * (jsp1-1) + 1 ; i1 = nsz(1) * jsp1 j = nsz(1) * (jsp2-1) + 1 ; j1 = nsz(1) * jsp2 `````` Daniel Wortmann committed Oct 07, 2016 256 257 `````` if (l_real) THEN sigma_xc(i:i1,j:j1) = & `````` Gregor Michalicek committed Jun 04, 2018 258 259 `````` & MATMUL ( TRANSPOSE(zmat(1)%data_r(:nbas,:)) ,& & MATMUL ( sigma_xc_apw, zmat(1)%data_r(:nbas,:) ) ) `````` Daniel Wortmann committed Oct 07, 2016 260 ``````else `````` Markus Betzinger committed Apr 26, 2016 261 `````` sigma_xc(i:i1,j:j1) = & `````` Gregor Michalicek committed Jun 04, 2018 262 263 `````` & MATMUL ( CONJG(TRANSPOSE(zmat(1)%data_c(:nbas,:))) ,& & MATMUL ( sigma_xc_apw, zmat(1)%data_c(:nbas,:) ) ) `````` Daniel Wortmann committed Oct 07, 2016 264 `````` endif `````` Markus Betzinger committed Apr 26, 2016 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 `````` hso(i:i1,j:j1) = hso(i:i1,j:j1) + CONJG(sigma_xc(i:i1,j:j1)) IF(jsp1.NE.jsp2) THEN sigma_xc(j:j1,i:i1) = TRANSPOSE(CONJG(sigma_xc(i:i1,j:j1))) hso(j:j1,i:i1) = hso(j:j1,i:i1)+CONJG(sigma_xc(j:j1,i:i1)) ENDIF ENDDO ENDDO DEALLOCATE ( sigma_xc_apw ) ENDIF ! ! diagonalize the hamiltonian using library-routines ! idim_c = 4*DIMENSION%neigd idim_r = 6*DIMENSION%neigd CALL timestart("alineso SOC: -diag") ALLOCATE ( cwork(idim_c),rwork(idim_r) ) IF (input%eonly) THEN vectors= 'N' ELSE vectors= 'V' ENDIF CALL CPP_LAPACK_cheev(vectors,'U',nsize,& & hso,2*DIMENSION%neigd,& & eig_so,& & cwork, idim_c, rwork, & & info) IF (info.NE.0) WRITE (6,FMT=8000) info 8000 FORMAT (' AFTER CPP_LAPACK_cheev: info=',i4) CALL timestop("alineso SOC: -diag") DEALLOCATE ( cwork,rwork ) IF (input%eonly) THEN IF(l_socvec) CALL juDFT_error& & ("EONLY set. Vectors not calculated.",calledby ="alineso") ELSE ALLOCATE ( zhelp2(DIMENSION%neigd,2*DIMENSION%neigd) ) ! ! proj. back to G - space: old eigenvector 'z' to new one 'Z' ! + ! s --- s | z(G,j,s) ... z(ig,i,jsp) ! Z (G) = > z (G) * C (i,j) | Z(G,j,s) ... zso(ig,j,jsp) ! j --- i | C(i,j) ... hso(i ,j) ! i + ! reorder new e.w. in 2x2 spin space : zhelp(,1),zhelp(,2) ! DO i1 = 1,2 jsp = i1 jsp2= i1 IF (input%jspins.EQ.1) jsp = 1 `````` Gregor Michalicek committed Nov 29, 2017 321 `````` IF (input%jspins.EQ.1 .AND..NOT.(input%l_wann.OR.l_socvec)) jsp2=1 `````` Markus Betzinger committed Apr 26, 2016 322 323 324 325 326 327 328 329 330 331 `````` IF (i1.EQ.1) nn = 0 IF (i1.EQ.2) nn = nsz(1) zhelp2(:,:) = 0.d0 DO j = 1,nsize DO i = 1,nsz(jsp) zhelp2(i,j) = CONJG(hso(i+nn,j)) ENDDO ENDDO ! j `````` Daniel Wortmann committed Feb 09, 2017 332 `````` if (l_real) THEN `````` Gregor Michalicek committed Jun 04, 2018 333 334 `````` CALL CPP_BLAS_cgemm("N","N",zmat(1)%matsize1,2*dimension%neigd,dimension%neigd,CMPLX(1.d0,0.d0),CMPLX(zmat(jsp)%data_r(:,:)),& zmat(1)%matsize1, zhelp2,DIMENSION%neigd,CMPLX(0.d0,0.d0), zso(1,1,jsp2),zmat(1)%matsize1) `````` Daniel Wortmann committed Feb 09, 2017 335 `````` else `````` Gregor Michalicek committed Jun 04, 2018 336 337 `````` CALL CPP_BLAS_cgemm("N","N",zmat(1)%matsize1,2*dimension%neigd,dimension%neigd, CMPLX(1.d0,0.d0),zmat(jsp)%data_c(:,:),& zmat(1)%matsize1, zhelp2,DIMENSION%neigd,CMPLX(0.d0,0.d0), zso(:,:,jsp2),zmat(1)%matsize1) `````` Daniel Wortmann committed Feb 09, 2017 338 `````` endif `````` Markus Betzinger committed Apr 26, 2016 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 `````` ENDDO !isp IF(l_socvec) THEN !RS: write SOC vectors to SOCVEC WRITE(4649) lapw%nmat,nsize,input%jspins,nsz,2*DIMENSION%neigd,CONJG(hso) !CF: write qsgw IF(l_qsgw) THEN nn = 2*nsz(1) sigma_xc = MATMUL ( TRANSPOSE(hso(:nn,:nn)) ,& & MATMUL ( sigma_xc , CONJG(hso(:nn,:nn)) ) ) WRITE(1014) nn WRITE(1014) ((sigma_xc(i,j),i=1,j),j=1,nn) DEALLOCATE ( sigma_xc ) ENDIF ENDIF DEALLOCATE ( zhelp2 ) ENDIF ! (.NOT.input%eonly) `````` Daniel Wortmann committed Feb 09, 2017 359 `````` DEALLOCATE ( hso ) `````` Gustav Bihlmayer committed Sep 20, 2018 360 `````` ENDIF ! (n_rank==0) `````` Markus Betzinger committed Apr 26, 2016 361 `````` ! `````` Daniel Wortmann committed May 02, 2016 362 `````` nmat=lapw%nmat `````` Markus Betzinger committed Apr 26, 2016 363 364 365 366 `````` RETURN END SUBROUTINE alineso END MODULE m_alineso``````