1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2019 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief g tensor calculation by dfpt 8!> Initialization of the epr_env, creation of the special neighbor lists 9!> Perturbation Hamiltonians by application of the p and rxp oprtators to psi0 10!> Write output 11!> Deallocate everything 12!> \note 13!> The psi0 should be localized 14!> the Sebastiani method works within the assumption that the orbitals are 15!> completely contained in the simulation box 16!> \par History 17!> created 07-2005 [MI] 18!> \author MI 19! ************************************************************************************************** 20MODULE qs_linres_epr_utils 21 USE atomic_kind_types, ONLY: atomic_kind_type 22 USE cell_types, ONLY: cell_type 23 USE cp_control_types, ONLY: dft_control_type 24 USE cp_log_handling, ONLY: cp_get_default_logger,& 25 cp_logger_type 26 USE cp_output_handling, ONLY: cp_print_key_finished_output,& 27 cp_print_key_unit_nr 28 USE input_section_types, ONLY: section_vals_get_subs_vals,& 29 section_vals_type 30 USE kinds, ONLY: dp 31 USE mathconstants, ONLY: fourpi,& 32 twopi 33 USE particle_types, ONLY: particle_type 34 USE physcon, ONLY: a_fine,& 35 e_gfactor 36 USE pw_env_types, ONLY: pw_env_get,& 37 pw_env_type 38 USE pw_pool_types, ONLY: pw_pool_create_pw,& 39 pw_pool_type 40 USE pw_types, ONLY: COMPLEXDATA1D,& 41 REALDATA3D,& 42 REALSPACE,& 43 RECIPROCALSPACE,& 44 pw_p_type 45 USE qs_environment_types, ONLY: get_qs_env,& 46 qs_environment_type 47 USE qs_kind_types, ONLY: qs_kind_type 48 USE qs_linres_types, ONLY: deallocate_nablavks_atom_set,& 49 epr_env_create,& 50 epr_env_type,& 51 init_nablavks_atom_set,& 52 linres_control_type,& 53 nablavks_atom_type,& 54 set_epr_env 55 USE qs_matrix_pools, ONLY: qs_matrix_pools_type 56 USE qs_mo_types, ONLY: mo_set_p_type 57 USE qs_rho_atom_types, ONLY: deallocate_rho_atom_set 58 USE qs_rho_types, ONLY: qs_rho_clear,& 59 qs_rho_create,& 60 qs_rho_set 61 USE scf_control_types, ONLY: scf_control_type 62#include "./base/base_uses.f90" 63 64 IMPLICIT NONE 65 66 PRIVATE 67 68 PUBLIC :: epr_env_cleanup, epr_env_init 69 70 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_linres_epr_utils' 71 72CONTAINS 73 74! ************************************************************************************************** 75!> \brief Initialize the epr environment 76!> \param epr_env ... 77!> \param qs_env ... 78!> \par History 79!> 07.2006 created [MI] 80!> \author MI 81! ************************************************************************************************** 82 SUBROUTINE epr_env_init(epr_env, qs_env) 83 ! 84 TYPE(epr_env_type) :: epr_env 85 TYPE(qs_environment_type), POINTER :: qs_env 86 87 CHARACTER(LEN=*), PARAMETER :: routineN = 'epr_env_init', routineP = moduleN//':'//routineN 88 89 INTEGER :: handle, i_B, idir, ispin, n_mo(2), nao, & 90 natom, nmoloc, nspins, output_unit 91 LOGICAL :: gapw 92 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set 93 TYPE(cell_type), POINTER :: cell 94 TYPE(cp_logger_type), POINTER :: logger 95 TYPE(dft_control_type), POINTER :: dft_control 96 TYPE(linres_control_type), POINTER :: linres_control 97 TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos 98 TYPE(nablavks_atom_type), DIMENSION(:), POINTER :: nablavks_atom_set 99 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set 100 TYPE(pw_env_type), POINTER :: pw_env 101 TYPE(pw_p_type), DIMENSION(:), POINTER :: rho_g, rho_r 102 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool 103 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set 104 TYPE(qs_matrix_pools_type), POINTER :: mpools 105 TYPE(scf_control_type), POINTER :: scf_control 106 TYPE(section_vals_type), POINTER :: lr_section 107 108 CALL timeset(routineN, handle) 109 110 NULLIFY (atomic_kind_set, qs_kind_set, cell, dft_control, linres_control, scf_control) 111 NULLIFY (logger, mos, mpools, particle_set) 112 NULLIFY (auxbas_pw_pool, pw_env) 113 NULLIFY (nablavks_atom_set) 114 115 n_mo(1:2) = 0 116 nao = 0 117 nmoloc = 0 118 119 logger => cp_get_default_logger() 120 !ionode = logger%para_env%mepos==logger%para_env%source 121 lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES") 122 123 output_unit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", & 124 extension=".linresLog") 125 126 IF (epr_env%ref_count /= 0) THEN 127 CALL epr_env_cleanup(epr_env) 128 END IF 129 130 IF (output_unit > 0) THEN 131 WRITE (output_unit, "(/,T20,A,/)") "*** Start EPR g tensor calculation ***" 132 WRITE (output_unit, "(T10,A,/)") "Initialization of the EPR environment" 133 ENDIF 134 135 CALL epr_env_create(epr_env) 136 137 CALL get_qs_env(qs_env=qs_env, & 138 atomic_kind_set=atomic_kind_set, & 139 qs_kind_set=qs_kind_set, & 140 cell=cell, & 141 dft_control=dft_control, & 142 linres_control=linres_control, & 143 mos=mos, & 144 mpools=mpools, & 145 particle_set=particle_set, & 146 pw_env=pw_env, & 147 scf_control=scf_control) 148 ! 149 ! Check if restat also psi0 should be restarted 150 !IF(epr_env%restart_epr .AND. scf_control%density_guess/=restart_guess)THEN 151 ! CPABORT("restart_epr requires density_guess=restart") 152 !ENDIF 153 ! 154 ! check that the psi0 are localized and you have all the centers 155 CPASSERT(linres_control%localized_psi0) 156 IF (output_unit > 0) THEN 157 WRITE (output_unit, '(A)') & 158 ' To get EPR parameters within PBC you need localized zero order orbitals ' 159 ENDIF 160 gapw = dft_control%qs_control%gapw 161 nspins = dft_control%nspins 162 natom = SIZE(particle_set, 1) 163 ! 164 ! Conversion factors 165 ! Magical constant twopi/cell%deth just like in NMR shift (basically undo scale_fac in qs_linres_nmr_current.F) 166 epr_env%g_free_factor = -1.0_dp*e_gfactor 167 epr_env%g_zke_factor = e_gfactor*(a_fine)**2 168 epr_env%g_so_factor = (a_fine)**2*(-1.0_dp*e_gfactor - 1.0_dp)/2.0_dp*twopi/cell%deth 169 epr_env%g_so_factor_gapw = (a_fine)**2*(-1.0_dp*e_gfactor - 1.0_dp)/2.0_dp 170 ! * 2 because B_ind = 2 * B_beta 171 epr_env%g_soo_factor = 2.0_dp*fourpi*(a_fine)**2*twopi/cell%deth 172 ! 2 * 2 * 1/4 * e^2 / m * a_0^2 * 2/3 * mu_0 / (omega * 1e-30 ) 173 epr_env%g_soo_chicorr_factor = 2.0/3.0_dp*fourpi*(a_fine)**2/cell%deth 174 ! 175 ! If the current density on the grid needs to be stored 176 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool) 177 ! 178 ! Initialize local current density if GAPW calculation 179 IF (gapw) THEN 180 CALL init_nablavks_atom_set(nablavks_atom_set, atomic_kind_set, qs_kind_set, nspins) 181 CALL set_epr_env(epr_env=epr_env, & 182 nablavks_atom_set=nablavks_atom_set) 183 ENDIF 184 ! 185 ! Bind 186 ALLOCATE (epr_env%bind_set(3, 3)) 187 DO i_B = 1, 3 188 DO idir = 1, 3 189 NULLIFY (epr_env%bind_set(idir, i_B)%rho, rho_r, rho_g) 190 CALL qs_rho_create(epr_env%bind_set(idir, i_B)%rho) 191 ALLOCATE (rho_r(1), rho_g(1)) 192 CALL pw_pool_create_pw(auxbas_pw_pool, rho_r(1)%pw, & 193 use_data=REALDATA3D, in_space=REALSPACE) 194 CALL pw_pool_create_pw(auxbas_pw_pool, rho_g(1)%pw, & 195 use_data=COMPLEXDATA1D, in_space=RECIPROCALSPACE) 196 CALL qs_rho_set(epr_env%bind_set(idir, i_B)%rho, rho_r=rho_r, rho_g=rho_g) 197 END DO 198 END DO 199 200 ! Nabla_V_ks 201 ALLOCATE (epr_env%nablavks_set(3, dft_control%nspins)) 202 DO idir = 1, 3 203 DO ispin = 1, nspins 204 NULLIFY (epr_env%nablavks_set(idir, ispin)%rho, rho_r, rho_g) 205 CALL qs_rho_create(epr_env%nablavks_set(idir, ispin)%rho) 206 ALLOCATE (rho_r(1), rho_g(1)) 207 CALL pw_pool_create_pw(auxbas_pw_pool, rho_r(1)%pw, & 208 use_data=REALDATA3D, in_space=REALSPACE) 209 CALL pw_pool_create_pw(auxbas_pw_pool, rho_g(1)%pw, & 210 use_data=COMPLEXDATA1D, in_space=RECIPROCALSPACE) 211 CALL qs_rho_set(epr_env%nablavks_set(idir, ispin)%rho, & 212 rho_r=rho_r, rho_g=rho_g) 213 END DO 214 END DO 215 216 ! Initialize the g tensor components 217 ALLOCATE (epr_env%g_total(3, 3)) 218 ALLOCATE (epr_env%g_so(3, 3)) 219 ALLOCATE (epr_env%g_soo(3, 3)) 220 epr_env%g_total = 0.0_dp 221 epr_env%g_zke = 0.0_dp 222 epr_env%g_so = 0.0_dp 223 epr_env%g_soo = 0.0_dp 224 225 CALL cp_print_key_finished_output(output_unit, logger, lr_section,& 226 & "PRINT%PROGRAM_RUN_INFO") 227 228 CALL timestop(handle) 229 230 END SUBROUTINE epr_env_init 231 232! ************************************************************************************************** 233!> \brief Deallocate the epr environment 234!> \param epr_env ... 235!> \par History 236!> 07.2005 created [MI] 237!> \author MI 238! ************************************************************************************************** 239 SUBROUTINE epr_env_cleanup(epr_env) 240 241 TYPE(epr_env_type) :: epr_env 242 243 CHARACTER(LEN=*), PARAMETER :: routineN = 'epr_env_cleanup', & 244 routineP = moduleN//':'//routineN 245 246 INTEGER :: i_B, idir, ispin 247 248 epr_env%ref_count = epr_env%ref_count - 1 249 IF (epr_env%ref_count == 0) THEN 250 ! nablavks_set 251 IF (ASSOCIATED(epr_env%nablavks_set)) THEN 252 DO ispin = 1, SIZE(epr_env%nablavks_set, 2) 253 DO idir = 1, SIZE(epr_env%nablavks_set, 1) 254 CALL qs_rho_clear(epr_env%nablavks_set(idir, ispin)%rho) 255 DEALLOCATE (epr_env%nablavks_set(idir, ispin)%rho) 256 ENDDO 257 ENDDO 258 DEALLOCATE (epr_env%nablavks_set) 259 END IF 260 ! nablavks_atom_set 261 IF (ASSOCIATED(epr_env%nablavks_atom_set)) THEN 262 CALL deallocate_nablavks_atom_set(epr_env%nablavks_atom_set) 263 END IF 264 ! vks_atom_set 265 IF (ASSOCIATED(epr_env%vks_atom_set)) THEN 266 CALL deallocate_rho_atom_set(epr_env%vks_atom_set) 267 END IF 268 ! bind_set 269 IF (ASSOCIATED(epr_env%bind_set)) THEN 270 DO i_B = 1, SIZE(epr_env%bind_set, 2) 271 DO idir = 1, SIZE(epr_env%bind_set, 1) 272 CALL qs_rho_clear(epr_env%bind_set(idir, i_B)%rho) 273 DEALLOCATE (epr_env%bind_set(idir, i_B)%rho) 274 ENDDO 275 ENDDO 276 DEALLOCATE (epr_env%bind_set) 277 END IF 278 ! bind_atom_set 279 IF (ASSOCIATED(epr_env%bind_atom_set)) THEN 280 DEALLOCATE (epr_env%bind_atom_set) 281 END IF 282 ! g_total 283 IF (ASSOCIATED(epr_env%g_total)) THEN 284 DEALLOCATE (epr_env%g_total) 285 END IF 286 ! g_so 287 IF (ASSOCIATED(epr_env%g_so)) THEN 288 DEALLOCATE (epr_env%g_so) 289 END IF 290 ! g_soo 291 IF (ASSOCIATED(epr_env%g_soo)) THEN 292 DEALLOCATE (epr_env%g_soo) 293 END IF 294 END IF ! ref count 295 296 END SUBROUTINE epr_env_cleanup 297 298END MODULE qs_linres_epr_utils 299