(git:ed6f26b)
Loading...
Searching...
No Matches
qs_linres_current_utils.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
8! **************************************************************************************************
9!> \brief Chemical shift calculation by dfpt
10!> Initialization of the nmr_env, creation of the special neighbor lists
11!> Perturbation Hamiltonians by application of the p and rxp oprtators to psi0
12!> Write output
13!> Deallocate everything
14!> \note
15!> The psi0 should be localized
16!> the Sebastiani method works within the assumption that the orbitals are
17!> completely contained in the simulation box
18!> \par History
19!> created 07-2005 [MI]
20!> \author MI
21! **************************************************************************************************
22MODULE qs_linres_current_utils
23 USE atomic_kind_types, ONLY: atomic_kind_type
24 USE cell_types, ONLY: cell_type,&
25 pbc
26 USE cp_array_utils, ONLY: cp_2d_i_p_type,&
27 cp_2d_r_p_type
28 USE cp_control_types, ONLY: dft_control_type
29 USE cp_dbcsr_api, ONLY: dbcsr_convert_offsets_to_sizes,&
30 dbcsr_copy,&
31 dbcsr_create,&
32 dbcsr_distribution_type,&
33 dbcsr_p_type,&
34 dbcsr_set,&
35 dbcsr_type_antisymmetric
36 USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
37 USE cp_dbcsr_operations, ONLY: cp_dbcsr_sm_fm_multiply,&
38 dbcsr_allocate_matrix_set,&
39 dbcsr_deallocate_matrix_set
40 USE cp_fm_basic_linalg, ONLY: cp_fm_column_scale,&
41 cp_fm_scale_and_add
42 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
43 cp_fm_struct_release,&
44 cp_fm_struct_type
45 USE cp_fm_types, ONLY: cp_fm_create,&
46 cp_fm_release,&
47 cp_fm_set_all,&
48 cp_fm_to_fm,&
49 cp_fm_type
50 USE cp_log_handling, ONLY: cp_get_default_logger,&
51 cp_logger_type,&
52 cp_to_string
53 USE cp_output_handling, ONLY: cp_p_file,&
54 cp_print_key_finished_output,&
55 cp_print_key_should_output,&
56 cp_print_key_unit_nr
57 USE input_constants, ONLY: current_gauge_atom,&
58 current_gauge_r,&
59 current_gauge_r_and_step_func,&
60 current_orb_center_atom,&
61 current_orb_center_box,&
62 current_orb_center_common,&
63 current_orb_center_wannier,&
64 ot_precond_full_all
65 USE input_section_types, ONLY: section_get_lval,&
66 section_vals_get_subs_vals,&
67 section_vals_type,&
68 section_vals_val_get
69 USE kinds, ONLY: default_path_length,&
70 dp
71 USE memory_utilities, ONLY: reallocate
72 USE message_passing, ONLY: mp_para_env_type
73 USE particle_methods, ONLY: get_particle_set
74 USE particle_types, ONLY: particle_type
75 USE pw_env_types, ONLY: pw_env_get,&
76 pw_env_type
77 USE pw_methods, ONLY: pw_zero
78 USE pw_pool_types, ONLY: pw_pool_type
79 USE pw_types, ONLY: pw_c1d_gs_type,&
80 pw_r3d_rs_type
81 USE qs_environment_types, ONLY: get_qs_env,&
82 qs_environment_type
83 USE qs_kind_types, ONLY: qs_kind_type
84 USE qs_linres_methods, ONLY: linres_read_restart,&
85 linres_solver,&
86 linres_write_restart
87 USE qs_linres_op, ONLY: set_vecp
88 USE qs_linres_types, ONLY: current_env_type,&
89 deallocate_jrho_atom_set,&
90 get_current_env,&
91 init_jrho_atom_set,&
92 jrho_atom_type,&
93 linres_control_type,&
94 set_current_env
95 USE qs_loc_methods, ONLY: qs_print_cubes
96 USE qs_loc_types, ONLY: get_qs_loc_env,&
97 localized_wfn_control_type,&
98 qs_loc_env_type
99 USE qs_matrix_pools, ONLY: qs_matrix_pools_type
100 USE qs_mo_types, ONLY: get_mo_set,&
101 mo_set_type
102 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
103 USE qs_operators_ao, ONLY: build_lin_mom_matrix
104 USE qs_p_env_types, ONLY: qs_p_env_type
105 USE qs_rho_types, ONLY: qs_rho_clear,&
106 qs_rho_create,&
107 qs_rho_set
108 USE realspace_grid_types, ONLY: rs_grid_release
109 USE scf_control_types, ONLY: scf_control_type
110#include "./base/base_uses.f90"
111
112 IMPLICIT NONE
113
114 PRIVATE
115 PUBLIC :: current_response, current_env_cleanup, current_env_init
116
117 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_linres_current_utils'
118
119CONTAINS
120
121! **************************************************************************************************
122!> \brief ...
123!> \param current_env ...
124!> \param p_env ...
125!> \param qs_env ...
126! **************************************************************************************************
127 SUBROUTINE current_response(current_env, p_env, qs_env)
128 !
129 TYPE(current_env_type) :: current_env
130 TYPE(qs_p_env_type) :: p_env
131 TYPE(qs_environment_type), POINTER :: qs_env
132
133 CHARACTER(LEN=*), PARAMETER :: routinen = 'current_response'
134
135 CHARACTER(LEN=default_path_length) :: my_pos
136 INTEGER :: first_center, handle, i, icenter, idir, ii, iii, ispin, ist_true, istate, j, &
137 jcenter, jstate, max_nbr_center, max_states, nao, natom, nbr_center(2), ncubes, nmo, &
138 nspins, nstates(2), output_unit
139 INTEGER, ALLOCATABLE, DIMENSION(:) :: first_sgf, last_sgf
140 INTEGER, DIMENSION(:), POINTER :: list_cubes, row_blk_sizes
141 INTEGER, DIMENSION(:, :, :), POINTER :: statetrueindex
142 LOGICAL :: append_cube, should_stop
143 REAL(dp) :: dk(3), dkl(3), dl(3)
144 REAL(dp), ALLOCATABLE, DIMENSION(:) :: dkl_vec_ii, dkl_vec_iii
145 REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: vecbuf_dklxp0
146 TYPE(cell_type), POINTER :: cell
147 TYPE(cp_2d_i_p_type), DIMENSION(:), POINTER :: center_list
148 TYPE(cp_2d_r_p_type), DIMENSION(:), POINTER :: centers_set
149 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
150 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:) :: fm_work_ii, fm_work_iii, h1_psi0, psi1
151 TYPE(cp_fm_type), DIMENSION(:), POINTER :: hpsi0_ptr, psi0_order
152 TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: p_psi0, psi1_d, psi1_p, psi1_rxp, &
153 rxp_psi0
154 TYPE(cp_fm_type), POINTER :: mo_coeff
155 TYPE(cp_logger_type), POINTER :: logger
156 TYPE(dbcsr_distribution_type), POINTER :: dbcsr_dist
157 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: op_p_ao
158 TYPE(dft_control_type), POINTER :: dft_control
159 TYPE(linres_control_type), POINTER :: linres_control
160 TYPE(mp_para_env_type), POINTER :: para_env
161 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
162 POINTER :: sab_orb
163 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
164 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
165 TYPE(qs_matrix_pools_type), POINTER :: mpools
166 TYPE(section_vals_type), POINTER :: current_section, lr_section, print_key
167
168 CALL timeset(routinen, handle)
169 !
170 NULLIFY (cell, dft_control, linres_control, lr_section, current_section, &
171 logger, mpools, mo_coeff, para_env, &
172 tmp_fm_struct, &
173 list_cubes, statetrueindex, centers_set, center_list, psi1_p, psi1_rxp, psi1_d, &
174 p_psi0, rxp_psi0, psi0_order, op_p_ao, sab_orb, particle_set)
175
176 logger => cp_get_default_logger()
177 lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
178 current_section => section_vals_get_subs_vals(qs_env%input, &
179 "PROPERTIES%LINRES%CURRENT")
180
181 output_unit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
182 extension=".linresLog")
183 IF (output_unit > 0) THEN
184 WRITE (unit=output_unit, fmt="(T10,A,/)") &
185 "*** Self consistent optimization of the response wavefunctions ***"
186 END IF
187
188 CALL get_qs_env(qs_env=qs_env, &
189 dft_control=dft_control, &
190 mpools=mpools, cell=cell, &
191 linres_control=linres_control, &
192 sab_orb=sab_orb, &
193 particle_set=particle_set, &
194 qs_kind_set=qs_kind_set, &
195 dbcsr_dist=dbcsr_dist, &
196 para_env=para_env)
197
198 nspins = dft_control%nspins
199
200 CALL get_current_env(current_env=current_env, &
201 nao=nao, &
202 nstates=nstates, &
203 centers_set=centers_set, &
204 nbr_center=nbr_center, &
205 center_list=center_list, &
206 list_cubes=list_cubes, &
207 statetrueindex=statetrueindex, &
208 psi1_p=psi1_p, &
209 psi1_rxp=psi1_rxp, &
210 psi1_d=psi1_d, &
211 p_psi0=p_psi0, &
212 rxp_psi0=rxp_psi0, &
213 psi0_order=psi0_order)
214 !
215 ! allocate the vectors
216 IF (current_env%full) THEN
217 ALLOCATE (psi1(nspins), h1_psi0(nspins))
218 DO ispin = 1, nspins
219 mo_coeff => psi0_order(ispin)
220 NULLIFY (tmp_fm_struct)
221 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
222 ncol_global=nstates(ispin), &
223 context=mo_coeff%matrix_struct%context)
224 CALL cp_fm_create(psi1(ispin), tmp_fm_struct)
225 CALL cp_fm_create(h1_psi0(ispin), tmp_fm_struct)
226 CALL cp_fm_struct_release(tmp_fm_struct)
227 END DO
228 !
229 ! prepare for allocation
230 natom = SIZE(particle_set, 1)
231 ALLOCATE (first_sgf(natom))
232 ALLOCATE (last_sgf(natom))
233 CALL get_particle_set(particle_set, qs_kind_set, &
234 first_sgf=first_sgf, &
235 last_sgf=last_sgf)
236 ALLOCATE (row_blk_sizes(natom))
237 CALL dbcsr_convert_offsets_to_sizes(first_sgf, row_blk_sizes, last_sgf)
238 DEALLOCATE (first_sgf, last_sgf)
239 !
240 ! rebuild the linear momentum matrices
241 CALL dbcsr_allocate_matrix_set(op_p_ao, 3)
242 ALLOCATE (op_p_ao(1)%matrix)
243 CALL dbcsr_create(matrix=op_p_ao(1)%matrix, &
244 name="OP_P", &
245 dist=dbcsr_dist, matrix_type=dbcsr_type_antisymmetric, &
246 row_blk_size=row_blk_sizes, col_blk_size=row_blk_sizes, &
247 mutable_work=.true.)
248 CALL cp_dbcsr_alloc_block_from_nbl(op_p_ao(1)%matrix, sab_orb)
249 !
250 !
251 DEALLOCATE (row_blk_sizes)
252 !
253 !
254 CALL dbcsr_set(op_p_ao(1)%matrix, 0.0_dp)
255 DO idir = 2, 3
256 ALLOCATE (op_p_ao(idir)%matrix)
257 CALL dbcsr_copy(op_p_ao(idir)%matrix, op_p_ao(1)%matrix, &
258 "current_env%op_p_ao"//"-"//trim(adjustl(cp_to_string(idir))))
259 CALL dbcsr_set(op_p_ao(idir)%matrix, 0.0_dp)
260 END DO
261 !
262 !CALL p_xyz_ao(op_p_ao,qs_env,minimum_image=.FALSE.)
263 CALL build_lin_mom_matrix(qs_env, op_p_ao)
264 !
265 END IF
266 !
267 ! set response to zero
268 !
269 DO idir = 1, 3
270 DO ispin = 1, nspins
271 CALL cp_fm_set_all(psi1_p(ispin, idir), 0.0_dp)
272 CALL cp_fm_set_all(psi1_rxp(ispin, idir), 0.0_dp)
273 IF (current_env%full) CALL cp_fm_set_all(psi1_d(ispin, idir), 0.0_dp)
274 END DO
275 END DO
276 !
277 !
278 !
279 ! load restart file
280 !
281 first_center = 0
282 IF (linres_control%linres_restart) THEN
283 DO idir = 1, 3
284 ! operator p
285 CALL linres_read_restart(qs_env, lr_section, psi1_p(:, idir), idir, "nmr_p")
286 ! operator rxp
287 CALL linres_read_restart(qs_env, lr_section, psi1_rxp(:, idir), idir, "nmr_rxp")
288 END DO
289 END IF
290 !
291 !
292 !
293 !
294 should_stop = .false.
295 current_env%all_pert_op_done = .false.
296 ! operator p
297 DO idir = 1, 3
298 IF (should_stop) EXIT
299 IF (output_unit > 0) THEN
300 WRITE (output_unit, "(T10,A)") "Response to the perturbation operator P_"//achar(idir + 119)
301 END IF
302 !
303 ! Initial guess for psi1
304 hpsi0_ptr => p_psi0(:, idir)
305 !
306 !
307 linres_control%converged = .false.
308 CALL linres_solver(p_env, qs_env, psi1_p(:, idir), hpsi0_ptr, psi0_order, output_unit, should_stop)
309 !
310 !
311 ! print response functions
312 IF (btest(cp_print_key_should_output(logger%iter_info, current_section,&
313 & "PRINT%RESPONSE_FUNCTION_CUBES"), cp_p_file)) THEN
314 ncubes = SIZE(list_cubes, 1)
315 print_key => section_vals_get_subs_vals(current_section, "PRINT%RESPONSE_FUNCTION_CUBES")
316 append_cube = section_get_lval(current_section, "PRINT%RESPONSE_FUNCTION_CUBES%APPEND")
317 my_pos = "REWIND"
318 IF (append_cube) THEN
319 my_pos = "APPEND"
320 END IF
321 !
322
323 DO ispin = 1, nspins
324 CALL qs_print_cubes(qs_env, psi1_p(ispin, idir), ncubes, list_cubes, &
325 centers_set(ispin)%array, print_key, 'psi1_p', &
326 idir=idir, ispin=ispin, file_position=my_pos)
327 END DO ! ispin
328 END IF ! print response functions
329 !
330 ! write restart file
331 CALL linres_write_restart(qs_env, lr_section, psi1_p(:, idir), idir, "nmr_p")
332 END DO ! idir
333 !
334 ! operator rxp
335 DO idir = 1, 3
336 IF (should_stop) EXIT
337 IF (output_unit > 0) THEN
338 WRITE (output_unit, "(T10,A)") "Response to the perturbation operator L_"//achar(idir + 119)
339 END IF
340 !
341 ! Initial guess for psi1
342 hpsi0_ptr => rxp_psi0(:, idir)
343 !
344 !
345 linres_control%converged = .false.
346 CALL linres_solver(p_env, qs_env, psi1_rxp(:, idir), hpsi0_ptr, psi0_order, output_unit, should_stop)
347 !
348 ! print response functions
349 IF (btest(cp_print_key_should_output(logger%iter_info, current_section,&
350 & "PRINT%RESPONSE_FUNCTION_CUBES"), cp_p_file)) THEN
351 ncubes = SIZE(list_cubes, 1)
352 print_key => section_vals_get_subs_vals(current_section, "PRINT%RESPONSE_FUNCTION_CUBES")
353
354 DO ispin = 1, nspins
355 CALL qs_print_cubes(qs_env, psi1_rxp(ispin, idir), ncubes, list_cubes, &
356 centers_set(ispin)%array, print_key, 'psi1_rxp', &
357 idir=idir, ispin=ispin, file_position=my_pos)
358 END DO ! ispin
359 END IF ! print response functions
360 !
361 ! write restart file
362 CALL linres_write_restart(qs_env, lr_section, psi1_rxp(:, idir), idir, "nmr_rxp")
363 END DO ! idir
364 IF (.NOT. should_stop) current_env%all_pert_op_done = .true.
365 !
366 ! operator D
367 IF (current_env%full) THEN
368 current_env%all_pert_op_done = .false.
369 !
370 !
371 DO ispin = 1, nspins
372 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
373 END DO
374 !
375 ! The correction is state depedent a loop over the states is necessary
376 max_nbr_center = maxval(nbr_center(1:nspins))
377 max_states = maxval(nstates(1:nspins))
378 !
379 ALLOCATE (vecbuf_dklxp0(1, nao), fm_work_ii(nspins), fm_work_iii(nspins), &
380 dkl_vec_ii(max_states), dkl_vec_iii(max_states))
381 vecbuf_dklxp0(1, nao) = 0.0_dp
382 !
383 DO ispin = 1, nspins
384 nmo = nstates(ispin)
385 mo_coeff => psi0_order(ispin)
386 NULLIFY (tmp_fm_struct)
387 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
388 ncol_global=nmo, para_env=para_env, &
389 context=mo_coeff%matrix_struct%context)
390
391 CALL cp_fm_create(fm_work_ii(ispin), tmp_fm_struct)
392 CALL cp_fm_set_all(fm_work_ii(ispin), 0.0_dp)
393 CALL cp_fm_create(fm_work_iii(ispin), tmp_fm_struct)
394 CALL cp_fm_set_all(fm_work_iii(ispin), 0.0_dp)
395 CALL cp_fm_struct_release(tmp_fm_struct)
396 END DO ! ispin
397 !
398 DO idir = 1, 3
399 IF (should_stop) EXIT
400 DO ispin = 1, nspins
401 CALL cp_fm_set_all(psi1_d(ispin, idir), 0.0_dp)
402 END DO
403 first_center = 0
404 IF (linres_control%linres_restart) THEN
405 CALL linres_read_restart(qs_env, lr_section, psi1_d(:, idir), idir, "nmr_dxp-", ind=first_center)
406 END IF
407 IF (first_center > 0) THEN
408 IF (output_unit > 0) THEN
409 WRITE (output_unit, "(T10,A,I4,A)")&
410 & "Response to the perturbation operators (dk-dl)xp up to state ", &
411 first_center, " have been read from restart"
412 END IF
413 END IF
414 !
415 ! here we run over the max number of states, then we need
416 ! to be careful with overflow while doing uks calculations.
417 DO icenter = 1, max_nbr_center
418 !
419 IF (should_stop) EXIT
420 IF (icenter > first_center) THEN
421 IF (output_unit > 0) THEN
422 WRITE (output_unit, "(T10,A,I4,A)")&
423 & "Response to the perturbation operator (dk-dl)xp for -state- ", &
424 icenter, " in dir. "//achar(idir + 119)
425 END IF
426 !
427 DO ispin = 1, nspins
428 nmo = nstates(ispin)
429 mo_coeff => psi0_order(ispin)
430 !
431 ! take care that no overflow can occur for uks
432 IF (icenter .GT. nbr_center(ispin)) THEN
433 !
434 ! set h1_psi0 and psi1 to zero to avoid problems in linres_scf
435 CALL cp_fm_set_all(h1_psi0(ispin), 0.0_dp)
436 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
437 cycle
438 END IF
439 !
440 dkl_vec_ii(:) = 0.0_dp
441 dkl_vec_iii(:) = 0.0_dp
442 !
443 ist_true = statetrueindex(idir, icenter, ispin)
444 !
445 ! the initial guess is the previous set of psi1, just optimized
446 CALL set_vecp(idir, ii, iii)
447 dk(1:3) = centers_set(ispin)%array(1:3, ist_true)
448 !
449 DO jcenter = 1, nbr_center(ispin)
450 dl(1:3) = centers_set(ispin)%array(1:3, jcenter)
451 dkl = pbc(dl, dk, cell)
452 DO j = center_list(ispin)%array(1, jcenter), center_list(ispin)%array(1, jcenter + 1) - 1
453 jstate = center_list(ispin)%array(2, j)
454 dkl_vec_ii(jstate) = dkl(ii)
455 dkl_vec_iii(jstate) = dkl(iii)
456 END DO
457 END DO
458 !
459 ! First term
460 ! Rescale the ground state orbitals by (dk-dl)_ii
461 CALL cp_fm_to_fm(mo_coeff, fm_work_ii(ispin))
462 CALL cp_fm_column_scale(fm_work_ii(ispin), dkl_vec_ii(1:nmo))
463 !
464 ! Apply the p_iii operator
465 ! fm_work_iii = -p_iii * (dk-dl)_ii * C0
466 CALL cp_dbcsr_sm_fm_multiply(op_p_ao(iii)%matrix, fm_work_ii(ispin), &
467 fm_work_iii(ispin), ncol=nmo, alpha=-1.0_dp)
468 !
469 ! Copy in h1_psi0
470 ! h1_psi0_i = fm_work_iii
471 CALL cp_fm_to_fm(fm_work_iii(ispin), h1_psi0(ispin))
472 !
473 ! Second term
474 ! Rescale the ground state orbitals by (dk-dl)_iii
475 CALL cp_fm_to_fm(mo_coeff, fm_work_iii(ispin))
476 CALL cp_fm_column_scale(fm_work_iii(ispin), dkl_vec_iii(1:nmo))
477 !
478 ! Apply the p_ii operator
479 ! fm_work_ii = -p_ii * (dk-dl)_iii * C0
480 CALL cp_dbcsr_sm_fm_multiply(op_p_ao(ii)%matrix, fm_work_iii(ispin), &
481 fm_work_ii(ispin), ncol=nmo, alpha=-1.0_dp)
482 !
483 ! Copy in h1_psi0
484 ! h1_psi0_i = fm_work_iii - fm_work_ii
485 CALL cp_fm_scale_and_add(1.0_dp, h1_psi0(ispin),&
486 & -1.0_dp, fm_work_ii(ispin))
487 END DO
488
489 !
490 ! Optimize the response wavefunctions
491 CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, output_unit, should_stop)
492 !
493 IF (output_unit > 0) THEN
494 WRITE (output_unit, "(T10,A,/)")&
495 & "Store the psi1 vector for the calculation of the response current density "
496 END IF
497 !
498 DO ispin = 1, nspins
499 !
500 ! take care that no overflow can occur for uks
501 IF (icenter .GT. nbr_center(ispin)) cycle
502 !
503 ! need to reset those guys
504 ist_true = statetrueindex(idir, icenter, ispin)
505 DO i = center_list(ispin)%array(1, ist_true), center_list(ispin)%array(1, ist_true + 1) - 1
506 istate = center_list(ispin)%array(2, i)
507 !
508 ! the optimized wfns are copied in the fm
509 CALL cp_fm_to_fm(psi1(ispin), psi1_d(ispin, idir), 1, istate, istate)
510 END DO
511 current_env%full_done(idir*ispin, icenter) = .true.
512 END DO ! ispin
513 !
514 ELSE
515 DO ispin = 1, nspins
516 current_env%full_done(idir*ispin, icenter) = .true.
517 END DO ! ispin
518 END IF
519 CALL linres_write_restart(qs_env, lr_section, psi1_d(:, idir), idir, "nmr_dxp-", ind=icenter)
520
521 END DO ! center
522 !
523 ! print response functions
524 IF (btest(cp_print_key_should_output(logger%iter_info, current_section,&
525 & "PRINT%RESPONSE_FUNCTION_CUBES"), cp_p_file)) THEN
526 ncubes = SIZE(list_cubes, 1)
527 print_key => section_vals_get_subs_vals(current_section, "PRINT%RESPONSE_FUNCTION_CUBES")
528 DO ispin = 1, nspins
529 CALL qs_print_cubes(qs_env, psi1_d(ispin, idir), &
530 ncubes, list_cubes, centers_set(ispin)%array, print_key, 'psi1_D', &
531 idir=idir, ispin=ispin, file_position=my_pos)
532 END DO
533 END IF ! print response functions
534 !
535 END DO ! idir
536 IF (.NOT. should_stop) current_env%all_pert_op_done = .true.
537 !
538 ! clean up
539 CALL cp_fm_release(fm_work_ii)
540 CALL cp_fm_release(fm_work_iii)
541 DEALLOCATE (dkl_vec_ii, dkl_vec_iii, vecbuf_dklxp0)
542
543 END IF
544 !
545 ! clean up
546 IF (current_env%full) THEN
547 CALL dbcsr_deallocate_matrix_set(op_p_ao)
548 CALL cp_fm_release(psi1)
549 CALL cp_fm_release(h1_psi0)
550 END IF
551 !
552 CALL cp_print_key_finished_output(output_unit, logger, lr_section,&
553 & "PRINT%PROGRAM_RUN_INFO")
554 !
555 CALL timestop(handle)
556 !
557 END SUBROUTINE current_response
558
559! **************************************************************************************************
560
561! **************************************************************************************************
562!> \brief ...
563!> \param current_env ...
564!> \param qs_env ...
565! **************************************************************************************************
566 SUBROUTINE current_env_init(current_env, qs_env)
567 !
568 TYPE(current_env_type) :: current_env
569 TYPE(qs_environment_type), POINTER :: qs_env
570
571 CHARACTER(LEN=*), PARAMETER :: routinen = 'current_env_init'
572
573 INTEGER :: handle, homo, i, iao, iatom, ibox, icenter, icount, idir, ii, ini, ir, is, ispin, &
574 istate, istate2, istate_next, ix, iy, iz, j, jstate, k, max_nbr_center, max_states, n, &
575 n_rep, nao, natom, nbr_box, ncubes, nmo, nspins, nstate, nstate_list(2), output_unit
576 INTEGER, ALLOCATABLE, DIMENSION(:) :: buff, first_sgf, last_sgf
577 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: state_list
578 INTEGER, DIMENSION(:), POINTER :: bounds, list, nbox, &
579 selected_states_on_atom_list
580 LOGICAL :: force_no_full, gapw, is0, &
581 uniform_occupation
582 LOGICAL, ALLOCATABLE, DIMENSION(:, :) :: state_done
583 REAL(dp) :: center(3), center2(3), dist, mdist, &
584 r(3), rab(3)
585 REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: rbuff
586 REAL(dp), DIMENSION(:), POINTER :: common_center
587 REAL(dp), DIMENSION(:, :), POINTER :: center_array
588 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
589 TYPE(cell_type), POINTER :: cell
590 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
591 TYPE(cp_fm_type), POINTER :: mo_coeff
592 TYPE(cp_logger_type), POINTER :: logger
593 TYPE(dft_control_type), POINTER :: dft_control
594 TYPE(jrho_atom_type), DIMENSION(:), POINTER :: jrho1_atom_set
595 TYPE(linres_control_type), POINTER :: linres_control
596 TYPE(localized_wfn_control_type), POINTER :: localized_wfn_control
597 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
598 TYPE(mp_para_env_type), POINTER :: para_env
599 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
600 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g
601 TYPE(pw_env_type), POINTER :: pw_env
602 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
603 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r
604 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
605 TYPE(qs_loc_env_type), POINTER :: qs_loc_env
606 TYPE(qs_matrix_pools_type), POINTER :: mpools
607 TYPE(scf_control_type), POINTER :: scf_control
608 TYPE(section_vals_type), POINTER :: current_section, lr_section
609
610 CALL timeset(routinen, handle)
611
612 NULLIFY (atomic_kind_set, qs_kind_set, cell, dft_control, linres_control, scf_control, &
613 logger, mos, mpools, current_section, particle_set, mo_coeff, &
614 auxbas_pw_pool, pw_env, jrho1_atom_set, common_center, tmp_fm_struct, &
615 para_env, qs_loc_env, localized_wfn_control, rho_g, rho_r)
616 !
617 !
618 CALL get_qs_env(qs_env=qs_env, &
619 atomic_kind_set=atomic_kind_set, &
620 qs_kind_set=qs_kind_set, &
621 cell=cell, &
622 dft_control=dft_control, &
623 linres_control=linres_control, &
624 mos=mos, &
625 mpools=mpools, &
626 particle_set=particle_set, &
627 pw_env=pw_env, &
628 scf_control=scf_control, &
629 para_env=para_env)
630 !
631 gapw = dft_control%qs_control%gapw
632 nspins = dft_control%nspins
633 natom = SIZE(particle_set, 1)
634 CALL get_mo_set(mo_set=mos(1), nao=nao)
635 !
636 max_states = 0
637 DO ispin = 1, nspins
638 CALL get_mo_set(mo_set=mos(ispin), nmo=nmo)
639 max_states = max(max_states, nmo)
640 END DO
641 !
642 !
643 !
644 ! some checks
645 DO ispin = 1, nspins
646 CALL get_mo_set(mo_set=mos(ispin), nmo=nmo, homo=homo, &
647 uniform_occupation=uniform_occupation)
648 !
649 ! check that homo=nmo
650 IF (nmo .NE. homo) cpabort("nmo != homo")
651 !
652 ! check that the nbr of localized states is equal to the nbr of states
653! nmoloc = SIZE(linres_control%localized_wfn_control%centers_set(ispin)%array,2)
654! IF (nmoloc.NE.nmo) CALL cp_abort(__LOCATION__,&
655! "The nbr of localized functions "//&
656! & "is not equal to the nbr of states.")
657 !
658 ! check that the occupation is uniform
659 IF (.NOT. uniform_occupation) cpabort("nmo != homo")
660 END DO
661 !
662 !
663 logger => cp_get_default_logger()
664 lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
665 output_unit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
666 extension=".linresLog")
667
668 IF (output_unit > 0) THEN
669 WRITE (output_unit, "(/,T20,A,/)") "*** Start current Calculation ***"
670 WRITE (output_unit, "(T10,A,/)") "Inizialization of the current environment"
671 END IF
672
673 CALL current_env_cleanup(current_env)
674 current_env%gauge_init = .false.
675 current_env%chi_tensor(:, :, :) = 0.0_dp
676 current_env%chi_tensor_loc(:, :, :) = 0.0_dp
677 current_env%nao = nao
678 current_env%full = .true.
679 current_env%do_selected_states = .false.
680 !
681 ! If current_density or full_nmr different allocations are required
682 current_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES%CURRENT")
683 !
684 ! Select the gauge
685 CALL section_vals_val_get(current_section, "GAUGE", i_val=current_env%gauge)
686 SELECT CASE (current_env%gauge)
687 CASE (current_gauge_r)
688 current_env%gauge_name = "R"
689 CASE (current_gauge_r_and_step_func)
690 current_env%gauge_name = "R_AND_STEP_FUNCTION"
691 CASE (current_gauge_atom)
692 current_env%gauge_name = "ATOM"
693 CASE DEFAULT
694 cpabort("Unknown gauge, try again...")
695 END SELECT
696 !
697 ! maximal radius to build the atom gauge
698 CALL section_vals_val_get(current_section, "GAUGE_ATOM_RADIUS", &
699 r_val=current_env%gauge_atom_radius)
700 ! use old gauge atom
701 CALL section_vals_val_get(current_section, "USE_OLD_GAUGE_ATOM", &
702 l_val=current_env%use_old_gauge_atom)
703 ! chi for pbc
704 CALL section_vals_val_get(current_section, "CHI_PBC", l_val=current_env%chi_pbc)
705 !
706 ! use old gauge atom
707 CALL section_vals_val_get(current_section, "USE_OLD_GAUGE_ATOM", &
708 l_val=current_env%use_old_gauge_atom)
709 !
710 ! chi for pbc
711 CALL section_vals_val_get(current_section, "CHI_PBC", l_val=current_env%chi_pbc)
712 !
713 ! which center for the orbitals shall we use
714 CALL section_vals_val_get(current_section, "ORBITAL_CENTER", i_val=current_env%orb_center)
715 SELECT CASE (current_env%orb_center)
716 CASE (current_orb_center_wannier)
717 !
718 current_env%orb_center_name = "WANNIER"
719 CASE (current_orb_center_common)
720 !
721 current_env%orb_center_name = "COMMON"
722 current_env%full = .false.
723 !
724 ! Is there a user specified common_center?
725 CALL section_vals_val_get(current_section, "COMMON_CENTER", r_vals=common_center)
726 CASE (current_orb_center_atom)
727 !
728 current_env%orb_center_name = "ATOM"
729 CASE (current_orb_center_box)
730 !
731 current_env%orb_center_name = "BOX"
732 !
733 ! Is there a user specified nbox?
734 CALL section_vals_val_get(current_section, "NBOX", i_vals=nbox)
735 CASE DEFAULT
736 cpabort("Unknown orbital center, try again...")
737 END SELECT
738
739 CALL section_vals_val_get(current_section, "FORCE_NO_FULL", l_val=force_no_full)
740 IF (force_no_full) current_env%full = .false.
741 !
742 ! Check if restat also psi0 should be restarted
743 !IF(current_env%restart_current .AND. scf_control%density_guess/=restart_guess)THEN
744 ! CPABORT("restart_nmr requires density_guess=restart")
745 !ENDIF
746 !
747 ! check that the psi0 are localized and you have all the centers
748 cpassert(linres_control%localized_psi0)
749 IF (output_unit > 0) THEN
750 WRITE (output_unit, '(A)') &
751 ' To get CURRENT parameters within PBC you need localized zero order orbitals '
752 END IF
753 qs_loc_env => linres_control%qs_loc_env
754 CALL get_qs_loc_env(qs_loc_env, localized_wfn_control=localized_wfn_control)
755 !
756 !
757 ALLOCATE (current_env%centers_set(nspins), current_env%center_list(nspins), &
758 state_list(max_states, nspins))
759 state_list(:, :) = huge(0)
760 nstate_list(:) = huge(0)
761 !
762 !
763 !
764 ! if qmmm is selected change the definition of the center of the box 0 -> L/2
765 IF (current_env%do_qmmm) THEN
766 DO ispin = 1, nspins
767 DO istate = 1, SIZE(localized_wfn_control%centers_set(ispin)%array, 2)
768 ! just to be sure...
769 r(:) = pbc(localized_wfn_control%centers_set(ispin)%array(:, istate), cell)
770 IF (r(1) .LT. 0.0_dp) THEN
771 localized_wfn_control%centers_set(ispin)%array(1, istate) = &
772 r(1) + cell%hmat(1, 1)
773 END IF
774 IF (r(2) .LT. 0.0_dp) THEN
775 localized_wfn_control%centers_set(ispin)%array(2, istate) = &
776 r(2) + cell%hmat(2, 2)
777 END IF
778 IF (r(3) .LT. 0.0_dp) THEN
779 localized_wfn_control%centers_set(ispin)%array(3, istate) = &
780 r(3) + cell%hmat(3, 3)
781 END IF
782 END DO
783 END DO
784 END IF
785 !
786 !
787 !
788 ! if the user has requested to compute the response for a subset of the states
789 ! we collect them here. it requies the states to be localized!
790 CALL section_vals_val_get(current_section, "SELECTED_STATES_ATOM_RADIUS", &
791 r_val=current_env%selected_states_atom_radius)
792 CALL section_vals_val_get(current_section, "SELECTED_STATES_ON_ATOM_LIST", n_rep_val=n_rep)
793 !
794 current_env%do_selected_states = n_rep .GT. 0
795 !
796 ! for the moment selected states doesnt work with the preconditioner FULL_ALL
797 IF (linres_control%preconditioner_type .EQ. ot_precond_full_all .AND. current_env%do_selected_states) &
798 cpabort("Selected states doesnt work with the preconditioner FULL_ALL")
799 !
800 !
801 NULLIFY (current_env%selected_states_on_atom_list)
802 n = 0
803 DO ir = 1, n_rep
804 NULLIFY (list)
805 CALL section_vals_val_get(current_section, "SELECTED_STATES_ON_ATOM_LIST", &
806 i_rep_val=ir, i_vals=list)
807 IF (ASSOCIATED(list)) THEN
808 CALL reallocate(current_env%selected_states_on_atom_list, 1, n + SIZE(list))
809 DO ini = 1, SIZE(list)
810 current_env%selected_states_on_atom_list(ini + n) = list(ini)
811 END DO
812 n = n + SIZE(list)
813 END IF
814 END DO
815 !
816 ! build the subset
817 IF (current_env%do_selected_states) THEN
818 selected_states_on_atom_list => current_env%selected_states_on_atom_list
819 DO ispin = 1, nspins
820 center_array => localized_wfn_control%centers_set(ispin)%array
821 nstate = 0
822 DO istate = 1, SIZE(center_array, 2)
823 DO i = 1, SIZE(selected_states_on_atom_list, 1)
824 iatom = selected_states_on_atom_list(i)
825 r(:) = pbc(center_array(1:3, istate) - particle_set(iatom)%r(:), cell)
826 ! SQRT(DOT_PRODUCT(r, r)) .LE. current_env%selected_states_atom_radius
827 IF ((dot_product(r, r)) .LE. (current_env%selected_states_atom_radius &
828 *current_env%selected_states_atom_radius)) &
829 THEN
830 !
831 ! add the state to the list
832 nstate = nstate + 1
833 state_list(nstate, ispin) = istate
834 EXIT
835 END IF
836 END DO
837 END DO
838 nstate_list(ispin) = nstate
839 END DO
840 ELSE
841 DO ispin = 1, nspins
842 center_array => localized_wfn_control%centers_set(ispin)%array
843 nstate = 0
844 DO istate = 1, SIZE(center_array, 2)
845 nstate = nstate + 1
846 state_list(nstate, ispin) = istate
847 END DO
848 nstate_list(ispin) = nstate
849 END DO
850 END IF
851 !
852 !
853 !
854 ! clustering the states
855 DO ispin = 1, nspins
856 nstate = nstate_list(ispin)
857 current_env%nstates(ispin) = nstate
858 !
859 ALLOCATE (current_env%center_list(ispin)%array(2, nstate + 1), &
860 current_env%centers_set(ispin)%array(3, nstate))
861 current_env%center_list(ispin)%array(:, :) = huge(0)
862 current_env%centers_set(ispin)%array(:, :) = huge(0.0_dp)
863 !
864 center_array => localized_wfn_control%centers_set(ispin)%array
865 !
866 ! point to the psi0 centers
867 SELECT CASE (current_env%orb_center)
868 CASE (current_orb_center_wannier)
869 !
870 ! use the wannier center as -center-
871 current_env%nbr_center(ispin) = nstate
872 DO is = 1, nstate
873 istate = state_list(is, ispin)
874 current_env%centers_set(ispin)%array(1:3, is) = center_array(1:3, istate)
875 current_env%center_list(ispin)%array(1, is) = is
876 current_env%center_list(ispin)%array(2, is) = istate
877 END DO
878 current_env%center_list(ispin)%array(1, nstate + 1) = nstate + 1
879 !
880 CASE (current_orb_center_common)
881 !
882 ! use a common -center-
883 current_env%centers_set(ispin)%array(:, 1) = common_center(:)
884 current_env%nbr_center(ispin) = 1
885 current_env%center_list(ispin)%array(1, 1) = 1
886 current_env%center_list(ispin)%array(1, 2) = nstate + 1
887 DO is = 1, nstate
888 istate = state_list(is, ispin)
889 current_env%center_list(ispin)%array(2, is) = istate
890 END DO
891 !
892 CASE (current_orb_center_atom)
893 !
894 ! use the atom as -center-
895 ALLOCATE (buff(nstate_list(ispin)))
896 buff(:) = 0
897 !
898 DO is = 1, nstate
899 istate = state_list(is, ispin)
900 mdist = huge(0.0_dp)
901 DO iatom = 1, natom
902 r = pbc(particle_set(iatom)%r(:), cell)
903 rab = pbc(r, center_array(1:3, istate), cell)
904 dist = sqrt(rab(1)**2 + rab(2)**2 + rab(3)**2)
905 IF (dist .LT. mdist) THEN
906 buff(is) = iatom
907 mdist = dist
908 END IF
909 END DO
910 END DO
911 !
912 i = 0
913 ii = 1
914 current_env%center_list(ispin)%array(1, 1) = 1
915 DO iatom = 1, natom
916 j = 0
917 is0 = .true.
918 DO is = 1, nstate
919 istate = state_list(is, ispin)
920 IF (buff(is) .EQ. iatom) THEN
921 j = j + 1
922 i = i + 1
923 is0 = .false.
924 current_env%center_list(ispin)%array(2, i) = istate
925 END IF
926 END DO
927 IF (.NOT. is0) THEN
928 IF (output_unit > 0) THEN
929 WRITE (output_unit, '(T2,A,I6,A,I6)') 'clustering ', j, ' center(s) on atom ', iatom
930 END IF
931 current_env%center_list(ispin)%array(1, ii + 1) = &
932 current_env%center_list(ispin)%array(1, ii) + j
933 current_env%centers_set(ispin)%array(:, ii) = &
934 pbc(particle_set(iatom)%r, cell)
935 ii = ii + 1
936 END IF
937 END DO
938 current_env%nbr_center(ispin) = ii - 1
939 !
940 DEALLOCATE (buff)
941 CASE (current_orb_center_box)
942 !
943 ! use boxes as -center-
944 nbr_box = nbox(1)*nbox(2)*nbox(3)
945 ALLOCATE (rbuff(3, nbr_box), buff(nstate))
946 rbuff(:, :) = huge(0.0_dp)
947 buff(:) = 0
948 !
949 ibox = 1
950 DO iz = 1, nbox(3)
951 DO iy = 1, nbox(2)
952 DO ix = 1, nbox(1)
953 rbuff(1, ibox) = cell%hmat(1, 1)*((real(ix, dp) - 0.5_dp)/real(nbox(1), dp) - 0.5_dp)
954 rbuff(2, ibox) = cell%hmat(2, 2)*((real(iy, dp) - 0.5_dp)/real(nbox(2), dp) - 0.5_dp)
955 rbuff(3, ibox) = cell%hmat(3, 3)*((real(iz, dp) - 0.5_dp)/real(nbox(3), dp) - 0.5_dp)
956 ibox = ibox + 1
957 END DO
958 END DO
959 END DO
960 !
961 DO is = 1, nstate
962 istate = state_list(is, ispin)
963 mdist = huge(0.0_dp)
964 DO ibox = 1, nbr_box
965 rab(:) = pbc(rbuff(:, ibox), center_array(1:3, istate), cell)
966 dist = sqrt(rab(1)**2 + rab(2)**2 + rab(3)**2)
967 IF (dist .LT. mdist) THEN
968 buff(is) = ibox
969 mdist = dist
970 END IF
971 END DO
972 END DO
973 !
974 i = 0
975 ii = 1
976 current_env%center_list(ispin)%array(1, 1) = 1
977 DO ibox = 1, nbr_box
978 j = 0
979 is0 = .true.
980 DO is = 1, nstate
981 istate = state_list(is, ispin)
982 IF (buff(is) .EQ. ibox) THEN
983 j = j + 1
984 i = i + 1
985 is0 = .false.
986 current_env%center_list(ispin)%array(2, i) = istate
987 END IF
988 END DO
989 IF (.NOT. is0) THEN
990 IF (output_unit > 0) THEN
991 WRITE (output_unit, '(T2,A,I6,A,I6)') 'clustering ', j, ' center(s) on box ', ibox
992 END IF
993 current_env%center_list(ispin)%array(1, ii + 1) = &
994 current_env%center_list(ispin)%array(1, ii) + j
995 current_env%centers_set(ispin)%array(:, ii) = rbuff(:, ibox)
996 ii = ii + 1
997 END IF
998 END DO
999 current_env%nbr_center(ispin) = ii - 1
1000 !
1001 DEALLOCATE (buff, rbuff)
1002 CASE DEFAULT
1003 cpabort("Unknown orbital center, try again...")
1004 END SELECT
1005 END DO
1006 !
1007 !
1008 !
1009 ! block the states for each center
1010 ALLOCATE (current_env%psi0_order(nspins))
1011 DO ispin = 1, nspins
1012 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
1013 NULLIFY (tmp_fm_struct)
1014 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
1015 ncol_global=nstate_list(ispin), para_env=para_env, &
1016 context=mo_coeff%matrix_struct%context)
1017 CALL cp_fm_create(current_env%psi0_order(ispin), tmp_fm_struct)
1018 CALL cp_fm_struct_release(tmp_fm_struct)
1019 CALL cp_fm_set_all(current_env%psi0_order(ispin), 0.0_dp)
1020 END DO
1021 !
1022 !
1023 DO ispin = 1, nspins
1024 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
1025 jstate = 0
1026 DO icenter = 1, current_env%nbr_center(ispin)
1027 DO i = current_env%center_list(ispin)%array(1, icenter), &
1028 current_env%center_list(ispin)%array(1, icenter + 1) - 1
1029 !
1030 istate = current_env%center_list(ispin)%array(2, i)
1031 jstate = jstate + 1
1032 !
1033 IF (current_env%do_selected_states) THEN ! this should be removed. always reorder the states
1034 ! the blocking works (so far) with all the precond except FULL_ALL
1035 ! the center_list(ispin)%array(2,i) array is obsolete then
1036 current_env%center_list(ispin)%array(2, i) = jstate
1037 CALL cp_fm_to_fm(mo_coeff, current_env%psi0_order(ispin), 1, istate, jstate)
1038 ELSE
1039 ! for the moment we just copy them
1040 CALL cp_fm_to_fm(mo_coeff, current_env%psi0_order(ispin), 1, istate, istate)
1041 END IF
1042 END DO
1043 END DO
1044 END DO
1045 !
1046 !
1047 !
1048 !
1049 IF (current_env%do_qmmm .AND.&
1050 & (current_env%orb_center .EQ. current_orb_center_wannier .OR. &
1051 current_env%orb_center .EQ. current_orb_center_atom .OR. &
1052 current_env%orb_center .EQ. current_orb_center_box)) THEN
1053 IF (output_unit > 0) THEN
1054 WRITE (output_unit, *) 'WARNING WARNING WARNING WARNING WARNING WARNING WARNING WARNING WARNING '
1055 WRITE (output_unit, *) 'orbital center shifted to match the '
1056 WRITE (output_unit, *) 'center of the box (L/2 L/2 L/2) (superdirty...)'
1057 WRITE (output_unit, *) 'WARNING WARNING WARNING WARNING WARNING WARNING WARNING WARNING WARNING '
1058 END IF
1059 DO ispin = 1, nspins
1060 DO istate = 1, current_env%nbr_center(ispin)
1061 IF (current_env%centers_set(ispin)%array(1, istate) .LE. 0.0_dp) THEN
1062 current_env%centers_set(ispin)%array(1, istate) = &
1063 current_env%centers_set(ispin)%array(1, istate) + cell%hmat(1, 1)
1064 END IF
1065 IF (current_env%centers_set(ispin)%array(2, istate) .LE. 0.0_dp) THEN
1066 current_env%centers_set(ispin)%array(2, istate) = &
1067 current_env%centers_set(ispin)%array(2, istate) + cell%hmat(2, 2)
1068 END IF
1069 IF (current_env%centers_set(ispin)%array(3, istate) .LE. 0.0_dp) THEN
1070 current_env%centers_set(ispin)%array(3, istate) = &
1071 current_env%centers_set(ispin)%array(3, istate) + cell%hmat(3, 3)
1072 END IF
1073 END DO
1074 END DO
1075 ! printing
1076 DO ispin = 1, nspins
1077 IF (output_unit > 0) THEN
1078 WRITE (output_unit, '(/,T2,A,I2)') "WANNIER CENTERS for spin ", ispin
1079 WRITE (output_unit, '(/,T18,A)') "--------------- Shifted Centers --------------- "
1080 DO istate = 1, current_env%nbr_center(ispin)
1081 WRITE (output_unit, '(T5,A,I6,3F12.6)') &
1082 'state ', istate, current_env%centers_set(ispin)%array(1:3, istate)
1083 END DO
1084 END IF
1085 END DO
1086 END IF
1087 !
1088 !
1089 !
1090 ! reorder the center dependent responses
1091 max_nbr_center = maxval(current_env%nbr_center(1:nspins))
1092 current_env%nao = nao
1093 ALLOCATE (current_env%statetrueindex(3, max_nbr_center, nspins))
1094 current_env%statetrueindex(:, :, :) = 0
1095 IF (.true.) THEN
1096 ALLOCATE (state_done(3, max_nbr_center))
1097 DO ispin = 1, nspins
1098 state_done(:, :) = .false.
1099 current_env%statetrueindex(1, 1, ispin) = 1
1100 center(1) = current_env%centers_set(ispin)%array(1, 1)
1101 center(2) = current_env%centers_set(ispin)%array(2, 1)
1102 center(3) = current_env%centers_set(ispin)%array(3, 1)
1103 state_done(1, 1) = .true.
1104 icount = 1
1105 CALL get_mo_set(mo_set=mos(ispin), nmo=nmo)
1106 !
1107 DO idir = 1, 3
1108 ini = 1
1109 IF (idir == 1) ini = 2
1110 DO istate = ini, current_env%nbr_center(ispin)
1111 mdist = huge(0.0_dp)
1112 !
1113 DO istate2 = 1, current_env%nbr_center(ispin)
1114 IF (.NOT. state_done(idir, istate2)) THEN
1115 center2(1) = current_env%centers_set(ispin)%array(1, istate2)
1116 center2(2) = current_env%centers_set(ispin)%array(2, istate2)
1117 center2(3) = current_env%centers_set(ispin)%array(3, istate2)
1118 !
1119 rab = pbc(center, center2, cell)
1120 CALL set_vecp(idir, j, k)
1121 dist = sqrt(rab(j)*rab(j) + rab(k)*rab(k))
1122 !
1123 IF (dist .LT. mdist) THEN
1124 mdist = dist
1125 istate_next = istate2
1126 END IF
1127 END IF
1128 END DO ! istate2
1129 !
1130 icount = icount + 1
1131 state_done(idir, istate_next) = .true.
1132 current_env%statetrueindex(idir, icount, ispin) = istate_next
1133 center(1) = current_env%centers_set(ispin)%array(1, istate_next)
1134 center(2) = current_env%centers_set(ispin)%array(2, istate_next)
1135 center(3) = current_env%centers_set(ispin)%array(3, istate_next)
1136 END DO ! istate
1137 icount = 0
1138 END DO ! idir
1139 END DO
1140 DEALLOCATE (state_done)
1141 ELSE
1142 DO ispin = 1, nspins
1143 DO icenter = 1, current_env%nbr_center(ispin)
1144 current_env%statetrueindex(:, icenter, ispin) = icenter
1145 END DO
1146 END DO
1147 END IF
1148 !
1149 !
1150 IF (output_unit > 0) THEN
1151 WRITE (output_unit, "(T2,A,T60,A)") "CURRENT| Gauge used", &
1152 repeat(' ', 20 - len_trim(current_env%gauge_name))//trim(current_env%gauge_name)
1153 WRITE (output_unit, "(T2,A,T79,L1)") "CURRENT| Use old gauge code", current_env%use_old_gauge_atom
1154 WRITE (output_unit, "(T2,A,T79,L1)") "CURRENT| Compute chi for PBC calculation", current_env%chi_pbc
1155 IF (current_env%gauge .EQ. current_gauge_atom) THEN
1156 WRITE (output_unit, "(T2,A,T68,E12.6)") "CURRENT| Radius for building the gauge", &
1157 current_env%gauge_atom_radius
1158 END IF
1159 WRITE (output_unit, "(T2,A,T60,A)") "CURRENT| Orbital center used", &
1160 repeat(' ', 20 - len_trim(current_env%orb_center_name))//trim(current_env%orb_center_name)
1161 IF (current_env%orb_center .EQ. current_orb_center_common) THEN
1162 WRITE (output_unit, "(T2,A,T50,3F10.6)") "CURRENT| Common center", common_center(1:3)
1163 ELSEIF (current_env%orb_center .EQ. current_orb_center_box) THEN
1164 WRITE (output_unit, "(T2,A,T60,3I5)") "CURRENT| Nbr boxes in each direction", nbox(1:3)
1165 END IF
1166 !
1167 DO ispin = 1, nspins
1168 CALL get_mo_set(mo_set=mos(ispin), nmo=nmo)
1169 WRITE (output_unit, '(T2,A,I6,A,I6,A,I2)') 'CURRENT| Compute', &
1170 nstate_list(ispin), ' selected response functions out of', &
1171 nmo, ' for spin ', ispin
1172 !
1173 WRITE (output_unit, '(T2,A,I6,A,I2)') 'CURRENT| There is a total of ', &
1174 current_env%nbr_center(ispin), ' (clustered) center(s) for spin ', ispin
1175 END DO
1176 END IF
1177
1178 IF (btest(cp_print_key_should_output(logger%iter_info, current_section,&
1179 & "PRINT%RESPONSE_FUNCTION_CUBES"), cp_p_file)) THEN
1180
1181 NULLIFY (bounds, list)
1182 ncubes = 0
1183 CALL section_vals_val_get(current_section,&
1184 & "PRINT%RESPONSE_FUNCTION_CUBES%CUBES_LU_BOUNDS", &
1185 i_vals=bounds)
1186 ncubes = bounds(2) - bounds(1) + 1
1187 IF (ncubes > 0) THEN
1188 ALLOCATE (current_env%list_cubes(ncubes))
1189 DO ir = 1, ncubes
1190 current_env%list_cubes(ir) = bounds(1) + (ir - 1)
1191 END DO
1192 END IF
1193 IF (.NOT. ASSOCIATED(current_env%list_cubes)) THEN
1194 CALL section_vals_val_get(current_section, "PRINT%RESPONSE_FUNCTION_CUBES%CUBES_LIST", &
1195 n_rep_val=n_rep)
1196 ncubes = 0
1197 DO ir = 1, n_rep
1198 NULLIFY (list)
1199 CALL section_vals_val_get(current_section, "PRINT%RESPONSE_FUNCTION_CUBES%CUBES_LIST", &
1200 i_rep_val=ir, i_vals=list)
1201 IF (ASSOCIATED(list)) THEN
1202 CALL reallocate(current_env%list_cubes, 1, ncubes + SIZE(list))
1203 DO ini = 1, SIZE(list)
1204 current_env%list_cubes(ini + ncubes) = list(ini)
1205 END DO
1206 ncubes = ncubes + SIZE(list)
1207 END IF
1208 END DO ! ir
1209 END IF
1210 IF (.NOT. ASSOCIATED(current_env%list_cubes)) THEN
1211 ALLOCATE (current_env%list_cubes(max_states))
1212 DO ir = 1, max_states
1213 current_env%list_cubes(ir) = ir
1214 END DO
1215 END IF
1216 END IF
1217 IF (btest(cp_print_key_should_output(logger%iter_info, current_section,&
1218 & "PRINT%CURRENT_CUBES"), cp_p_file)) THEN
1219 END IF
1220
1221 ! for the chemical shift we need 6 psi1, i.e. 6 optimization procedures
1222 ! They become 9 if full nmr is calculated, i.e. with the correction term too
1223 ! All of them are required at the end of the optimization procedure
1224 ! if the current density and the induced field have to be calculated
1225 ! If instead only the shift is needed, only one psi1 should be enough, providing
1226 ! that after every optimization the corresponding shift contribution is calculated
1227 ! prepare the psi1
1228 !
1229 ! for the rxp we cannot calculate it a priori because it is in facts (r-dk)xp
1230 ! where dk is the center of the orbital it is applied to. We would need nstate operators
1231 ! What we can store is (r-dk)xp|psi0k> for each k, which is a full matrix only
1232 ! Therefore we prepare here the full matrix p_psi0 and rxp_psi0
1233 ! We also need a temporary sparse matrix where to store the integrals during the calculation
1234 ALLOCATE (current_env%p_psi0(nspins, 3), current_env%rxp_psi0(nspins, 3), &
1235 current_env%psi1_p(nspins, 3), current_env%psi1_rxp(nspins, 3))
1236 IF (current_env%full) THEN
1237 ALLOCATE (current_env%psi1_D(nspins, 3))
1238 END IF
1239 DO ispin = 1, nspins
1240 mo_coeff => current_env%psi0_order(ispin)
1241 NULLIFY (tmp_fm_struct)
1242 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
1243 ncol_global=current_env%nstates(ispin), &
1244 context=mo_coeff%matrix_struct%context)
1245 DO idir = 1, 3
1246 CALL cp_fm_create(current_env%psi1_p(ispin, idir), tmp_fm_struct)
1247 CALL cp_fm_create(current_env%psi1_rxp(ispin, idir), tmp_fm_struct)
1248 CALL cp_fm_create(current_env%p_psi0(ispin, idir), tmp_fm_struct)
1249 CALL cp_fm_create(current_env%rxp_psi0(ispin, idir), tmp_fm_struct)
1250 IF (current_env%full) THEN
1251 CALL cp_fm_create(current_env%psi1_D(ispin, idir), tmp_fm_struct)
1252 END IF
1253 END DO
1254 CALL cp_fm_struct_release(tmp_fm_struct)
1255 END DO
1256 !
1257 ! If the current density on the grid needs to be stored
1258 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
1259 ALLOCATE (current_env%jrho1_set(3))
1260 DO idir = 1, 3
1261 NULLIFY (rho_r, rho_g)
1262 ALLOCATE (rho_r(nspins), rho_g(nspins))
1263 DO ispin = 1, nspins
1264 CALL auxbas_pw_pool%create_pw(rho_r(ispin))
1265 CALL pw_zero(rho_r(ispin))
1266 CALL auxbas_pw_pool%create_pw(rho_g(ispin))
1267 CALL pw_zero(rho_g(ispin))
1268 END DO
1269 NULLIFY (current_env%jrho1_set(idir)%rho)
1270 ALLOCATE (current_env%jrho1_set(idir)%rho)
1271 CALL qs_rho_create(current_env%jrho1_set(idir)%rho)
1272 CALL qs_rho_set(current_env%jrho1_set(idir)%rho, &
1273 rho_r=rho_r, rho_g=rho_g)
1274 END DO
1275 !
1276 ! Initialize local current density if GAPW calculation
1277 IF (gapw) THEN
1278 CALL init_jrho_atom_set(jrho1_atom_set, atomic_kind_set, nspins)
1279 CALL set_current_env(current_env=current_env, jrho1_atom_set=jrho1_atom_set)
1280 END IF
1281 !
1282 ALLOCATE (current_env%basisfun_center(3, current_env%nao), &
1283 first_sgf(natom), last_sgf(natom))
1284 current_env%basisfun_center = 0.0_dp
1285 CALL get_particle_set(particle_set, qs_kind_set, &
1286 first_sgf=first_sgf, &
1287 last_sgf=last_sgf)
1288 !Build 3 arrays where for each contracted basis function
1289 !the x y and z coordinates of the center are given
1290 DO iatom = 1, natom
1291 DO iao = first_sgf(iatom), last_sgf(iatom)
1292 DO idir = 1, 3
1293 current_env%basisfun_center(idir, iao) = particle_set(iatom)%r(idir)
1294 END DO
1295 END DO
1296 END DO
1297
1298 DEALLOCATE (first_sgf, last_sgf, state_list)
1299
1300 current_env%simple_done(1:6) = .false.
1301
1302 ALLOCATE (current_env%full_done(3*nspins, max_nbr_center))
1303 current_env%full_done = .false.
1304 !
1305 ! allocate pointer for the gauge
1306 ALLOCATE (current_env%rs_gauge(3, pw_env%gridlevel_info%ngrid_levels))
1307
1308 CALL cp_print_key_finished_output(output_unit, logger, lr_section, "PRINT%PROGRAM_RUN_INFO")
1309 CALL timestop(handle)
1310
1311 END SUBROUTINE current_env_init
1312
1313! **************************************************************************************************
1314!> \brief ...
1315!> \param current_env ...
1316! **************************************************************************************************
1317 SUBROUTINE current_env_cleanup(current_env)
1318
1319 TYPE(current_env_type) :: current_env
1320
1321 INTEGER :: i, idir, ispin
1322
1323 CALL cp_fm_release(current_env%psi0_order)
1324 !
1325 !psi1_p
1326 CALL cp_fm_release(current_env%psi1_p)
1327 !
1328 !psi1_rxp
1329 CALL cp_fm_release(current_env%psi1_rxp)
1330 !
1331 !psi1_D
1332 CALL cp_fm_release(current_env%psi1_D)
1333 !
1334 !p_psi0
1335 CALL cp_fm_release(current_env%p_psi0)
1336 !
1337 !rxp_psi0
1338 CALL cp_fm_release(current_env%rxp_psi0)
1339 !
1340 IF (ASSOCIATED(current_env%centers_set)) THEN
1341 DO ispin = 1, SIZE(current_env%centers_set, 1)
1342 DEALLOCATE (current_env%centers_set(ispin)%array)
1343 END DO
1344 DEALLOCATE (current_env%centers_set)
1345 END IF
1346
1347 IF (ASSOCIATED(current_env%center_list)) THEN
1348 DO ispin = 1, SIZE(current_env%center_list, 1)
1349 DEALLOCATE (current_env%center_list(ispin)%array)
1350 END DO
1351 DEALLOCATE (current_env%center_list)
1352 END IF
1353
1354 IF (ASSOCIATED(current_env%list_cubes)) THEN
1355 DEALLOCATE (current_env%list_cubes)
1356 END IF
1357 ! Current density on the grid
1358 IF (ASSOCIATED(current_env%jrho1_set)) THEN
1359 DO idir = 1, 3
1360 CALL qs_rho_clear(current_env%jrho1_set(idir)%rho)
1361 DEALLOCATE (current_env%jrho1_set(idir)%rho)
1362 END DO
1363 DEALLOCATE (current_env%jrho1_set)
1364 END IF
1365 !
1366 ! Local current density, atom by atom (only gapw)
1367 IF (ASSOCIATED(current_env%jrho1_atom_set)) THEN
1368 CALL deallocate_jrho_atom_set(current_env%jrho1_atom_set)
1369 END IF
1370
1371 !fullnmr_done
1372 IF (ASSOCIATED(current_env%full_done)) THEN
1373 DEALLOCATE (current_env%full_done)
1374 END IF
1375 IF (ASSOCIATED(current_env%basisfun_center)) THEN
1376 DEALLOCATE (current_env%basisfun_center)
1377 END IF
1378 IF (ASSOCIATED(current_env%statetrueindex)) THEN
1379 DEALLOCATE (current_env%statetrueindex)
1380 END IF
1381 IF (ASSOCIATED(current_env%selected_states_on_atom_list)) THEN
1382 DEALLOCATE (current_env%selected_states_on_atom_list)
1383 END IF
1384 ! give back the gauge
1385 IF (ASSOCIATED(current_env%rs_gauge)) THEN
1386 DO idir = 1, 3
1387 DO i = 1, SIZE(current_env%rs_gauge, 2)
1388 CALL rs_grid_release(current_env%rs_gauge(idir, i))
1389 END DO
1390 END DO
1391 DEALLOCATE (current_env%rs_gauge)
1392 END IF
1393
1394 ! give back the buf
1395 IF (ASSOCIATED(current_env%rs_buf)) THEN
1396 DO i = 1, SIZE(current_env%rs_buf)
1397 CALL rs_grid_release(current_env%rs_buf(i))
1398 END DO
1399 DEALLOCATE (current_env%rs_buf)
1400 END IF
1401
1402 END SUBROUTINE current_env_cleanup
1403
1404END MODULE qs_linres_current_utils
cell_types::pbc
Definition cell_types.F:92
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:192
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition cp_dbcsr_operations.F:77
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:85
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:87
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:82
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
memory_utilities::reallocate
Definition memory_utilities.F:27
pw_methods::pw_zero
Definition pw_methods.F:82
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cp_array_utils
various utilities that regard array of different kinds: output, allocation,... maybe it is not a good...
Definition cp_array_utils.F:18
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_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_copy
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
Definition cp_dbcsr_api.F:368
cp_dbcsr_api::dbcsr_set
subroutine, public dbcsr_set(matrix, alpha)
...
Definition cp_dbcsr_api.F:1181
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_dbcsr_operations::cp_dbcsr_sm_fm_multiply
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
Definition cp_dbcsr_operations.F:552
cp_fm_basic_linalg
basic linear algebra operations for full matrices
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_column_scale
subroutine, public cp_fm_column_scale(matrixa, scaling)
scales column i of matrix a with scaling(i)
Definition cp_fm_basic_linalg.F:1885
cp_fm_basic_linalg::cp_fm_scale_and_add
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....
Definition cp_fm_basic_linalg.F:167
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:351
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_set_all
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
Definition cp_fm_types.F:528
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Definition cp_fm_types.F:164
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
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
cp_output_handling
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
Definition cp_output_handling.F:25
cp_output_handling::cp_print_key_unit_nr
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
Definition cp_output_handling.F:853
cp_output_handling::cp_print_key_finished_output
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
Definition cp_output_handling.F:1063
cp_output_handling::cp_p_file
integer, parameter, public cp_p_file
Definition cp_output_handling.F:103
cp_output_handling::cp_print_key_should_output
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
Definition cp_output_handling.F:345
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::current_orb_center_wannier
integer, parameter, public current_orb_center_wannier
Definition input_constants.F:502
input_constants::current_gauge_atom
integer, parameter, public current_gauge_atom
Definition input_constants.F:502
input_constants::current_gauge_r
integer, parameter, public current_gauge_r
Definition input_constants.F:502
input_constants::current_gauge_r_and_step_func
integer, parameter, public current_gauge_r_and_step_func
Definition input_constants.F:502
input_constants::current_orb_center_box
integer, parameter, public current_orb_center_box
Definition input_constants.F:502
input_constants::current_orb_center_common
integer, parameter, public current_orb_center_common
Definition input_constants.F:502
input_constants::current_orb_center_atom
integer, parameter, public current_orb_center_atom
Definition input_constants.F:502
input_constants::ot_precond_full_all
integer, parameter, public ot_precond_full_all
Definition input_constants.F:514
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_get_subs_vals
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
Definition input_section_types.F:731
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1047
input_section_types::section_get_lval
logical function, public section_get_lval(section_vals, keyword_name)
...
Definition input_section_types.F:1012
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_path_length
integer, parameter, public default_path_length
Definition kinds.F:58
list
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
Definition list.F:24
memory_utilities
Utility routines for the memory handling.
Definition memory_utilities.F:14
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
particle_methods
Define methods related to particle_type.
Definition particle_methods.F:14
particle_methods::get_particle_set
subroutine, public get_particle_set(particle_set, qs_kind_set, first_sgf, last_sgf, nsgf, nmao, basis)
Get the components of a particle set.
Definition particle_methods.F:98
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
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_environment_types
Definition qs_environment_types.F:14
qs_environment_types::get_qs_env
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, 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, 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)
Get the QUICKSTEP environment.
Definition qs_environment_types.F:514
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_linres_current_utils
Chemical shift calculation by dfpt Initialization of the nmr_env, creation of the special neighbor li...
Definition qs_linres_current_utils.F:22
qs_linres_current_utils::current_response
subroutine, public current_response(current_env, p_env, qs_env)
...
Definition qs_linres_current_utils.F:128
qs_linres_current_utils::current_env_init
subroutine, public current_env_init(current_env, qs_env)
...
Definition qs_linres_current_utils.F:567
qs_linres_current_utils::current_env_cleanup
subroutine, public current_env_cleanup(current_env)
...
Definition qs_linres_current_utils.F:1318
qs_linres_methods
localize wavefunctions linear response scf
Definition qs_linres_methods.F:15
qs_linres_methods::linres_write_restart
subroutine, public linres_write_restart(qs_env, linres_section, vec, ivec, tag, ind)
...
Definition qs_linres_methods.F:811
qs_linres_methods::linres_read_restart
subroutine, public linres_read_restart(qs_env, linres_section, vec, ivec, tag, ind)
...
Definition qs_linres_methods.F:925
qs_linres_methods::linres_solver
subroutine, public linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, iounit, should_stop)
scf loop to optimize the first order wavefunctions (psi1) given a perturbation as an operator applied...
Definition qs_linres_methods.F:216
qs_linres_op
Calculate the operators p rxp and D needed in the optimization of the different contribution of the f...
Definition qs_linres_op.F:18
qs_linres_op::set_vecp
subroutine, public set_vecp(i1, i2, i3)
...
Definition qs_linres_op.F:1320
qs_linres_types
Type definitiona for linear response calculations.
Definition qs_linres_types.F:12
qs_linres_types::init_jrho_atom_set
subroutine, public init_jrho_atom_set(jrho1_atom_set, atomic_kind_set, nspins)
...
Definition qs_linres_types.F:1091
qs_linres_types::get_current_env
subroutine, public get_current_env(current_env, simple_done, simple_converged, full_done, nao, nstates, gauge, list_cubes, statetrueindex, gauge_name, basisfun_center, nbr_center, center_list, centers_set, psi1_p, psi1_rxp, psi1_d, p_psi0, rxp_psi0, jrho1_atom_set, jrho1_set, chi_tensor, chi_tensor_loc, gauge_atom_radius, rs_gauge, use_old_gauge_atom, chi_pbc, psi0_order)
...
Definition qs_linres_types.F:429
qs_linres_types::set_current_env
subroutine, public set_current_env(current_env, jrho1_atom_set, jrho1_set)
...
Definition qs_linres_types.F:612
qs_linres_types::deallocate_jrho_atom_set
subroutine, public deallocate_jrho_atom_set(jrho_atom_set)
...
Definition qs_linres_types.F:786
qs_loc_methods
Driver for the localization that should be general for all the methods available and all the definiti...
Definition qs_loc_methods.F:21
qs_loc_methods::qs_print_cubes
subroutine, public qs_print_cubes(qs_env, mo_coeff, nstates, state_list, centers, print_key, root, ispin, idir, state0, file_position)
write the cube files for a set of selected states
Definition qs_loc_methods.F:786
qs_loc_types
New version of the module for the localization of the molecular orbitals This should be able to use d...
Definition qs_loc_types.F:25
qs_loc_types::get_qs_loc_env
subroutine, public get_qs_loc_env(qs_loc_env, cell, local_molecules, localized_wfn_control, moloc_coeff, op_sm_set, op_fm_set, para_env, particle_set, weights, dim_op)
...
Definition qs_loc_types.F:317
qs_matrix_pools
wrapper for the pools of matrixes
Definition qs_matrix_pools.F:14
qs_mo_types
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
qs_mo_types::get_mo_set
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Definition qs_mo_types.F:397
qs_neighbor_list_types
Define the neighbor list data types and the corresponding functionality.
Definition qs_neighbor_list_types.F:17
qs_operators_ao
Definition qs_operators_ao.F:13
qs_operators_ao::build_lin_mom_matrix
subroutine, public build_lin_mom_matrix(qs_env, matrix)
Calculation of the linear momentum matrix <mu|∂|nu> over Cartesian Gaussian functions.
Definition qs_operators_ao.F:74
qs_p_env_types
basis types for the calculation of the perturbation of density theory.
Definition qs_p_env_types.F:14
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_create
subroutine, public qs_rho_create(rho)
Allocates a new instance of rho.
Definition qs_rho_types.F:101
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:127
realspace_grid_types
Definition realspace_grid_types.F:18
realspace_grid_types::rs_grid_release
subroutine, public rs_grid_release(rs_grid)
releases the given rs grid (see doc/ReferenceCounting.html)
Definition realspace_grid_types.F:2012
scf_control_types
parameters that control an scf iteration
Definition scf_control_types.F:17
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
cp_array_utils::cp_2d_i_p_type
represent a pointer to a 2d array
Definition cp_array_utils.F:231
cp_array_utils::cp_2d_r_p_type
represent a pointer to a 2d array
Definition cp_array_utils.F:109
cp_control_types::dft_control_type
Definition cp_control_types.F:570
cp_dbcsr_api::dbcsr_distribution_type
Definition cp_dbcsr_api.F:180
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:170
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:113
cp_log_handling::cp_logger_type
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Definition cp_log_handling.F:140
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:763
particle_types::particle_type
Definition particle_types.F:35
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:217
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:171
qs_linres_types::current_env_type
Definition qs_linres_types.F:94
qs_linres_types::jrho_atom_type
Definition qs_linres_types.F:219
qs_linres_types::linres_control_type
General settings for linear response calculations.
Definition qs_linres_types.F:53
qs_loc_types::localized_wfn_control_type
A type that holds controlling information for the calculation of the spread of wfn and the optimizati...
Definition qs_loc_types.F:121
qs_loc_types::qs_loc_env_type
contains all the info needed by quickstep to calculate the spread of a selected set of orbitals and i...
Definition qs_loc_types.F:80
qs_matrix_pools::qs_matrix_pools_type
container for the pools of matrixes used by qs
Definition qs_matrix_pools.F:66
qs_mo_types::mo_set_type
Definition qs_mo_types.F:46
qs_neighbor_list_types::neighbor_list_set_p_type
Definition qs_neighbor_list_types.F:67
qs_p_env_types::qs_p_env_type
Represent a qs system that is perturbed. Can calculate the linear operator and the rhs of the system ...
Definition qs_p_env_types.F:58
scf_control_types::scf_control_type
Definition scf_control_types.F:124