1!------------------------------------------------------------------------ 2PROGRAM open_grid 3 !------------------------------------------------------------------------ 4 ! 5 USE kinds, ONLY : DP 6 USE io_global, ONLY : stdout, ionode, ionode_id 7 USE mp_global, ONLY : mp_startup 8 USE mp_images, ONLY : intra_image_comm 9 USE mp_pools, ONLY : npool 10 USE mp, ONLY : mp_bcast 11 USE cell_base, ONLY : at, bg, tpiba2, alat 12 USE klist, ONLY : nks, nkstot, xk, wk, igk_k, ngk, qnorm 13 USE io_files, ONLY : prefix, tmp_dir, nwordwfc, iunwfc, diropn 14 USE noncollin_module, ONLY : noncolin, m_loc, angle1, angle2, nspin_lsda 15 USE spin_orb, ONLY : domag 16 USE control_flags, ONLY : gamma_only 17 USE environment, ONLY : environment_start, environment_end 18 USE ions_base, ONLY : nat, tau, ityp 19 USE symm_base, ONLY : nrot, nsym, s, t_rev, fft_fact, find_sym 20 USE parameters, ONLY : npk 21 USE exx_base, ONLY : nq1,nq2,nq3, xkq_collect, & 22 nkqs, exx_mp_init, index_xk, exx_grid_init 23 USE exx, ONLY : exxbuff, exxinit, use_ace, ecutfock 24 USE gvecw, ONLY : ecutwfc, gcutw 25 USE gvect, ONLY : g, ngm 26 USE funct, ONLY : dft_force_hybrid 27 USE wvfct, ONLY : nbnd, npwx, g2kin, et, wg 28 USE wavefunctions, ONLY : evc 29 USE buffers, ONLY : save_buffer, open_buffer, close_buffer 30 USE scf, ONLY : rho 31 USE lsda_mod, ONLY : nspin, isk, lsda, starting_magnetization 32 USE io_rho_xml, ONLY : write_scf 33 USE noncollin_module, ONLY : nspin_mag, npol 34 USE fft_interfaces, ONLY : fwfft 35 ! 36 USE fft_base, ONLY : dffts 37 USE control_flags, ONLY : gamma_only, io_level 38 USE start_k, ONLY : init_start_k 39 USE extfield, ONLY : gate 40 ! 41 IMPLICIT NONE 42 ! 43 CHARACTER(LEN=256), EXTERNAL :: trimcheck 44 ! 45 INTEGER :: ios, ik, ibnd, ik_idx, ik_idx_kpt, ik_idx_exx, is, na 46 CHARACTER(len=4) :: spin_component 47 CHARACTER(len=256) :: outdir 48 INTEGER :: k1, k2, k3 49 LOGICAL :: exst, opnd, exst_mem, magnetic_sym 50 REAL(DP),ALLOCATABLE :: et0(:,:), wg0(:,:), yk(:,:), wk0(:) 51 INTEGER, EXTERNAL :: n_plane_waves 52 COMPLEX(DP),ALLOCATABLE :: psic(:), evx(:,:) 53 ! 54 LOGICAL :: use_ace_back, exx_status_back 55 REAL(DP) :: ecutfock_back 56 INTEGER :: nq_back(3) 57 58 ! these are in wannier module.....-> integer :: ispinw, ikstart, ikstop, iknum 59 NAMELIST / inputpp / outdir, prefix !, nq 60 ! 61 ! initialise environment 62 ! 63#if defined(__MPI) 64 CALL mp_startup ( ) 65#endif 66 !! not sure if this should be called also in 'library' mode or not !! 67 CALL environment_start ( 'OPEN_GRID' ) 68 ! 69 ios = 0 70 IF(ionode) THEN 71 ! 72 ! Check to see if we are reading from a file 73 ! 74 CALL input_from_file() 75 ! 76 ! set default values for variables in namelist 77 ! 78 CALL get_environment_variable( 'ESPRESSO_TMPDIR', outdir ) 79 IF ( trim( outdir ) == ' ' ) outdir = './' 80 prefix = ' ' 81 !nq = 0 82 ! 83 READ (5, inputpp, iostat=ios) 84 tmp_dir = trimcheck(outdir) 85 ENDIF 86 ! 87 ! 88 CALL mp_bcast(outdir,ionode_id, intra_image_comm) 89 CALL mp_bcast(tmp_dir,ionode_id, intra_image_comm) 90 CALL mp_bcast(prefix,ionode_id, intra_image_comm) 91 ! 92 WRITE(stdout,*) 93 WRITE(stdout,*) ' Reading nscf_save data' 94 CALL read_file 95 !print*, "initial", nks, nkstot 96 CALL open_buffer(iunwfc, 'wfc', nwordwfc, io_level, exst_mem, exst) 97 ! 98 WRITE(stdout,*) 99 IF ( npool > 1 .and. nspin_mag>1) CALL errore( 'open_grid', & 100 'pools+spin not implemented, run this tool with only one pool', npool ) 101 ! 102 IF(gamma_only) CALL errore("open_grid", & 103 "not implemented, and pointless, for gamma-only",1) 104 ! 105 ! Store some variables related to exx to be put back before writing 106 use_ace_back = use_ace 107 ecutfock_back = ecutfock 108 nq_back = (/ nq1, nq2, nq3 /) 109 exx_status_back = .true. 110 CALL dft_force_hybrid(exx_status_back) 111 112 magnetic_sym = noncolin .AND. domag 113 ALLOCATE(m_loc(3,nat)) 114 IF (noncolin.and.domag) THEN 115 DO na = 1, nat 116 m_loc(1,na) = starting_magnetization(ityp(na)) * & 117 SIN( angle1(ityp(na)) ) * COS( angle2(ityp(na)) ) 118 m_loc(2,na) = starting_magnetization(ityp(na)) * & 119 SIN( angle1(ityp(na)) ) * SIN( angle2(ityp(na)) ) 120 m_loc(3,na) = starting_magnetization(ityp(na)) * & 121 COS( angle1(ityp(na)) ) 122 ENDDO 123 ENDIF 124 CALL find_sym ( nat, tau, ityp, magnetic_sym, m_loc, gate ) 125 126 nq1 = -1 127 nq2 = -1 128 nq3 = -1 129 ecutfock = 4*ecutwfc 130 use_ace = .false. 131 132 CALL exx_grid_init() 133 CALL exx_mp_init() 134 ! 135 CALL exxinit(.false.) 136 qnorm = 0._dp 137 ! 138 CALL close_buffer(iunwfc,'keep') 139 ! 140 ALLOCATE(et0(nbnd,nks), wg0(nbnd,nks), wk0(nks)) 141 et0 = et 142 wg0 = wg 143 wk0 = wk(1:nks) 144 ! 145 nks = nkqs 146 nkstot = nks 147 ! Temporary order of xk points, they will be resorted later (if nspin>1) 148 xk(:,1:nks) = xkq_collect(:,1:nks) 149 wk(1:nks) = 1._dp/DBLE(nks) !/DBLE(nspin_mag) 150 npwx = n_plane_waves(gcutw, nks, xk, g, ngm) 151 IF (nspin==2) THEN 152 isk(1:nks/2) = 1 153 isk(nks/2+1:nks) = 2 154 ENDIF 155 ! 156 DEALLOCATE(igk_k, ngk, et, wg, g2kin) 157 ALLOCATE(igk_k(npwx,nks), ngk(nks), g2kin(npwx)) 158 ALLOCATE(et(nbnd,nks), wg(nbnd,nks)) 159 ! 160 DEALLOCATE(evc) 161 ALLOCATE(evc(npwx*npol,nbnd)) 162 ! 163 prefix = TRIM(prefix)//"_open" 164 nwordwfc = nbnd * npwx * npol 165 CALL open_buffer(iunwfc, 'wfc', nwordwfc, +1, exst_mem, exst) 166 ! 167 ! Write everything again with the new prefix 168 CALL write_scf(rho, nspin) 169 ! 170 ALLOCATE(psic(dffts%nnr), evx(npol*npwx, nbnd)) 171 172 DO ik = 1, nks 173 ik_idx_kpt = ik !+ (is-1)*(nks/nspin_mag) !(ik-1)*nspin_mag + is 174 ik_idx_exx = ik !+ (is-1)*(nks/nspin_mag) 175 xk(:,ik_idx_kpt) = xkq_collect(:,ik_idx_exx) 176! wk(ik_idx_kpt) = 1._dp/nkqs*DBLE(nspin_mag) 177 ! 178 CALL gk_sort (xk(:,ik_idx_kpt), ngm, g, ecutwfc / tpiba2, & 179 ngk(ik_idx_kpt), igk_k(:,ik_idx_kpt), g2kin) 180 evx(:,:) = 0._dp 181 DO ibnd = 1, nbnd 182 psic(1:dffts%nnr) = exxbuff(1:dffts%nnr,ibnd,ik_idx_exx) 183 CALL fwfft('Wave', psic, dffts) 184 evx(1:ngk(ik_idx_kpt),ibnd) = psic(dffts%nl(igk_k(1:ngk(ik_idx_kpt),ik_idx_kpt))) 185 ! 186 IF(noncolin)THEN 187 psic(1:dffts%nnr) = exxbuff(dffts%nnr+1:2*dffts%nnr,ibnd,ik_idx_exx) 188 CALL fwfft('Wave', psic, dffts) 189 evx(npwx+1:npwx+ngk(ik_idx_kpt),ibnd) = & 190 psic(dffts%nl(igk_k(1:ngk(ik_idx_kpt),ik_idx_kpt))) 191 ENDIF 192 ! 193 ENDDO 194 ! 195 CALL save_buffer( evx, nwordwfc, iunwfc, ik_idx_kpt ) 196 et(:,ik_idx_kpt) = et0(:,index_xk(ik_idx_kpt)) 197 wg(:,ik_idx_kpt) = wg0(:,index_xk(ik_idx_kpt))/wk0(index_xk(ik_idx_kpt))*wk(ik_idx_kpt) 198 !ENDDO 199 ENDDO 200 DEALLOCATE(psic, et0, wg0) 201 ! 202 ! reconstruct input variables 203 k1 = 0 204 k2 = 0 205 k3 = 0 206 CALL init_start_k(nq1,nq2,nq3, k1, k2, k3, "automatic",nks/nspin_lsda, xk, wk) 207 ! 208 ! Restore EXX variables 209 use_ace = use_ace_back 210 ecutfock = ecutfock_back 211 nq1 = nq_back(1) 212 nq2 = nq_back(2) 213 nq3 = nq_back(3) 214 CALL dft_force_hybrid(exx_status_back) 215 ! 216 CALL punch('all') 217 ! 218 ALLOCATE(yk(3,nks)) 219 yk = xk 220 CALL cryst_to_cart(nks, yk, at, -1) 221 WRITE(stdout,'(5x,a)') "Grid of q-points" 222 WRITE(stdout,'(5x,a,3i4)') "Dimensions:", nq1, nq2, nq3 223 WRITE(stdout,'(5x,a,3i4)') "Shift: ", k1,k2,k3 224 WRITE(stdout,'(5x,a)') "List to be put in the .win file of wannier90: & 225 &(already in crystal/fractionary coordinates):" 226 227 DO ik = 1, nks/nspin_lsda 228 WRITE(stdout,'(3f21.15,3x,f13.10)') yk(:,ik), wk(ik) 229 ENDDO 230 DEALLOCATE(yk) 231 ! 232 CALL environment_end ( 'OPEN_GRID' ) 233 WRITE( stdout, * ) 234 CALL stop_pp 235 ! 236END PROGRAM open_grid 237! 238