(git:b17b328)
Loading...
Searching...
No Matches
qs_tddfpt2_operators.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
9 USE admm_types, ONLY: admm_type
10 USE cell_types, ONLY: cell_type,&
11 pbc
12 USE cp_dbcsr_api, ONLY: &
15 dbcsr_release, dbcsr_set, dbcsr_type, dbcsr_type_no_symmetry
22 USE cp_fm_types, ONLY: cp_fm_create,&
30 USE hfx_types, ONLY: hfx_type
34 USE kinds, ONLY: dp
38 USE pw_env_types, ONLY: pw_env_get,&
40 USE pw_methods, ONLY: pw_axpy,&
42 pw_scale,&
47 USE pw_pool_types, ONLY: pw_pool_p_type,&
49 USE pw_types, ONLY: pw_c1d_gs_type,&
56 USE qs_rho_types, ONLY: qs_rho_get,&
68#include "./base/base_uses.f90"
69
70 IMPLICIT NONE
71
72 PRIVATE
73
74 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_tddfpt2_operators'
75
76 LOGICAL, PARAMETER, PRIVATE :: debug_this_module = .false.
77 ! number of first derivative components (3: d/dx, d/dy, d/dz)
78 INTEGER, PARAMETER, PRIVATE :: nderivs = 3
79 INTEGER, PARAMETER, PRIVATE :: maxspins = 2
80
83
84! **************************************************************************************************
85
86CONTAINS
87
88! **************************************************************************************************
89!> \brief Apply orbital energy difference term:
90!> Aop_evects(spin,state) += KS(spin) * evects(spin,state) -
91!> S * evects(spin,state) * diag(evals_occ(spin))
92!> \param Aop_evects action of TDDFPT operator on trial vectors (modified on exit)
93!> \param evects trial vectors C_{1,i}
94!> \param S_evects S * C_{1,i}
95!> \param gs_mos molecular orbitals optimised for the ground state (only occupied orbital
96!> energies [component %evals_occ] are needed)
97!> \param matrix_ks Kohn-Sham matrix
98!> \par History
99!> * 05.2016 initialise all matrix elements in one go [Sergey Chulkov]
100!> * 03.2017 renamed from tddfpt_init_energy_diff(), altered prototype [Sergey Chulkov]
101!> \note Based on the subroutine p_op_l1() which was originally created by
102!> Thomas Chassaing on 08.2002.
103! **************************************************************************************************
104 SUBROUTINE tddfpt_apply_energy_diff(Aop_evects, evects, S_evects, gs_mos, matrix_ks)
105 TYPE(cp_fm_type), DIMENSION(:, :), INTENT(INOUT) :: aop_evects
106 TYPE(cp_fm_type), DIMENSION(:, :), INTENT(IN) :: evects, s_evects
107 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
108 INTENT(in) :: gs_mos
109 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(in) :: matrix_ks
110
111 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_apply_energy_diff'
112
113 INTEGER :: handle, ispin, ivect, nactive, nao, &
114 nspins, nvects
115 TYPE(cp_fm_struct_type), POINTER :: matrix_struct
116 TYPE(cp_fm_type) :: hevec
117
118 CALL timeset(routinen, handle)
119
120 nspins = SIZE(evects, 1)
121 nvects = SIZE(evects, 2)
122
123 DO ispin = 1, nspins
124 CALL cp_fm_get_info(matrix=evects(ispin, 1), matrix_struct=matrix_struct, &
125 nrow_global=nao, ncol_global=nactive)
126 CALL cp_fm_create(hevec, matrix_struct)
127
128 DO ivect = 1, nvects
129 CALL cp_dbcsr_sm_fm_multiply(matrix_ks(ispin)%matrix, evects(ispin, ivect), &
130 aop_evects(ispin, ivect), ncol=nactive, &
131 alpha=1.0_dp, beta=1.0_dp)
132
133 IF (ASSOCIATED(gs_mos(ispin)%evals_occ_matrix)) THEN
134 ! orbital energy correction: evals_occ_matrix is not a diagonal matrix
135 CALL parallel_gemm('N', 'N', nao, nactive, nactive, 1.0_dp, &
136 s_evects(ispin, ivect), gs_mos(ispin)%evals_occ_matrix, &
137 0.0_dp, hevec)
138 ELSE
139 CALL cp_fm_to_fm(s_evects(ispin, ivect), hevec)
140 CALL cp_fm_column_scale(hevec, gs_mos(ispin)%evals_occ)
141 END IF
142
143 ! KS * C1 - S * C1 * occupied_orbital_energies
144 CALL cp_fm_scale_and_add(1.0_dp, aop_evects(ispin, ivect), -1.0_dp, hevec)
145 END DO
146 CALL cp_fm_release(hevec)
147 END DO
148
149 CALL timestop(handle)
150
151 END SUBROUTINE tddfpt_apply_energy_diff
152
153! **************************************************************************************************
154!> \brief Update v_rspace by adding coulomb term.
155!> \param A_ia_rspace action of TDDFPT operator on the trial vector expressed in a plane wave
156!> representation (modified on exit)
157!> \param rho_ia_g response density in reciprocal space for the given trial vector
158!> \param local_rho_set ...
159!> \param hartree_local ...
160!> \param qs_env ...
161!> \param sub_env the full sub_environment needed for calculation
162!> \param gapw Flag indicating GAPW cacluation
163!> \param work_v_gspace work reciprocal-space grid to store Coulomb potential (modified on exit)
164!> \param work_v_rspace work real-space grid to store Coulomb potential (modified on exit)
165!> \param tddfpt_mgrid ...
166!> \par History
167!> * 05.2016 compute all coulomb terms in one go [Sergey Chulkov]
168!> * 03.2017 proceed excited states sequentially; minimise the number of conversions between
169!> DBCSR and FM matrices [Sergey Chulkov]
170!> * 06.2018 return the action expressed in the plane wave representation instead of the one
171!> in the atomic basis set representation
172!> \note Based on the subroutine kpp1_calc_k_p_p1() which was originally created by
173!> Mohamed Fawzi on 10.2002.
174! **************************************************************************************************
175 SUBROUTINE tddfpt_apply_coulomb(A_ia_rspace, rho_ia_g, local_rho_set, hartree_local, &
176 qs_env, sub_env, gapw, work_v_gspace, work_v_rspace, tddfpt_mgrid)
177 TYPE(pw_r3d_rs_type), DIMENSION(:), INTENT(INOUT) :: a_ia_rspace
178 TYPE(pw_c1d_gs_type), INTENT(INOUT) :: rho_ia_g
179 TYPE(local_rho_type), POINTER :: local_rho_set
180 TYPE(hartree_local_type), POINTER :: hartree_local
181 TYPE(qs_environment_type), POINTER :: qs_env
182 TYPE(tddfpt_subgroup_env_type), INTENT(in) :: sub_env
183 LOGICAL, INTENT(IN) :: gapw
184 TYPE(pw_c1d_gs_type), INTENT(INOUT) :: work_v_gspace
185 TYPE(pw_r3d_rs_type), INTENT(INOUT) :: work_v_rspace
186 LOGICAL, INTENT(IN) :: tddfpt_mgrid
187
188 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_apply_coulomb'
189
190 INTEGER :: handle, ispin, nspins
191 REAL(kind=dp) :: alpha, pair_energy
192 TYPE(pw_env_type), POINTER :: pw_env
193 TYPE(pw_poisson_type), POINTER :: poisson_env
194 TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: my_pools
195 TYPE(realspace_grid_desc_p_type), DIMENSION(:), &
196 POINTER :: my_rs_descs
197 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: my_rs_grids
198
199 CALL timeset(routinen, handle)
200
201 nspins = SIZE(a_ia_rspace)
202 pw_env => sub_env%pw_env
203 IF (tddfpt_mgrid) THEN
204 CALL pw_env_get(pw_env, poisson_env=poisson_env, rs_grids=my_rs_grids, &
205 rs_descs=my_rs_descs, pw_pools=my_pools)
206 ELSE
207 CALL pw_env_get(pw_env, poisson_env=poisson_env)
208 END IF
209
210 IF (nspins > 1) THEN
211 alpha = 1.0_dp
212 ELSE
213 ! spin-restricted case: alpha == 2 due to singlet state.
214 ! In case of triplet states alpha == 0, so we should not call this subroutine at all.
215 alpha = 2.0_dp
216 END IF
217
218 IF (gapw) THEN
219 cpassert(ASSOCIATED(local_rho_set))
220 CALL pw_axpy(local_rho_set%rho0_mpole%rho0_s_gs, rho_ia_g)
221 IF (ASSOCIATED(local_rho_set%rho0_mpole%rhoz_cneo_s_gs)) THEN
222 CALL pw_axpy(local_rho_set%rho0_mpole%rhoz_cneo_s_gs, rho_ia_g)
223 END IF
224 END IF
225
226 CALL pw_poisson_solve(poisson_env, rho_ia_g, pair_energy, work_v_gspace)
227 CALL pw_transfer(work_v_gspace, work_v_rspace)
228
229 ! (i a || j b) = ( i_alpha a_alpha + i_beta a_beta || j_alpha b_alpha + j_beta b_beta) =
230 ! tr (Cj_alpha^T * [J_i{alpha}a{alpha}_munu + J_i{beta}a{beta}_munu] * Cb_alpha) +
231 ! tr (Cj_beta^T * [J_i{alpha}a{alpha}_munu + J_i{beta}a{beta}_munu] * Cb_beta)
232 DO ispin = 1, nspins
233 CALL pw_axpy(work_v_rspace, a_ia_rspace(ispin), alpha)
234 END DO
235
236 IF (gapw) THEN
237 CALL vh_1c_gg_integrals(qs_env, pair_energy, &
238 hartree_local%ecoul_1c, &
239 local_rho_set, &
240 sub_env%para_env, tddft=.true., core_2nd=.true.)
241 CALL pw_scale(work_v_rspace, work_v_rspace%pw_grid%dvol)
242 IF (tddfpt_mgrid) THEN
243 CALL integrate_vhg0_rspace(qs_env, work_v_rspace, sub_env%para_env, &
244 calculate_forces=.false., &
245 local_rho_set=local_rho_set, my_pools=my_pools, &
246 my_rs_descs=my_rs_descs)
247 ELSE
248 CALL integrate_vhg0_rspace(qs_env, work_v_rspace, sub_env%para_env, &
249 calculate_forces=.false., &
250 local_rho_set=local_rho_set)
251 END IF
252 END IF
253
254 CALL timestop(handle)
255
256 END SUBROUTINE tddfpt_apply_coulomb
257
258! **************************************************************************************************
259!> \brief Driver routine for applying fxc (analyic vs. finite difference for testing
260!> \param A_ia_rspace action of TDDFPT operator on trial vectors expressed in a plane wave
261!> representation (modified on exit)
262!> \param kernel_env kernel environment
263!> \param rho_ia_struct response density for the given trial vector
264!> \param is_rks_triplets indicates that the triplet excited states calculation using
265!> spin-unpolarised molecular orbitals has been requested
266!> \param pw_env plain wave environment
267!> \param work_v_xc work real-space grid to store the gradient of the exchange-correlation
268!> potential with respect to the response density (modified on exit)
269!> \param work_v_xc_tau ...
270! **************************************************************************************************
271 SUBROUTINE tddfpt_apply_xc(A_ia_rspace, kernel_env, rho_ia_struct, is_rks_triplets, &
272 pw_env, work_v_xc, work_v_xc_tau)
273
274 TYPE(pw_r3d_rs_type), DIMENSION(:), INTENT(INOUT) :: a_ia_rspace
275 TYPE(full_kernel_env_type), INTENT(IN) :: kernel_env
276 TYPE(qs_rho_type), POINTER :: rho_ia_struct
277 LOGICAL, INTENT(in) :: is_rks_triplets
278 TYPE(pw_env_type), POINTER :: pw_env
279 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: work_v_xc, work_v_xc_tau
280
281 INTEGER :: ispin, nspins
282
283 nspins = SIZE(a_ia_rspace)
284
285 IF (kernel_env%deriv2_analytic) THEN
286 CALL tddfpt_apply_xc_analytic(kernel_env, rho_ia_struct, is_rks_triplets, nspins, &
287 pw_env, work_v_xc, work_v_xc_tau)
288 ELSE
289 CALL tddfpt_apply_xc_fd(kernel_env, rho_ia_struct, is_rks_triplets, nspins, &
290 pw_env, work_v_xc, work_v_xc_tau)
291 END IF
292
293 DO ispin = 1, nspins
294 ! pw2 = pw2 + alpha * pw1
295 CALL pw_axpy(work_v_xc(ispin), a_ia_rspace(ispin), kernel_env%alpha)
296 END DO
297
298 END SUBROUTINE tddfpt_apply_xc
299
300! **************************************************************************************************
301!> \brief Routine for applying fxc potential
302!> \param A_ia_rspace action of TDDFPT operator on trial vectors expressed in a plane wave
303!> representation (modified on exit)
304!> \param fxc_rspace ...
305!> \param rho_ia_struct response density for the given trial vector
306!> \param is_rks_triplets ...
307! **************************************************************************************************
308 SUBROUTINE tddfpt_apply_xc_potential(A_ia_rspace, fxc_rspace, rho_ia_struct, is_rks_triplets)
309
310 TYPE(pw_r3d_rs_type), DIMENSION(:), INTENT(INOUT) :: a_ia_rspace
311 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: fxc_rspace
312 TYPE(qs_rho_type), POINTER :: rho_ia_struct
313 LOGICAL, INTENT(in) :: is_rks_triplets
314
315 INTEGER :: nspins
316 REAL(kind=dp) :: alpha
317 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho1_r
318
319 nspins = SIZE(a_ia_rspace)
320
321 alpha = 1.0_dp
322
323 CALL qs_rho_get(rho_ia_struct, rho_r=rho1_r)
324
325 IF (nspins == 2) THEN
326 CALL pw_multiply(a_ia_rspace(1), fxc_rspace(1), rho1_r(1), alpha)
327 CALL pw_multiply(a_ia_rspace(1), fxc_rspace(2), rho1_r(2), alpha)
328 CALL pw_multiply(a_ia_rspace(2), fxc_rspace(3), rho1_r(2), alpha)
329 CALL pw_multiply(a_ia_rspace(2), fxc_rspace(2), rho1_r(1), alpha)
330 ELSE IF (is_rks_triplets) THEN
331 CALL pw_multiply(a_ia_rspace(1), fxc_rspace(1), rho1_r(1), alpha)
332 CALL pw_multiply(a_ia_rspace(1), fxc_rspace(2), rho1_r(1), -alpha)
333 ELSE
334 CALL pw_multiply(a_ia_rspace(1), fxc_rspace(1), rho1_r(1), alpha)
335 CALL pw_multiply(a_ia_rspace(1), fxc_rspace(2), rho1_r(1), alpha)
336 END IF
337
338 END SUBROUTINE tddfpt_apply_xc_potential
339
340! **************************************************************************************************
341!> \brief Update A_ia_munu by adding exchange-correlation term.
342!> \param kernel_env kernel environment
343!> \param rho_ia_struct response density for the given trial vector
344!> \param is_rks_triplets indicates that the triplet excited states calculation using
345!> spin-unpolarised molecular orbitals has been requested
346!> \param nspins ...
347!> \param pw_env plain wave environment
348!> \param work_v_xc work real-space grid to store the gradient of the exchange-correlation
349!> potential with respect to the response density (modified on exit)
350!> \param work_v_xc_tau ...
351!> \par History
352!> * 05.2016 compute all kernel terms in one go [Sergey Chulkov]
353!> * 03.2017 proceed excited states sequentially; minimise the number of conversions between
354!> DBCSR and FM matrices [Sergey Chulkov]
355!> * 06.2018 return the action expressed in the plane wave representation instead of the one
356!> in the atomic basis set representation
357!> \note Based on the subroutine kpp1_calc_k_p_p1() which was originally created by
358!> Mohamed Fawzi on 10.2002.
359! **************************************************************************************************
360 SUBROUTINE tddfpt_apply_xc_analytic(kernel_env, rho_ia_struct, is_rks_triplets, nspins, &
361 pw_env, work_v_xc, work_v_xc_tau)
362 TYPE(full_kernel_env_type), INTENT(in) :: kernel_env
363 TYPE(qs_rho_type), POINTER :: rho_ia_struct
364 LOGICAL, INTENT(in) :: is_rks_triplets
365 INTEGER, INTENT(in) :: nspins
366 TYPE(pw_env_type), POINTER :: pw_env
367 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: work_v_xc, work_v_xc_tau
368
369 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_apply_xc_analytic'
370
371 INTEGER :: handle, ispin
372 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_ia_g, rho_ia_g2
373 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
374 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_ia_r, rho_ia_r2, tau_ia_r, tau_ia_r2
375
376 CALL timeset(routinen, handle)
377
378 CALL qs_rho_get(rho_ia_struct, rho_g=rho_ia_g, rho_r=rho_ia_r, tau_r=tau_ia_r)
379 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
380
381 IF (debug_this_module) THEN
382 cpassert(SIZE(rho_ia_g) == nspins)
383 cpassert(SIZE(rho_ia_r) == nspins)
384 cpassert((.NOT. ASSOCIATED(tau_ia_r)) .OR. SIZE(tau_ia_r) == nspins)
385 cpassert((.NOT. is_rks_triplets) .OR. nspins == 1)
386 END IF
387
388 NULLIFY (tau_ia_r2)
389 IF (is_rks_triplets) THEN
390 ALLOCATE (rho_ia_r2(2))
391 ALLOCATE (rho_ia_g2(2))
392 rho_ia_r2(1) = rho_ia_r(1)
393 rho_ia_r2(2) = rho_ia_r(1)
394 rho_ia_g2(1) = rho_ia_g(1)
395 rho_ia_g2(2) = rho_ia_g(1)
396
397 IF (ASSOCIATED(tau_ia_r)) THEN
398 ALLOCATE (tau_ia_r2(2))
399 tau_ia_r2(1) = tau_ia_r(1)
400 tau_ia_r2(2) = tau_ia_r(1)
401 END IF
402 ELSE
403 rho_ia_r2 => rho_ia_r
404 rho_ia_g2 => rho_ia_g
405
406 tau_ia_r2 => tau_ia_r
407 END IF
408
409 DO ispin = 1, nspins
410 CALL pw_zero(work_v_xc(ispin))
411 IF (ASSOCIATED(work_v_xc_tau)) CALL pw_zero(work_v_xc_tau(ispin))
412 END DO
413
414 CALL xc_rho_set_update(rho_set=kernel_env%xc_rho1_set, rho_r=rho_ia_r2, rho_g=rho_ia_g2, tau=tau_ia_r2, &
415 needs=kernel_env%xc_rho1_cflags, xc_deriv_method_id=kernel_env%deriv_method_id, &
416 xc_rho_smooth_id=kernel_env%rho_smooth_id, pw_pool=auxbas_pw_pool)
417
418 CALL xc_calc_2nd_deriv_analytical(v_xc=work_v_xc, v_xc_tau=work_v_xc_tau, deriv_set=kernel_env%xc_deriv_set, &
419 rho_set=kernel_env%xc_rho_set, &
420 rho1_set=kernel_env%xc_rho1_set, pw_pool=auxbas_pw_pool, &
421 xc_section=kernel_env%xc_section, gapw=.false., tddfpt_fac=kernel_env%beta)
422
423 IF (is_rks_triplets) THEN
424 DEALLOCATE (rho_ia_r2)
425 DEALLOCATE (rho_ia_g2)
426 IF (ASSOCIATED(tau_ia_r2)) DEALLOCATE (tau_ia_r2)
427 END IF
428
429 CALL timestop(handle)
430
431 END SUBROUTINE tddfpt_apply_xc_analytic
432
433! **************************************************************************************************
434!> \brief Update A_ia_munu by adding exchange-correlation term using finite difference methods.
435!> \param kernel_env kernel environment
436!> \param rho_ia_struct response density for the given trial vector
437!> \param is_rks_triplets indicates that the triplet excited states calculation using
438!> spin-unpolarised molecular orbitals has been requested
439!> \param nspins ...
440!> \param pw_env plain wave environment
441!> \param work_v_xc work real-space grid to store the gradient of the exchange-correlation
442!> potential with respect to the response density (modified on exit)
443!> \param work_v_xc_tau ...
444! **************************************************************************************************
445 SUBROUTINE tddfpt_apply_xc_fd(kernel_env, rho_ia_struct, is_rks_triplets, nspins, &
446 pw_env, work_v_xc, work_v_xc_tau)
447 TYPE(full_kernel_env_type), INTENT(in) :: kernel_env
448 TYPE(qs_rho_type), POINTER :: rho_ia_struct
449 LOGICAL, INTENT(in) :: is_rks_triplets
450 INTEGER, INTENT(in) :: nspins
451 TYPE(pw_env_type), POINTER :: pw_env
452 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: work_v_xc, work_v_xc_tau
453
454 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_apply_xc_fd'
455
456 INTEGER :: handle, ispin
457 LOGICAL :: lsd, singlet, triplet
458 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho1_g
459 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
460 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho1_r, tau1_r
461 TYPE(xc_rho_set_type), POINTER :: rho_set
462
463 CALL timeset(routinen, handle)
464
465 CALL qs_rho_get(rho_ia_struct, rho_r=rho1_r, rho_g=rho1_g, tau_r=tau1_r)
466 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
467 DO ispin = 1, nspins
468 CALL pw_zero(work_v_xc(ispin))
469 END DO
470 rho_set => kernel_env%xc_rho_set
471
472 singlet = .false.
473 triplet = .false.
474 lsd = .false.
475 IF (nspins == 1 .AND. .NOT. is_rks_triplets) THEN
476 singlet = .true.
477 ELSE IF (nspins == 1 .AND. is_rks_triplets) THEN
478 triplet = .true.
479 ELSE IF (nspins == 2) THEN
480 lsd = .true.
481 ELSE
482 cpabort("illegal options")
483 END IF
484
485 IF (ASSOCIATED(tau1_r)) THEN
486 DO ispin = 1, nspins
487 CALL pw_zero(work_v_xc_tau(ispin))
488 END DO
489 END IF
490
491 CALL xc_calc_2nd_deriv_numerical(work_v_xc, work_v_xc_tau, rho_set, rho1_r, rho1_g, tau1_r, &
492 auxbas_pw_pool, kernel_env%xc_section, &
493 is_rks_triplets)
494
495 CALL timestop(handle)
496
497 END SUBROUTINE tddfpt_apply_xc_fd
498
499! **************************************************************************************************
500!> \brief Update action of TDDFPT operator on trial vectors by adding exact-exchange term.
501!> \param Aop_evects action of TDDFPT operator on trial vectors (modified on exit)
502!> \param evects trial vectors
503!> \param gs_mos molecular orbitals optimised for the ground state (only occupied
504!> molecular orbitals [component %mos_occ] are needed)
505!> \param do_admm perform auxiliary density matrix method calculations
506!> \param qs_env Quickstep environment
507!> \param work_rho_ia_ao_symm ...
508!> \param work_hmat_symm ...
509!> \param work_rho_ia_ao_asymm ...
510!> \param work_hmat_asymm ...
511!> \param wfm_rho_orb ...
512!> \par History
513!> * 05.2016 compute all exact-exchange terms in one go [Sergey Chulkov]
514!> * 03.2017 code related to ADMM correction is now moved to tddfpt_apply_admm_correction()
515!> in order to compute this correction within parallel groups [Sergey Chulkov]
516!> \note Based on the subroutine kpp1_calc_k_p_p1() which was originally created by
517!> Mohamed Fawzi on 10.2002.
518! **************************************************************************************************
519 SUBROUTINE tddfpt_apply_hfx(Aop_evects, evects, gs_mos, do_admm, qs_env, &
520 work_rho_ia_ao_symm, work_hmat_symm, work_rho_ia_ao_asymm, &
521 work_hmat_asymm, wfm_rho_orb)
522 TYPE(cp_fm_type), DIMENSION(:, :), INTENT(INOUT) :: aop_evects
523 TYPE(cp_fm_type), DIMENSION(:, :), INTENT(IN) :: evects
524 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
525 INTENT(in) :: gs_mos
526 LOGICAL, INTENT(in) :: do_admm
527 TYPE(qs_environment_type), POINTER :: qs_env
528 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT) :: work_rho_ia_ao_symm
529 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
530 TARGET :: work_hmat_symm
531 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT) :: work_rho_ia_ao_asymm
532 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
533 TARGET :: work_hmat_asymm
534 TYPE(cp_fm_type), INTENT(IN) :: wfm_rho_orb
535
536 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_apply_hfx'
537
538 INTEGER :: handle, ispin, ivect, nao, nao_aux, &
539 nspins, nvects
540 INTEGER, DIMENSION(maxspins) :: nactive
541 LOGICAL :: do_hfx
542 REAL(kind=dp) :: alpha
543 TYPE(admm_type), POINTER :: admm_env
544 TYPE(section_vals_type), POINTER :: hfx_section, input
545
546 CALL timeset(routinen, handle)
547
548 ! Check for hfx section
549 CALL get_qs_env(qs_env, input=input)
550 hfx_section => section_vals_get_subs_vals(input, "DFT%XC%HF")
551 CALL section_vals_get(hfx_section, explicit=do_hfx)
552
553 IF (do_hfx) THEN
554 nspins = SIZE(evects, 1)
555 nvects = SIZE(evects, 2)
556
557 IF (nspins > 1) THEN
558 alpha = 1.0_dp
559 ELSE
560 alpha = 2.0_dp
561 END IF
562
563 CALL cp_fm_get_info(gs_mos(1)%mos_occ, nrow_global=nao)
564 DO ispin = 1, nspins
565 CALL cp_fm_get_info(evects(ispin, 1), ncol_global=nactive(ispin))
566 END DO
567
568 IF (do_admm) THEN
569 CALL get_qs_env(qs_env, admm_env=admm_env)
570 CALL cp_fm_get_info(admm_env%A, nrow_global=nao_aux)
571 END IF
572
573 !Note: the symmetrized transition density matrix is P = 0.5*(C*evect^T + evect*C^T)
574 ! in the end, we only want evect*C^T for consistency with the MO formulation of TDDFT
575 ! therefore, we go in 2 steps: with the symmetric 0.5*(C*evect^T + evect*C^T) and
576 ! the antisymemtric 0.5*(C*evect^T - evect*C^T)
577
578 ! some stuff from qs_ks_build_kohn_sham_matrix
579 ! TO DO: add SIC support
580 DO ivect = 1, nvects
581 DO ispin = 1, nspins
582
583 !The symmetric density matrix
584 CALL parallel_gemm('N', 'T', nao, nao, nactive(ispin), 0.5_dp, evects(ispin, ivect), &
585 gs_mos(ispin)%mos_occ, 0.0_dp, wfm_rho_orb)
586 CALL parallel_gemm('N', 'T', nao, nao, nactive(ispin), 0.5_dp, gs_mos(ispin)%mos_occ, &
587 evects(ispin, ivect), 1.0_dp, wfm_rho_orb)
588
589 CALL dbcsr_set(work_hmat_symm(ispin)%matrix, 0.0_dp)
590 IF (do_admm) THEN
591 CALL parallel_gemm('N', 'N', nao_aux, nao, nao, 1.0_dp, admm_env%A, &
592 wfm_rho_orb, 0.0_dp, admm_env%work_aux_orb)
593 CALL parallel_gemm('N', 'T', nao_aux, nao_aux, nao, 1.0_dp, admm_env%work_aux_orb, admm_env%A, &
594 0.0_dp, admm_env%work_aux_aux)
595 CALL copy_fm_to_dbcsr(admm_env%work_aux_aux, work_rho_ia_ao_symm(ispin)%matrix, keep_sparsity=.true.)
596 ELSE
597 CALL copy_fm_to_dbcsr(wfm_rho_orb, work_rho_ia_ao_symm(ispin)%matrix, keep_sparsity=.true.)
598 END IF
599 END DO
600
601 CALL tddft_hfx_matrix(work_hmat_symm, work_rho_ia_ao_symm, qs_env)
602
603 IF (do_admm) THEN
604 DO ispin = 1, nspins
605 CALL cp_dbcsr_sm_fm_multiply(work_hmat_symm(ispin)%matrix, admm_env%A, admm_env%work_aux_orb, &
606 ncol=nao, alpha=1.0_dp, beta=0.0_dp)
607
608 CALL parallel_gemm('T', 'N', nao, nao, nao_aux, 1.0_dp, admm_env%A, &
609 admm_env%work_aux_orb, 0.0_dp, wfm_rho_orb)
610
611 CALL parallel_gemm('N', 'N', nao, nactive(ispin), nao, alpha, wfm_rho_orb, &
612 gs_mos(ispin)%mos_occ, 1.0_dp, aop_evects(ispin, ivect))
613 END DO
614 ELSE
615 DO ispin = 1, nspins
616 CALL cp_dbcsr_sm_fm_multiply(work_hmat_symm(ispin)%matrix, gs_mos(ispin)%mos_occ, &
617 aop_evects(ispin, ivect), ncol=nactive(ispin), &
618 alpha=alpha, beta=1.0_dp)
619 END DO
620 END IF
621
622 !The anti-symmetric density matrix
623 DO ispin = 1, nspins
624
625 !The symmetric density matrix
626 CALL parallel_gemm('N', 'T', nao, nao, nactive(ispin), 0.5_dp, evects(ispin, ivect), &
627 gs_mos(ispin)%mos_occ, 0.0_dp, wfm_rho_orb)
628 CALL parallel_gemm('N', 'T', nao, nao, nactive(ispin), -0.5_dp, gs_mos(ispin)%mos_occ, &
629 evects(ispin, ivect), 1.0_dp, wfm_rho_orb)
630
631 CALL dbcsr_set(work_hmat_asymm(ispin)%matrix, 0.0_dp)
632 IF (do_admm) THEN
633 CALL parallel_gemm('N', 'N', nao_aux, nao, nao, 1.0_dp, admm_env%A, &
634 wfm_rho_orb, 0.0_dp, admm_env%work_aux_orb)
635 CALL parallel_gemm('N', 'T', nao_aux, nao_aux, nao, 1.0_dp, admm_env%work_aux_orb, admm_env%A, &
636 0.0_dp, admm_env%work_aux_aux)
637 CALL copy_fm_to_dbcsr(admm_env%work_aux_aux, work_rho_ia_ao_asymm(ispin)%matrix, keep_sparsity=.true.)
638 ELSE
639 CALL copy_fm_to_dbcsr(wfm_rho_orb, work_rho_ia_ao_asymm(ispin)%matrix, keep_sparsity=.true.)
640 END IF
641 END DO
642
643 CALL tddft_hfx_matrix(work_hmat_asymm, work_rho_ia_ao_asymm, qs_env)
644
645 IF (do_admm) THEN
646 DO ispin = 1, nspins
647 CALL cp_dbcsr_sm_fm_multiply(work_hmat_asymm(ispin)%matrix, admm_env%A, admm_env%work_aux_orb, &
648 ncol=nao, alpha=1.0_dp, beta=0.0_dp)
649
650 CALL parallel_gemm('T', 'N', nao, nao, nao_aux, 1.0_dp, admm_env%A, &
651 admm_env%work_aux_orb, 0.0_dp, wfm_rho_orb)
652
653 CALL parallel_gemm('N', 'N', nao, nactive(ispin), nao, alpha, wfm_rho_orb, &
654 gs_mos(ispin)%mos_occ, 1.0_dp, aop_evects(ispin, ivect))
655 END DO
656 ELSE
657 DO ispin = 1, nspins
658 CALL cp_dbcsr_sm_fm_multiply(work_hmat_asymm(ispin)%matrix, gs_mos(ispin)%mos_occ, &
659 aop_evects(ispin, ivect), ncol=nactive(ispin), &
660 alpha=alpha, beta=1.0_dp)
661 END DO
662 END IF
663 END DO
664 END IF
665
666 CALL timestop(handle)
667
668 END SUBROUTINE tddfpt_apply_hfx
669
670! **************************************************************************************************
671!> \brief Update action of TDDFPT operator on trial vectors by adding exact-exchange term.
672!> \param Aop_evects action of TDDFPT operator on trial vectors (modified on exit)
673!> \param evects trial vectors
674!> \param gs_mos molecular orbitals optimised for the ground state (only occupied
675!> molecular orbitals [component %mos_occ] are needed)
676!> \param qs_env Quickstep environment
677!> \param admm_env ...
678!> \param hfx_section ...
679!> \param x_data ...
680!> \param symmetry ...
681!> \param recalc_integrals ...
682!> \param work_rho_ia_ao ...
683!> \param work_hmat ...
684!> \param wfm_rho_orb ...
685! **************************************************************************************************
686 SUBROUTINE tddfpt_apply_hfxsr_kernel(Aop_evects, evects, gs_mos, qs_env, admm_env, &
687 hfx_section, x_data, symmetry, recalc_integrals, &
688 work_rho_ia_ao, work_hmat, wfm_rho_orb)
689 TYPE(cp_fm_type), DIMENSION(:, :), INTENT(in) :: aop_evects, evects
690 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
691 INTENT(in) :: gs_mos
692 TYPE(qs_environment_type), POINTER :: qs_env
693 TYPE(admm_type), POINTER :: admm_env
694 TYPE(section_vals_type), POINTER :: hfx_section
695 TYPE(hfx_type), DIMENSION(:, :), POINTER :: x_data
696 INTEGER, INTENT(IN) :: symmetry
697 LOGICAL, INTENT(IN) :: recalc_integrals
698 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT) :: work_rho_ia_ao
699 TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
700 TARGET :: work_hmat
701 TYPE(cp_fm_type), INTENT(IN) :: wfm_rho_orb
702
703 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_apply_hfxsr_kernel'
704
705 INTEGER :: handle, ispin, ivect, nao, nao_aux, &
706 nspins, nvects
707 INTEGER, DIMENSION(maxspins) :: nactive
708 LOGICAL :: reint
709 REAL(kind=dp) :: alpha
710
711 CALL timeset(routinen, handle)
712
713 nspins = SIZE(evects, 1)
714 nvects = SIZE(evects, 2)
715
716 alpha = 2.0_dp
717 IF (nspins > 1) alpha = 1.0_dp
718
719 CALL cp_fm_get_info(gs_mos(1)%mos_occ, nrow_global=nao)
720 CALL cp_fm_get_info(admm_env%A, nrow_global=nao_aux)
721 DO ispin = 1, nspins
722 CALL cp_fm_get_info(evects(ispin, 1), ncol_global=nactive(ispin))
723 END DO
724
725 reint = recalc_integrals
726
727 DO ivect = 1, nvects
728 DO ispin = 1, nspins
729 CALL parallel_gemm('N', 'T', nao, nao, nactive(ispin), 0.5_dp, evects(ispin, ivect), &
730 gs_mos(ispin)%mos_occ, 0.0_dp, wfm_rho_orb)
731 CALL parallel_gemm('N', 'T', nao, nao, nactive(ispin), 0.5_dp*symmetry, gs_mos(ispin)%mos_occ, &
732 evects(ispin, ivect), 1.0_dp, wfm_rho_orb)
733 CALL dbcsr_set(work_hmat(ispin)%matrix, 0.0_dp)
734 CALL parallel_gemm('N', 'N', nao_aux, nao, nao, 1.0_dp, admm_env%A, &
735 wfm_rho_orb, 0.0_dp, admm_env%work_aux_orb)
736 CALL parallel_gemm('N', 'T', nao_aux, nao_aux, nao, 1.0_dp, admm_env%work_aux_orb, admm_env%A, &
737 0.0_dp, admm_env%work_aux_aux)
738 CALL copy_fm_to_dbcsr(admm_env%work_aux_aux, work_rho_ia_ao(ispin)%matrix, keep_sparsity=.true.)
739 END DO
740
741 CALL tddft_hfx_matrix(work_hmat, work_rho_ia_ao, qs_env, .false., reint, hfx_section, x_data)
742 reint = .false.
743
744 DO ispin = 1, nspins
745 CALL cp_dbcsr_sm_fm_multiply(work_hmat(ispin)%matrix, admm_env%A, admm_env%work_aux_orb, &
746 ncol=nao, alpha=1.0_dp, beta=0.0_dp)
747 CALL parallel_gemm('T', 'N', nao, nao, nao_aux, 1.0_dp, admm_env%A, &
748 admm_env%work_aux_orb, 0.0_dp, wfm_rho_orb)
749 CALL parallel_gemm('N', 'N', nao, nactive(ispin), nao, alpha, wfm_rho_orb, &
750 gs_mos(ispin)%mos_occ, 1.0_dp, aop_evects(ispin, ivect))
751 END DO
752 END DO
753
754 CALL timestop(handle)
755
756 END SUBROUTINE tddfpt_apply_hfxsr_kernel
757
758! **************************************************************************************************
759!> \brief ...Calculate the HFXLR kernel contribution by contracting the Lowdin MO coefficients --
760!> transition charges with the exchange-type integrals using the sTDA approximation
761!> \param qs_env ...
762!> \param sub_env ...
763!> \param rcut ...
764!> \param hfx_scale ...
765!> \param work ...
766!> \param X ...
767!> \param res ... vector AX with A being the sTDA matrix and X the Davidson trial vector of the
768!> eigenvalue problem A*X = omega*X
769! **************************************************************************************************
770 SUBROUTINE tddfpt_apply_hfxlr_kernel(qs_env, sub_env, rcut, hfx_scale, work, X, res)
771
772 TYPE(qs_environment_type), POINTER :: qs_env
773 TYPE(tddfpt_subgroup_env_type) :: sub_env
774 REAL(kind=dp), INTENT(IN) :: rcut, hfx_scale
775 TYPE(tddfpt_work_matrices) :: work
776 TYPE(cp_fm_type), DIMENSION(:), INTENT(IN) :: x
777 TYPE(cp_fm_type), DIMENSION(:), INTENT(INOUT) :: res
778
779 CHARACTER(len=*), PARAMETER :: routinen = 'tddfpt_apply_hfxlr_kernel'
780
781 INTEGER :: handle, iatom, ispin, jatom, natom, &
782 nsgf, nspins
783 INTEGER, DIMENSION(2) :: nactive
784 REAL(kind=dp) :: dr, eps_filter, fcut, gabr
785 REAL(kind=dp), DIMENSION(3) :: rij
786 REAL(kind=dp), DIMENSION(:, :), POINTER :: pblock
787 TYPE(cell_type), POINTER :: cell
788 TYPE(cp_fm_struct_type), POINTER :: fmstruct
789 TYPE(cp_fm_type) :: cvec
790 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:) :: xtransformed
791 TYPE(cp_fm_type), POINTER :: ct
792 TYPE(dbcsr_iterator_type) :: iter
793 TYPE(dbcsr_type) :: pdens
794 TYPE(dbcsr_type), POINTER :: tempmat
795 TYPE(mp_para_env_type), POINTER :: para_env
796 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
797
798 CALL timeset(routinen, handle)
799
800 ! parameters
801 eps_filter = 1.e-08_dp
802
803 nspins = SIZE(x)
804 DO ispin = 1, nspins
805 CALL cp_fm_get_info(x(ispin), ncol_global=nactive(ispin))
806 END DO
807
808 para_env => sub_env%para_env
809
810 CALL get_qs_env(qs_env, natom=natom, cell=cell, particle_set=particle_set)
811
812 ! calculate Loewdin transformed Davidson trial vector tilde(X)=S^1/2*X
813 ! and tilde(tilde(X))=S^1/2_A*tilde(X)_A
814 ALLOCATE (xtransformed(nspins))
815 DO ispin = 1, nspins
816 NULLIFY (fmstruct)
817 ct => work%ctransformed(ispin)
818 CALL cp_fm_get_info(ct, matrix_struct=fmstruct)
819 CALL cp_fm_create(matrix=xtransformed(ispin), matrix_struct=fmstruct, name="XTRANSFORMED")
820 END DO
821 CALL get_lowdin_x(work%shalf, x, xtransformed)
822
823 DO ispin = 1, nspins
824 ct => work%ctransformed(ispin)
825 CALL cp_fm_get_info(ct, matrix_struct=fmstruct, nrow_global=nsgf)
826 CALL cp_fm_create(cvec, fmstruct)
827 !
828 tempmat => work%shalf
829 CALL dbcsr_create(pdens, template=tempmat, matrix_type=dbcsr_type_no_symmetry)
830 ! P(nu,mu) = SUM_j XT(nu,j)*CT(mu,j)
831 ct => work%ctransformed(ispin)
832 CALL dbcsr_set(pdens, 0.0_dp)
833 CALL cp_dbcsr_plus_fm_fm_t(pdens, xtransformed(ispin), ct, nactive(ispin), &
834 1.0_dp, keep_sparsity=.false.)
835 CALL dbcsr_filter(pdens, eps_filter)
836 ! Apply PP*gab -> PP; gab = gamma_coulomb
837 ! P(nu,mu) = P(nu,mu)*g(a of nu,b of mu)
838 CALL dbcsr_iterator_start(iter, pdens)
839 DO WHILE (dbcsr_iterator_blocks_left(iter))
840 CALL dbcsr_iterator_next_block(iter, iatom, jatom, pblock)
841 rij = particle_set(iatom)%r - particle_set(jatom)%r
842 rij = pbc(rij, cell)
843 dr = sqrt(sum(rij(:)**2))
844 gabr = 1._dp/rcut
845 IF (dr < 1.e-6) THEN
846 gabr = 2._dp*gabr/sqrt(3.1415926_dp)
847 ELSE
848 gabr = erf(gabr*dr)/dr
849 fcut = exp(dr - 4._dp*rcut)
850 fcut = fcut/(fcut + 1._dp)
851 END IF
852 pblock = hfx_scale*gabr*pblock
853 END DO
854 CALL dbcsr_iterator_stop(iter)
855 ! CV(mu,i) = P(nu,mu)*CT(mu,i)
856 CALL cp_dbcsr_sm_fm_multiply(pdens, ct, cvec, nactive(ispin), 1.0_dp, 0.0_dp)
857 ! rho(nu,i) = rho(nu,i) + ShalfP(nu,mu)*CV(mu,i)
858 CALL cp_dbcsr_sm_fm_multiply(work%shalf, cvec, res(ispin), nactive(ispin), &
859 -1.0_dp, 1.0_dp)
860 !
861 CALL dbcsr_release(pdens)
862 !
863 CALL cp_fm_release(cvec)
864 END DO
865
866 CALL cp_fm_release(xtransformed)
867
868 CALL timestop(handle)
869
870 END SUBROUTINE tddfpt_apply_hfxlr_kernel
871
872! **************************************************************************************************
873
874END MODULE qs_tddfpt2_operators
Types and set/get functions for auxiliary density matrix methods.
Definition admm_types.F:15
Handles all functions related to the CELL.
Definition cell_types.F:15
subroutine, public dbcsr_iterator_next_block(iterator, row, column, block, block_number_argument_has_been_removed, row_size, col_size, row_offset, col_offset)
...
logical function, public dbcsr_iterator_blocks_left(iterator)
...
subroutine, public dbcsr_iterator_stop(iterator)
...
subroutine, public dbcsr_filter(matrix, eps)
...
subroutine, public dbcsr_iterator_start(iterator, matrix, shared, dynamic, dynamic_byrows)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_release(matrix)
...
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
subroutine, public cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
performs the multiplication sparse_matrix+dense_mat*dens_mat^T if matrix_g is not explicitly given,...
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
Basic linear algebra operations for full matrices.
subroutine, public cp_fm_column_scale(matrixa, scaling)
scales column i of matrix a with scaling(i)
subroutine, public cp_fm_scale_and_add(alpha, matrix_a, beta, matrix_b)
calc A <- alpha*A + beta*B optimized for alpha == 1.0 (just add beta*B) and beta == 0....
represent the structure of a full matrix
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
subroutine, public vh_1c_gg_integrals(qs_env, energy_hartree_1c, ecoul_1c, local_rho_set, para_env, tddft, local_rho_set_2nd, core_2nd)
Calculates one center GAPW Hartree energies and matrix elements Hartree potentials are input Takes po...
Utilities for hfx and admm methods.
subroutine, public tddft_hfx_matrix(matrix_ks, rho_ao, qs_env, update_energy, recalc_integrals, external_hfx_sections, external_x_data, external_para_env)
Add the hfx contributions to the Hamiltonian.
Types and set/get functions for HFX.
Definition hfx_types.F:15
objects that represent the structure of input sections and the data contained in an input section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Interface to the message passing library MPI.
basic linear algebra operations for full matrixes
Define the data structure for the particle information.
container for various plainwaves related things
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
functions related to the poisson solver on regular grids
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
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_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, 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, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_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, harris_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, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
subroutine, public integrate_vhg0_rspace(qs_env, v_rspace, para_env, calculate_forces, local_rho_set, local_rho_set_2nd, atener, kforce, my_pools, my_rs_descs)
...
superstucture that hold various representations of the density and keeps track of which ones are vali...
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...
subroutine, public tddfpt_apply_hfx(aop_evects, evects, gs_mos, do_admm, qs_env, work_rho_ia_ao_symm, work_hmat_symm, work_rho_ia_ao_asymm, work_hmat_asymm, wfm_rho_orb)
Update action of TDDFPT operator on trial vectors by adding exact-exchange term.
subroutine, public tddfpt_apply_coulomb(a_ia_rspace, rho_ia_g, local_rho_set, hartree_local, qs_env, sub_env, gapw, work_v_gspace, work_v_rspace, tddfpt_mgrid)
Update v_rspace by adding coulomb term.
subroutine, public tddfpt_apply_hfxlr_kernel(qs_env, sub_env, rcut, hfx_scale, work, x, res)
...Calculate the HFXLR kernel contribution by contracting the Lowdin MO coefficients – transition cha...
subroutine, public tddfpt_apply_energy_diff(aop_evects, evects, s_evects, gs_mos, matrix_ks)
Apply orbital energy difference term: Aop_evects(spin,state) += KS(spin) * evects(spin,...
subroutine, public tddfpt_apply_xc(a_ia_rspace, kernel_env, rho_ia_struct, is_rks_triplets, pw_env, work_v_xc, work_v_xc_tau)
Driver routine for applying fxc (analyic vs. finite difference for testing.
subroutine, public tddfpt_apply_hfxsr_kernel(aop_evects, evects, gs_mos, qs_env, admm_env, hfx_section, x_data, symmetry, recalc_integrals, work_rho_ia_ao, work_hmat, wfm_rho_orb)
Update action of TDDFPT operator on trial vectors by adding exact-exchange term.
subroutine, public tddfpt_apply_xc_potential(a_ia_rspace, fxc_rspace, rho_ia_struct, is_rks_triplets)
Routine for applying fxc potential.
Simplified Tamm Dancoff approach (sTDA).
subroutine, public get_lowdin_x(shalf, xvec, xt)
Calculate Lowdin transformed Davidson trial vector X shalf (dbcsr), xvec, xt (fm) are defined in the ...
contains the structure
subroutine, public xc_rho_set_update(rho_set, rho_r, rho_g, tau, needs, xc_deriv_method_id, xc_rho_smooth_id, pw_pool)
updates the given rho set with the density given by rho_r (and rho_g). The rho set will contain the c...
Exchange and Correlation functional calculations.
Definition xc.F:17
subroutine, public xc_calc_2nd_deriv_analytical(v_xc, v_xc_tau, deriv_set, rho_set, rho1_set, pw_pool, xc_section, gapw, vxg, tddfpt_fac, compute_virial, virial_xc)
Calculates the second derivative of E_xc at rho in the direction rho1 (if you see the second derivati...
Definition xc.F:2640
subroutine, public xc_calc_2nd_deriv_numerical(v_xc, v_tau, rho_set, rho1_r, rho1_g, tau1_r, pw_pool, xc_section, do_triplet, calc_virial, virial_xc, deriv_set)
calculates 2nd derivative numerically
Definition xc.F:1653
stores some data used in wavefunction fitting
Definition admm_types.F:120
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
keeps the information about the structure of a full matrix
represent a full matrix
stores some data used in construction of Kohn-Sham matrix
Definition hfx_types.F:510
stores all the informations relevant to an mpi environment
contained for different pw related things
environment for the poisson solver
to create arrays of pools
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Collection of variables required to evaluate adiabatic TDDFPT kernel.
keeps the density in various representations, keeping track of which ones are valid.
Ground state molecular orbitals.
Set of temporary ("work") matrices.
represent a density, with all the representation and data needed to perform a functional evaluation