1! 2! Copyright (C) 2001-2020 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8! 9!---------------------------------------------------------------------------- 10MODULE io_rho_xml 11 !---------------------------------------------------------------------------- 12 ! 13 USE kinds, ONLY : DP 14 USE io_files, ONLY : create_directory 15 USE io_base, ONLY : write_rhog, read_rhog 16 ! 17 PRIVATE 18 PUBLIC :: write_scf, read_scf 19 ! 20 ! {read|write}_rho: read or write the charge density 21 ! {read|write}_scf: as above, plus ldaU ns, PAW becsum, meta-GGA 22 ! 23 CONTAINS 24 25 SUBROUTINE write_scf ( rho, nspin ) 26 ! 27 USE paw_variables, ONLY : okpaw 28 USE ldaU, ONLY : lda_plus_u, hub_back, lda_plus_u_kind, nsg 29 USE funct, ONLY : dft_is_meta 30 USE noncollin_module, ONLY : noncolin 31 USE spin_orb, ONLY : domag 32 USE scf, ONLY : scf_type 33 ! 34 USE cell_base, ONLY : bg, tpiba 35 USE gvect, ONLY : ig_l2g, mill 36 USE control_flags, ONLY : gamma_only 37 USE io_files, ONLY : restart_dir 38 USE io_global, ONLY : ionode, ionode_id, stdout 39 USE mp_pools, ONLY : my_pool_id 40 USE mp_bands, ONLY : my_bgrp_id, root_bgrp_id, & 41 root_bgrp, intra_bgrp_comm 42 USE mp_images, ONLY : intra_image_comm 43 USE mp, ONLY : mp_bcast 44 45 ! 46 IMPLICIT NONE 47 TYPE(scf_type), INTENT(IN) :: rho 48 INTEGER, INTENT(IN) :: nspin 49 ! 50 CHARACTER (LEN=256) :: dirname 51 INTEGER :: nspin_, iunocc, iunpaw, ierr 52 INTEGER, EXTERNAL :: find_free_unit 53 54 dirname = restart_dir ( ) 55 CALL create_directory( dirname ) 56 ! in the following case do not read or write polarization 57 IF ( noncolin .AND. .NOT.domag ) THEN 58 nspin_ = 1 59 ELSE 60 nspin_ = nspin 61 ENDIF 62 ! Write G-space density 63 IF ( my_pool_id == 0 .AND. my_bgrp_id == root_bgrp_id ) & 64 CALL write_rhog( TRIM(dirname) // "charge-density", & 65 root_bgrp, intra_bgrp_comm, & 66 bg(:,1)*tpiba, bg(:,2)*tpiba, bg(:,3)*tpiba, & 67 gamma_only, mill, ig_l2g, rho%of_g(:,1:nspin_) ) 68 ! 69 ! Write kinetic energy density density 70 IF ( dft_is_meta() ) THEN 71 IF ( my_pool_id == 0 .AND. my_bgrp_id == root_bgrp_id ) & 72 CALL write_rhog( TRIM(dirname) // "ekin-density", & 73 root_bgrp, intra_bgrp_comm, & 74 bg(:,1)*tpiba, bg(:,2)*tpiba, bg(:,3)*tpiba, & 75 gamma_only, mill, ig_l2g, rho%kin_g(:,1:nspin_) ) 76 WRITE(stdout,'(5x,"Writing meta-gga kinetic term")') 77 ENDIF 78 79 ! Then write the other terms to separate files 80 81 IF ( lda_plus_u ) THEN 82 ! 83 iunocc = find_free_unit () 84 IF ( ionode ) THEN 85 OPEN ( UNIT=iunocc, FILE = TRIM(dirname) // 'occup.txt', & 86 FORM='formatted', STATUS='unknown' ) 87 IF (lda_plus_u_kind.EQ.0) THEN 88 WRITE( iunocc, * , iostat = ierr) rho%ns 89 IF (hub_back) WRITE( iunocc, * , iostat = ierr) rho%nsb 90 ELSEIF (lda_plus_u_kind.EQ.1) THEN 91 IF (noncolin) THEN 92 WRITE( iunocc, * , iostat = ierr) rho%ns_nc 93 ELSE 94 WRITE( iunocc, * , iostat = ierr) rho%ns 95 ENDIF 96 ELSEIF (lda_plus_u_kind.EQ.2) THEN 97 WRITE( iunocc, * , iostat = ierr) nsg 98 ! Write Hubbard_V to file 99 CALL write_V 100 ENDIF 101 ENDIF 102 CALL mp_bcast( ierr, ionode_id, intra_image_comm ) 103 IF ( ierr/=0 ) CALL errore('write_scf', 'Writing ldaU ns', 1) 104 IF ( ionode ) THEN 105 CLOSE( UNIT = iunocc, STATUS = 'KEEP' ) 106 ENDIF 107 ! 108 END IF 109 ! 110 IF ( okpaw ) THEN 111 ! 112 iunpaw = find_free_unit () 113 IF ( ionode ) THEN 114 OPEN ( UNIT=iunpaw, FILE = TRIM(dirname) // 'paw.txt', & 115 FORM='formatted', STATUS='unknown' ) 116 WRITE( iunpaw, * , iostat = ierr) rho%bec 117 END IF 118 CALL mp_bcast( ierr, ionode_id, intra_image_comm ) 119 IF ( ierr/=0 ) CALL errore('write_scf', 'Writing PAW becsum',1) 120 IF ( ionode ) THEN 121 CLOSE( UNIT = iunpaw, STATUS = 'KEEP' ) 122 ENDIF 123 ! 124 END IF 125 126 RETURN 127 END SUBROUTINE write_scf 128 129 SUBROUTINE read_scf ( rho, nspin, gamma_only ) 130 ! 131 USE scf, ONLY : scf_type 132 USE paw_variables, ONLY : okpaw 133 USE ldaU, ONLY : lda_plus_u, starting_ns, hub_back, & 134 lda_plus_u_kind, nsg 135 USE noncollin_module, ONLY : noncolin 136 USE spin_orb, ONLY : domag 137 USE gvect, ONLY : ig_l2g 138 USE funct, ONLY : dft_is_meta 139 USE io_files, ONLY : restart_dir 140 USE io_global, ONLY : ionode, ionode_id, stdout 141 USE mp_bands, ONLY : root_bgrp, intra_bgrp_comm 142 USE mp_images, ONLY : intra_image_comm 143 USE mp, ONLY : mp_bcast, mp_sum 144 ! 145 IMPLICIT NONE 146 TYPE(scf_type), INTENT(INOUT) :: rho 147 INTEGER, INTENT(IN) :: nspin 148 LOGICAL, OPTIONAL,INTENT(IN) :: gamma_only 149 ! 150 CHARACTER(LEN=256) :: dirname 151 LOGICAL :: lexist 152 INTEGER :: nspin_, iunocc, iunpaw, ierr 153 INTEGER, EXTERNAL :: find_free_unit 154 155 dirname = restart_dir ( ) 156 ! in the following case do not read or write polarization 157 IF ( noncolin .AND. .NOT.domag ) THEN 158 nspin_=1 159 ELSE 160 nspin_=nspin 161 ENDIF 162 ! read charge density 163 CALL read_rhog( TRIM(dirname) // "charge-density", & 164 root_bgrp, intra_bgrp_comm, & 165 ig_l2g, nspin_, rho%of_g, gamma_only ) 166 IF ( nspin > nspin_) rho%of_g(:,nspin_+1:nspin) = (0.0_dp, 0.0_dp) 167 ! 168 ! read kinetic energy density 169 IF ( dft_is_meta() ) THEN 170 CALL read_rhog( TRIM(dirname) // "ekin-density", & 171 root_bgrp, intra_bgrp_comm, & 172 ig_l2g, nspin_, rho%kin_g, gamma_only, ierr ) 173 IF ( ierr == 0 ) THEN 174 WRITE(stdout,'(5x,"Reading meta-gga kinetic term")') 175 ELSE 176 rho%kin_g(:,:) = (0.0_dp, 0.0_dp) 177 WRITE(stdout,'(5x,"BEWARE: kinetic-energy density file not found,",& 178 & " Kinetic-energy density set to 0")') 179 ENDIF 180 END IF 181 182 IF ( lda_plus_u ) THEN 183 ! 184 ! The occupations ns also need to be read in order to build up 185 ! the potential 186 ! 187 iunocc = find_free_unit () 188 IF ( ionode ) THEN 189 OPEN ( UNIT=iunocc, FILE = TRIM(dirname) // 'occup.txt', & 190 FORM='formatted', STATUS='old', IOSTAT=ierr ) 191 IF (lda_plus_u_kind.EQ.0) THEN 192 READ( UNIT = iunocc, FMT = *, iostat = ierr ) rho%ns 193 IF (hub_back) READ( UNIT = iunocc, FMT = * , iostat = ierr) rho%nsb 194 ELSEIF (lda_plus_u_kind.EQ.1) THEN 195 IF (noncolin) THEN 196 READ( UNIT = iunocc, FMT = *, iostat = ierr ) rho%ns_nc 197 ELSE 198 READ( UNIT = iunocc, FMT = *, iostat = ierr ) rho%ns 199 ENDIF 200 ELSEIF (lda_plus_u_kind.EQ.2) THEN 201 READ( UNIT = iunocc, FMT = * , iostat = ierr) nsg 202 ENDIF 203 ENDIF 204 ! 205 CALL mp_bcast( ierr, ionode_id, intra_image_comm ) 206 IF ( ierr/=0 ) CALL errore('read_scf', 'Reading ldaU ns', 1) 207 ! 208 IF ( ionode ) THEN 209 CLOSE( UNIT = iunocc, STATUS = 'KEEP') 210 ELSE 211 IF (lda_plus_u_kind.EQ.0) THEN 212 rho%ns(:,:,:,:) = 0.D0 213 IF (hub_back) rho%nsb(:,:,:,:) = 0.D0 214 ELSEIF (lda_plus_u_kind.EQ.1) THEN 215 IF (noncolin) THEN 216 rho%ns_nc(:,:,:,:) = 0.D0 217 ELSE 218 rho%ns(:,:,:,:) = 0.D0 219 ENDIF 220 ELSEIF (lda_plus_u_kind.EQ.2) THEN 221 nsg(:,:,:,:,:) = (0.d0, 0.d0) 222 ENDIF 223 ENDIF 224 ! 225 IF (lda_plus_u_kind.EQ.0) THEN 226 CALL mp_sum(rho%ns, intra_image_comm) 227 IF (hub_back) CALL mp_sum(rho%nsb, intra_image_comm) 228 ELSEIF (lda_plus_u_kind.EQ.1) THEN 229 IF (noncolin) THEN 230 CALL mp_sum(rho%ns_nc, intra_image_comm) 231 ELSE 232 CALL mp_sum(rho%ns, intra_image_comm) 233 ENDIF 234 ELSEIF (lda_plus_u_kind.EQ.2) THEN 235 CALL mp_sum(nsg, intra_image_comm) 236 ENDIF 237 ! 238 ! If projections on Hubbard manifold are read from file, there is no 239 ! need to set starting values: reset them to prevent any problem 240 starting_ns = -1.0_dp 241 ! 242 ENDIF 243 ! 244 IF ( okpaw ) THEN 245 ! 246 ! Also the PAW coefficients are needed: 247 ! 248 iunpaw = find_free_unit () 249 IF ( ionode ) THEN 250 OPEN ( UNIT=iunpaw, FILE = TRIM(dirname) // 'paw.txt', & 251 FORM='formatted', STATUS='old', IOSTAT=ierr ) 252 READ( UNIT = iunpaw, FMT = *, iostat=ierr ) rho%bec 253 END IF 254 CALL mp_bcast( ierr, ionode_id, intra_image_comm ) 255 IF ( ierr/=0 ) CALL errore('read_scf', 'Reading PAW becsum',1) 256 IF ( ionode ) THEN 257 CLOSE( UNIT = iunpaw, STATUS = 'KEEP') 258 ELSE 259 rho%bec(:,:,:) = 0.D0 260 END IF 261 CALL mp_sum(rho%bec, intra_image_comm) 262 ! 263 END IF 264 ! 265 RETURN 266 END SUBROUTINE read_scf 267 ! 268END MODULE io_rho_xml 269