Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Utility subroutine for qs energy calculation
10 : !> \par History
11 : !> none
12 : !> \author MK (29.10.2002)
13 : ! **************************************************************************************************
14 : MODULE qs_energy_utils
15 : USE atomic_kind_types, ONLY: atomic_kind_type
16 : USE atprop_types, ONLY: atprop_array_add,&
17 : atprop_array_init,&
18 : atprop_type
19 : USE bse_bs_env_link, ONLY: bse_from_bs_env
20 : USE cp_control_types, ONLY: dft_control_type
21 : USE cp_control_utils, ONLY: read_ddapc_section
22 : USE cp_dbcsr_api, ONLY: dbcsr_copy,&
23 : dbcsr_create,&
24 : dbcsr_p_type,&
25 : dbcsr_release,&
26 : dbcsr_set
27 : USE et_coupling, ONLY: calc_et_coupling
28 : USE et_coupling_proj, ONLY: calc_et_coupling_proj
29 : USE hartree_local_methods, ONLY: Vh_1c_gg_integrals
30 : USE hartree_local_types, ONLY: ecoul_1center_type
31 : USE input_section_types, ONLY: section_vals_get,&
32 : section_vals_get_subs_vals,&
33 : section_vals_type
34 : USE kinds, ONLY: dp
35 : USE message_passing, ONLY: mp_para_env_type
36 : USE mulliken, ONLY: atom_trace
37 : USE post_scf_bandstructure_methods, ONLY: post_scf_bandstructure
38 : USE pw_env_types, ONLY: pw_env_get,&
39 : pw_env_type
40 : USE pw_methods, ONLY: pw_axpy,&
41 : pw_scale
42 : USE pw_pool_types, ONLY: pw_pool_type
43 : USE pw_types, ONLY: pw_r3d_rs_type
44 : USE qs_core_matrices, ONLY: core_matrices
45 : USE qs_dispersion_types, ONLY: qs_dispersion_type
46 : USE qs_energy_types, ONLY: qs_energy_type
47 : USE qs_environment_types, ONLY: get_qs_env,&
48 : qs_environment_type
49 : USE qs_integrate_potential, ONLY: integrate_v_core_rspace,&
50 : integrate_v_rspace
51 : USE qs_kind_types, ONLY: qs_kind_type
52 : USE qs_ks_atom, ONLY: update_ks_atom
53 : USE qs_ks_methods, ONLY: qs_ks_update_qs_env
54 : USE qs_ks_types, ONLY: qs_ks_env_type
55 : USE qs_linres_module, ONLY: linres_calculation_low
56 : USE qs_local_rho_types, ONLY: local_rho_type
57 : USE qs_rho0_ggrid, ONLY: integrate_vhg0_rspace
58 : USE qs_rho_atom_types, ONLY: rho_atom_type,&
59 : zero_rho_atom_integrals
60 : USE qs_rho_types, ONLY: qs_rho_get,&
61 : qs_rho_type
62 : USE qs_scf, ONLY: scf
63 : USE qs_tddfpt2_methods, ONLY: tddfpt
64 : USE qs_vxc, ONLY: qs_xc_density
65 : USE qs_vxc_atom, ONLY: calculate_vxc_atom
66 : USE rixs_methods, ONLY: rixs
67 : USE tip_scan_methods, ONLY: tip_scanning
68 : USE xas_methods, ONLY: xas
69 : USE xas_tdp_methods, ONLY: xas_tdp
70 : USE xc_derivatives, ONLY: xc_functionals_get_needs
71 : USE xc_rho_cflags_types, ONLY: xc_rho_cflags_type
72 : #include "./base/base_uses.f90"
73 :
74 : IMPLICIT NONE
75 :
76 : PRIVATE
77 :
78 : ! *** Global parameters ***
79 :
80 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_energy_utils'
81 :
82 : PUBLIC :: qs_energies_properties
83 :
84 : CONTAINS
85 :
86 : ! **************************************************************************************************
87 : !> \brief Refactoring of qs_energies_scf. Moves computation of properties
88 : !> into separate subroutine. Calculates energy and if calc_forces then the response vector
89 : !> and some contributions to the force
90 : !> \param qs_env ...
91 : !> \param calc_forces ...
92 : !> \par History
93 : !> 05.2013 created [Florian Schiffmann]
94 : ! **************************************************************************************************
95 :
96 144755 : SUBROUTINE qs_energies_properties(qs_env, calc_forces)
97 : TYPE(qs_environment_type), POINTER :: qs_env
98 : LOGICAL, INTENT(IN) :: calc_forces
99 :
100 : CHARACTER(len=*), PARAMETER :: routineN = 'qs_energies_properties'
101 :
102 : INTEGER :: handle, natom
103 : LOGICAL :: do_bse_props, do_et, do_et_proj, &
104 : do_post_scf_bandstructure, do_tip_scan
105 : REAL(KIND=dp) :: ekts
106 : TYPE(atprop_type), POINTER :: atprop
107 : TYPE(dft_control_type), POINTER :: dft_control
108 : TYPE(mp_para_env_type), POINTER :: para_env
109 : TYPE(pw_env_type), POINTER :: pw_env
110 : TYPE(pw_r3d_rs_type), POINTER :: v_hartree_rspace
111 : TYPE(qs_energy_type), POINTER :: energy
112 : TYPE(section_vals_type), POINTER :: bse_props_section, input, &
113 : post_scf_bands_section, proj_section, &
114 : rest_b_section, tip_scan_section
115 :
116 28951 : NULLIFY (atprop, energy, pw_env)
117 28951 : CALL timeset(routineN, handle)
118 :
119 : ! atomic energies using Mulliken partition
120 : CALL get_qs_env(qs_env, &
121 : dft_control=dft_control, &
122 : input=input, &
123 : atprop=atprop, &
124 : energy=energy, &
125 : v_hartree_rspace=v_hartree_rspace, &
126 : para_env=para_env, &
127 28951 : pw_env=pw_env)
128 28951 : IF (atprop%energy) THEN
129 140 : CALL qs_energies_mulliken(qs_env)
130 140 : CALL get_qs_env(qs_env, natom=natom)
131 : IF (.NOT. dft_control%qs_control%semi_empirical .AND. &
132 140 : .NOT. dft_control%qs_control%xtb .AND. &
133 : .NOT. dft_control%qs_control%dftb) THEN
134 : ! Nuclear charge correction
135 36 : CALL integrate_v_core_rspace(v_hartree_rspace, qs_env)
136 36 : IF (.NOT. ASSOCIATED(atprop%ateb)) THEN
137 8 : CALL atprop_array_init(atprop%ateb, natom)
138 : END IF
139 : ! Kohn-Sham Functional corrections
140 36 : CALL ks_xc_correction(qs_env)
141 : END IF
142 140 : CALL atprop_array_add(atprop%atener, atprop%ateb)
143 140 : CALL atprop_array_add(atprop%atener, atprop%ateself)
144 140 : CALL atprop_array_add(atprop%atener, atprop%atexc)
145 140 : CALL atprop_array_add(atprop%atener, atprop%atecoul)
146 140 : CALL atprop_array_add(atprop%atener, atprop%atevdw)
147 140 : CALL atprop_array_add(atprop%atener, atprop%ategcp)
148 140 : CALL atprop_array_add(atprop%atener, atprop%atecc)
149 140 : CALL atprop_array_add(atprop%atener, atprop%ate1c)
150 : ! entropic energy
151 140 : ekts = energy%kts/REAL(natom, KIND=dp)/REAL(para_env%num_pe, KIND=dp)
152 7504 : atprop%atener(:) = atprop%atener(:) + ekts
153 : END IF
154 :
155 : ! ET coupling - projection-operator approach
156 28951 : NULLIFY (proj_section)
157 : proj_section => &
158 28951 : section_vals_get_subs_vals(input, "PROPERTIES%ET_COUPLING%PROJECTION")
159 28951 : CALL section_vals_get(proj_section, explicit=do_et_proj)
160 28951 : IF (do_et_proj) THEN
161 10 : CALL calc_et_coupling_proj(qs_env)
162 : END IF
163 :
164 : ! ********** Calculate the electron transfer coupling elements********
165 28951 : do_et = .FALSE.
166 28951 : do_et = dft_control%qs_control%et_coupling_calc
167 28951 : IF (do_et) THEN
168 0 : qs_env%et_coupling%energy = energy%total
169 0 : qs_env%et_coupling%keep_matrix = .TRUE.
170 0 : qs_env%et_coupling%first_run = .TRUE.
171 0 : CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE., just_energy=.TRUE.)
172 0 : qs_env%et_coupling%first_run = .FALSE.
173 0 : IF (dft_control%qs_control%ddapc_restraint) THEN
174 0 : rest_b_section => section_vals_get_subs_vals(input, "PROPERTIES%ET_COUPLING%DDAPC_RESTRAINT_B")
175 : CALL read_ddapc_section(qs_control=dft_control%qs_control, &
176 0 : ddapc_restraint_section=rest_b_section)
177 : END IF
178 0 : CALL scf(qs_env=qs_env)
179 0 : qs_env%et_coupling%keep_matrix = .TRUE.
180 :
181 0 : CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE., just_energy=.TRUE.)
182 0 : CALL calc_et_coupling(qs_env)
183 : END IF
184 :
185 : !Properties
186 :
187 28951 : IF (dft_control%do_xas_calculation) THEN
188 42 : CALL xas(qs_env, dft_control)
189 : END IF
190 :
191 28951 : IF (dft_control%do_xas_tdp_calculation) THEN
192 52 : CALL xas_tdp(qs_env)
193 : END IF
194 :
195 : ! Compute Linear Response properties as post-scf
196 : ! Initialization of TDDFPT calculation
197 28951 : IF (.NOT. qs_env%linres_run) THEN
198 28281 : CALL linres_calculation_low(qs_env)
199 : END IF
200 :
201 : ! Calculates the energy, response vector and some contributions to the force
202 28951 : IF (dft_control%tddfpt2_control%enabled) THEN
203 1428 : CALL tddfpt(qs_env, calc_forces)
204 : END IF
205 :
206 : ! stand-alone RIXS, does not depend on previous xas_tdp and/or tddfpt2 calculations
207 28951 : IF (qs_env%do_rixs) CALL rixs(qs_env)
208 :
209 : ! post-SCF bandstructure calculation from higher level methods
210 28951 : NULLIFY (post_scf_bands_section)
211 28951 : post_scf_bands_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%BANDSTRUCTURE")
212 28951 : CALL section_vals_get(post_scf_bands_section, explicit=do_post_scf_bandstructure)
213 28951 : IF (do_post_scf_bandstructure) THEN
214 140 : CALL post_scf_bandstructure(qs_env, post_scf_bands_section)
215 : END IF
216 :
217 : ! Linear-response Bethe-Salpeter on the GW of the bandstructure section
218 28951 : NULLIFY (bse_props_section)
219 28951 : bse_props_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%BSE")
220 28951 : CALL section_vals_get(bse_props_section, explicit=do_bse_props)
221 28951 : IF (do_bse_props) CALL bse_from_bs_env(qs_env, bse_props_section)
222 :
223 : ! tip scan
224 28951 : NULLIFY (tip_scan_section)
225 28951 : tip_scan_section => section_vals_get_subs_vals(input, "PROPERTIES%TIP_SCAN")
226 28951 : CALL section_vals_get(tip_scan_section, explicit=do_tip_scan)
227 28951 : IF (do_tip_scan) THEN
228 0 : CALL tip_scanning(qs_env, tip_scan_section)
229 : END IF
230 :
231 28951 : CALL timestop(handle)
232 :
233 28951 : END SUBROUTINE qs_energies_properties
234 :
235 : ! **************************************************************************************************
236 : !> \brief Use a simple Mulliken-like energy decomposition
237 : !> \param qs_env ...
238 : !> \date 07.2011
239 : !> \author JHU
240 : !> \version 1.0
241 : ! **************************************************************************************************
242 140 : SUBROUTINE qs_energies_mulliken(qs_env)
243 :
244 : TYPE(qs_environment_type), POINTER :: qs_env
245 :
246 : INTEGER :: ispin, natom, nspin
247 140 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: atcore
248 : TYPE(atprop_type), POINTER :: atprop
249 : TYPE(dbcsr_p_type), ALLOCATABLE, DIMENSION(:), &
250 140 : TARGET :: core_mat
251 140 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_h, matrix_ks, rho_ao
252 140 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: math, matp
253 : TYPE(dft_control_type), POINTER :: dft_control
254 : TYPE(qs_rho_type), POINTER :: rho
255 :
256 140 : CALL get_qs_env(qs_env=qs_env, atprop=atprop)
257 140 : IF (atprop%energy) THEN
258 140 : CALL get_qs_env(qs_env=qs_env, matrix_ks=matrix_ks, matrix_h=matrix_h, rho=rho)
259 140 : CALL qs_rho_get(rho, rho_ao=rho_ao)
260 : ! E = 0.5*Tr(H*P+F*P)
261 7504 : atprop%atener = 0._dp
262 140 : nspin = SIZE(rho_ao)
263 284 : DO ispin = 1, nspin
264 : CALL atom_trace(matrix_h(1)%matrix, rho_ao(ispin)%matrix, &
265 144 : 0.5_dp, atprop%atener)
266 : CALL atom_trace(matrix_ks(ispin)%matrix, rho_ao(ispin)%matrix, &
267 284 : 0.5_dp, atprop%atener)
268 : END DO
269 : !
270 140 : CALL get_qs_env(qs_env=qs_env, dft_control=dft_control)
271 : IF (.NOT. dft_control%qs_control%semi_empirical .AND. &
272 140 : .NOT. dft_control%qs_control%xtb .AND. &
273 : .NOT. dft_control%qs_control%dftb) THEN
274 36 : CALL get_qs_env(qs_env=qs_env, natom=natom)
275 108 : ALLOCATE (atcore(natom))
276 36 : atcore = 0.0_dp
277 108 : ALLOCATE (core_mat(1))
278 36 : ALLOCATE (core_mat(1)%matrix)
279 36 : CALL dbcsr_create(core_mat(1)%matrix, template=matrix_h(1)%matrix)
280 36 : CALL dbcsr_copy(core_mat(1)%matrix, matrix_h(1)%matrix)
281 36 : CALL dbcsr_set(core_mat(1)%matrix, 0.0_dp)
282 36 : math(1:1, 1:1) => core_mat(1:1)
283 36 : matp(1:nspin, 1:1) => rho_ao(1:nspin)
284 36 : CALL core_matrices(qs_env, math, matp, .FALSE., 0, atcore=atcore)
285 160 : atprop%atener = atprop%atener + 0.5_dp*atcore
286 76 : DO ispin = 1, nspin
287 : CALL atom_trace(core_mat(1)%matrix, rho_ao(ispin)%matrix, &
288 76 : -0.5_dp, atprop%atener)
289 : END DO
290 36 : DEALLOCATE (atcore)
291 36 : CALL dbcsr_release(core_mat(1)%matrix)
292 36 : DEALLOCATE (core_mat(1)%matrix)
293 36 : DEALLOCATE (core_mat)
294 : END IF
295 : END IF
296 :
297 280 : END SUBROUTINE qs_energies_mulliken
298 :
299 : ! **************************************************************************************************
300 : !> \brief ...
301 : !> \param qs_env ...
302 : ! **************************************************************************************************
303 36 : SUBROUTINE ks_xc_correction(qs_env)
304 : TYPE(qs_environment_type), POINTER :: qs_env
305 :
306 : CHARACTER(len=*), PARAMETER :: routineN = 'ks_xc_correction'
307 :
308 : INTEGER :: handle, iatom, ispin, natom, nspins
309 : LOGICAL :: gapw, gapw_xc
310 : REAL(KIND=dp) :: eh1, exc1
311 36 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
312 : TYPE(atprop_type), POINTER :: atprop
313 36 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s, rho_ao, xcmat
314 36 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_p
315 : TYPE(dft_control_type), POINTER :: dft_control
316 36 : TYPE(ecoul_1center_type), DIMENSION(:), POINTER :: ecoul_1c
317 : TYPE(local_rho_type), POINTER :: local_rho_set
318 : TYPE(mp_para_env_type), POINTER :: para_env
319 : TYPE(pw_env_type), POINTER :: pw_env
320 : TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
321 : TYPE(pw_r3d_rs_type) :: xc_den
322 36 : TYPE(pw_r3d_rs_type), ALLOCATABLE, DIMENSION(:) :: vtau, vxc
323 : TYPE(pw_r3d_rs_type), POINTER :: v_hartree_rspace
324 : TYPE(qs_dispersion_type), POINTER :: dispersion_env
325 36 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
326 : TYPE(qs_ks_env_type), POINTER :: ks_env
327 : TYPE(qs_rho_type), POINTER :: rho_struct
328 36 : TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom_set
329 : TYPE(section_vals_type), POINTER :: xc_fun_section, xc_section
330 : TYPE(xc_rho_cflags_type) :: needs
331 :
332 36 : CALL timeset(routineN, handle)
333 :
334 36 : CALL get_qs_env(qs_env, ks_env=ks_env, dft_control=dft_control, pw_env=pw_env, atprop=atprop)
335 :
336 36 : IF (atprop%energy) THEN
337 :
338 36 : nspins = dft_control%nspins
339 36 : xc_section => section_vals_get_subs_vals(qs_env%input, "DFT%XC")
340 36 : xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
341 36 : needs = xc_functionals_get_needs(xc_fun_section, (nspins == 2), .TRUE.)
342 36 : gapw = dft_control%qs_control%gapw
343 36 : gapw_xc = dft_control%qs_control%gapw_xc
344 :
345 : ! Nuclear charge correction
346 36 : CALL get_qs_env(qs_env, v_hartree_rspace=v_hartree_rspace)
347 36 : IF (gapw .OR. gapw_xc) THEN
348 : CALL get_qs_env(qs_env=qs_env, local_rho_set=local_rho_set, &
349 : rho_atom_set=rho_atom_set, ecoul_1c=ecoul_1c, &
350 12 : natom=natom, para_env=para_env)
351 12 : CALL zero_rho_atom_integrals(rho_atom_set)
352 12 : CALL calculate_vxc_atom(qs_env, .FALSE., exc1)
353 12 : IF (gapw) THEN
354 8 : CALL Vh_1c_gg_integrals(qs_env, eh1, ecoul_1c, local_rho_set, para_env, tddft=.FALSE.)
355 8 : CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set, qs_kind_set=qs_kind_set)
356 : CALL integrate_vhg0_rspace(qs_env, v_hartree_rspace, para_env, calculate_forces=.FALSE., &
357 8 : local_rho_set=local_rho_set, atener=atprop%ateb)
358 : END IF
359 : END IF
360 :
361 36 : CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
362 36 : CALL auxbas_pw_pool%create_pw(xc_den)
363 148 : ALLOCATE (vxc(nspins))
364 76 : DO ispin = 1, nspins
365 76 : CALL auxbas_pw_pool%create_pw(vxc(ispin))
366 : END DO
367 36 : IF (needs%tau .OR. needs%tau_spin) THEN
368 36 : ALLOCATE (vtau(nspins))
369 24 : DO ispin = 1, nspins
370 48 : CALL auxbas_pw_pool%create_pw(vtau(ispin))
371 : END DO
372 : END IF
373 :
374 36 : IF (gapw_xc) THEN
375 4 : CALL get_qs_env(qs_env, rho_xc=rho_struct, dispersion_env=dispersion_env)
376 : ELSE
377 32 : CALL get_qs_env(qs_env, rho=rho_struct, dispersion_env=dispersion_env)
378 : END IF
379 36 : IF (needs%tau .OR. needs%tau_spin) THEN
380 : CALL qs_xc_density(ks_env, rho_struct, xc_section, dispersion_env=dispersion_env, &
381 12 : xc_den=xc_den, vxc=vxc, vtau=vtau)
382 : ELSE
383 : CALL qs_xc_density(ks_env, rho_struct, xc_section, dispersion_env=dispersion_env, &
384 24 : xc_den=xc_den, vxc=vxc)
385 : END IF
386 36 : CALL get_qs_env(qs_env, rho=rho_struct)
387 36 : CALL qs_rho_get(rho_struct, rho_ao=rho_ao)
388 36 : CALL get_qs_env(qs_env, natom=natom, matrix_s=matrix_s)
389 36 : CALL atprop_array_init(atprop%atexc, natom)
390 148 : ALLOCATE (xcmat(nspins))
391 76 : DO ispin = 1, nspins
392 40 : ALLOCATE (xcmat(ispin)%matrix)
393 40 : CALL dbcsr_create(xcmat(ispin)%matrix, template=matrix_s(1)%matrix)
394 40 : CALL dbcsr_copy(xcmat(ispin)%matrix, matrix_s(1)%matrix)
395 40 : CALL dbcsr_set(xcmat(ispin)%matrix, 0.0_dp)
396 40 : CALL pw_scale(vxc(ispin), -0.5_dp)
397 40 : CALL pw_axpy(xc_den, vxc(ispin))
398 40 : CALL pw_scale(vxc(ispin), vxc(ispin)%pw_grid%dvol)
399 : CALL integrate_v_rspace(qs_env=qs_env, v_rspace=vxc(ispin), hmat=xcmat(ispin), &
400 40 : calculate_forces=.FALSE., gapw=(gapw .OR. gapw_xc))
401 76 : IF (needs%tau .OR. needs%tau_spin) THEN
402 12 : CALL pw_scale(vtau(ispin), -0.5_dp*vtau(ispin)%pw_grid%dvol)
403 : CALL integrate_v_rspace(qs_env=qs_env, v_rspace=vtau(ispin), &
404 : hmat=xcmat(ispin), calculate_forces=.FALSE., &
405 12 : gapw=(gapw .OR. gapw_xc), compute_tau=.TRUE.)
406 : END IF
407 : END DO
408 36 : IF (gapw .OR. gapw_xc) THEN
409 : ! remove one-center potential matrix part
410 12 : CALL qs_rho_get(rho_struct, rho_ao_kp=matrix_p)
411 12 : CALL update_ks_atom(qs_env, xcmat, matrix_p, forces=.FALSE., kscale=-0.5_dp)
412 12 : CALL get_qs_env(qs_env=qs_env, rho_atom_set=rho_atom_set)
413 12 : CALL atprop_array_init(atprop%ate1c, natom)
414 48 : atprop%ate1c = 0.0_dp
415 48 : DO iatom = 1, natom
416 : atprop%ate1c(iatom) = atprop%ate1c(iatom) + &
417 48 : rho_atom_set(iatom)%exc_h - rho_atom_set(iatom)%exc_s
418 : END DO
419 12 : IF (gapw) THEN
420 8 : CALL get_qs_env(qs_env=qs_env, ecoul_1c=ecoul_1c)
421 32 : DO iatom = 1, natom
422 : atprop%ate1c(iatom) = atprop%ate1c(iatom) + &
423 : ecoul_1c(iatom)%ecoul_1_h - ecoul_1c(iatom)%ecoul_1_s + &
424 32 : ecoul_1c(iatom)%ecoul_1_z - ecoul_1c(iatom)%ecoul_1_0
425 : END DO
426 : END IF
427 : END IF
428 76 : DO ispin = 1, nspins
429 40 : CALL atom_trace(xcmat(ispin)%matrix, rho_ao(ispin)%matrix, 1.0_dp, atprop%atexc)
430 40 : CALL dbcsr_release(xcmat(ispin)%matrix)
431 76 : DEALLOCATE (xcmat(ispin)%matrix)
432 : END DO
433 36 : DEALLOCATE (xcmat)
434 :
435 36 : CALL auxbas_pw_pool%give_back_pw(xc_den)
436 76 : DO ispin = 1, nspins
437 76 : CALL auxbas_pw_pool%give_back_pw(vxc(ispin))
438 : END DO
439 36 : IF (needs%tau .OR. needs%tau_spin) THEN
440 24 : DO ispin = 1, nspins
441 48 : CALL auxbas_pw_pool%give_back_pw(vtau(ispin))
442 : END DO
443 : END IF
444 :
445 : END IF
446 :
447 36 : CALL timestop(handle)
448 :
449 72 : END SUBROUTINE ks_xc_correction
450 :
451 : END MODULE qs_energy_utils
|
