Line data Source code
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 : MODULE qs_tddfpt2_densities
9 : USE admm_types, ONLY: admm_type,&
10 : get_admm_env
11 : USE cp_control_types, ONLY: dft_control_type
12 : USE cp_dbcsr_api, ONLY: dbcsr_p_type,&
13 : dbcsr_scale
14 : USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
15 : copy_fm_to_dbcsr
16 : USE cp_fm_types, ONLY: cp_fm_get_info,&
17 : cp_fm_type
18 : USE kinds, ONLY: default_string_length,&
19 : dp
20 : USE parallel_gemm_api, ONLY: parallel_gemm
21 : USE pw_env_types, ONLY: pw_env_get
22 : USE pw_pool_types, ONLY: pw_pool_type
23 : USE pw_types, ONLY: pw_c1d_gs_type,&
24 : pw_r3d_rs_type
25 : USE qs_collocate_density, ONLY: calculate_rho_elec
26 : USE qs_environment_types, ONLY: get_qs_env,&
27 : qs_environment_type
28 : USE qs_gapw_densities, ONLY: prepare_gapw_den
29 : USE qs_ks_types, ONLY: qs_ks_env_type
30 : USE qs_local_rho_types, ONLY: local_rho_type
31 : USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
32 : USE qs_rho_atom_methods, ONLY: calculate_rho_atom_coeff
33 : USE qs_rho_methods, ONLY: qs_rho_copy,&
34 : qs_rho_update_rho
35 : USE qs_rho_types, ONLY: qs_rho_get,&
36 : qs_rho_type
37 : USE qs_tddfpt2_subgroups, ONLY: tddfpt_subgroup_env_type
38 : USE task_list_types, ONLY: task_list_type
39 : #include "./base/base_uses.f90"
40 :
41 : IMPLICIT NONE
42 :
43 : PRIVATE
44 :
45 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_tddfpt2_densities'
46 :
47 : INTEGER, PARAMETER, PRIVATE :: maxspins = 2
48 :
49 : PUBLIC :: tddfpt_construct_ground_state_orb_density, tddfpt_construct_aux_fit_density
50 :
51 : ! **************************************************************************************************
52 :
53 : CONTAINS
54 :
55 : ! **************************************************************************************************
56 : !> \brief Compute the ground-state charge density expressed in primary basis set.
57 : !> \param rho_orb_struct ground-state density in primary basis set
58 : !> \param rho_xc_struct ground-state density in primary basis set for GAPW_XC
59 : !> \param is_rks_triplets indicates that the triplet excited states calculation using
60 : !> spin-unpolarised molecular orbitals has been requested
61 : !> \param qs_env Quickstep environment
62 : !> \param sub_env parallel (sub)group environment
63 : !> \param wfm_rho_orb work dense matrix with shape [nao x nao] distributed among
64 : !> processors of the given parallel group (modified on exit)
65 : !> \par History
66 : !> * 06.2018 created by splitting the subroutine tddfpt_apply_admm_correction() in two
67 : !> subroutines tddfpt_construct_ground_state_orb_density() and
68 : !> tddfpt_construct_aux_fit_density [Sergey Chulkov]
69 : ! **************************************************************************************************
70 566 : SUBROUTINE tddfpt_construct_ground_state_orb_density(rho_orb_struct, rho_xc_struct, is_rks_triplets, &
71 : qs_env, sub_env, wfm_rho_orb)
72 : TYPE(qs_rho_type), POINTER :: rho_orb_struct, rho_xc_struct
73 : LOGICAL, INTENT(in) :: is_rks_triplets
74 : TYPE(qs_environment_type), POINTER :: qs_env
75 : TYPE(tddfpt_subgroup_env_type), INTENT(in) :: sub_env
76 : TYPE(cp_fm_type), INTENT(IN) :: wfm_rho_orb
77 :
78 : CHARACTER(LEN=*), PARAMETER :: routineN = 'tddfpt_construct_ground_state_orb_density'
79 :
80 : INTEGER :: handle, ispin, nao, nspins
81 : INTEGER, DIMENSION(maxspins) :: nmo_occ
82 566 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ij_ao
83 : TYPE(dft_control_type), POINTER :: dft_control
84 : TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
85 :
86 566 : CALL timeset(routineN, handle)
87 :
88 566 : nspins = SIZE(sub_env%mos_occ)
89 1224 : DO ispin = 1, nspins
90 1224 : CALL cp_fm_get_info(sub_env%mos_occ(ispin), nrow_global=nao, ncol_global=nmo_occ(ispin))
91 : END DO
92 :
93 566 : CALL qs_rho_get(rho_orb_struct, rho_ao=rho_ij_ao)
94 1224 : DO ispin = 1, nspins
95 : CALL parallel_gemm('N', 'T', nao, nao, nmo_occ(ispin), 1.0_dp, &
96 : sub_env%mos_occ(ispin), sub_env%mos_occ(ispin), &
97 658 : 0.0_dp, wfm_rho_orb)
98 :
99 1224 : CALL copy_fm_to_dbcsr(wfm_rho_orb, rho_ij_ao(ispin)%matrix, keep_sparsity=.TRUE.)
100 : END DO
101 :
102 : ! take into account that all MOs are doubly occupied in spin-restricted case
103 566 : IF (nspins == 1 .AND. (.NOT. is_rks_triplets)) &
104 394 : CALL dbcsr_scale(rho_ij_ao(1)%matrix, 2.0_dp)
105 :
106 566 : CALL get_qs_env(qs_env, dft_control=dft_control)
107 :
108 566 : IF (dft_control%qs_control%gapw) THEN
109 : CALL qs_rho_update_rho(rho_orb_struct, qs_env, &
110 : local_rho_set=sub_env%local_rho_set, &
111 : pw_env_external=sub_env%pw_env, &
112 : task_list_external=sub_env%task_list_orb_soft, &
113 160 : para_env_external=sub_env%para_env)
114 160 : CALL prepare_gapw_den(qs_env, local_rho_set=sub_env%local_rho_set)
115 406 : ELSEIF (dft_control%qs_control%gapw_xc) THEN
116 : CALL qs_rho_update_rho(rho_orb_struct, qs_env, &
117 : rho_xc_external=rho_xc_struct, &
118 : local_rho_set=sub_env%local_rho_set, &
119 : pw_env_external=sub_env%pw_env, &
120 : task_list_external=sub_env%task_list_orb, &
121 : task_list_external_soft=sub_env%task_list_orb_soft, &
122 32 : para_env_external=sub_env%para_env)
123 32 : CALL pw_env_get(sub_env%pw_env, auxbas_pw_pool=auxbas_pw_pool)
124 32 : CALL qs_rho_copy(rho_xc_struct, rho_orb_struct, auxbas_pw_pool, nspins)
125 32 : CALL prepare_gapw_den(qs_env, local_rho_set=sub_env%local_rho_set, do_rho0=.FALSE.)
126 : ELSE
127 : CALL qs_rho_update_rho(rho_orb_struct, qs_env, &
128 : pw_env_external=sub_env%pw_env, &
129 : task_list_external=sub_env%task_list_orb, &
130 374 : para_env_external=sub_env%para_env)
131 : END IF
132 :
133 566 : CALL timestop(handle)
134 :
135 566 : END SUBROUTINE tddfpt_construct_ground_state_orb_density
136 :
137 : ! **************************************************************************************************
138 : !> \brief Project a charge density expressed in primary basis set into the auxiliary basis set.
139 : !> \param rho_orb_struct response density in primary basis set
140 : !> \param rho_aux_fit_struct response density in auxiliary basis set (modified on exit)
141 : !> \param local_rho_set GAPW density of auxiliary basis set density
142 : !> \param qs_env Quickstep environment
143 : !> \param sub_env parallel (sub)group environment
144 : !> \param wfm_rho_orb work dense matrix with shape [nao x nao] distributed among
145 : !> processors of the given parallel group (modified on exit)
146 : !> \param wfm_rho_aux_fit work dense matrix with shape [nao_aux x nao_aux] distributed among
147 : !> processors of the given parallel group (modified on exit)
148 : !> \param wfm_aux_orb work dense matrix with shape [nao_aux x nao] distributed among
149 : !> processors of the given parallel group (modified on exit)
150 : !> \par History
151 : !> * 03.2017 the subroutine tddfpt_apply_admm_correction() was originally created by splitting
152 : !> the subroutine tddfpt_apply_hfx() in two parts [Sergey Chulkov]
153 : !> * 06.2018 created by splitting the subroutine tddfpt_apply_admm_correction() in two subroutines
154 : !> tddfpt_construct_ground_state_orb_density() and tddfpt_construct_aux_fit_density()
155 : !> in order to avoid code duplication [Sergey Chulkov]
156 : ! **************************************************************************************************
157 2028 : SUBROUTINE tddfpt_construct_aux_fit_density(rho_orb_struct, rho_aux_fit_struct, local_rho_set, &
158 : qs_env, sub_env, &
159 : wfm_rho_orb, wfm_rho_aux_fit, wfm_aux_orb)
160 : TYPE(qs_rho_type), POINTER :: rho_orb_struct, rho_aux_fit_struct
161 : TYPE(local_rho_type), POINTER :: local_rho_set
162 : TYPE(qs_environment_type), POINTER :: qs_env
163 : TYPE(tddfpt_subgroup_env_type), INTENT(in) :: sub_env
164 : TYPE(cp_fm_type), INTENT(IN) :: wfm_rho_orb, wfm_rho_aux_fit, wfm_aux_orb
165 :
166 : CHARACTER(LEN=*), PARAMETER :: routineN = 'tddfpt_construct_aux_fit_density'
167 :
168 : CHARACTER(LEN=default_string_length) :: basis_type
169 : INTEGER :: handle, ispin, nao, nao_aux, nspins
170 1014 : REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_aux_fit_r
171 : TYPE(admm_type), POINTER :: admm_env
172 1014 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao_aux_fit, rho_ao_orb
173 : TYPE(neighbor_list_set_p_type), DIMENSION(:), &
174 1014 : POINTER :: sab_aux_fit
175 1014 : TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_aux_fit_g
176 1014 : TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_aux_fit_r
177 : TYPE(qs_ks_env_type), POINTER :: ks_env
178 : TYPE(task_list_type), POINTER :: task_list
179 :
180 1014 : CALL timeset(routineN, handle)
181 :
182 1014 : CPASSERT(ASSOCIATED(sub_env%admm_A))
183 :
184 1014 : CALL get_qs_env(qs_env, ks_env=ks_env, admm_env=admm_env)
185 1014 : CALL qs_rho_get(rho_orb_struct, rho_ao=rho_ao_orb)
186 : CALL qs_rho_get(rho_aux_fit_struct, rho_ao=rho_ao_aux_fit, rho_g=rho_aux_fit_g, &
187 1014 : rho_r=rho_aux_fit_r, tot_rho_r=tot_rho_aux_fit_r)
188 :
189 1014 : nspins = SIZE(rho_ao_orb)
190 :
191 1014 : IF (admm_env%do_gapw) THEN
192 150 : basis_type = "AUX_FIT_SOFT"
193 150 : task_list => sub_env%task_list_aux_fit_soft
194 : ELSE
195 864 : basis_type = "AUX_FIT"
196 864 : task_list => sub_env%task_list_aux_fit
197 : END IF
198 :
199 1014 : CALL cp_fm_get_info(sub_env%admm_A, nrow_global=nao_aux, ncol_global=nao)
200 2032 : DO ispin = 1, nspins
201 : ! TO DO: consider sub_env%admm_A to be a DBCSR matrix
202 1018 : CALL copy_dbcsr_to_fm(rho_ao_orb(ispin)%matrix, wfm_rho_orb)
203 : CALL parallel_gemm('N', 'N', nao_aux, nao, nao, 1.0_dp, sub_env%admm_A, &
204 1018 : wfm_rho_orb, 0.0_dp, wfm_aux_orb)
205 : CALL parallel_gemm('N', 'T', nao_aux, nao_aux, nao, 1.0_dp, sub_env%admm_A, wfm_aux_orb, &
206 1018 : 0.0_dp, wfm_rho_aux_fit)
207 1018 : CALL copy_fm_to_dbcsr(wfm_rho_aux_fit, rho_ao_aux_fit(ispin)%matrix, keep_sparsity=.TRUE.)
208 :
209 : CALL calculate_rho_elec(matrix_p=rho_ao_aux_fit(ispin)%matrix, &
210 : rho=rho_aux_fit_r(ispin), rho_gspace=rho_aux_fit_g(ispin), &
211 : total_rho=tot_rho_aux_fit_r(ispin), ks_env=ks_env, &
212 : soft_valid=.FALSE., basis_type=basis_type, &
213 2032 : pw_env_external=sub_env%pw_env, task_list_external=task_list)
214 : END DO
215 1014 : IF (admm_env%do_gapw) THEN
216 150 : CALL get_admm_env(qs_env%admm_env, sab_aux_fit=sab_aux_fit)
217 : CALL calculate_rho_atom_coeff(qs_env, rho_ao_aux_fit, &
218 : rho_atom_set=local_rho_set%rho_atom_set, &
219 : qs_kind_set=admm_env%admm_gapw_env%admm_kind_set, &
220 150 : oce=admm_env%admm_gapw_env%oce, sab=sab_aux_fit, para_env=sub_env%para_env)
221 : CALL prepare_gapw_den(qs_env, local_rho_set=local_rho_set, &
222 150 : do_rho0=.FALSE., kind_set_external=admm_env%admm_gapw_env%admm_kind_set)
223 : END IF
224 :
225 1014 : CALL timestop(handle)
226 :
227 1014 : END SUBROUTINE tddfpt_construct_aux_fit_density
228 :
229 : END MODULE qs_tddfpt2_densities
|
