(git:374b731)
Loading...
Searching...
No Matches
qs_rho_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 methods of the rho structure (defined in qs_rho_types)
10!> \par History
11!> 08.2002 created [fawzi]
12!> 08.2014 kpoints [JGH]
13!> \author Fawzi Mohamed
14! **************************************************************************************************
15MODULE qs_rho_methods
16 USE admm_types, ONLY: get_admm_env
17 USE atomic_kind_types, ONLY: atomic_kind_type
18 USE cp_control_types, ONLY: dft_control_type
19 USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
20 USE cp_dbcsr_operations, ONLY: dbcsr_allocate_matrix_set,&
21 dbcsr_deallocate_matrix_set
22 USE cp_log_handling, ONLY: cp_to_string
23 USE dbcsr_api, ONLY: &
24 dbcsr_add, dbcsr_copy, dbcsr_create, dbcsr_p_type, dbcsr_scale, dbcsr_set, dbcsr_type, &
25 dbcsr_type_antisymmetric, dbcsr_type_symmetric
26 USE kinds, ONLY: default_string_length,&
27 dp
28 USE kpoint_types, ONLY: get_kpoint_info,&
29 kpoint_type
30 USE lri_environment_methods, ONLY: calculate_lri_densities
31 USE lri_environment_types, ONLY: lri_density_type,&
32 lri_environment_type
33 USE message_passing, ONLY: mp_para_env_type
34 USE pw_env_types, ONLY: pw_env_get,&
35 pw_env_type
36 USE pw_methods, ONLY: pw_axpy,&
37 pw_copy,&
38 pw_scale,&
39 pw_zero
40 USE pw_pool_types, ONLY: pw_pool_type
41 USE pw_types, ONLY: pw_c1d_gs_type,&
42 pw_r3d_rs_type
43 USE qs_collocate_density, ONLY: calculate_drho_elec,&
44 calculate_rho_elec
45 USE qs_environment_types, ONLY: get_qs_env,&
46 qs_environment_type,&
47 set_qs_env
48 USE qs_kind_types, ONLY: qs_kind_type
49 USE qs_ks_types, ONLY: get_ks_env,&
50 qs_ks_env_type
51 USE qs_local_rho_types, ONLY: local_rho_type
52 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
53 USE qs_oce_types, ONLY: oce_matrix_type
54 USE qs_rho_atom_methods, ONLY: calculate_rho_atom_coeff
55 USE qs_rho_atom_types, ONLY: rho_atom_type
56 USE qs_rho_types, ONLY: qs_rho_clear,&
57 qs_rho_get,&
58 qs_rho_set,&
59 qs_rho_type
60 USE ri_environment_methods, ONLY: calculate_ri_densities
61 USE task_list_types, ONLY: task_list_type
62#include "./base/base_uses.f90"
63
64 IMPLICIT NONE
65 PRIVATE
66
67 LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .true.
68 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_rho_methods'
69
70 PUBLIC :: qs_rho_update_rho, qs_rho_update_tddfpt, &
71 qs_rho_rebuild, qs_rho_copy, qs_rho_scale_and_add
72 PUBLIC :: duplicate_rho_type, allocate_rho_ao_imag_from_real
73
74CONTAINS
75
76! **************************************************************************************************
77!> \brief rebuilds rho (if necessary allocating and initializing it)
78!> \param rho the rho type to rebuild (defaults to qs_env%rho)
79!> \param qs_env the environment to which rho belongs
80!> \param rebuild_ao if it is necessary to rebuild rho_ao. Defaults to true.
81!> \param rebuild_grids if it in necessary to rebuild rho_r and rho_g.
82!> Defaults to false.
83!> \param admm (use aux_fit basis)
84!> \param pw_env_external external plane wave environment
85!> \par History
86!> 11.2002 created replacing qs_rho_create and qs_env_rebuild_rho[fawzi]
87!> \author Fawzi Mohamed
88!> \note
89!> needs updated pw pools, s, s_mstruct and h in qs_env.
90!> The use of p to keep the structure of h (needed for the forces)
91!> is ugly and should be removed.
92!> Change so that it does not allocate a subcomponent if it is not
93!> associated and not requested?
94! **************************************************************************************************
95 SUBROUTINE qs_rho_rebuild(rho, qs_env, rebuild_ao, rebuild_grids, admm, pw_env_external)
96 TYPE(qs_rho_type), INTENT(INOUT) :: rho
97 TYPE(qs_environment_type), POINTER :: qs_env
98 LOGICAL, INTENT(in), OPTIONAL :: rebuild_ao, rebuild_grids, admm
99 TYPE(pw_env_type), OPTIONAL, POINTER :: pw_env_external
100
101 CHARACTER(LEN=*), PARAMETER :: routinen = 'qs_rho_rebuild'
102
103 CHARACTER(LEN=default_string_length) :: headline
104 INTEGER :: handle, i, ic, j, nimg, nspins
105 LOGICAL :: do_kpoints, my_admm, my_rebuild_ao, &
106 my_rebuild_grids, rho_ao_is_complex
107 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_r
108 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_kp, rho_ao_im_kp, rho_ao_kp
109 TYPE(dbcsr_type), POINTER :: refmatrix, tmatrix
110 TYPE(dft_control_type), POINTER :: dft_control
111 TYPE(kpoint_type), POINTER :: kpoints
112 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
113 POINTER :: sab_orb
114 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g, tau_g
115 TYPE(pw_c1d_gs_type), DIMENSION(:, :), POINTER :: drho_g
116 TYPE(pw_env_type), POINTER :: pw_env
117 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
118 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r, tau_r
119 TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_r
120 TYPE(pw_r3d_rs_type), POINTER :: rho_r_sccs
121
122 CALL timeset(routinen, handle)
123
124 NULLIFY (pw_env, auxbas_pw_pool, matrix_s_kp, dft_control)
125 NULLIFY (tot_rho_r, rho_ao_kp, rho_r, rho_g, drho_r, drho_g, tau_r, tau_g, rho_ao_im_kp)
126 NULLIFY (rho_r_sccs)
127 NULLIFY (sab_orb)
128 my_rebuild_ao = .true.
129 my_rebuild_grids = .true.
130 my_admm = .false.
131 IF (PRESENT(rebuild_ao)) my_rebuild_ao = rebuild_ao
132 IF (PRESENT(rebuild_grids)) my_rebuild_grids = rebuild_grids
133 IF (PRESENT(admm)) my_admm = admm
134
135 CALL get_qs_env(qs_env, &
136 kpoints=kpoints, &
137 do_kpoints=do_kpoints, &
138 pw_env=pw_env, &
139 dft_control=dft_control)
140 IF (PRESENT(pw_env_external)) &
141 pw_env => pw_env_external
142
143 nimg = dft_control%nimages
144
145 IF (my_admm) THEN
146 CALL get_admm_env(qs_env%admm_env, sab_aux_fit=sab_orb, matrix_s_aux_fit_kp=matrix_s_kp)
147 ELSE
148 CALL get_qs_env(qs_env, matrix_s_kp=matrix_s_kp)
149
150 IF (do_kpoints) THEN
151 CALL get_kpoint_info(kpoints, sab_nl=sab_orb)
152 ELSE
153 CALL get_qs_env(qs_env, sab_orb=sab_orb)
154 END IF
155 END IF
156 refmatrix => matrix_s_kp(1, 1)%matrix
157
158 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
159 nspins = dft_control%nspins
160
161 CALL qs_rho_get(rho, &
162 tot_rho_r=tot_rho_r, &
163 rho_ao_kp=rho_ao_kp, &
164 rho_ao_im_kp=rho_ao_im_kp, &
165 rho_r=rho_r, &
166 rho_g=rho_g, &
167 drho_r=drho_r, &
168 drho_g=drho_g, &
169 tau_r=tau_r, &
170 tau_g=tau_g, &
171 rho_r_sccs=rho_r_sccs, &
172 complex_rho_ao=rho_ao_is_complex)
173
174 IF (.NOT. ASSOCIATED(tot_rho_r)) THEN
175 ALLOCATE (tot_rho_r(nspins))
176 tot_rho_r = 0.0_dp
177 CALL qs_rho_set(rho, tot_rho_r=tot_rho_r)
178 END IF
179
180 ! rho_ao
181 IF (my_rebuild_ao .OR. (.NOT. ASSOCIATED(rho_ao_kp))) THEN
182 IF (ASSOCIATED(rho_ao_kp)) &
183 CALL dbcsr_deallocate_matrix_set(rho_ao_kp)
184 ! Create a new density matrix set
185 CALL dbcsr_allocate_matrix_set(rho_ao_kp, nspins, nimg)
186 CALL qs_rho_set(rho, rho_ao_kp=rho_ao_kp)
187 DO i = 1, nspins
188 DO ic = 1, nimg
189 IF (nspins > 1) THEN
190 IF (i == 1) THEN
191 headline = "DENSITY MATRIX FOR ALPHA SPIN"
192 ELSE
193 headline = "DENSITY MATRIX FOR BETA SPIN"
194 END IF
195 ELSE
196 headline = "DENSITY MATRIX"
197 END IF
198 ALLOCATE (rho_ao_kp(i, ic)%matrix)
199 tmatrix => rho_ao_kp(i, ic)%matrix
200 CALL dbcsr_create(matrix=tmatrix, template=refmatrix, name=trim(headline), &
201 matrix_type=dbcsr_type_symmetric, nze=0)
202 CALL cp_dbcsr_alloc_block_from_nbl(tmatrix, sab_orb)
203 CALL dbcsr_set(tmatrix, 0.0_dp)
204 END DO
205 END DO
206 IF (rho_ao_is_complex) THEN
207 IF (ASSOCIATED(rho_ao_im_kp)) THEN
208 CALL dbcsr_deallocate_matrix_set(rho_ao_im_kp)
209 END IF
210 CALL dbcsr_allocate_matrix_set(rho_ao_im_kp, nspins, nimg)
211 CALL qs_rho_set(rho, rho_ao_im_kp=rho_ao_im_kp)
212 DO i = 1, nspins
213 DO ic = 1, nimg
214 IF (nspins > 1) THEN
215 IF (i == 1) THEN
216 headline = "IMAGINARY PART OF DENSITY MATRIX FOR ALPHA SPIN"
217 ELSE
218 headline = "IMAGINARY PART OF DENSITY MATRIX FOR BETA SPIN"
219 END IF
220 ELSE
221 headline = "IMAGINARY PART OF DENSITY MATRIX"
222 END IF
223 ALLOCATE (rho_ao_im_kp(i, ic)%matrix)
224 tmatrix => rho_ao_im_kp(i, ic)%matrix
225 CALL dbcsr_create(matrix=tmatrix, template=refmatrix, name=trim(headline), &
226 matrix_type=dbcsr_type_antisymmetric, nze=0)
227 CALL cp_dbcsr_alloc_block_from_nbl(tmatrix, sab_orb)
228 CALL dbcsr_set(tmatrix, 0.0_dp)
229 END DO
230 END DO
231 END IF
232 END IF
233
234 ! rho_r
235 IF (my_rebuild_grids .OR. .NOT. ASSOCIATED(rho_r)) THEN
236 IF (ASSOCIATED(rho_r)) THEN
237 DO i = 1, SIZE(rho_r)
238 CALL rho_r(i)%release()
239 END DO
240 DEALLOCATE (rho_r)
241 END IF
242 ALLOCATE (rho_r(nspins))
243 CALL qs_rho_set(rho, rho_r=rho_r)
244 DO i = 1, nspins
245 CALL auxbas_pw_pool%create_pw(rho_r(i))
246 END DO
247 END IF
248
249 ! rho_g
250 IF (my_rebuild_grids .OR. .NOT. ASSOCIATED(rho_g)) THEN
251 IF (ASSOCIATED(rho_g)) THEN
252 DO i = 1, SIZE(rho_g)
253 CALL rho_g(i)%release()
254 END DO
255 DEALLOCATE (rho_g)
256 END IF
257 ALLOCATE (rho_g(nspins))
258 CALL qs_rho_set(rho, rho_g=rho_g)
259 DO i = 1, nspins
260 CALL auxbas_pw_pool%create_pw(rho_g(i))
261 END DO
262 END IF
263
264 ! SCCS
265 IF (dft_control%do_sccs) THEN
266 IF (my_rebuild_grids .OR. (.NOT. ASSOCIATED(rho_r_sccs))) THEN
267 IF (ASSOCIATED(rho_r_sccs)) THEN
268 CALL rho_r_sccs%release()
269 DEALLOCATE (rho_r_sccs)
270 END IF
271 ALLOCATE (rho_r_sccs)
272 CALL qs_rho_set(rho, rho_r_sccs=rho_r_sccs)
273 CALL auxbas_pw_pool%create_pw(rho_r_sccs)
274 CALL pw_zero(rho_r_sccs)
275 END IF
276 END IF
277
278 ! allocate drho_r and drho_g if xc_deriv_collocate
279 IF (dft_control%drho_by_collocation) THEN
280 ! drho_r
281 IF (my_rebuild_grids .OR. .NOT. ASSOCIATED(drho_r)) THEN
282 IF (ASSOCIATED(drho_r)) THEN
283 DO j = 1, SIZE(drho_r, 2)
284 DO i = 1, SIZE(drho_r, 1)
285 CALL drho_r(i, j)%release()
286 END DO
287 END DO
288 DEALLOCATE (drho_r)
289 END IF
290 ALLOCATE (drho_r(3, nspins))
291 CALL qs_rho_set(rho, drho_r=drho_r)
292 DO j = 1, nspins
293 DO i = 1, 3
294 CALL auxbas_pw_pool%create_pw(drho_r(i, j))
295 END DO
296 END DO
297 END IF
298 ! drho_g
299 IF (my_rebuild_grids .OR. .NOT. ASSOCIATED(drho_g)) THEN
300 IF (ASSOCIATED(drho_g)) THEN
301 DO j = 1, SIZE(drho_g, 2)
302 DO i = 1, SIZE(drho_r, 1)
303 CALL drho_g(i, j)%release()
304 END DO
305 END DO
306 DEALLOCATE (drho_g)
307 END IF
308 ALLOCATE (drho_g(3, nspins))
309 CALL qs_rho_set(rho, drho_g=drho_g)
310 DO j = 1, nspins
311 DO i = 1, 3
312 CALL auxbas_pw_pool%create_pw(drho_g(i, j))
313 END DO
314 END DO
315 END IF
316 END IF
317
318 ! allocate tau_r and tau_g if use_kinetic_energy_density
319 IF (dft_control%use_kinetic_energy_density) THEN
320 ! tau_r
321 IF (my_rebuild_grids .OR. .NOT. ASSOCIATED(tau_r)) THEN
322 IF (ASSOCIATED(tau_r)) THEN
323 DO i = 1, SIZE(tau_r)
324 CALL tau_r(i)%release()
325 END DO
326 DEALLOCATE (tau_r)
327 END IF
328 ALLOCATE (tau_r(nspins))
329 CALL qs_rho_set(rho, tau_r=tau_r)
330 DO i = 1, nspins
331 CALL auxbas_pw_pool%create_pw(tau_r(i))
332 END DO
333 END IF
334
335 ! tau_g
336 IF (my_rebuild_grids .OR. .NOT. ASSOCIATED(tau_g)) THEN
337 IF (ASSOCIATED(tau_g)) THEN
338 DO i = 1, SIZE(tau_g)
339 CALL tau_g(i)%release()
340 END DO
341 DEALLOCATE (tau_g)
342 END IF
343 ALLOCATE (tau_g(nspins))
344 CALL qs_rho_set(rho, tau_g=tau_g)
345 DO i = 1, nspins
346 CALL auxbas_pw_pool%create_pw(tau_g(i))
347 END DO
348 END IF
349 END IF ! use_kinetic_energy_density
350
351 CALL timestop(handle)
352
353 END SUBROUTINE qs_rho_rebuild
354
355! **************************************************************************************************
356!> \brief updates rho_r and rho_g to the rho%rho_ao.
357!> if use_kinetic_energy_density also computes tau_r and tau_g
358!> this works for all ground state and ground state response methods
359!> \param rho_struct the rho structure that should be updated
360!> \param qs_env the qs_env rho_struct refers to
361!> the integrated charge in r space
362!> \param rho_xc_external ...
363!> \param local_rho_set ...
364!> \param task_list_external external task list
365!> \param task_list_external_soft external task list (soft_version)
366!> \param pw_env_external external plane wave environment
367!> \param para_env_external external MPI environment
368!> \par History
369!> 08.2002 created [fawzi]
370!> \author Fawzi Mohamed
371! **************************************************************************************************
372 SUBROUTINE qs_rho_update_rho(rho_struct, qs_env, rho_xc_external, local_rho_set, &
373 task_list_external, task_list_external_soft, &
374 pw_env_external, para_env_external)
375 TYPE(qs_rho_type), INTENT(INOUT) :: rho_struct
376 TYPE(qs_environment_type), POINTER :: qs_env
377 TYPE(qs_rho_type), OPTIONAL, POINTER :: rho_xc_external
378 TYPE(local_rho_type), OPTIONAL, POINTER :: local_rho_set
379 TYPE(task_list_type), OPTIONAL, POINTER :: task_list_external, &
380 task_list_external_soft
381 TYPE(pw_env_type), OPTIONAL, POINTER :: pw_env_external
382 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env_external
383
384 INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
385 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
386 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
387 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_kp
388 TYPE(dft_control_type), POINTER :: dft_control
389 TYPE(kpoint_type), POINTER :: kpoints
390 TYPE(lri_density_type), POINTER :: lri_density
391 TYPE(lri_environment_type), POINTER :: lri_env
392 TYPE(mp_para_env_type), POINTER :: para_env
393 TYPE(qs_ks_env_type), POINTER :: ks_env
394
395 CALL get_qs_env(qs_env, dft_control=dft_control, &
396 atomic_kind_set=atomic_kind_set, &
397 para_env=para_env)
398 IF (PRESENT(para_env_external)) para_env => para_env_external
399
400 IF (dft_control%qs_control%semi_empirical .OR. &
401 dft_control%qs_control%dftb .OR. dft_control%qs_control%xtb) THEN
402
403 CALL qs_rho_set(rho_struct, rho_r_valid=.false., rho_g_valid=.false.)
404
405 ELSEIF (dft_control%qs_control%lrigpw) THEN
406 cpassert(.NOT. dft_control%use_kinetic_energy_density)
407 cpassert(.NOT. dft_control%drho_by_collocation)
408 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
409 CALL get_qs_env(qs_env, ks_env=ks_env)
410 CALL get_ks_env(ks_env=ks_env, kpoints=kpoints)
411 CALL get_kpoint_info(kpoint=kpoints, cell_to_index=cell_to_index)
412 CALL get_qs_env(qs_env, lri_env=lri_env, lri_density=lri_density)
413 CALL calculate_lri_densities(lri_env, lri_density, qs_env, rho_ao_kp, cell_to_index, &
414 lri_rho_struct=rho_struct, &
415 atomic_kind_set=atomic_kind_set, &
416 para_env=para_env, &
417 response_density=.false.)
418 CALL set_qs_env(qs_env, lri_density=lri_density)
419 CALL qs_rho_set(rho_struct, rho_r_valid=.true., rho_g_valid=.true.)
420
421 ELSEIF (dft_control%qs_control%rigpw) THEN
422 cpassert(.NOT. dft_control%use_kinetic_energy_density)
423 cpassert(.NOT. dft_control%drho_by_collocation)
424 CALL get_qs_env(qs_env, lri_env=lri_env)
425 CALL qs_rho_get(rho_struct, rho_ao=rho_ao)
426 CALL calculate_ri_densities(lri_env, qs_env, rho_ao, &
427 lri_rho_struct=rho_struct, &
428 atomic_kind_set=atomic_kind_set, &
429 para_env=para_env)
430 CALL qs_rho_set(rho_struct, rho_r_valid=.true., rho_g_valid=.true.)
431
432 ELSE
433 CALL qs_rho_update_rho_low(rho_struct=rho_struct, qs_env=qs_env, &
434 rho_xc_external=rho_xc_external, &
435 local_rho_set=local_rho_set, &
436 task_list_external=task_list_external, &
437 task_list_external_soft=task_list_external_soft, &
438 pw_env_external=pw_env_external, &
439 para_env_external=para_env_external)
440
441 END IF
442
443 END SUBROUTINE qs_rho_update_rho
444
445! **************************************************************************************************
446!> \brief updates rho_r and rho_g to the rho%rho_ao.
447!> if use_kinetic_energy_density also computes tau_r and tau_g
448!> \param rho_struct the rho structure that should be updated
449!> \param qs_env the qs_env rho_struct refers to
450!> the integrated charge in r space
451!> \param rho_xc_external rho structure for GAPW_XC
452!> \param local_rho_set ...
453!> \param pw_env_external external plane wave environment
454!> \param task_list_external external task list (use for default and GAPW)
455!> \param task_list_external_soft external task list (soft density for GAPW_XC)
456!> \param para_env_external ...
457!> \par History
458!> 08.2002 created [fawzi]
459!> \author Fawzi Mohamed
460! **************************************************************************************************
461 SUBROUTINE qs_rho_update_rho_low(rho_struct, qs_env, rho_xc_external, &
462 local_rho_set, pw_env_external, &
463 task_list_external, task_list_external_soft, &
464 para_env_external)
465 TYPE(qs_rho_type), INTENT(INOUT) :: rho_struct
466 TYPE(qs_environment_type), POINTER :: qs_env
467 TYPE(qs_rho_type), OPTIONAL, POINTER :: rho_xc_external
468 TYPE(local_rho_type), OPTIONAL, POINTER :: local_rho_set
469 TYPE(pw_env_type), OPTIONAL, POINTER :: pw_env_external
470 TYPE(task_list_type), OPTIONAL, POINTER :: task_list_external, &
471 task_list_external_soft
472 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env_external
473
474 CHARACTER(len=*), PARAMETER :: routinen = 'qs_rho_update_rho_low'
475
476 INTEGER :: handle, img, ispin, nimg, nspins
477 LOGICAL :: gapw, gapw_xc
478 REAL(kind=dp) :: dum
479 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_r, tot_rho_r_xc
480 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
481 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
482 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_kp, rho_xc_ao
483 TYPE(dft_control_type), POINTER :: dft_control
484 TYPE(mp_para_env_type), POINTER :: para_env
485 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
486 POINTER :: sab
487 TYPE(oce_matrix_type), POINTER :: oce
488 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g, rho_xc_g, tau_g, tau_xc_g
489 TYPE(pw_c1d_gs_type), DIMENSION(:, :), POINTER :: drho_g, drho_xc_g
490 TYPE(pw_env_type), POINTER :: pw_env
491 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r, rho_xc_r, tau_r, tau_xc_r
492 TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_r, drho_xc_r
493 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
494 TYPE(qs_ks_env_type), POINTER :: ks_env
495 TYPE(qs_rho_type), POINTER :: rho_xc
496 TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom_set
497 TYPE(task_list_type), POINTER :: task_list
498
499 CALL timeset(routinen, handle)
500
501 NULLIFY (dft_control, rho_xc, ks_env, rho_ao, rho_r, rho_g, drho_r, drho_g, tau_r, tau_g)
502 NULLIFY (rho_xc_ao, rho_xc_g, rho_xc_r, drho_xc_g, tau_xc_r, tau_xc_g, tot_rho_r, tot_rho_r_xc)
503 NULLIFY (para_env, pw_env, atomic_kind_set)
504
505 CALL get_qs_env(qs_env, &
506 ks_env=ks_env, &
507 dft_control=dft_control, &
508 atomic_kind_set=atomic_kind_set)
509
510 CALL qs_rho_get(rho_struct, &
511 rho_r=rho_r, &
512 rho_g=rho_g, &
513 tot_rho_r=tot_rho_r, &
514 drho_r=drho_r, &
515 drho_g=drho_g, &
516 tau_r=tau_r, &
517 tau_g=tau_g)
518
519 CALL get_qs_env(qs_env, task_list=task_list, &
520 para_env=para_env, pw_env=pw_env)
521 IF (PRESENT(pw_env_external)) pw_env => pw_env_external
522 IF (PRESENT(task_list_external)) task_list => task_list_external
523 IF (PRESENT(para_env_external)) para_env => para_env_external
524
525 nspins = dft_control%nspins
526 nimg = dft_control%nimages
527 gapw = dft_control%qs_control%gapw
528 gapw_xc = dft_control%qs_control%gapw_xc
529
530 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
531 DO ispin = 1, nspins
532 rho_ao => rho_ao_kp(ispin, :)
533 CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
534 rho=rho_r(ispin), &
535 rho_gspace=rho_g(ispin), &
536 total_rho=tot_rho_r(ispin), &
537 ks_env=ks_env, soft_valid=gapw, &
538 task_list_external=task_list_external, &
539 pw_env_external=pw_env_external)
540 END DO
541 CALL qs_rho_set(rho_struct, rho_r_valid=.true., rho_g_valid=.true.)
542
543 IF (gapw_xc) THEN
544 IF (PRESENT(rho_xc_external)) THEN
545 rho_xc => rho_xc_external
546 ELSE
547 CALL get_qs_env(qs_env=qs_env, rho_xc=rho_xc)
548 END IF
549 CALL qs_rho_get(rho_xc, &
550 rho_ao_kp=rho_xc_ao, &
551 rho_r=rho_xc_r, &
552 rho_g=rho_xc_g, &
553 tot_rho_r=tot_rho_r_xc)
554 ! copy rho_ao into rho_xc_ao
555 DO ispin = 1, nspins
556 DO img = 1, nimg
557 CALL dbcsr_copy(rho_xc_ao(ispin, img)%matrix, rho_ao_kp(ispin, img)%matrix)
558 END DO
559 END DO
560 DO ispin = 1, nspins
561 rho_ao => rho_xc_ao(ispin, :)
562 CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
563 rho=rho_xc_r(ispin), &
564 rho_gspace=rho_xc_g(ispin), &
565 total_rho=tot_rho_r_xc(ispin), &
566 ks_env=ks_env, soft_valid=gapw_xc, &
567 task_list_external=task_list_external_soft, &
568 pw_env_external=pw_env_external)
569 END DO
570 CALL qs_rho_set(rho_xc, rho_r_valid=.true., rho_g_valid=.true.)
571 END IF
572
573 ! GAPW o GAPW_XC require the calculation of hard and soft local densities
574 IF (gapw .OR. gapw_xc) THEN
575 CALL get_qs_env(qs_env=qs_env, &
576 rho_atom_set=rho_atom_set, &
577 qs_kind_set=qs_kind_set, &
578 oce=oce, sab_orb=sab)
579 IF (PRESENT(local_rho_set)) rho_atom_set => local_rho_set%rho_atom_set
580 cpassert(ASSOCIATED(rho_atom_set))
581 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
582 CALL calculate_rho_atom_coeff(qs_env, rho_ao_kp, rho_atom_set, qs_kind_set, oce, sab, para_env)
583 END IF
584
585 IF (.NOT. gapw_xc) THEN
586 ! if needed compute also the gradient of the density
587 IF (dft_control%drho_by_collocation) THEN
588 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
589 cpassert(.NOT. PRESENT(task_list_external))
590 DO ispin = 1, nspins
591 rho_ao => rho_ao_kp(ispin, :)
592 CALL calculate_drho_elec(matrix_p_kp=rho_ao, &
593 drho=drho_r(:, ispin), &
594 drho_gspace=drho_g(:, ispin), &
595 qs_env=qs_env, soft_valid=gapw)
596 END DO
597 CALL qs_rho_set(rho_struct, drho_r_valid=.true., drho_g_valid=.true.)
598 END IF
599 ! if needed compute also the kinetic energy density
600 IF (dft_control%use_kinetic_energy_density) THEN
601 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
602 DO ispin = 1, nspins
603 rho_ao => rho_ao_kp(ispin, :)
604 CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
605 rho=tau_r(ispin), &
606 rho_gspace=tau_g(ispin), &
607 total_rho=dum, & ! presumably not meaningful
608 ks_env=ks_env, soft_valid=gapw, &
609 compute_tau=.true., &
610 task_list_external=task_list_external, &
611 pw_env_external=pw_env_external)
612 END DO
613 CALL qs_rho_set(rho_struct, tau_r_valid=.true., tau_g_valid=.true.)
614 END IF
615 ELSE
616 CALL qs_rho_get(rho_xc, &
617 drho_r=drho_xc_r, &
618 drho_g=drho_xc_g, &
619 tau_r=tau_xc_r, &
620 tau_g=tau_xc_g)
621 ! if needed compute also the gradient of the density
622 IF (dft_control%drho_by_collocation) THEN
623 cpassert(.NOT. PRESENT(task_list_external))
624 DO ispin = 1, nspins
625 rho_ao => rho_xc_ao(ispin, :)
626 CALL calculate_drho_elec(matrix_p_kp=rho_ao, &
627 drho=drho_xc_r(:, ispin), &
628 drho_gspace=drho_xc_g(:, ispin), &
629 qs_env=qs_env, soft_valid=gapw_xc)
630 END DO
631 CALL qs_rho_set(rho_xc, drho_r_valid=.true., drho_g_valid=.true.)
632 END IF
633 ! if needed compute also the kinetic energy density
634 IF (dft_control%use_kinetic_energy_density) THEN
635 DO ispin = 1, nspins
636 rho_ao => rho_xc_ao(ispin, :)
637 CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
638 rho=tau_xc_r(ispin), &
639 rho_gspace=tau_xc_g(ispin), &
640 ks_env=ks_env, soft_valid=gapw_xc, &
641 compute_tau=.true., &
642 task_list_external=task_list_external_soft, &
643 pw_env_external=pw_env_external)
644 END DO
645 CALL qs_rho_set(rho_xc, tau_r_valid=.true., tau_g_valid=.true.)
646 END IF
647 END IF
648
649 CALL timestop(handle)
650
651 END SUBROUTINE qs_rho_update_rho_low
652
653! **************************************************************************************************
654!> \brief updates rho_r and rho_g to the rho%rho_ao.
655!> if use_kinetic_energy_density also computes tau_r and tau_g
656!> \param rho_struct the rho structure that should be updated
657!> \param qs_env the qs_env rho_struct refers to
658!> the integrated charge in r space
659!> \param pw_env_external external plane wave environment
660!> \param task_list_external external task list
661!> \param para_env_external ...
662!> \param tddfpt_lri_env ...
663!> \param tddfpt_lri_density ...
664!> \par History
665!> 08.2002 created [fawzi]
666!> \author Fawzi Mohamed
667! **************************************************************************************************
668 SUBROUTINE qs_rho_update_tddfpt(rho_struct, qs_env, pw_env_external, task_list_external, &
669 para_env_external, tddfpt_lri_env, tddfpt_lri_density)
670 TYPE(qs_rho_type), INTENT(INOUT) :: rho_struct
671 TYPE(qs_environment_type), POINTER :: qs_env
672 TYPE(pw_env_type), OPTIONAL, POINTER :: pw_env_external
673 TYPE(task_list_type), OPTIONAL, POINTER :: task_list_external
674 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env_external
675 TYPE(lri_environment_type), OPTIONAL, POINTER :: tddfpt_lri_env
676 TYPE(lri_density_type), OPTIONAL, POINTER :: tddfpt_lri_density
677
678 CHARACTER(len=*), PARAMETER :: routinen = 'qs_rho_update_tddfpt'
679
680 INTEGER :: handle, ispin, nspins
681 INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
682 LOGICAL :: lri_response
683 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_r
684 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
685 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
686 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_kp
687 TYPE(dft_control_type), POINTER :: dft_control
688 TYPE(kpoint_type), POINTER :: kpoints
689 TYPE(mp_para_env_type), POINTER :: para_env
690 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g
691 TYPE(pw_env_type), POINTER :: pw_env
692 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r
693 TYPE(qs_ks_env_type), POINTER :: ks_env
694 TYPE(task_list_type), POINTER :: task_list
695
696 CALL timeset(routinen, handle)
697
698 CALL get_qs_env(qs_env, &
699 ks_env=ks_env, &
700 dft_control=dft_control, &
701 atomic_kind_set=atomic_kind_set, &
702 task_list=task_list, &
703 para_env=para_env, &
704 pw_env=pw_env)
705 IF (PRESENT(pw_env_external)) pw_env => pw_env_external
706 IF (PRESENT(task_list_external)) task_list => task_list_external
707 IF (PRESENT(para_env_external)) para_env => para_env_external
708
709 CALL qs_rho_get(rho_struct, &
710 rho_r=rho_r, &
711 rho_g=rho_g, &
712 tot_rho_r=tot_rho_r)
713
714 nspins = dft_control%nspins
715
716 lri_response = PRESENT(tddfpt_lri_env)
717 IF (lri_response) THEN
718 cpassert(PRESENT(tddfpt_lri_density))
719 END IF
720
721 cpassert(.NOT. dft_control%drho_by_collocation)
722 cpassert(.NOT. dft_control%use_kinetic_energy_density)
723 cpassert(.NOT. dft_control%qs_control%gapw)
724 cpassert(.NOT. dft_control%qs_control%gapw_xc)
725
726 IF (lri_response) THEN
727 CALL get_ks_env(ks_env=ks_env, kpoints=kpoints)
728 CALL get_kpoint_info(kpoint=kpoints, cell_to_index=cell_to_index)
729 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
730 CALL calculate_lri_densities(tddfpt_lri_env, tddfpt_lri_density, qs_env, rho_ao_kp, cell_to_index, &
731 lri_rho_struct=rho_struct, &
732 atomic_kind_set=atomic_kind_set, &
733 para_env=para_env, &
734 response_density=lri_response)
735 CALL qs_rho_set(rho_struct, rho_r_valid=.true., rho_g_valid=.true.)
736 ELSE
737 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
738 DO ispin = 1, nspins
739 rho_ao => rho_ao_kp(ispin, :)
740 CALL calculate_rho_elec(matrix_p_kp=rho_ao, &
741 rho=rho_r(ispin), &
742 rho_gspace=rho_g(ispin), &
743 total_rho=tot_rho_r(ispin), &
744 ks_env=ks_env, &
745 task_list_external=task_list_external, &
746 pw_env_external=pw_env_external)
747 END DO
748 CALL qs_rho_set(rho_struct, rho_r_valid=.true., rho_g_valid=.true.)
749 END IF
750
751 CALL timestop(handle)
752
753 END SUBROUTINE qs_rho_update_tddfpt
754
755! **************************************************************************************************
756!> \brief Allocate a density structure and fill it with data from an input structure
757!> SIZE(rho_input) == mspin == 1 direct copy
758!> SIZE(rho_input) == mspin == 2 direct copy of alpha and beta spin
759!> SIZE(rho_input) == 1 AND mspin == 2 copy rho/2 into alpha and beta spin
760!> \param rho_input ...
761!> \param rho_output ...
762!> \param auxbas_pw_pool ...
763!> \param mspin ...
764! **************************************************************************************************
765 SUBROUTINE qs_rho_copy(rho_input, rho_output, auxbas_pw_pool, mspin)
766
767 TYPE(qs_rho_type), INTENT(IN) :: rho_input
768 TYPE(qs_rho_type), INTENT(INOUT) :: rho_output
769 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
770 INTEGER, INTENT(IN) :: mspin
771
772 CHARACTER(len=*), PARAMETER :: routinen = 'qs_rho_copy'
773
774 INTEGER :: handle, i, j, nspins
775 LOGICAL :: complex_rho_ao, drho_g_valid_in, drho_r_valid_in, rho_g_valid_in, rho_r_valid_in, &
776 soft_valid_in, tau_g_valid_in, tau_r_valid_in
777 REAL(kind=dp) :: ospin
778 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_g_in, tot_rho_g_out, &
779 tot_rho_r_in, tot_rho_r_out
780 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao_im_in, rho_ao_im_out, rho_ao_in, &
781 rho_ao_out
782 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_kp_in
783 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_in, rho_g_out, tau_g_in, tau_g_out
784 TYPE(pw_c1d_gs_type), DIMENSION(:, :), POINTER :: drho_g_in, drho_g_out
785 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_in, rho_r_out, tau_r_in, tau_r_out
786 TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_r_in, drho_r_out
787 TYPE(pw_r3d_rs_type), POINTER :: rho_r_sccs_in, rho_r_sccs_out
788
789 CALL timeset(routinen, handle)
790
791 cpassert(mspin == 1 .OR. mspin == 2)
792 ospin = 1._dp/real(mspin, kind=dp)
793
794 CALL qs_rho_clear(rho_output)
795
796 NULLIFY (rho_ao_in, rho_ao_kp_in, rho_ao_im_in, rho_r_in, rho_g_in, drho_r_in, &
797 drho_g_in, tau_r_in, tau_g_in, tot_rho_r_in, tot_rho_g_in, rho_r_sccs_in)
798
799 CALL qs_rho_get(rho_input, &
800 rho_ao=rho_ao_in, &
801 rho_ao_kp=rho_ao_kp_in, &
802 rho_ao_im=rho_ao_im_in, &
803 rho_r=rho_r_in, &
804 rho_g=rho_g_in, &
805 drho_r=drho_r_in, &
806 drho_g=drho_g_in, &
807 tau_r=tau_r_in, &
808 tau_g=tau_g_in, &
809 tot_rho_r=tot_rho_r_in, &
810 tot_rho_g=tot_rho_g_in, &
811 rho_g_valid=rho_g_valid_in, &
812 rho_r_valid=rho_r_valid_in, &
813 drho_g_valid=drho_g_valid_in, &
814 drho_r_valid=drho_r_valid_in, &
815 tau_r_valid=tau_r_valid_in, &
816 tau_g_valid=tau_g_valid_in, &
817 rho_r_sccs=rho_r_sccs_in, &
818 soft_valid=soft_valid_in, &
819 complex_rho_ao=complex_rho_ao)
820
821 NULLIFY (rho_ao_out, rho_ao_im_out, rho_r_out, rho_g_out, drho_r_out, &
822 drho_g_out, tau_r_out, tau_g_out, tot_rho_r_out, tot_rho_g_out, rho_r_sccs_out)
823 ! rho_ao
824 IF (ASSOCIATED(rho_ao_in)) THEN
825 nspins = SIZE(rho_ao_in)
826 cpassert(mspin >= nspins)
827 CALL dbcsr_allocate_matrix_set(rho_ao_out, mspin)
828 CALL qs_rho_set(rho_output, rho_ao=rho_ao_out)
829 IF (mspin > nspins) THEN
830 DO i = 1, mspin
831 ALLOCATE (rho_ao_out(i)%matrix)
832 CALL dbcsr_copy(rho_ao_out(i)%matrix, rho_ao_in(1)%matrix, name="RHO copy")
833 CALL dbcsr_scale(rho_ao_out(i)%matrix, ospin)
834 END DO
835 ELSE
836 DO i = 1, nspins
837 ALLOCATE (rho_ao_out(i)%matrix)
838 CALL dbcsr_copy(rho_ao_out(i)%matrix, rho_ao_in(i)%matrix, name="RHO copy")
839 END DO
840 END IF
841 END IF
842
843 ! rho_ao_kp
844 ! only for non-kp, we could probably just copy this pointer, should work also for non-kp?
845 !IF (ASSOCIATED(rho_ao_kp_in)) THEN
846 ! CPABORT("Copy not available")
847 !END IF
848
849 ! rho_ao_im
850 IF (ASSOCIATED(rho_ao_im_in)) THEN
851 nspins = SIZE(rho_ao_im_in)
852 cpassert(mspin >= nspins)
853 CALL dbcsr_allocate_matrix_set(rho_ao_im_out, mspin)
854 CALL qs_rho_set(rho_output, rho_ao_im=rho_ao_im_out)
855 IF (mspin > nspins) THEN
856 DO i = 1, mspin
857 ALLOCATE (rho_ao_im_out(i)%matrix)
858 CALL dbcsr_copy(rho_ao_im_out(i)%matrix, rho_ao_im_in(1)%matrix, name="RHO copy")
859 CALL dbcsr_scale(rho_ao_im_out(i)%matrix, ospin)
860 END DO
861 ELSE
862 DO i = 1, nspins
863 ALLOCATE (rho_ao_im_out(i)%matrix)
864 CALL dbcsr_copy(rho_ao_im_out(i)%matrix, rho_ao_im_in(i)%matrix, name="RHO copy")
865 END DO
866 END IF
867 END IF
868
869 ! rho_r
870 IF (ASSOCIATED(rho_r_in)) THEN
871 nspins = SIZE(rho_r_in)
872 cpassert(mspin >= nspins)
873 ALLOCATE (rho_r_out(mspin))
874 CALL qs_rho_set(rho_output, rho_r=rho_r_out)
875 IF (mspin > nspins) THEN
876 DO i = 1, mspin
877 CALL auxbas_pw_pool%create_pw(rho_r_out(i))
878 CALL pw_copy(rho_r_in(1), rho_r_out(i))
879 CALL pw_scale(rho_r_out(i), ospin)
880 END DO
881 ELSE
882 DO i = 1, nspins
883 CALL auxbas_pw_pool%create_pw(rho_r_out(i))
884 CALL pw_copy(rho_r_in(i), rho_r_out(i))
885 END DO
886 END IF
887 END IF
888
889 ! rho_g
890 IF (ASSOCIATED(rho_g_in)) THEN
891 nspins = SIZE(rho_g_in)
892 cpassert(mspin >= nspins)
893 ALLOCATE (rho_g_out(mspin))
894 CALL qs_rho_set(rho_output, rho_g=rho_g_out)
895 IF (mspin > nspins) THEN
896 DO i = 1, mspin
897 CALL auxbas_pw_pool%create_pw(rho_g_out(i))
898 CALL pw_copy(rho_g_in(1), rho_g_out(i))
899 CALL pw_scale(rho_g_out(i), ospin)
900 END DO
901 ELSE
902 DO i = 1, nspins
903 CALL auxbas_pw_pool%create_pw(rho_g_out(i))
904 CALL pw_copy(rho_g_in(i), rho_g_out(i))
905 END DO
906 END IF
907 END IF
908
909 ! SCCS
910 IF (ASSOCIATED(rho_r_sccs_in)) THEN
911 CALL qs_rho_set(rho_output, rho_r_sccs=rho_r_sccs_out)
912 CALL auxbas_pw_pool%create_pw(rho_r_sccs_out)
913 CALL pw_copy(rho_r_sccs_in, rho_r_sccs_out)
914 END IF
915
916 ! drho_r
917 IF (ASSOCIATED(drho_r_in)) THEN
918 nspins = SIZE(drho_r_in)
919 cpassert(mspin >= nspins)
920 ALLOCATE (drho_r_out(3, mspin))
921 CALL qs_rho_set(rho_output, drho_r=drho_r_out)
922 IF (mspin > nspins) THEN
923 DO j = 1, mspin
924 DO i = 1, 3
925 CALL auxbas_pw_pool%create_pw(drho_r_out(i, j))
926 CALL pw_copy(drho_r_in(i, 1), drho_r_out(i, j))
927 CALL pw_scale(drho_r_out(i, j), ospin)
928 END DO
929 END DO
930 ELSE
931 DO j = 1, nspins
932 DO i = 1, 3
933 CALL auxbas_pw_pool%create_pw(drho_r_out(i, j))
934 CALL pw_copy(drho_r_in(i, j), drho_r_out(i, j))
935 END DO
936 END DO
937 END IF
938 END IF
939
940 ! drho_g
941 IF (ASSOCIATED(drho_g_in)) THEN
942 nspins = SIZE(drho_g_in)
943 cpassert(mspin >= nspins)
944 ALLOCATE (drho_g_out(3, mspin))
945 CALL qs_rho_set(rho_output, drho_g=drho_g_out)
946 IF (mspin > nspins) THEN
947 DO j = 1, mspin
948 DO i = 1, 3
949 CALL auxbas_pw_pool%create_pw(drho_g_out(i, j))
950 CALL pw_copy(drho_g_in(i, 1), drho_g_out(i, j))
951 CALL pw_scale(drho_g_out(i, j), ospin)
952 END DO
953 END DO
954 ELSE
955 DO j = 1, nspins
956 DO i = 1, 3
957 CALL auxbas_pw_pool%create_pw(drho_g_out(i, j))
958 CALL pw_copy(drho_g_in(i, j), drho_g_out(i, j))
959 END DO
960 END DO
961 END IF
962 END IF
963
964 ! tau_r
965 IF (ASSOCIATED(tau_r_in)) THEN
966 nspins = SIZE(tau_r_in)
967 cpassert(mspin >= nspins)
968 ALLOCATE (tau_r_out(mspin))
969 CALL qs_rho_set(rho_output, tau_r=tau_r_out)
970 IF (mspin > nspins) THEN
971 DO i = 1, mspin
972 CALL auxbas_pw_pool%create_pw(tau_r_out(i))
973 CALL pw_copy(tau_r_in(1), tau_r_out(i))
974 CALL pw_scale(tau_r_out(i), ospin)
975 END DO
976 ELSE
977 DO i = 1, nspins
978 CALL auxbas_pw_pool%create_pw(tau_r_out(i))
979 CALL pw_copy(tau_r_in(i), tau_r_out(i))
980 END DO
981 END IF
982 END IF
983
984 ! tau_g
985 IF (ASSOCIATED(tau_g_in)) THEN
986 nspins = SIZE(tau_g_in)
987 cpassert(mspin >= nspins)
988 ALLOCATE (tau_g_out(mspin))
989 CALL qs_rho_set(rho_output, tau_g=tau_g_out)
990 IF (mspin > nspins) THEN
991 DO i = 1, mspin
992 CALL auxbas_pw_pool%create_pw(tau_g_out(i))
993 CALL pw_copy(tau_g_in(1), tau_g_out(i))
994 CALL pw_scale(tau_g_out(i), ospin)
995 END DO
996 ELSE
997 DO i = 1, nspins
998 CALL auxbas_pw_pool%create_pw(tau_g_out(i))
999 CALL pw_copy(tau_g_in(i), tau_g_out(i))
1000 END DO
1001 END IF
1002 END IF
1003
1004 ! tot_rho_r
1005 IF (ASSOCIATED(tot_rho_r_in)) THEN
1006 nspins = SIZE(tot_rho_r_in)
1007 cpassert(mspin >= nspins)
1008 ALLOCATE (tot_rho_r_out(mspin))
1009 CALL qs_rho_set(rho_output, tot_rho_r=tot_rho_r_out)
1010 IF (mspin > nspins) THEN
1011 DO i = 1, mspin
1012 tot_rho_r_out(i) = tot_rho_r_in(1)*ospin
1013 END DO
1014 ELSE
1015 DO i = 1, nspins
1016 tot_rho_r_out(i) = tot_rho_r_in(i)
1017 END DO
1018 END IF
1019 END IF
1020
1021 ! tot_rho_g
1022 IF (ASSOCIATED(tot_rho_g_in)) THEN
1023 nspins = SIZE(tot_rho_g_in)
1024 cpassert(mspin >= nspins)
1025 ALLOCATE (tot_rho_g_out(mspin))
1026 CALL qs_rho_set(rho_output, tot_rho_g=tot_rho_g_out)
1027 IF (mspin > nspins) THEN
1028 DO i = 1, mspin
1029 tot_rho_g_out(i) = tot_rho_g_in(1)*ospin
1030 END DO
1031 ELSE
1032 DO i = 1, nspins
1033 tot_rho_g_out(i) = tot_rho_g_in(i)
1034 END DO
1035 END IF
1036 END IF
1037
1038 CALL qs_rho_set(rho_output, &
1039 rho_g_valid=rho_g_valid_in, &
1040 rho_r_valid=rho_r_valid_in, &
1041 drho_g_valid=drho_g_valid_in, &
1042 drho_r_valid=drho_r_valid_in, &
1043 tau_r_valid=tau_r_valid_in, &
1044 tau_g_valid=tau_g_valid_in, &
1045 soft_valid=soft_valid_in, &
1046 complex_rho_ao=complex_rho_ao)
1047
1048 CALL timestop(handle)
1049
1050 END SUBROUTINE qs_rho_copy
1051
1052! **************************************************************************************************
1053!> \brief rhoa = alpha*rhoa+beta*rhob
1054!> \param rhoa ...
1055!> \param rhob ...
1056!> \param alpha ...
1057!> \param beta ...
1058! **************************************************************************************************
1059 SUBROUTINE qs_rho_scale_and_add(rhoa, rhob, alpha, beta)
1060
1061 TYPE(qs_rho_type), INTENT(IN) :: rhoa, rhob
1062 REAL(kind=dp), INTENT(IN) :: alpha, beta
1063
1064 CHARACTER(len=*), PARAMETER :: routinen = 'qs_rho_scale_and_add'
1065
1066 INTEGER :: handle, i, j, nspina, nspinb, nspins
1067 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_g_a, tot_rho_g_b, tot_rho_r_a, &
1068 tot_rho_r_b
1069 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao_a, rho_ao_b, rho_ao_im_a, &
1070 rho_ao_im_b
1071 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_a, rho_g_b, tau_g_a, tau_g_b
1072 TYPE(pw_c1d_gs_type), DIMENSION(:, :), POINTER :: drho_g_a, drho_g_b
1073 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_a, rho_r_b, tau_r_a, tau_r_b
1074 TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_r_a, drho_r_b
1075 TYPE(pw_r3d_rs_type), POINTER :: rho_r_sccs_a, rho_r_sccs_b
1076
1077 CALL timeset(routinen, handle)
1078
1079 NULLIFY (rho_ao_a, rho_ao_im_a, rho_r_a, rho_g_a, drho_r_a, &
1080 drho_g_a, tau_r_a, tau_g_a, tot_rho_r_a, tot_rho_g_a, rho_r_sccs_a)
1081
1082 CALL qs_rho_get(rhoa, &
1083 rho_ao=rho_ao_a, &
1084 rho_ao_im=rho_ao_im_a, &
1085 rho_r=rho_r_a, &
1086 rho_g=rho_g_a, &
1087 drho_r=drho_r_a, &
1088 drho_g=drho_g_a, &
1089 tau_r=tau_r_a, &
1090 tau_g=tau_g_a, &
1091 tot_rho_r=tot_rho_r_a, &
1092 tot_rho_g=tot_rho_g_a, &
1093 rho_r_sccs=rho_r_sccs_a)
1094
1095 NULLIFY (rho_ao_b, rho_ao_im_b, rho_r_b, rho_g_b, drho_r_b, &
1096 drho_g_b, tau_r_b, tau_g_b, tot_rho_r_b, tot_rho_g_b, rho_r_sccs_b)
1097
1098 CALL qs_rho_get(rhob, &
1099 rho_ao=rho_ao_b, &
1100 rho_ao_im=rho_ao_im_b, &
1101 rho_r=rho_r_b, &
1102 rho_g=rho_g_b, &
1103 drho_r=drho_r_b, &
1104 drho_g=drho_g_b, &
1105 tau_r=tau_r_b, &
1106 tau_g=tau_g_b, &
1107 tot_rho_r=tot_rho_r_b, &
1108 tot_rho_g=tot_rho_g_b, &
1109 rho_r_sccs=rho_r_sccs_b)
1110 ! rho_ao
1111 IF (ASSOCIATED(rho_ao_a) .AND. ASSOCIATED(rho_ao_b)) THEN
1112 nspina = SIZE(rho_ao_a)
1113 nspinb = SIZE(rho_ao_b)
1114 nspins = min(nspina, nspinb)
1115 DO i = 1, nspins
1116 CALL dbcsr_add(rho_ao_a(i)%matrix, rho_ao_b(i)%matrix, alpha, beta)
1117 END DO
1118 END IF
1119
1120 ! rho_ao_im
1121 IF (ASSOCIATED(rho_ao_im_a) .AND. ASSOCIATED(rho_ao_im_b)) THEN
1122 nspina = SIZE(rho_ao_im_a)
1123 nspinb = SIZE(rho_ao_im_b)
1124 nspins = min(nspina, nspinb)
1125 DO i = 1, nspins
1126 CALL dbcsr_add(rho_ao_im_a(i)%matrix, rho_ao_im_b(i)%matrix, alpha, beta)
1127 END DO
1128 END IF
1129
1130 ! rho_r
1131 IF (ASSOCIATED(rho_r_a) .AND. ASSOCIATED(rho_r_b)) THEN
1132 nspina = SIZE(rho_ao_a)
1133 nspinb = SIZE(rho_ao_b)
1134 nspins = min(nspina, nspinb)
1135 DO i = 1, nspins
1136 CALL pw_axpy(rho_r_b(i), rho_r_a(i), beta, alpha)
1137 END DO
1138 END IF
1139
1140 ! rho_g
1141 IF (ASSOCIATED(rho_g_a) .AND. ASSOCIATED(rho_g_b)) THEN
1142 nspina = SIZE(rho_ao_a)
1143 nspinb = SIZE(rho_ao_b)
1144 nspins = min(nspina, nspinb)
1145 DO i = 1, nspins
1146 CALL pw_axpy(rho_g_b(i), rho_g_a(i), beta, alpha)
1147 END DO
1148 END IF
1149
1150 ! SCCS
1151 IF (ASSOCIATED(rho_r_sccs_a) .AND. ASSOCIATED(rho_r_sccs_b)) THEN
1152 CALL pw_axpy(rho_r_sccs_b, rho_r_sccs_a, beta, alpha)
1153 END IF
1154
1155 ! drho_r
1156 IF (ASSOCIATED(drho_r_a) .AND. ASSOCIATED(drho_r_b)) THEN
1157 cpassert(all(shape(drho_r_a) == shape(drho_r_b))) ! not implemented
1158 DO j = 1, SIZE(drho_r_a, 2)
1159 DO i = 1, SIZE(drho_r_a, 1)
1160 CALL pw_axpy(drho_r_b(i, j), drho_r_a(i, j), beta, alpha)
1161 END DO
1162 END DO
1163 END IF
1164
1165 ! drho_g
1166 IF (ASSOCIATED(drho_g_a) .AND. ASSOCIATED(drho_g_b)) THEN
1167 cpassert(all(shape(drho_g_a) == shape(drho_g_b))) ! not implemented
1168 DO j = 1, SIZE(drho_g_a, 2)
1169 DO i = 1, SIZE(drho_g_a, 1)
1170 CALL pw_axpy(drho_g_b(i, j), drho_g_a(i, j), beta, alpha)
1171 END DO
1172 END DO
1173 END IF
1174
1175 ! tau_r
1176 IF (ASSOCIATED(tau_r_a) .AND. ASSOCIATED(tau_r_b)) THEN
1177 nspina = SIZE(rho_ao_a)
1178 nspinb = SIZE(rho_ao_b)
1179 nspins = min(nspina, nspinb)
1180 DO i = 1, nspins
1181 CALL pw_axpy(tau_r_b(i), tau_r_a(i), beta, alpha)
1182 END DO
1183 END IF
1184
1185 ! tau_g
1186 IF (ASSOCIATED(tau_g_a) .AND. ASSOCIATED(tau_g_b)) THEN
1187 nspina = SIZE(rho_ao_a)
1188 nspinb = SIZE(rho_ao_b)
1189 nspins = min(nspina, nspinb)
1190 DO i = 1, nspins
1191 CALL pw_axpy(tau_g_b(i), tau_g_a(i), beta, alpha)
1192 END DO
1193 END IF
1194
1195 ! tot_rho_r
1196 IF (ASSOCIATED(tot_rho_r_a) .AND. ASSOCIATED(tot_rho_r_b)) THEN
1197 nspina = SIZE(rho_ao_a)
1198 nspinb = SIZE(rho_ao_b)
1199 nspins = min(nspina, nspinb)
1200 DO i = 1, nspins
1201 tot_rho_r_a(i) = alpha*tot_rho_r_a(i) + beta*tot_rho_r_b(i)
1202 END DO
1203 END IF
1204
1205 ! tot_rho_g
1206 IF (ASSOCIATED(tot_rho_g_a) .AND. ASSOCIATED(tot_rho_g_b)) THEN
1207 nspina = SIZE(rho_ao_a)
1208 nspinb = SIZE(rho_ao_b)
1209 nspins = min(nspina, nspinb)
1210 DO i = 1, nspins
1211 tot_rho_g_a(i) = alpha*tot_rho_g_a(i) + beta*tot_rho_g_b(i)
1212 END DO
1213 END IF
1214
1215 CALL timestop(handle)
1216
1217 END SUBROUTINE qs_rho_scale_and_add
1218
1219! **************************************************************************************************
1220!> \brief Duplicates a pointer physically
1221!> \param rho_input The rho structure to be duplicated
1222!> \param rho_output The duplicate rho structure
1223!> \param qs_env The QS environment from which the auxiliary PW basis-set
1224!> pool is taken
1225!> \par History
1226!> 07.2005 initial create [tdk]
1227!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
1228!> \note
1229!> Associated pointers are deallocated, nullified pointers are NOT accepted!
1230! **************************************************************************************************
1231 SUBROUTINE duplicate_rho_type(rho_input, rho_output, qs_env)
1232
1233 TYPE(qs_rho_type), INTENT(INOUT) :: rho_input, rho_output
1234 TYPE(qs_environment_type), POINTER :: qs_env
1235
1236 CHARACTER(len=*), PARAMETER :: routinen = 'duplicate_rho_type'
1237
1238 INTEGER :: handle, i, j, nspins
1239 LOGICAL :: complex_rho_ao_in, drho_g_valid_in, drho_r_valid_in, rho_g_valid_in, &
1240 rho_r_valid_in, soft_valid_in, tau_g_valid_in, tau_r_valid_in
1241 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_g_in, tot_rho_g_out, &
1242 tot_rho_r_in, tot_rho_r_out
1243 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao_im_in, rho_ao_im_out, rho_ao_in, &
1244 rho_ao_out
1245 TYPE(dft_control_type), POINTER :: dft_control
1246 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_in, rho_g_out, tau_g_in, tau_g_out
1247 TYPE(pw_c1d_gs_type), DIMENSION(:, :), POINTER :: drho_g_in, drho_g_out
1248 TYPE(pw_env_type), POINTER :: pw_env
1249 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
1250 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_in, rho_r_out, tau_r_in, tau_r_out
1251 TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_r_in, drho_r_out
1252 TYPE(pw_r3d_rs_type), POINTER :: rho_r_sccs_in, rho_r_sccs_out
1253
1254 CALL timeset(routinen, handle)
1255
1256 NULLIFY (dft_control, pw_env, auxbas_pw_pool)
1257 NULLIFY (rho_ao_in, rho_ao_out, rho_ao_im_in, rho_ao_im_out)
1258 NULLIFY (rho_r_in, rho_r_out, rho_g_in, rho_g_out, drho_r_in, drho_r_out)
1259 NULLIFY (drho_g_in, drho_g_out, tau_r_in, tau_r_out, tau_g_in, tau_g_out)
1260 NULLIFY (tot_rho_r_in, tot_rho_r_out, tot_rho_g_in, tot_rho_g_out)
1261 NULLIFY (rho_r_sccs_in, rho_r_sccs_out)
1262
1263 cpassert(ASSOCIATED(qs_env))
1264
1265 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env, dft_control=dft_control)
1266 CALL pw_env_get(pw_env=pw_env, auxbas_pw_pool=auxbas_pw_pool)
1267 nspins = dft_control%nspins
1268
1269 CALL qs_rho_clear(rho_output)
1270
1271 CALL qs_rho_get(rho_input, &
1272 rho_ao=rho_ao_in, &
1273 rho_ao_im=rho_ao_im_in, &
1274 rho_r=rho_r_in, &
1275 rho_g=rho_g_in, &
1276 drho_r=drho_r_in, &
1277 drho_g=drho_g_in, &
1278 tau_r=tau_r_in, &
1279 tau_g=tau_g_in, &
1280 tot_rho_r=tot_rho_r_in, &
1281 tot_rho_g=tot_rho_g_in, &
1282 rho_g_valid=rho_g_valid_in, &
1283 rho_r_valid=rho_r_valid_in, &
1284 drho_g_valid=drho_g_valid_in, &
1285 drho_r_valid=drho_r_valid_in, &
1286 tau_r_valid=tau_r_valid_in, &
1287 tau_g_valid=tau_g_valid_in, &
1288 rho_r_sccs=rho_r_sccs_in, &
1289 soft_valid=soft_valid_in, &
1290 complex_rho_ao=complex_rho_ao_in)
1291
1292 ! rho_ao
1293 IF (ASSOCIATED(rho_ao_in)) THEN
1294 CALL dbcsr_allocate_matrix_set(rho_ao_out, nspins)
1295 CALL qs_rho_set(rho_output, rho_ao=rho_ao_out)
1296 DO i = 1, nspins
1297 ALLOCATE (rho_ao_out(i)%matrix)
1298 CALL dbcsr_copy(rho_ao_out(i)%matrix, rho_ao_in(i)%matrix, &
1299 name="myDensityMatrix_for_Spin_"//trim(adjustl(cp_to_string(i))))
1300 CALL dbcsr_set(rho_ao_out(i)%matrix, 0.0_dp)
1301 END DO
1302 END IF
1303
1304 ! rho_ao_im
1305 IF (ASSOCIATED(rho_ao_im_in)) THEN
1306 CALL dbcsr_allocate_matrix_set(rho_ao_im_out, nspins)
1307 CALL qs_rho_set(rho_output, rho_ao=rho_ao_im_out)
1308 DO i = 1, nspins
1309 ALLOCATE (rho_ao_im_out(i)%matrix)
1310 CALL dbcsr_copy(rho_ao_im_out(i)%matrix, rho_ao_im_in(i)%matrix, &
1311 name="myImagDensityMatrix_for_Spin_"//trim(adjustl(cp_to_string(i))))
1312 CALL dbcsr_set(rho_ao_im_out(i)%matrix, 0.0_dp)
1313 END DO
1314 END IF
1315
1316 ! rho_r
1317 IF (ASSOCIATED(rho_r_in)) THEN
1318 ALLOCATE (rho_r_out(nspins))
1319 CALL qs_rho_set(rho_output, rho_r=rho_r_out)
1320 DO i = 1, nspins
1321 CALL auxbas_pw_pool%create_pw(rho_r_out(i))
1322 CALL pw_copy(rho_r_in(i), rho_r_out(i))
1323 END DO
1324 END IF
1325
1326 ! rho_g
1327 IF (ASSOCIATED(rho_g_in)) THEN
1328 ALLOCATE (rho_g_out(nspins))
1329 CALL qs_rho_set(rho_output, rho_g=rho_g_out)
1330 DO i = 1, nspins
1331 CALL auxbas_pw_pool%create_pw(rho_g_out(i))
1332 CALL pw_copy(rho_g_in(i), rho_g_out(i))
1333 END DO
1334 END IF
1335
1336 ! SCCS
1337 IF (ASSOCIATED(rho_r_sccs_in)) THEN
1338 CALL qs_rho_set(rho_output, rho_r_sccs=rho_r_sccs_out)
1339 CALL auxbas_pw_pool%create_pw(rho_r_sccs_out)
1340 CALL pw_copy(rho_r_sccs_in, rho_r_sccs_out)
1341 END IF
1342
1343 ! drho_r and drho_g are only needed if calculated by collocation
1344 IF (dft_control%drho_by_collocation) THEN
1345 ! drho_r
1346 IF (ASSOCIATED(drho_r_in)) THEN
1347 ALLOCATE (drho_r_out(3, nspins))
1348 CALL qs_rho_set(rho_output, drho_r=drho_r_out)
1349 DO j = 1, nspins
1350 DO i = 1, 3
1351 CALL auxbas_pw_pool%create_pw(drho_r_out(i, j))
1352 CALL pw_copy(drho_r_in(i, j), drho_r_out(i, j))
1353 END DO
1354 END DO
1355 END IF
1356
1357 ! drho_g
1358 IF (ASSOCIATED(drho_g_in)) THEN
1359 ALLOCATE (drho_g_out(3, nspins))
1360 CALL qs_rho_set(rho_output, drho_g=drho_g_out)
1361 DO j = 1, nspins
1362 DO i = 1, 3
1363 CALL auxbas_pw_pool%create_pw(drho_g_out(i, j))
1364 CALL pw_copy(drho_g_in(i, j), drho_g_out(i, j))
1365 END DO
1366 END DO
1367 END IF
1368 END IF
1369
1370 ! tau_r and tau_g are only needed in the case of Meta-GGA XC-functionals
1371 ! are used. Therefore they are only allocated if
1372 ! dft_control%use_kinetic_energy_density is true
1373 IF (dft_control%use_kinetic_energy_density) THEN
1374 ! tau_r
1375 IF (ASSOCIATED(tau_r_in)) THEN
1376 ALLOCATE (tau_r_out(nspins))
1377 CALL qs_rho_set(rho_output, tau_r=tau_r_out)
1378 DO i = 1, nspins
1379 CALL auxbas_pw_pool%create_pw(tau_r_out(i))
1380 CALL pw_copy(tau_r_in(i), tau_r_out(i))
1381 END DO
1382 END IF
1383
1384 ! tau_g
1385 IF (ASSOCIATED(tau_g_in)) THEN
1386 ALLOCATE (tau_g_out(nspins))
1387 CALL qs_rho_set(rho_output, tau_g=tau_g_out)
1388 DO i = 1, nspins
1389 CALL auxbas_pw_pool%create_pw(tau_g_out(i))
1390 CALL pw_copy(tau_g_in(i), tau_g_out(i))
1391 END DO
1392 END IF
1393 END IF
1394
1395 CALL qs_rho_set(rho_output, &
1396 rho_g_valid=rho_g_valid_in, &
1397 rho_r_valid=rho_r_valid_in, &
1398 drho_g_valid=drho_g_valid_in, &
1399 drho_r_valid=drho_r_valid_in, &
1400 tau_r_valid=tau_r_valid_in, &
1401 tau_g_valid=tau_g_valid_in, &
1402 soft_valid=soft_valid_in, &
1403 complex_rho_ao=complex_rho_ao_in)
1404
1405 ! tot_rho_r
1406 IF (ASSOCIATED(tot_rho_r_in)) THEN
1407 ALLOCATE (tot_rho_r_out(nspins))
1408 CALL qs_rho_set(rho_output, tot_rho_r=tot_rho_r_out)
1409 DO i = 1, nspins
1410 tot_rho_r_out(i) = tot_rho_r_in(i)
1411 END DO
1412 END IF
1413
1414 ! tot_rho_g
1415 IF (ASSOCIATED(tot_rho_g_in)) THEN
1416 ALLOCATE (tot_rho_g_out(nspins))
1417 CALL qs_rho_set(rho_output, tot_rho_g=tot_rho_g_out)
1418 DO i = 1, nspins
1419 tot_rho_g_out(i) = tot_rho_g_in(i)
1420 END DO
1421
1422 END IF
1423
1424 CALL timestop(handle)
1425
1426 END SUBROUTINE duplicate_rho_type
1427
1428! **************************************************************************************************
1429!> \brief (Re-)allocates rho_ao_im from real part rho_ao
1430!> \param rho ...
1431!> \param qs_env ...
1432! **************************************************************************************************
1433 SUBROUTINE allocate_rho_ao_imag_from_real(rho, qs_env)
1434 TYPE(qs_rho_type), POINTER :: rho
1435 TYPE(qs_environment_type), POINTER :: qs_env
1436
1437 CHARACTER(LEN=default_string_length) :: headline
1438 INTEGER :: i, ic, nimages, nspins
1439 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_im_kp, rho_ao_kp
1440 TYPE(dbcsr_type), POINTER :: template
1441 TYPE(dft_control_type), POINTER :: dft_control
1442 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
1443 POINTER :: sab_orb
1444
1445 NULLIFY (rho_ao_im_kp, rho_ao_kp, dft_control, template, sab_orb)
1446
1447 CALL get_qs_env(qs_env, &
1448 dft_control=dft_control, &
1449 sab_orb=sab_orb)
1450
1451 CALL qs_rho_get(rho, rho_ao_im_kp=rho_ao_im_kp, rho_ao_kp=rho_ao_kp)
1452
1453 nspins = dft_control%nspins
1454 nimages = dft_control%nimages
1455
1456 cpassert(nspins .EQ. SIZE(rho_ao_kp, 1))
1457 cpassert(nimages .EQ. SIZE(rho_ao_kp, 2))
1458
1459 CALL dbcsr_allocate_matrix_set(rho_ao_im_kp, nspins, nimages)
1460 CALL qs_rho_set(rho, rho_ao_im_kp=rho_ao_im_kp)
1461 DO i = 1, nspins
1462 DO ic = 1, nimages
1463 IF (nspins > 1) THEN
1464 IF (i == 1) THEN
1465 headline = "IMAGINARY PART OF DENSITY MATRIX FOR ALPHA SPIN"
1466 ELSE
1467 headline = "IMAGINARY PART OF DENSITY MATRIX FOR BETA SPIN"
1468 END IF
1469 ELSE
1470 headline = "IMAGINARY PART OF DENSITY MATRIX"
1471 END IF
1472 ALLOCATE (rho_ao_im_kp(i, ic)%matrix)
1473 template => rho_ao_kp(i, ic)%matrix ! base on real part, but anti-symmetric
1474 CALL dbcsr_create(matrix=rho_ao_im_kp(i, ic)%matrix, template=template, &
1475 name=trim(headline), matrix_type=dbcsr_type_antisymmetric, nze=0)
1476 CALL cp_dbcsr_alloc_block_from_nbl(rho_ao_im_kp(i, ic)%matrix, sab_orb)
1477 CALL dbcsr_set(rho_ao_im_kp(i, ic)%matrix, 0.0_dp)
1478 END DO
1479 END DO
1480
1481 END SUBROUTINE allocate_rho_ao_imag_from_real
1482
1483END MODULE qs_rho_methods
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition cp_dbcsr_operations.F:81
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:89
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
pw_methods::pw_axpy
Definition pw_methods.F:164
pw_methods::pw_copy
Definition pw_methods.F:145
pw_methods::pw_scale
Definition pw_methods.F:93
pw_methods::pw_zero
Definition pw_methods.F:82
admm_types
Types and set/get functions for auxiliary density matrix methods.
Definition admm_types.F:15
admm_types::get_admm_env
subroutine, public get_admm_env(admm_env, mo_derivs_aux_fit, mos_aux_fit, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, matrix_s_aux_fit, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp, task_list_aux_fit, matrix_ks_aux_fit, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp, rho_aux_fit, rho_aux_fit_buffer, admm_dm)
Get routine for the ADMM env.
Definition admm_types.F:593
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
cp_control_types
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Definition cp_control_types.F:12
cp_dbcsr_cp2k_link
Routines that link DBCSR and CP2K concepts together.
Definition cp_dbcsr_cp2k_link.F:14
cp_dbcsr_cp2k_link::cp_dbcsr_alloc_block_from_nbl
subroutine, public cp_dbcsr_alloc_block_from_nbl(matrix, sab_orb, desymmetrize)
allocate the blocks of a dbcsr based on the neighbor list
Definition cp_dbcsr_cp2k_link.F:445
cp_dbcsr_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_string_length
integer, parameter, public default_string_length
Definition kinds.F:57
kpoint_types
Types and basic routines needed for a kpoint calculation.
Definition kpoint_types.F:15
kpoint_types::get_kpoint_info
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym)
Retrieve information from a kpoint environment.
Definition kpoint_types.F:333
lri_environment_methods
Calculates integral matrices for LRIGPW method lri : local resolution of the identity.
Definition lri_environment_methods.F:18
lri_environment_methods::calculate_lri_densities
subroutine, public calculate_lri_densities(lri_env, lri_density, qs_env, pmatrix, cell_to_index, lri_rho_struct, atomic_kind_set, para_env, response_density)
performs the fitting of the density and distributes the fitted density on the grid
Definition lri_environment_methods.F:569
lri_environment_types
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
Definition lri_environment_types.F:18
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
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_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_drho_elec
subroutine, public calculate_drho_elec(matrix_p, matrix_p_kp, drho, drho_gspace, qs_env, soft_valid, basis_type)
computes the gradient of the density corresponding to a given density matrix on the grid
Definition qs_collocate_density.F:1876
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_environment_types
Definition qs_environment_types.F:14
qs_environment_types::set_qs_env
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, wanniercentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Set the QUICKSTEP environment.
Definition qs_environment_types.F:1035
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_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_ks_types
Definition qs_ks_types.F:15
qs_ks_types::get_ks_env
subroutine, public get_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, rho, rho_xc, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, task_list, task_list_soft, kpoints, do_kpoints, atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, particle_set, energy, force, local_particles, local_molecules, molecule_kind_set, molecule_set, subsys, cp_subsys, virial, results, atprop, nkind, natom, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env, nelectron_total, nelectron_spin)
...
Definition qs_ks_types.F:330
qs_local_rho_types
Definition qs_local_rho_types.F:8
qs_neighbor_list_types
Define the neighbor list data types and the corresponding functionality.
Definition qs_neighbor_list_types.F:17
qs_oce_types
Definition qs_oce_types.F:8
qs_rho_atom_methods
Definition qs_rho_atom_methods.F:7
qs_rho_atom_methods::calculate_rho_atom_coeff
subroutine, public calculate_rho_atom_coeff(qs_env, rho_ao, rho_atom_set, qs_kind_set, oce, sab, para_env)
...
Definition qs_rho_atom_methods.F:420
qs_rho_atom_types
Definition qs_rho_atom_types.F:8
qs_rho_methods
methods of the rho structure (defined in qs_rho_types)
Definition qs_rho_methods.F:15
qs_rho_methods::qs_rho_update_tddfpt
subroutine, public qs_rho_update_tddfpt(rho_struct, qs_env, pw_env_external, task_list_external, para_env_external, tddfpt_lri_env, tddfpt_lri_density)
updates rho_r and rho_g to the rhorho_ao. if use_kinetic_energy_density also computes tau_r and tau_g
Definition qs_rho_methods.F:670
qs_rho_methods::allocate_rho_ao_imag_from_real
subroutine, public allocate_rho_ao_imag_from_real(rho, qs_env)
(Re-)allocates rho_ao_im from real part rho_ao
Definition qs_rho_methods.F:1434
qs_rho_methods::qs_rho_update_rho
subroutine, public qs_rho_update_rho(rho_struct, qs_env, rho_xc_external, local_rho_set, task_list_external, task_list_external_soft, pw_env_external, para_env_external)
updates rho_r and rho_g to the rhorho_ao. if use_kinetic_energy_density also computes tau_r and tau_g...
Definition qs_rho_methods.F:375
qs_rho_methods::duplicate_rho_type
subroutine, public duplicate_rho_type(rho_input, rho_output, qs_env)
Duplicates a pointer physically.
Definition qs_rho_methods.F:1232
qs_rho_methods::qs_rho_copy
subroutine, public qs_rho_copy(rho_input, rho_output, auxbas_pw_pool, mspin)
Allocate a density structure and fill it with data from an input structure SIZE(rho_input) == mspin =...
Definition qs_rho_methods.F:766
qs_rho_methods::qs_rho_scale_and_add
subroutine, public qs_rho_scale_and_add(rhoa, rhob, alpha, beta)
rhoa = alpha*rhoa+beta*rhob
Definition qs_rho_methods.F:1060
qs_rho_methods::qs_rho_rebuild
subroutine, public qs_rho_rebuild(rho, qs_env, rebuild_ao, rebuild_grids, admm, pw_env_external)
rebuilds rho (if necessary allocating and initializing it)
Definition qs_rho_methods.F:96
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_set
subroutine, public qs_rho_set(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)
...
Definition qs_rho_types.F:308
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
qs_rho_types::qs_rho_clear
subroutine, public qs_rho_clear(rho_struct)
Deallocates all components, without deallocating rho_struct itself.
Definition qs_rho_types.F:125
ri_environment_methods
Calculates integral matrices for RIGPW method.
Definition ri_environment_methods.F:18
ri_environment_methods::calculate_ri_densities
subroutine, public calculate_ri_densities(lri_env, qs_env, pmatrix, lri_rho_struct, atomic_kind_set, para_env)
performs the fitting of the density and distributes the fitted density on the grid
Definition ri_environment_methods.F:358
task_list_types
types for task lists
Definition task_list_types.F:14
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
cp_control_types::dft_control_type
Definition cp_control_types.F:559
kpoint_types::kpoint_type
Contains information about kpoints.
Definition kpoint_types.F:138
lri_environment_types::lri_density_type
Definition lri_environment_types.F:357
lri_environment_types::lri_environment_type
Definition lri_environment_types.F:269
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:763
pw_env_types::pw_env_type
contained for different pw related things
Definition pw_env_types.F:70
pw_pool_types::pw_pool_type
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:83
pw_types::pw_c1d_gs_type
Definition pw_types.F:127
pw_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:215
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:171
qs_ks_types::qs_ks_env_type
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
Definition qs_ks_types.F:135
qs_local_rho_types::local_rho_type
Definition qs_local_rho_types.F:43
qs_neighbor_list_types::neighbor_list_set_p_type
Definition qs_neighbor_list_types.F:67
qs_oce_types::oce_matrix_type
Definition qs_oce_types.F:37
qs_rho_atom_types::rho_atom_type
Definition qs_rho_atom_types.F:23
qs_rho_types::qs_rho_type
keeps the density in various representations, keeping track of which ones are valid.
Definition qs_rho_types.F:62
task_list_types::task_list_type
Definition task_list_types.F:58