(git:374b731)
Loading...
Searching...
No Matches
qs_linres_issc_utils.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 Chemical shift calculation by dfpt
10!> Initialization of the issc_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! **************************************************************************************************
19MODULE qs_linres_issc_utils
20 USE atomic_kind_types, ONLY: atomic_kind_type,&
21 get_atomic_kind
22 USE cell_types, ONLY: cell_type,&
23 pbc
24 USE cp_control_types, ONLY: dft_control_type
25 USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
26 USE cp_dbcsr_operations, ONLY: cp_dbcsr_sm_fm_multiply,&
27 dbcsr_allocate_matrix_set,&
28 dbcsr_deallocate_matrix_set
29 USE cp_fm_basic_linalg, ONLY: cp_fm_frobenius_norm,&
30 cp_fm_trace
31 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
32 cp_fm_struct_release,&
33 cp_fm_struct_type
34 USE cp_fm_types, ONLY: cp_fm_create,&
35 cp_fm_get_info,&
36 cp_fm_release,&
37 cp_fm_set_all,&
38 cp_fm_to_fm,&
39 cp_fm_type
40 USE cp_log_handling, ONLY: cp_get_default_logger,&
41 cp_logger_get_default_io_unit,&
42 cp_logger_type
43 USE cp_output_handling, ONLY: cp_p_file,&
44 cp_print_key_finished_output,&
45 cp_print_key_should_output,&
46 cp_print_key_unit_nr
47 USE dbcsr_api, ONLY: dbcsr_convert_offsets_to_sizes,&
48 dbcsr_copy,&
49 dbcsr_create,&
50 dbcsr_distribution_type,&
51 dbcsr_p_type,&
52 dbcsr_set,&
53 dbcsr_type_antisymmetric,&
54 dbcsr_type_symmetric
55 USE input_section_types, ONLY: section_vals_get_subs_vals,&
56 section_vals_type,&
57 section_vals_val_get
58 USE kinds, ONLY: default_string_length,&
59 dp
60 USE mathlib, ONLY: diamat_all
61 USE memory_utilities, ONLY: reallocate
62 USE message_passing, ONLY: mp_para_env_type
63 USE particle_methods, ONLY: get_particle_set
64 USE particle_types, ONLY: particle_type
65 USE physcon, ONLY: a_fine,&
66 e_mass,&
67 hertz,&
68 p_mass
69 USE qs_elec_field, ONLY: build_efg_matrix
70 USE qs_environment_types, ONLY: get_qs_env,&
71 qs_environment_type
72 USE qs_fermi_contact, ONLY: build_fermi_contact_matrix
73 USE qs_kind_types, ONLY: qs_kind_type
74 USE qs_linres_methods, ONLY: linres_solver
75 USE qs_linres_types, ONLY: get_issc_env,&
76 issc_env_type,&
77 linres_control_type
78 USE qs_matrix_pools, ONLY: qs_matrix_pools_type
79 USE qs_mo_types, ONLY: get_mo_set,&
80 mo_set_type
81 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
82 USE qs_p_env_types, ONLY: qs_p_env_type
83 USE qs_spin_orbit, ONLY: build_pso_matrix
84#include "./base/base_uses.f90"
85
86 IMPLICIT NONE
87
88 PRIVATE
89 PUBLIC :: issc_env_cleanup, issc_env_init, issc_response, issc_issc, issc_print
90
91 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_linres_issc_utils'
92
93CONTAINS
94
95! **************************************************************************************************
96!> \brief Initialize the issc environment
97!> \param issc_env ...
98!> \param p_env ...
99!> \param qs_env ...
100! **************************************************************************************************
101 SUBROUTINE issc_response(issc_env, p_env, qs_env)
102 !
103 TYPE(issc_env_type) :: issc_env
104 TYPE(qs_p_env_type) :: p_env
105 TYPE(qs_environment_type), POINTER :: qs_env
106
107 CHARACTER(LEN=*), PARAMETER :: routinen = 'issc_response'
108
109 INTEGER :: handle, idir, ijdir, ispin, jdir, nao, &
110 nmo, nspins, output_unit
111 LOGICAL :: do_dso, do_fc, do_pso, do_sd, should_stop
112 REAL(dp) :: chk, fro
113 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
114 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:) :: h1_psi0, psi0_order, psi1
115 TYPE(cp_fm_type), DIMENSION(:), POINTER :: fc_psi0, psi1_fc
116 TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: dso_psi0, efg_psi0, psi1_dso, psi1_efg, &
117 psi1_pso, pso_psi0
118 TYPE(cp_fm_type), POINTER :: mo_coeff
119 TYPE(cp_logger_type), POINTER :: logger
120 TYPE(dft_control_type), POINTER :: dft_control
121 TYPE(linres_control_type), POINTER :: linres_control
122 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
123 TYPE(mp_para_env_type), POINTER :: para_env
124 TYPE(qs_matrix_pools_type), POINTER :: mpools
125 TYPE(section_vals_type), POINTER :: issc_section, lr_section
126
127 CALL timeset(routinen, handle)
128 !
129 NULLIFY (dft_control, linres_control, lr_section, issc_section)
130 NULLIFY (logger, mpools, mo_coeff, para_env)
131 NULLIFY (tmp_fm_struct, psi1_fc, psi1_efg, psi1_pso, pso_psi0, fc_psi0, efg_psi0)
132
133 logger => cp_get_default_logger()
134 lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
135 issc_section => section_vals_get_subs_vals(qs_env%input, &
136 "PROPERTIES%LINRES%SPINSPIN")
137
138 output_unit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
139 extension=".linresLog")
140 IF (output_unit > 0) THEN
141 WRITE (unit=output_unit, fmt="(T10,A,/)") &
142 "*** Self consistent optimization of the response wavefunctions ***"
143 END IF
144
145 CALL get_qs_env(qs_env=qs_env, &
146 dft_control=dft_control, &
147 mpools=mpools, &
148 linres_control=linres_control, &
149 mos=mos, &
150 para_env=para_env)
151
152 nspins = dft_control%nspins
153
154 CALL get_issc_env(issc_env=issc_env, &
155 !list_cubes=list_cubes, &
156 psi1_efg=psi1_efg, &
157 psi1_pso=psi1_pso, &
158 psi1_dso=psi1_dso, &
159 psi1_fc=psi1_fc, &
160 efg_psi0=efg_psi0, &
161 pso_psi0=pso_psi0, &
162 dso_psi0=dso_psi0, &
163 fc_psi0=fc_psi0, &
164 do_fc=do_fc, &
165 do_sd=do_sd, &
166 do_pso=do_pso, &
167 do_dso=do_dso)
168 !
169 ! allocate the vectors
170 ALLOCATE (psi0_order(nspins))
171 ALLOCATE (psi1(nspins), h1_psi0(nspins))
172 DO ispin = 1, nspins
173 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
174 psi0_order(ispin) = mo_coeff
175 CALL cp_fm_get_info(mo_coeff, ncol_global=nmo, nrow_global=nao)
176 NULLIFY (tmp_fm_struct)
177 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
178 ncol_global=nmo, &
179 context=mo_coeff%matrix_struct%context)
180 CALL cp_fm_create(psi1(ispin), tmp_fm_struct)
181 CALL cp_fm_create(h1_psi0(ispin), tmp_fm_struct)
182 CALL cp_fm_struct_release(tmp_fm_struct)
183 END DO
184 chk = 0.0_dp
185 should_stop = .false.
186 !
187 ! operator efg
188 IF (do_sd) THEN
189 ijdir = 0
190 DO idir = 1, 3
191 DO jdir = idir, 3
192 ijdir = ijdir + 1
193 DO ispin = 1, nspins
194 CALL cp_fm_set_all(psi1_efg(ispin, ijdir), 0.0_dp)
195 END DO
196 IF (output_unit > 0) THEN
197 WRITE (output_unit, "(T10,A)") "Response to the perturbation operator efg_"//achar(idir + 119)//achar(jdir + 119)
198 END IF
199 !
200 !Initial guess for psi1
201 DO ispin = 1, nspins
202 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
203 !CALL cp_fm_to_fm(p_psi0(ispin,ijdir)%matrix, psi1(ispin))
204 !CALL cp_fm_scale(-1.0_dp,psi1(ispin))
205 END DO
206 !
207 !DO scf cycle to optimize psi1
208 DO ispin = 1, nspins
209 CALL cp_fm_to_fm(efg_psi0(ispin, ijdir), h1_psi0(ispin))
210 END DO
211 !
212 !
213 linres_control%lr_triplet = .false.
214 linres_control%do_kernel = .false.
215 linres_control%converged = .false.
216 CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, output_unit, should_stop)
217 !
218 !
219 ! copy the response
220 DO ispin = 1, nspins
221 CALL cp_fm_to_fm(psi1(ispin), psi1_efg(ispin, ijdir))
222 fro = cp_fm_frobenius_norm(psi1(ispin))
223 chk = chk + fro
224 END DO
225 !
226 ! print response functions
227 !IF(BTEST(cp_print_key_should_output(logger%iter_info,issc_section,&
228 ! & "PRINT%RESPONSE_FUNCTION_CUBES"),cp_p_file)) THEN
229 ! ncubes = SIZE(list_cubes,1)
230 ! print_key => section_vals_get_subs_vals(issc_section,"PRINT%RESPONSE_FUNCTION_CUBES")
231 ! DO ispin = 1,nspins
232 ! CALL qs_print_cubes(qs_env,psi1(ispin),ncubes,list_cubes,&
233 ! centers_set(ispin)%array,print_key,'psi1_efg',&
234 ! idir=ijdir,ispin=ispin)
235 ! ENDDO ! ispin
236 !ENDIF ! print response functions
237 !
238 !
239 IF (output_unit > 0) THEN
240 WRITE (output_unit, "(T10,A)") "Write the resulting psi1 in restart file... not implemented yet"
241 END IF
242 !
243 ! Write the result in the restart file
244 END DO ! jdir
245 END DO ! idir
246 END IF
247 !
248 ! operator pso
249 IF (do_pso) THEN
250 DO idir = 1, 3
251 DO ispin = 1, nspins
252 CALL cp_fm_set_all(psi1_pso(ispin, idir), 0.0_dp)
253 END DO
254 IF (output_unit > 0) THEN
255 WRITE (output_unit, "(T10,A)") "Response to the perturbation operator pso_"//achar(idir + 119)
256 END IF
257 !
258 !Initial guess for psi1
259 DO ispin = 1, nspins
260 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
261 !CALL cp_fm_to_fm(rxp_psi0(ispin,idir)%matrix, psi1(ispin))
262 !CALL cp_fm_scale(-1.0_dp,psi1(ispin))
263 END DO
264 !
265 !DO scf cycle to optimize psi1
266 DO ispin = 1, nspins
267 CALL cp_fm_to_fm(pso_psi0(ispin, idir), h1_psi0(ispin))
268 END DO
269 !
270 !
271 linres_control%lr_triplet = .false. ! we do singlet response
272 linres_control%do_kernel = .false. ! we do uncoupled response
273 linres_control%converged = .false.
274 CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, output_unit, should_stop)
275 !
276 !
277 ! copy the response
278 DO ispin = 1, nspins
279 CALL cp_fm_to_fm(psi1(ispin), psi1_pso(ispin, idir))
280 fro = cp_fm_frobenius_norm(psi1(ispin))
281 chk = chk + fro
282 END DO
283 !
284 ! print response functions
285 !IF(BTEST(cp_print_key_should_output(logger%iter_info,issc_section,&
286 ! & "PRINT%RESPONSE_FUNCTION_CUBES"),cp_p_file)) THEN
287 ! ncubes = SIZE(list_cubes,1)
288 ! print_key => section_vals_get_subs_vals(issc_section,"PRINT%RESPONSE_FUNCTION_CUBES")
289 ! DO ispin = 1,nspins
290 ! CALL qs_print_cubes(qs_env,psi1(ispin),ncubes,list_cubes,&
291 ! centers_set(ispin)%array,print_key,'psi1_pso',&
292 ! idir=idir,ispin=ispin)
293 ! ENDDO ! ispin
294 !ENDIF ! print response functions
295 !
296 !
297 IF (output_unit > 0) THEN
298 WRITE (output_unit, "(T10,A)") "Write the resulting psi1 in restart file... not implemented yet"
299 END IF
300 !
301 ! Write the result in the restart file
302 END DO ! idir
303 END IF
304 !
305 ! operator fc
306 IF (do_fc) THEN
307 DO ispin = 1, nspins
308 CALL cp_fm_set_all(psi1_fc(ispin), 0.0_dp)
309 END DO
310 IF (output_unit > 0) THEN
311 WRITE (output_unit, "(T10,A)") "Response to the perturbation operator fc"
312 END IF
313 !
314 !Initial guess for psi1
315 DO ispin = 1, nspins
316 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
317 !CALL cp_fm_to_fm(rxp_psi0(ispin,idir)%matrix, psi1(ispin))
318 !CALL cp_fm_scale(-1.0_dp,psi1(ispin))
319 END DO
320 !
321 !DO scf cycle to optimize psi1
322 DO ispin = 1, nspins
323 CALL cp_fm_to_fm(fc_psi0(ispin), h1_psi0(ispin))
324 END DO
325 !
326 !
327 linres_control%lr_triplet = .true. ! we do triplet response
328 linres_control%do_kernel = .true. ! we do coupled response
329 linres_control%converged = .false.
330 CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, output_unit, should_stop)
331 !
332 !
333 ! copy the response
334 DO ispin = 1, nspins
335 CALL cp_fm_to_fm(psi1(ispin), psi1_fc(ispin))
336 fro = cp_fm_frobenius_norm(psi1(ispin))
337 chk = chk + fro
338 END DO
339 !
340 ! print response functions
341 !IF(BTEST(cp_print_key_should_output(logger%iter_info,issc_section,&
342 ! & "PRINT%RESPONSE_FUNCTION_CUBES"),cp_p_file)) THEN
343 ! ncubes = SIZE(list_cubes,1)
344 ! print_key => section_vals_get_subs_vals(issc_section,"PRINT%RESPONSE_FUNCTION_CUBES")
345 ! DO ispin = 1,nspins
346 ! CALL qs_print_cubes(qs_env,psi1(ispin),ncubes,list_cubes,&
347 ! centers_set(ispin)%array,print_key,'psi1_pso',&
348 ! idir=idir,ispin=ispin)
349 ! ENDDO ! ispin
350 !ENDIF ! print response functions
351 !
352 !
353 IF (output_unit > 0) THEN
354 WRITE (output_unit, "(T10,A)") "Write the resulting psi1 in restart file... not implemented yet"
355 END IF
356 !
357 ! Write the result in the restart file
358 END IF
359
360 !>>>> debugging only
361 !
362 ! here we have the operator r and compute the polarizability for debugging the kernel only
363 IF (do_dso) THEN
364 DO idir = 1, 3
365 DO ispin = 1, nspins
366 CALL cp_fm_set_all(psi1_dso(ispin, idir), 0.0_dp)
367 END DO
368 IF (output_unit > 0) THEN
369 WRITE (output_unit, "(T10,A)") "Response to the perturbation operator r_"//achar(idir + 119)
370 END IF
371 !
372 !Initial guess for psi1
373 DO ispin = 1, nspins
374 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
375 !CALL cp_fm_to_fm(rxp_psi0(ispin,idir)%matrix, psi1(ispin))
376 !CALL cp_fm_scale(-1.0_dp,psi1(ispin))
377 END DO
378 !
379 !DO scf cycle to optimize psi1
380 DO ispin = 1, nspins
381 CALL cp_fm_to_fm(dso_psi0(ispin, idir), h1_psi0(ispin))
382 END DO
383 !
384 !
385 linres_control%lr_triplet = .false. ! we do singlet response
386 linres_control%do_kernel = .true. ! we do uncoupled response
387 linres_control%converged = .false.
388 CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, output_unit, should_stop)
389 !
390 !
391 ! copy the response
392 DO ispin = 1, nspins
393 CALL cp_fm_to_fm(psi1(ispin), psi1_dso(ispin, idir))
394 fro = cp_fm_frobenius_norm(psi1(ispin))
395 chk = chk + fro
396 END DO
397 !
398 ! print response functions
399 !IF(BTEST(cp_print_key_should_output(logger%iter_info,issc_section,&
400 ! & "PRINT%RESPONSE_FUNCTION_CUBES"),cp_p_file)) THEN
401 ! ncubes = SIZE(list_cubes,1)
402 ! print_key => section_vals_get_subs_vals(issc_section,"PRINT%RESPONSE_FUNCTION_CUBES")
403 ! DO ispin = 1,nspins
404 ! CALL qs_print_cubes(qs_env,psi1(ispin),ncubes,list_cubes,&
405 ! centers_set(ispin)%array,print_key,'psi1_pso',&
406 ! idir=idir,ispin=ispin)
407 ! ENDDO ! ispin
408 !ENDIF ! print response functions
409 !
410 !
411 IF (output_unit > 0) THEN
412 WRITE (output_unit, "(T10,A)") "Write the resulting psi1 in restart file... not implemented yet"
413 END IF
414 !
415 ! Write the result in the restart file
416 END DO ! idir
417 END IF
418 !<<<< debugging only
419
420 !
421 !
422 ! print the checksum
423 IF (output_unit > 0) THEN
424 WRITE (output_unit, '(T2,A,E23.16)') 'ISSC| response: CheckSum =', chk
425 END IF
426 !
427 !
428 ! clean up
429 CALL cp_fm_release(psi1)
430 CALL cp_fm_release(h1_psi0)
431 DEALLOCATE (psi0_order)
432 !
433 CALL cp_print_key_finished_output(output_unit, logger, lr_section,&
434 & "PRINT%PROGRAM_RUN_INFO")
435 !
436 CALL timestop(handle)
437 !
438 END SUBROUTINE issc_response
439
440! **************************************************************************************************
441!> \brief ...
442!> \param issc_env ...
443!> \param qs_env ...
444!> \param iatom ...
445! **************************************************************************************************
446 SUBROUTINE issc_issc(issc_env, qs_env, iatom)
447
448 TYPE(issc_env_type) :: issc_env
449 TYPE(qs_environment_type), POINTER :: qs_env
450 INTEGER, INTENT(IN) :: iatom
451
452 CHARACTER(LEN=*), PARAMETER :: routinen = 'issc_issc'
453
454 INTEGER :: handle, ispin, ixyz, jatom, jxyz, natom, &
455 nmo, nspins
456 LOGICAL :: do_dso, do_fc, do_pso, do_sd, gapw
457 REAL(dp) :: buf, facdso, facfc, facpso, facsd, g, &
458 issc_dso, issc_fc, issc_pso, issc_sd, &
459 maxocc
460 REAL(dp), DIMENSION(3) :: r_i, r_j
461 REAL(dp), DIMENSION(:, :, :, :, :), POINTER :: issc
462 TYPE(cell_type), POINTER :: cell
463 TYPE(cp_fm_type), DIMENSION(:), POINTER :: fc_psi0, psi1_fc
464 TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: psi1_dso, psi1_efg, psi1_pso
465 TYPE(cp_fm_type), POINTER :: mo_coeff
466 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_dso, matrix_efg, matrix_fc, &
467 matrix_pso
468 TYPE(dft_control_type), POINTER :: dft_control
469 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
470 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
471 TYPE(section_vals_type), POINTER :: issc_section
472
473 CALL timeset(routinen, handle)
474
475 NULLIFY (cell, dft_control, particle_set, issc, psi1_fc, psi1_efg, psi1_pso)
476 NULLIFY (matrix_efg, matrix_fc, matrix_pso, mos, mo_coeff, fc_psi0)
477
478 CALL get_qs_env(qs_env=qs_env, &
479 cell=cell, &
480 dft_control=dft_control, &
481 particle_set=particle_set, &
482 mos=mos)
483
484 gapw = dft_control%qs_control%gapw
485 natom = SIZE(particle_set, 1)
486 nspins = dft_control%nspins
487
488 CALL get_issc_env(issc_env=issc_env, &
489 matrix_efg=matrix_efg, &
490 matrix_pso=matrix_pso, &
491 matrix_fc=matrix_fc, &
492 matrix_dso=matrix_dso, &
493 psi1_fc=psi1_fc, &
494 psi1_efg=psi1_efg, &
495 psi1_pso=psi1_pso, &
496 psi1_dso=psi1_dso, &
497 fc_psi0=fc_psi0, &
498 issc=issc, &
499 do_fc=do_fc, &
500 do_sd=do_sd, &
501 do_pso=do_pso, &
502 do_dso=do_dso)
503
504 g = e_mass/(2.0_dp*p_mass)
505 facfc = hertz*g**2*a_fine**4
506 facpso = hertz*g**2*a_fine**4
507 facsd = hertz*g**2*a_fine**4
508 facdso = hertz*g**2*a_fine**4
509
510 !
511 !
512 issc_section => section_vals_get_subs_vals(qs_env%input, &
513 & "PROPERTIES%LINRES%SPINSPIN")
514 !
515 ! Initialize
516 DO ispin = 1, nspins
517 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff, maxocc=maxocc)
518 CALL cp_fm_get_info(mo_coeff, ncol_global=nmo)
519
520 DO jatom = 1, natom
521 r_i = particle_set(iatom)%r
522 r_j = particle_set(jatom)%r
523 r_j = pbc(r_i, r_j, cell) + r_i
524 !
525 !
526 !
527 !write(*,*) 'iatom =',iatom,' r_i=',r_i
528 !write(*,*) 'jatom =',jatom,' r_j=',r_j
529 !
530 ! FC term
531 !
532 IF (do_fc .AND. iatom .NE. jatom) THEN
533 !
534 ! build the integral for the jatom
535 CALL dbcsr_set(matrix_fc(1)%matrix, 0.0_dp)
536 CALL build_fermi_contact_matrix(qs_env, matrix_fc, r_j)
537 CALL cp_dbcsr_sm_fm_multiply(matrix_fc(1)%matrix, mo_coeff, &
538 fc_psi0(ispin), ncol=nmo,& ! fc_psi0 a buffer
539 & alpha=1.0_dp)
540
541 CALL cp_fm_trace(fc_psi0(ispin), mo_coeff, buf)
542 WRITE (*, *) ' jatom', jatom, 'tr(P*fc)=', buf
543
544 CALL cp_fm_trace(fc_psi0(ispin), psi1_fc(ispin), buf)
545 issc_fc = 2.0_dp*2.0_dp*maxocc*facfc*buf
546 issc(1, 1, iatom, jatom, 1) = issc(1, 1, iatom, jatom, 1) + issc_fc
547 issc(2, 2, iatom, jatom, 1) = issc(2, 2, iatom, jatom, 1) + issc_fc
548 issc(3, 3, iatom, jatom, 1) = issc(3, 3, iatom, jatom, 1) + issc_fc
549 END IF
550 !
551 ! SD term
552 !
553 IF (do_sd .AND. iatom .NE. jatom) THEN
554 !
555 ! build the integral for the jatom
556 CALL dbcsr_set(matrix_efg(1)%matrix, 0.0_dp)
557 CALL dbcsr_set(matrix_efg(2)%matrix, 0.0_dp)
558 CALL dbcsr_set(matrix_efg(3)%matrix, 0.0_dp)
559 CALL dbcsr_set(matrix_efg(4)%matrix, 0.0_dp)
560 CALL dbcsr_set(matrix_efg(5)%matrix, 0.0_dp)
561 CALL dbcsr_set(matrix_efg(6)%matrix, 0.0_dp)
562 CALL build_efg_matrix(qs_env, matrix_efg, r_j)
563 DO ixyz = 1, 6
564 CALL cp_dbcsr_sm_fm_multiply(matrix_efg(ixyz)%matrix, mo_coeff, &
565 fc_psi0(ispin), ncol=nmo,& ! fc_psi0 a buffer
566 & alpha=1.0_dp, beta=0.0_dp)
567 CALL cp_fm_trace(fc_psi0(ispin), mo_coeff, buf)
568 WRITE (*, *) ' jatom', jatom, ixyz, 'tr(P*efg)=', buf
569 DO jxyz = 1, 6
570 CALL cp_fm_trace(fc_psi0(ispin), psi1_efg(ispin, jxyz), buf)
571 issc_sd = 2.0_dp*maxocc*facsd*buf
572 !issc(ixyz,jxyz,iatom,jatom) = issc_sd
573 !write(*,*) 'pso_',ixyz,jxyz,' iatom',iatom,'jatom',jatom,issc_pso
574 END DO
575 END DO
576 END IF
577 !
578 ! PSO term
579 !
580 IF (do_pso .AND. iatom .NE. jatom) THEN
581 !
582 ! build the integral for the jatom
583 CALL dbcsr_set(matrix_pso(1)%matrix, 0.0_dp)
584 CALL dbcsr_set(matrix_pso(2)%matrix, 0.0_dp)
585 CALL dbcsr_set(matrix_pso(3)%matrix, 0.0_dp)
586 CALL build_pso_matrix(qs_env, matrix_pso, r_j)
587 DO ixyz = 1, 3
588 CALL cp_dbcsr_sm_fm_multiply(matrix_pso(ixyz)%matrix, mo_coeff, &
589 fc_psi0(ispin), ncol=nmo,& ! fc_psi0 a buffer
590 & alpha=1.0_dp, beta=0.0_dp)
591 DO jxyz = 1, 3
592 CALL cp_fm_trace(fc_psi0(ispin), psi1_pso(ispin, jxyz), buf)
593 issc_pso = -2.0_dp*maxocc*facpso*buf
594 issc(ixyz, jxyz, iatom, jatom, 3) = issc(ixyz, jxyz, iatom, jatom, 3) + issc_pso
595 END DO
596 END DO
597 END IF
598 !
599 ! DSO term
600 !
601 !>>>>> for debugging we compute here the polarizability and NOT the DSO term!
602 IF (do_dso .AND. iatom .EQ. natom .AND. jatom .EQ. natom) THEN
603 !
604 ! build the integral for the jatom
605 !CALL dbcsr_set(matrix_dso(1)%matrix,0.0_dp)
606 !CALL dbcsr_set(matrix_dso(2)%matrix,0.0_dp)
607 !CALL dbcsr_set(matrix_dso(3)%matrix,0.0_dp)
608 !CALL dbcsr_set(matrix_dso(4)%matrix,0.0_dp)
609 !CALL dbcsr_set(matrix_dso(5)%matrix,0.0_dp)
610 !CALL dbcsr_set(matrix_dso(6)%matrix,0.0_dp)
611 !CALL build_dso_matrix(qs_env,matrix_dso,r_i,r_j)
612 !DO ixyz = 1,6
613 ! CALL cp_sm_fm_multiply(matrix_dso(ixyz)%matrix,mo_coeff,&
614 ! & fc_psi0(ispin),ncol=nmo,& ! fc_psi0 a buffer
615 ! & alpha=1.0_dp,beta=0.0_dp)
616 ! CALL cp_fm_trace(fc_psi0(ispin),mo_coeff,buf)
617 ! issc_dso = 2.0_dp * maxocc * facdso * buf
618 ! issc(ixyz,jxyz,iatom,jatom,4) = issc_dso
619 !ENDDO
620 !CALL dbcsr_set(matrix_dso(1)%matrix,0.0_dp)
621 !CALL dbcsr_set(matrix_dso(2)%matrix,0.0_dp)
622 !CALL dbcsr_set(matrix_dso(3)%matrix,0.0_dp)
623 !CALL rRc_xyz_ao(matrix_dso,qs_env,(/0.0_dp,0.0_dp,0.0_dp/),1)
624 DO ixyz = 1, 3
625 CALL cp_dbcsr_sm_fm_multiply(matrix_dso(ixyz)%matrix, mo_coeff, &
626 fc_psi0(ispin), ncol=nmo,& ! fc_psi0 a buffer
627 & alpha=1.0_dp, beta=0.0_dp)
628 DO jxyz = 1, 3
629 CALL cp_fm_trace(psi1_dso(ispin, jxyz), fc_psi0(ispin), buf)
630 ! we save the polarizability for a checksum later on !
631 issc_dso = 2.0_dp*maxocc*buf
632 issc(ixyz, jxyz, iatom, jatom, 4) = issc(ixyz, jxyz, iatom, jatom, 4) + issc_dso
633 END DO
634 END DO
635
636 END IF
637 !
638 END DO ! jatom
639 END DO ! ispin
640 !
641 !
642 ! Finalize
643 CALL timestop(handle)
644 !
645 END SUBROUTINE issc_issc
646
647! **************************************************************************************************
648!> \brief ...
649!> \param issc_env ...
650!> \param qs_env ...
651! **************************************************************************************************
652 SUBROUTINE issc_print(issc_env, qs_env)
653 TYPE(issc_env_type) :: issc_env
654 TYPE(qs_environment_type), POINTER :: qs_env
655
656 CHARACTER(LEN=2) :: element_symbol_i, element_symbol_j
657 CHARACTER(LEN=default_string_length) :: name_i, name_j, title
658 INTEGER :: iatom, jatom, natom, output_unit, &
659 unit_atoms
660 LOGICAL :: do_dso, do_fc, do_pso, do_sd, gapw
661 REAL(dp) :: eig(3), issc_iso_dso, issc_iso_fc, &
662 issc_iso_pso, issc_iso_sd, &
663 issc_iso_tot, issc_tmp(3, 3)
664 REAL(dp), DIMENSION(:, :, :, :, :), POINTER :: issc
665 REAL(dp), EXTERNAL :: ddot
666 TYPE(atomic_kind_type), POINTER :: atom_kind_i, atom_kind_j
667 TYPE(cp_logger_type), POINTER :: logger
668 TYPE(dft_control_type), POINTER :: dft_control
669 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
670 TYPE(section_vals_type), POINTER :: issc_section
671
672 NULLIFY (logger, particle_set, atom_kind_i, atom_kind_j, dft_control)
673
674 logger => cp_get_default_logger()
675 output_unit = cp_logger_get_default_io_unit(logger)
676
677 issc_section => section_vals_get_subs_vals(qs_env%input, &
678 "PROPERTIES%LINRES%SPINSPIN")
679
680 CALL get_issc_env(issc_env=issc_env, &
681 issc=issc, &
682 do_fc=do_fc, &
683 do_sd=do_sd, &
684 do_pso=do_pso, &
685 do_dso=do_dso)
686 !
687 CALL get_qs_env(qs_env=qs_env, &
688 dft_control=dft_control, &
689 particle_set=particle_set)
690
691 natom = SIZE(particle_set, 1)
692 gapw = dft_control%qs_control%gapw
693
694 !
695 IF (output_unit > 0) THEN
696 WRITE (output_unit, '(T2,A,E14.6)') 'ISSC| CheckSum K =', &
697 sqrt(ddot(SIZE(issc), issc, 1, issc, 1))
698 END IF
699 !
700 IF (btest(cp_print_key_should_output(logger%iter_info, issc_section, &
701 "PRINT%K_MATRIX"), cp_p_file)) THEN
702
703 unit_atoms = cp_print_key_unit_nr(logger, issc_section, "PRINT%K_MATRIX", &
704 extension=".data", middle_name="K", log_filename=.false.)
705
706 IF (unit_atoms > 0) THEN
707 WRITE (unit_atoms, *)
708 WRITE (unit_atoms, *)
709 WRITE (title, '(A)') "Indirect spin-spin coupling matrix"
710 WRITE (unit_atoms, '(T2,A)') title
711 DO iatom = 1, natom
712 atom_kind_i => particle_set(iatom)%atomic_kind
713 CALL get_atomic_kind(atom_kind_i, name=name_i, element_symbol=element_symbol_i)
714 DO jatom = 1, natom
715 atom_kind_j => particle_set(jatom)%atomic_kind
716 CALL get_atomic_kind(atom_kind_j, name=name_j, element_symbol=element_symbol_j)
717 !
718 IF (iatom .EQ. jatom .AND. .NOT. do_dso) cycle
719 !
720 !
721 ! FC
722 issc_tmp(:, :) = issc(:, :, iatom, jatom, 1)
723 issc_tmp(:, :) = 0.5_dp*(issc_tmp(:, :) + transpose(issc_tmp(:, :)))
724 CALL diamat_all(issc_tmp, eig)
725 issc_iso_fc = (eig(1) + eig(2) + eig(3))/3.0_dp
726 !
727 ! SD
728 issc_tmp(:, :) = issc(:, :, iatom, jatom, 2)
729 issc_tmp(:, :) = 0.5_dp*(issc_tmp(:, :) + transpose(issc_tmp(:, :)))
730 CALL diamat_all(issc_tmp, eig)
731 issc_iso_sd = (eig(1) + eig(2) + eig(3))/3.0_dp
732 !
733 ! PSO
734 issc_tmp(:, :) = issc(:, :, iatom, jatom, 3)
735 issc_tmp(:, :) = 0.5_dp*(issc_tmp(:, :) + transpose(issc_tmp(:, :)))
736 CALL diamat_all(issc_tmp, eig)
737 issc_iso_pso = (eig(1) + eig(2) + eig(3))/3.0_dp
738 !
739 ! DSO
740 issc_tmp(:, :) = issc(:, :, iatom, jatom, 4)
741 issc_tmp(:, :) = 0.5_dp*(issc_tmp(:, :) + transpose(issc_tmp(:, :)))
742 CALL diamat_all(issc_tmp, eig)
743 issc_iso_dso = (eig(1) + eig(2) + eig(3))/3.0_dp
744 !
745 ! TOT
746 issc_iso_tot = issc_iso_fc + issc_iso_sd + issc_iso_dso + issc_iso_pso
747 !
748 !
749 WRITE (unit_atoms, *)
750 WRITE (unit_atoms, '(T2,2(A,I5,A,2X,A2))') 'Indirect spin-spin coupling between ', &
751 iatom, trim(name_i), element_symbol_i, ' and ', &
752 jatom, trim(name_j), element_symbol_j
753 !
754 IF (do_fc) WRITE (unit_atoms, '(T1,A,f12.4,A)') ' Isotropic FC contribution = ', issc_iso_fc, ' Hz'
755 IF (do_sd) WRITE (unit_atoms, '(T1,A,f12.4,A)') ' Isotropic SD contribution = ', issc_iso_sd, ' Hz'
756 IF (do_pso) WRITE (unit_atoms, '(T1,A,f12.4,A)') ' Isotropic PSO contribution = ', issc_iso_pso, ' Hz'
757 !IF(do_dso) WRITE(unit_atoms,'(T1,A,f12.4,A)') ' Isotropic DSO contribution = ',issc_iso_dso,' Hz'
758 IF (do_dso) WRITE (unit_atoms, '(T1,A,f12.4,A)') ' !!! POLARIZABILITY (for the moment) = ', issc_iso_dso, ' Hz'
759 IF (.NOT. do_dso) WRITE (unit_atoms, '(T1,A,f12.4,A)') ' Isotropic coupling = ', issc_iso_tot, ' Hz'
760 END DO
761 END DO
762 END IF
763 CALL cp_print_key_finished_output(unit_atoms, logger, issc_section,&
764 & "PRINT%K_MATRIX")
765 END IF
766 !
767 !
768 END SUBROUTINE issc_print
769
770! **************************************************************************************************
771!> \brief Initialize the issc environment
772!> \param issc_env ...
773!> \param qs_env ...
774! **************************************************************************************************
775 SUBROUTINE issc_env_init(issc_env, qs_env)
776 !
777 TYPE(issc_env_type) :: issc_env
778 TYPE(qs_environment_type), POINTER :: qs_env
779
780 CHARACTER(LEN=*), PARAMETER :: routinen = 'issc_env_init'
781
782 INTEGER :: handle, iatom, idir, ini, ir, ispin, &
783 istat, m, n, n_rep, nao, natom, &
784 nspins, output_unit
785 INTEGER, ALLOCATABLE, DIMENSION(:) :: first_sgf, last_sgf
786 INTEGER, DIMENSION(:), POINTER :: list, row_blk_sizes
787 LOGICAL :: gapw
788 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
789 TYPE(cp_fm_type), POINTER :: mo_coeff
790 TYPE(cp_logger_type), POINTER :: logger
791 TYPE(dbcsr_distribution_type), POINTER :: dbcsr_dist
792 TYPE(dft_control_type), POINTER :: dft_control
793 TYPE(linres_control_type), POINTER :: linres_control
794 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
795 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
796 POINTER :: sab_orb
797 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
798 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
799 TYPE(section_vals_type), POINTER :: issc_section, lr_section
800
801!
802
803 CALL timeset(routinen, handle)
804
805 NULLIFY (linres_control)
806 NULLIFY (logger, issc_section)
807 NULLIFY (tmp_fm_struct)
808 NULLIFY (particle_set, qs_kind_set)
809 NULLIFY (sab_orb)
810
811 logger => cp_get_default_logger()
812 lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
813
814 output_unit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
815 extension=".linresLog")
816
817 CALL issc_env_cleanup(issc_env)
818
819 IF (output_unit > 0) THEN
820 WRITE (output_unit, "(/,T20,A,/)") "*** Start indirect spin-spin coupling Calculation ***"
821 WRITE (output_unit, "(T10,A,/)") "Inizialization of the ISSC environment"
822 END IF
823
824 issc_section => section_vals_get_subs_vals(qs_env%input, &
825 & "PROPERTIES%LINRES%SPINSPIN")
826 !CALL section_vals_val_get(nmr_section,"INTERPOLATE_SHIFT",l_val=nmr_env%interpolate_shift)
827 !CALL section_vals_val_get(nmr_section,"SHIFT_GAPW_RADIUS",r_val=nmr_env%shift_gapw_radius)
828
829 CALL get_qs_env(qs_env=qs_env, &
830 dft_control=dft_control, &
831 linres_control=linres_control, &
832 mos=mos, &
833 sab_orb=sab_orb, &
834 particle_set=particle_set, &
835 qs_kind_set=qs_kind_set, &
836 dbcsr_dist=dbcsr_dist)
837 !
838 !
839 gapw = dft_control%qs_control%gapw
840 nspins = dft_control%nspins
841 natom = SIZE(particle_set, 1)
842 !
843 ! check that the psi0 are localized and you have all the centers
844 IF (.NOT. linres_control%localized_psi0) &
845 CALL cp_warn(__location__, 'To get indirect spin-spin coupling parameters within '// &
846 'PBC you need to localize zero order orbitals')
847 !
848 !
849 ! read terms need to be calculated
850 ! FC
851 CALL section_vals_val_get(issc_section, "DO_FC", l_val=issc_env%do_fc)
852 ! SD
853 CALL section_vals_val_get(issc_section, "DO_SD", l_val=issc_env%do_sd)
854 ! PSO
855 CALL section_vals_val_get(issc_section, "DO_PSO", l_val=issc_env%do_pso)
856 ! DSO
857 CALL section_vals_val_get(issc_section, "DO_DSO", l_val=issc_env%do_dso)
858 !
859 !
860 ! read the list of atoms on which the issc need to be calculated
861 CALL section_vals_val_get(issc_section, "ISSC_ON_ATOM_LIST", n_rep_val=n_rep)
862 !
863 !
864 NULLIFY (issc_env%issc_on_atom_list)
865 n = 0
866 DO ir = 1, n_rep
867 NULLIFY (list)
868 CALL section_vals_val_get(issc_section, "ISSC_ON_ATOM_LIST", i_rep_val=ir, i_vals=list)
869 IF (ASSOCIATED(list)) THEN
870 CALL reallocate(issc_env%issc_on_atom_list, 1, n + SIZE(list))
871 DO ini = 1, SIZE(list)
872 issc_env%issc_on_atom_list(ini + n) = list(ini)
873 END DO
874 n = n + SIZE(list)
875 END IF
876 END DO
877 !
878 IF (.NOT. ASSOCIATED(issc_env%issc_on_atom_list)) THEN
879 ALLOCATE (issc_env%issc_on_atom_list(natom), stat=istat)
880 cpassert(istat .EQ. 0)
881 DO iatom = 1, natom
882 issc_env%issc_on_atom_list(iatom) = iatom
883 END DO
884 END IF
885 issc_env%issc_natms = SIZE(issc_env%issc_on_atom_list)
886 !
887 !
888 ! Initialize the issc tensor
889 ALLOCATE (issc_env%issc(3, 3, issc_env%issc_natms, issc_env%issc_natms, 4), &
890 issc_env%issc_loc(3, 3, issc_env%issc_natms, issc_env%issc_natms, 4), &
891 stat=istat)
892 cpassert(istat == 0)
893 issc_env%issc(:, :, :, :, :) = 0.0_dp
894 issc_env%issc_loc(:, :, :, :, :) = 0.0_dp
895 !
896 ! allocation
897 ALLOCATE (issc_env%efg_psi0(nspins, 6), issc_env%pso_psi0(nspins, 3), issc_env%fc_psi0(nspins), &
898 issc_env%psi1_efg(nspins, 6), issc_env%psi1_pso(nspins, 3), issc_env%psi1_fc(nspins), &
899 issc_env%dso_psi0(nspins, 3), issc_env%psi1_dso(nspins, 3), &
900 stat=istat)
901 cpassert(istat == 0)
902 DO ispin = 1, nspins
903 !mo_coeff => current_env%psi0_order(ispin)%matrix
904 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
905 CALL cp_fm_get_info(mo_coeff, ncol_global=m, nrow_global=nao)
906
907 NULLIFY (tmp_fm_struct)
908 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
909 ncol_global=m, &
910 context=mo_coeff%matrix_struct%context)
911 DO idir = 1, 6
912 CALL cp_fm_create(issc_env%psi1_efg(ispin, idir), tmp_fm_struct)
913 CALL cp_fm_create(issc_env%efg_psi0(ispin, idir), tmp_fm_struct)
914 END DO
915 DO idir = 1, 3
916 CALL cp_fm_create(issc_env%psi1_pso(ispin, idir), tmp_fm_struct)
917 CALL cp_fm_create(issc_env%pso_psi0(ispin, idir), tmp_fm_struct)
918 CALL cp_fm_create(issc_env%psi1_dso(ispin, idir), tmp_fm_struct)
919 CALL cp_fm_create(issc_env%dso_psi0(ispin, idir), tmp_fm_struct)
920 END DO
921 CALL cp_fm_create(issc_env%psi1_fc(ispin), tmp_fm_struct)
922 CALL cp_fm_create(issc_env%fc_psi0(ispin), tmp_fm_struct)
923 CALL cp_fm_struct_release(tmp_fm_struct)
924 END DO
925 !
926 ! prepare for allocation
927 ALLOCATE (first_sgf(natom))
928 ALLOCATE (last_sgf(natom))
929 CALL get_particle_set(particle_set, qs_kind_set, &
930 first_sgf=first_sgf, &
931 last_sgf=last_sgf)
932 ALLOCATE (row_blk_sizes(natom))
933 CALL dbcsr_convert_offsets_to_sizes(first_sgf, row_blk_sizes, last_sgf)
934 DEALLOCATE (first_sgf)
935 DEALLOCATE (last_sgf)
936
937 !
938 ! efg, pso and fc operators
939 CALL dbcsr_allocate_matrix_set(issc_env%matrix_efg, 6)
940 ALLOCATE (issc_env%matrix_efg(1)%matrix)
941 CALL dbcsr_create(matrix=issc_env%matrix_efg(1)%matrix, &
942 name="efg (3xx-rr)/3", &
943 dist=dbcsr_dist, matrix_type=dbcsr_type_symmetric, &
944 row_blk_size=row_blk_sizes, col_blk_size=row_blk_sizes, &
945 nze=0, mutable_work=.true.)
946 CALL cp_dbcsr_alloc_block_from_nbl(issc_env%matrix_efg(1)%matrix, sab_orb)
947
948 ALLOCATE (issc_env%matrix_efg(2)%matrix, &
949 issc_env%matrix_efg(3)%matrix, issc_env%matrix_efg(4)%matrix, &
950 issc_env%matrix_efg(5)%matrix, issc_env%matrix_efg(6)%matrix)
951 CALL dbcsr_copy(issc_env%matrix_efg(2)%matrix, issc_env%matrix_efg(1)%matrix, &
952 'efg xy')
953 CALL dbcsr_copy(issc_env%matrix_efg(3)%matrix, issc_env%matrix_efg(1)%matrix, &
954 'efg xz')
955 CALL dbcsr_copy(issc_env%matrix_efg(4)%matrix, issc_env%matrix_efg(1)%matrix, &
956 'efg (3yy-rr)/3')
957 CALL dbcsr_copy(issc_env%matrix_efg(5)%matrix, issc_env%matrix_efg(1)%matrix, &
958 'efg yz')
959 CALL dbcsr_copy(issc_env%matrix_efg(6)%matrix, issc_env%matrix_efg(1)%matrix, &
960 'efg (3zz-rr)/3')
961
962 CALL dbcsr_set(issc_env%matrix_efg(1)%matrix, 0.0_dp)
963 CALL dbcsr_set(issc_env%matrix_efg(2)%matrix, 0.0_dp)
964 CALL dbcsr_set(issc_env%matrix_efg(3)%matrix, 0.0_dp)
965 CALL dbcsr_set(issc_env%matrix_efg(4)%matrix, 0.0_dp)
966 CALL dbcsr_set(issc_env%matrix_efg(5)%matrix, 0.0_dp)
967 CALL dbcsr_set(issc_env%matrix_efg(6)%matrix, 0.0_dp)
968 !
969 ! PSO
970 CALL dbcsr_allocate_matrix_set(issc_env%matrix_pso, 3)
971 ALLOCATE (issc_env%matrix_pso(1)%matrix)
972 CALL dbcsr_create(matrix=issc_env%matrix_pso(1)%matrix, &
973 name="pso x", &
974 dist=dbcsr_dist, matrix_type=dbcsr_type_antisymmetric, &
975 row_blk_size=row_blk_sizes, col_blk_size=row_blk_sizes, &
976 nze=0, mutable_work=.true.)
977 CALL cp_dbcsr_alloc_block_from_nbl(issc_env%matrix_pso(1)%matrix, sab_orb)
978
979 ALLOCATE (issc_env%matrix_pso(2)%matrix, issc_env%matrix_pso(3)%matrix)
980 CALL dbcsr_copy(issc_env%matrix_pso(2)%matrix, issc_env%matrix_pso(1)%matrix, &
981 'pso y')
982 CALL dbcsr_copy(issc_env%matrix_pso(3)%matrix, issc_env%matrix_pso(1)%matrix, &
983 'pso z')
984 CALL dbcsr_set(issc_env%matrix_pso(1)%matrix, 0.0_dp)
985 CALL dbcsr_set(issc_env%matrix_pso(2)%matrix, 0.0_dp)
986 CALL dbcsr_set(issc_env%matrix_pso(3)%matrix, 0.0_dp)
987 !
988 ! DSO
989 CALL dbcsr_allocate_matrix_set(issc_env%matrix_dso, 3)
990 ALLOCATE (issc_env%matrix_dso(1)%matrix, issc_env%matrix_dso(2)%matrix, issc_env%matrix_dso(3)%matrix)
991 CALL dbcsr_copy(issc_env%matrix_dso(1)%matrix, issc_env%matrix_efg(1)%matrix, &
992 'dso x')
993 CALL dbcsr_copy(issc_env%matrix_dso(2)%matrix, issc_env%matrix_efg(1)%matrix, &
994 'dso y')
995 CALL dbcsr_copy(issc_env%matrix_dso(3)%matrix, issc_env%matrix_efg(1)%matrix, &
996 'dso z')
997 CALL dbcsr_set(issc_env%matrix_dso(1)%matrix, 0.0_dp)
998 CALL dbcsr_set(issc_env%matrix_dso(2)%matrix, 0.0_dp)
999 CALL dbcsr_set(issc_env%matrix_dso(3)%matrix, 0.0_dp)
1000 !
1001 ! FC
1002 CALL dbcsr_allocate_matrix_set(issc_env%matrix_fc, 1)
1003 ALLOCATE (issc_env%matrix_fc(1)%matrix)
1004 CALL dbcsr_copy(issc_env%matrix_fc(1)%matrix, issc_env%matrix_efg(1)%matrix, &
1005 'fc')
1006 CALL dbcsr_set(issc_env%matrix_fc(1)%matrix, 0.0_dp)
1007
1008 DEALLOCATE (row_blk_sizes)
1009 !
1010 ! Conversion factors
1011 IF (output_unit > 0) THEN
1012 WRITE (output_unit, "(T2,A,T60,I4,A)")&
1013 & "ISSC| spin-spin coupling computed for ", issc_env%issc_natms, ' atoms'
1014 END IF
1015
1016 CALL cp_print_key_finished_output(output_unit, logger, lr_section,&
1017 & "PRINT%PROGRAM_RUN_INFO")
1018
1019 CALL timestop(handle)
1020
1021 END SUBROUTINE issc_env_init
1022
1023! **************************************************************************************************
1024!> \brief Deallocate the issc environment
1025!> \param issc_env ...
1026!> \par History
1027! **************************************************************************************************
1028 SUBROUTINE issc_env_cleanup(issc_env)
1029
1030 TYPE(issc_env_type), INTENT(INOUT) :: issc_env
1031
1032 IF (ASSOCIATED(issc_env%issc_on_atom_list)) THEN
1033 DEALLOCATE (issc_env%issc_on_atom_list)
1034 END IF
1035 IF (ASSOCIATED(issc_env%issc)) THEN
1036 DEALLOCATE (issc_env%issc)
1037 END IF
1038 IF (ASSOCIATED(issc_env%issc_loc)) THEN
1039 DEALLOCATE (issc_env%issc_loc)
1040 END IF
1041 !
1042 !efg_psi0
1043 CALL cp_fm_release(issc_env%efg_psi0)
1044 !
1045 !pso_psi0
1046 CALL cp_fm_release(issc_env%pso_psi0)
1047 !
1048 !dso_psi0
1049 CALL cp_fm_release(issc_env%dso_psi0)
1050 !
1051 !fc_psi0
1052 CALL cp_fm_release(issc_env%fc_psi0)
1053 !
1054 !psi1_efg
1055 CALL cp_fm_release(issc_env%psi1_efg)
1056 !
1057 !psi1_pso
1058 CALL cp_fm_release(issc_env%psi1_pso)
1059 !
1060 !psi1_dso
1061 CALL cp_fm_release(issc_env%psi1_dso)
1062 !
1063 !psi1_fc
1064 CALL cp_fm_release(issc_env%psi1_fc)
1065 !
1066 !matrix_efg
1067 IF (ASSOCIATED(issc_env%matrix_efg)) THEN
1068 CALL dbcsr_deallocate_matrix_set(issc_env%matrix_efg)
1069 END IF
1070 !
1071 !matrix_pso
1072 IF (ASSOCIATED(issc_env%matrix_pso)) THEN
1073 CALL dbcsr_deallocate_matrix_set(issc_env%matrix_pso)
1074 END IF
1075 !
1076 !matrix_dso
1077 IF (ASSOCIATED(issc_env%matrix_dso)) THEN
1078 CALL dbcsr_deallocate_matrix_set(issc_env%matrix_dso)
1079 END IF
1080 !
1081 !matrix_fc
1082 IF (ASSOCIATED(issc_env%matrix_fc)) THEN
1083 CALL dbcsr_deallocate_matrix_set(issc_env%matrix_fc)
1084 END IF
1085
1086 END SUBROUTINE issc_env_cleanup
1087
1088END MODULE qs_linres_issc_utils
cell_types::pbc
Definition cell_types.F:92
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_fm_basic_linalg::cp_fm_trace
Definition cp_fm_basic_linalg.F:74
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:90
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:85
memory_utilities::reallocate
Definition memory_utilities.F:27
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Definition atomic_kind_types.F:138
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
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_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:744
cp_fm_basic_linalg
basic linear algebra operations for full matrices
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_frobenius_norm
real(kind=dp) function, public cp_fm_frobenius_norm(matrix_a)
computes the Frobenius norm of matrix_a
Definition cp_fm_basic_linalg.F:629
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:125
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:320
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1016
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:535
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:167
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_logger_get_default_io_unit
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Definition cp_log_handling.F:515
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:803
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:1013
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_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:724
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:1040
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
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
mathlib
Collection of simple mathematical functions and subroutines.
Definition mathlib.F:15
mathlib::diamat_all
subroutine, public diamat_all(a, eigval, dac)
Diagonalize the symmetric n by n matrix a using the LAPACK library. Only the upper triangle of matrix...
Definition mathlib.F:372
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
physcon
Definition of physical constants:
Definition physcon.F:68
physcon::a_fine
real(kind=dp), parameter, public a_fine
Definition physcon.F:119
physcon::hertz
real(kind=dp), parameter, public hertz
Definition physcon.F:186
physcon::p_mass
real(kind=dp), parameter, public p_mass
Definition physcon.F:112
physcon::e_mass
real(kind=dp), parameter, public e_mass
Definition physcon.F:109
qs_elec_field
Distribution of the electric field gradient integral matrix.
Definition qs_elec_field.F:13
qs_elec_field::build_efg_matrix
subroutine, public build_efg_matrix(qs_env, matrix_efg, rc)
Calculation of the electric field gradient matrix over Cartesian Gaussian functions.
Definition qs_elec_field.F:75
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_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_fermi_contact
Distribution of the Fermi contact integral matrix.
Definition qs_fermi_contact.F:13
qs_fermi_contact::build_fermi_contact_matrix
subroutine, public build_fermi_contact_matrix(qs_env, matrix_fc, rc)
Calculation of the Fermi contact matrix over Cartesian Gaussian functions.
Definition qs_fermi_contact.F:75
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_linres_issc_utils
Chemical shift calculation by dfpt Initialization of the issc_env, creation of the special neighbor l...
Definition qs_linres_issc_utils.F:19
qs_linres_issc_utils::issc_issc
subroutine, public issc_issc(issc_env, qs_env, iatom)
...
Definition qs_linres_issc_utils.F:447
qs_linres_issc_utils::issc_response
subroutine, public issc_response(issc_env, p_env, qs_env)
Initialize the issc environment.
Definition qs_linres_issc_utils.F:102
qs_linres_issc_utils::issc_env_cleanup
subroutine, public issc_env_cleanup(issc_env)
Deallocate the issc environment.
Definition qs_linres_issc_utils.F:1029
qs_linres_issc_utils::issc_env_init
subroutine, public issc_env_init(issc_env, qs_env)
Initialize the issc environment.
Definition qs_linres_issc_utils.F:776
qs_linres_issc_utils::issc_print
subroutine, public issc_print(issc_env, qs_env)
...
Definition qs_linres_issc_utils.F:653
qs_linres_methods
localize wavefunctions linear response scf
Definition qs_linres_methods.F:15
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_types
Type definitiona for linear response calculations.
Definition qs_linres_types.F:12
qs_linres_types::get_issc_env
subroutine, public get_issc_env(issc_env, issc_on_atom_list, issc_gapw_radius, issc_loc, do_fc, do_sd, do_pso, do_dso, issc, interpolate_issc, psi1_efg, psi1_pso, psi1_dso, psi1_fc, efg_psi0, pso_psi0, dso_psi0, fc_psi0, matrix_efg, matrix_pso, matrix_dso, matrix_fc)
...
Definition qs_linres_types.F:558
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_p_env_types
basis types for the calculation of the perturbation of density theory.
Definition qs_p_env_types.F:14
qs_spin_orbit
Distribution of the spin orbit integral matrix.
Definition qs_spin_orbit.F:13
qs_spin_orbit::build_pso_matrix
subroutine, public build_pso_matrix(qs_env, matrix_so, rc)
Calculation of the paramagnetic spin orbit matrix over Cartesian Gaussian functions.
Definition qs_spin_orbit.F:74
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_control_types::dft_control_type
Definition cp_control_types.F:559
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:116
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:126
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
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_linres_types::issc_env_type
Definition qs_linres_types.F:146
qs_linres_types::linres_control_type
General settings for linear response calculations.
Definition qs_linres_types.F:53
qs_matrix_pools::qs_matrix_pools_type
container for the pools of matrixes used by qs
Definition qs_matrix_pools.F:68
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