(git:1f285aa)
qs_elf_methods.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief Does all kind of post scf calculations for GPW/GAPW
10 !> \par History
11 !> Taken out from qs_scf_post_gpw
12 !> \author JGH
13 ! **************************************************************************************************
14 MODULE qs_elf_methods
15  USE dbcsr_api, ONLY: dbcsr_p_type
16  USE kinds, ONLY: dp
17  USE mathconstants, ONLY: pi
18  USE pw_env_types, ONLY: pw_env_get,&
19  pw_env_type
20  USE pw_methods, ONLY: pw_derive,&
21  pw_transfer,&
22  pw_zero
23  USE pw_pool_types, ONLY: pw_pool_p_type,&
24  pw_pool_type
25  USE pw_types, ONLY: pw_c1d_gs_type,&
26  pw_r3d_rs_type
27  USE qs_collocate_density, ONLY: calculate_rho_elec
28  USE qs_environment_types, ONLY: get_qs_env,&
29  qs_environment_type
30  USE qs_ks_types, ONLY: qs_ks_env_type
31  USE qs_rho_types, ONLY: qs_rho_get,&
32  qs_rho_type
33 #include "./base/base_uses.f90"
34 
35  IMPLICIT NONE
36  PRIVATE
37 
38  ! Global parameters
39  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_elf_methods'
40 
41  PUBLIC :: qs_elf_calc
42 
43 ! **************************************************************************************************
44 
45 CONTAINS
46 
47 ! **************************************************************************************************
48 !> \brief ...
49 !> \param qs_env ...
50 !> \param elf_r ...
51 !> \param rho_cutoff ...
52 ! **************************************************************************************************
53  SUBROUTINE qs_elf_calc(qs_env, elf_r, rho_cutoff)
54 
55  TYPE(qs_environment_type), POINTER :: qs_env
56  TYPE(pw_r3d_rs_type), DIMENSION(:), INTENT(IN) :: elf_r
57  REAL(kind=dp), INTENT(IN) :: rho_cutoff
58 
59  CHARACTER(len=*), PARAMETER :: routinen = 'qs_elf_calc'
60  INTEGER, DIMENSION(3, 3), PARAMETER :: nd = reshape((/1, 0, 0, 0, 1, 0, 0, 0, 1/), (/3, 3/))
61  REAL(kind=dp), PARAMETER :: elfcut = 0.0001_dp, &
62  f18 = (1.0_dp/8.0_dp), &
63  f23 = (2.0_dp/3.0_dp), &
64  f53 = (5.0_dp/3.0_dp)
65 
66  INTEGER :: handle, i, idir, ispin, j, k, nspin
67  INTEGER, DIMENSION(2, 3) :: bo
68  LOGICAL :: deriv_pw, drho_r_valid, tau_r_valid
69  REAL(kind=dp) :: cfermi, elf_kernel, norm_drho, rho_53, &
70  udvol
71  TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
72  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_struct_ao
73  TYPE(pw_c1d_gs_type) :: tmp_g
74  TYPE(pw_env_type), POINTER :: pw_env
75  TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: pw_pools
76  TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
77  TYPE(pw_r3d_rs_type), DIMENSION(3) :: drho_r
78  TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_struct_r, tau_struct_r
79  TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_struct_r
80  TYPE(pw_r3d_rs_type), POINTER :: rho_r, tau_r
81  TYPE(qs_ks_env_type), POINTER :: ks_env
82  TYPE(qs_rho_type), POINTER :: rho_struct
83 
84  CALL timeset(routinen, handle)
85 
86  NULLIFY (rho_struct, rho_r, tau_r, pw_env, auxbas_pw_pool, pw_pools, ks_env)
87  NULLIFY (rho_struct_ao, rho_struct_r, tau_struct_r, drho_struct_r)
88 
89  CALL get_qs_env(qs_env, ks_env=ks_env, pw_env=pw_env, rho=rho_struct)
90 
91  CALL qs_rho_get(rho_struct, &
92  rho_ao_kp=rho_struct_ao, &
93  rho_r=rho_struct_r, &
94  tau_r=tau_struct_r, &
95  drho_r=drho_struct_r, &
96  tau_r_valid=tau_r_valid, &
97  drho_r_valid=drho_r_valid)
98 
99  CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
100  pw_pools=pw_pools)
101  nspin = SIZE(rho_struct_r)
102  bo = rho_struct_r(1)%pw_grid%bounds_local
103  cfermi = (3.0_dp/10.0_dp)*(pi*pi*3.0_dp)**f23
104 
105  ! In this case, we need a work matrix containing tau in g space
106  ! We will not have further use for it, so we will need only one
107  IF (.NOT. tau_r_valid) THEN
108  ALLOCATE (tau_r)
109  CALL auxbas_pw_pool%create_pw(tau_r)
110  END IF
111  IF (.NOT. tau_r_valid .OR. .NOT. drho_r_valid) THEN
112  CALL auxbas_pw_pool%create_pw(tmp_g)
113  END IF
114  IF (.NOT. drho_r_valid) THEN
115  DO idir = 1, 3
116  CALL auxbas_pw_pool%create_pw(drho_r(idir))
117  END DO
118  END IF
119 
120  DO ispin = 1, nspin
121  rho_r => rho_struct_r(ispin)
122  IF (tau_r_valid) THEN
123  tau_r => tau_struct_r(ispin)
124  ELSE
125  rho_ao => rho_struct_ao(ispin, :)
126  CALL pw_zero(tau_r)
127  CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
128  rho=tau_r, &
129  rho_gspace=tmp_g, &
130  ks_env=ks_env, soft_valid=.false., &
131  compute_tau=.true.)
132  END IF
133 
134  IF (drho_r_valid) THEN
135  drho_r(:) = drho_struct_r(:, ispin)
136  ELSE
137  deriv_pw = .false.
138  IF (deriv_pw) THEN
139  udvol = 1.0_dp/rho_r%pw_grid%dvol
140  DO idir = 1, 3
141  CALL pw_transfer(rho_r, tmp_g)
142  CALL pw_derive(tmp_g, nd(:, idir))
143  CALL pw_transfer(tmp_g, drho_r(idir))
144  END DO
145 
146  ELSE
147  DO idir = 1, 3
148  rho_ao => rho_struct_ao(ispin, :)
149  CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
150  rho=drho_r(idir), &
151  rho_gspace=tmp_g, &
152  ks_env=ks_env, soft_valid=.false., &
153  compute_tau=.false., compute_grad=.true., idir=idir)
154 
155  END DO
156  END IF
157  END IF
158 
159  ! Calculate elf_r
160 !$OMP PARALLEL DO DEFAULT(NONE) SHARED(bo,elf_r, ispin, drho_r,rho_r, tau_r, cfermi, rho_cutoff)&
161 !$OMP PRIVATE(k,j,i, norm_drho, rho_53, elf_kernel)
162  DO k = bo(1, 3), bo(2, 3)
163  DO j = bo(1, 2), bo(2, 2)
164  DO i = bo(1, 1), bo(2, 1)
165  norm_drho = drho_r(1)%array(i, j, k)**2 + &
166  drho_r(2)%array(i, j, k)**2 + &
167  drho_r(3)%array(i, j, k)**2
168  norm_drho = norm_drho/max(rho_r%array(i, j, k), rho_cutoff)
169  rho_53 = cfermi*max(rho_r%array(i, j, k), rho_cutoff)**f53
170  elf_kernel = (tau_r%array(i, j, k) - f18*norm_drho) + 2.87e-5_dp
171  elf_kernel = (elf_kernel/rho_53)**2
172  elf_r(ispin)%array(i, j, k) = 1.0_dp/(1.0_dp + elf_kernel)
173  IF (elf_r(ispin)%array(i, j, k) < elfcut) elf_r(ispin)%array(i, j, k) = 0.0_dp
174  END DO
175  END DO
176  END DO
177  END DO
178 
179  IF (.NOT. drho_r_valid) THEN
180  DO idir = 1, 3
181  CALL auxbas_pw_pool%give_back_pw(drho_r(idir))
182  END DO
183  END IF
184  IF (.NOT. tau_r_valid) THEN
185  CALL auxbas_pw_pool%give_back_pw(tau_r)
186  DEALLOCATE (tau_r)
187  END IF
188  IF (.NOT. tau_r_valid .OR. .NOT. drho_r_valid) THEN
189  CALL auxbas_pw_pool%give_back_pw(tmp_g)
190  END IF
191 
192  CALL timestop(handle)
193 
194  END SUBROUTINE qs_elf_calc
195 
196 END MODULE qs_elf_methods
kinds
Defines the basic variable types.
Definition: kinds.F:23
kinds::dp
integer, parameter, public dp
Definition: kinds.F:34
mathconstants
Definition of mathematical constants and functions.
Definition: mathconstants.F:16
mathconstants::pi
real(kind=dp), parameter, public pi
Definition: mathconstants.F:110
pw_env_types
container for various plainwaves related things
Definition: pw_env_types.F:14
pw_env_types::pw_env_get
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Definition: pw_env_types.F:113
pw_methods
Definition: pw_methods.F:25
pw_methods::pw_derive
subroutine, public pw_derive(pw, n)
Calculate the derivative of a plane wave vector.
Definition: pw_methods.F:10106
pw_pool_types
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition: pw_pool_types.F:24
pw_types
Definition: pw_types.F:24
qs_collocate_density
Calculate the plane wave density by collocating the primitive Gaussian functions (pgf).
Definition: qs_collocate_density.F:38
qs_collocate_density::calculate_rho_elec
subroutine, public calculate_rho_elec(matrix_p, matrix_p_kp, rho, rho_gspace, total_rho, ks_env, soft_valid, compute_tau, compute_grad, basis_type, der_type, idir, task_list_external, pw_env_external)
computes the density corresponding to a given density matrix on the grid
Definition: qs_collocate_density.F:1710
qs_elf_methods
Does all kind of post scf calculations for GPW/GAPW.
Definition: qs_elf_methods.F:14
qs_elf_methods::qs_elf_calc
subroutine, public qs_elf_calc(qs_env, elf_r, rho_cutoff)
...
Definition: qs_elf_methods.F:54
qs_environment_types
Definition: qs_environment_types.F:14
qs_environment_types::get_qs_env
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_RI_aux_kp, matrix_s, matrix_s_RI_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, WannierCentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Definition: qs_environment_types.F:502
qs_ks_types
Definition: qs_ks_types.F:15
qs_rho_types
superstucture that hold various representations of the density and keeps track of which ones are vali...
Definition: qs_rho_types.F:18
qs_rho_types::qs_rho_get
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
Definition: qs_rho_types.F:229