(git:e7e05ae)
qs_linres_polar_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 Polarizability calculation by dfpt
10 !> Initialization of the polar_env,
11 !> Perturbation Hamiltonian by application of the Berry phase operator to psi0
12 !> Write output
13 !> Deallocate everything
14 !> periodic Raman SL February 2013
15 !> \note
16 ! **************************************************************************************************
17 MODULE qs_linres_polar_utils
18  USE bibliography, ONLY: luber2014,&
19  cite_reference
20  USE cell_types, ONLY: cell_type
21  USE cp_control_types, ONLY: dft_control_type
22  USE cp_fm_basic_linalg, ONLY: cp_fm_trace
23  USE cp_fm_struct, ONLY: cp_fm_struct_create,&
24  cp_fm_struct_release,&
25  cp_fm_struct_type
26  USE cp_fm_types, ONLY: cp_fm_create,&
27  cp_fm_get_info,&
28  cp_fm_release,&
29  cp_fm_set_all,&
30  cp_fm_to_fm,&
31  cp_fm_type
32  USE cp_log_handling, ONLY: cp_get_default_logger,&
33  cp_logger_get_default_io_unit,&
34  cp_logger_type
35  USE cp_output_handling, ONLY: cp_p_file,&
36  cp_print_key_finished_output,&
37  cp_print_key_should_output,&
38  cp_print_key_unit_nr
39  USE cp_result_methods, ONLY: cp_results_erase,&
40  put_results
41  USE cp_result_types, ONLY: cp_result_type
42  USE dbcsr_api, ONLY: dbcsr_p_type
43  USE force_env_types, ONLY: force_env_get,&
44  force_env_type
45  USE input_section_types, ONLY: section_vals_get_subs_vals,&
46  section_vals_type,&
47  section_vals_val_get
48  USE kinds, ONLY: default_string_length,&
49  dp
50  USE machine, ONLY: m_flush
51  USE mathconstants, ONLY: twopi
52  USE message_passing, ONLY: mp_para_env_type
53  USE physcon, ONLY: angstrom
54  USE qs_environment_types, ONLY: get_qs_env,&
55  qs_environment_type,&
56  set_qs_env
57  USE qs_linres_methods, ONLY: linres_read_restart,&
58  linres_solver,&
59  linres_write_restart
60  USE qs_linres_types, ONLY: get_polar_env,&
61  linres_control_type,&
62  polar_env_type,&
63  set_polar_env
64  USE qs_matrix_pools, ONLY: qs_matrix_pools_type
65  USE qs_mo_types, ONLY: get_mo_set,&
66  mo_set_type
67  USE qs_p_env_types, ONLY: qs_p_env_type
68 #include "./base/base_uses.f90"
69 
70  IMPLICIT NONE
71 
72  PRIVATE
73 
74  PUBLIC :: polar_env_init, polar_polar, polar_print, polar_response, write_polarisability_tensor
75 
76  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_linres_polar_utils'
77 
78 CONTAINS
79 
80 ! **************************************************************************************************
81 !> \brief Initialize the polar environment
82 !> \param qs_env ...
83 !> \par History
84 !> 06.2018 polar_env integrated into qs_env (MK)
85 ! **************************************************************************************************
86  SUBROUTINE polar_env_init(qs_env)
87 
88  TYPE(qs_environment_type), POINTER :: qs_env
89 
90  CHARACTER(LEN=*), PARAMETER :: routinen = 'polar_env_init'
91 
92  INTEGER :: handle, idir, iounit, ispin, m, nao, &
93  nmo, nspins
94  TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
95  TYPE(cp_fm_type), POINTER :: mo_coeff
96  TYPE(cp_logger_type), POINTER :: logger
97  TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
98  TYPE(dft_control_type), POINTER :: dft_control
99  TYPE(linres_control_type), POINTER :: linres_control
100  TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
101  TYPE(polar_env_type), POINTER :: polar_env
102  TYPE(section_vals_type), POINTER :: lr_section, polar_section
103 
104  CALL timeset(routinen, handle)
105 
106  NULLIFY (dft_control)
107  NULLIFY (linres_control)
108  NULLIFY (logger)
109  NULLIFY (matrix_s)
110  NULLIFY (mos)
111  NULLIFY (polar_env)
112  NULLIFY (lr_section, polar_section)
113 
114  logger => cp_get_default_logger()
115  lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
116 
117  iounit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
118  extension=".linresLog")
119 
120  IF (iounit > 0) THEN
121  WRITE (iounit, "(/,(T2,A))") "POLAR| Starting polarizability calculation", &
122  "POLAR| Initialization of the polar environment"
123  END IF
124 
125  polar_section => section_vals_get_subs_vals(qs_env%input, &
126  "PROPERTIES%LINRES%POLAR")
127 
128  CALL get_qs_env(qs_env=qs_env, &
129  polar_env=polar_env, &
130  dft_control=dft_control, &
131  matrix_s=matrix_s, &
132  linres_control=linres_control, &
133  mos=mos)
134 
135  ! Create polar environment if needed
136  IF (.NOT. ASSOCIATED(polar_env)) THEN
137  ALLOCATE (polar_env)
138  CALL set_qs_env(qs_env=qs_env, polar_env=polar_env)
139  END IF
140 
141  nspins = dft_control%nspins
142 
143  CALL section_vals_val_get(polar_section, "DO_RAMAN", l_val=polar_env%do_raman)
144  CALL section_vals_val_get(polar_section, "PERIODIC_DIPOLE_OPERATOR", l_val=polar_env%do_periodic)
145 
146  ! Allocate components of the polar environment if needed
147  IF (.NOT. ASSOCIATED(polar_env%polar)) THEN
148  ALLOCATE (polar_env%polar(3, 3))
149  polar_env%polar(:, :) = 0.0_dp
150  END IF
151  IF (.NOT. ASSOCIATED(polar_env%dBerry_psi0)) THEN
152  ALLOCATE (polar_env%dBerry_psi0(3, nspins))
153  ELSE
154  ! Remove previous matrices
155  DO ispin = 1, nspins
156  DO idir = 1, 3
157  CALL cp_fm_release(polar_env%dBerry_psi0(idir, ispin))
158  END DO
159  END DO
160  END IF
161  IF (.NOT. ASSOCIATED(polar_env%psi1_dBerry)) THEN
162  ALLOCATE (polar_env%psi1_dBerry(3, nspins))
163  ELSE
164  ! Remove previous matrices
165  DO ispin = 1, nspins
166  DO idir = 1, 3
167  CALL cp_fm_release(polar_env%psi1_dBerry(idir, ispin))
168  END DO
169  END DO
170  END IF
171  DO ispin = 1, nspins
172  CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff, nao=nao, nmo=nmo)
173  CALL cp_fm_get_info(mo_coeff, ncol_global=m, nrow_global=nao)
174  NULLIFY (tmp_fm_struct)
175  CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
176  ncol_global=m, &
177  context=mo_coeff%matrix_struct%context)
178  DO idir = 1, 3
179  CALL cp_fm_create(polar_env%dBerry_psi0(idir, ispin), tmp_fm_struct)
180  CALL cp_fm_create(polar_env%psi1_dBerry(idir, ispin), tmp_fm_struct)
181  END DO
182  CALL cp_fm_struct_release(tmp_fm_struct)
183  END DO
184 
185  CALL cp_print_key_finished_output(iounit, logger, lr_section, &
186  "PRINT%PROGRAM_RUN_INFO")
187 
188  CALL timestop(handle)
189 
190  END SUBROUTINE polar_env_init
191 
192 ! **************************************************************************************************
193 !> \brief ...
194 !> \param qs_env ...
195 !> \par History
196 !> 06.2018 polar_env integrated into qs_env (MK)
197 ! **************************************************************************************************
198  SUBROUTINE polar_polar(qs_env)
199 
200  TYPE(qs_environment_type), POINTER :: qs_env
201 
202  CHARACTER(LEN=*), PARAMETER :: routinen = 'polar_polar'
203 
204  INTEGER :: handle, i, iounit, ispin, nspins, z
205  LOGICAL :: do_periodic, do_raman, run_stopped
206  REAL(dp) :: ptmp
207  REAL(dp), DIMENSION(:, :), POINTER :: polar, polar_tmp
208  TYPE(cell_type), POINTER :: cell
209  TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: dberry_psi0, psi1_dberry
210  TYPE(cp_logger_type), POINTER :: logger
211  TYPE(dft_control_type), POINTER :: dft_control
212  TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
213  TYPE(polar_env_type), POINTER :: polar_env
214 
215  CALL timeset(routinen, handle)
216 
217  NULLIFY (cell, dft_control, polar, psi1_dberry, logger)
218  NULLIFY (mos, dberry_psi0)
219  logger => cp_get_default_logger()
220  iounit = cp_logger_get_default_io_unit(logger)
221 
222  CALL get_qs_env(qs_env=qs_env, &
223  cell=cell, &
224  dft_control=dft_control, &
225  mos=mos, &
226  polar_env=polar_env)
227 
228  nspins = dft_control%nspins
229 
230  CALL get_polar_env(polar_env=polar_env, &
231  do_raman=do_raman, &
232  run_stopped=run_stopped)
233 
234  IF (.NOT. run_stopped .AND. do_raman) THEN
235 
236  CALL cite_reference(luber2014)
237 
238  CALL get_polar_env(polar_env=polar_env, &
239  do_periodic=do_periodic, &
240  dberry_psi0=dberry_psi0, &
241  polar=polar, &
242  psi1_dberry=psi1_dberry)
243 
244  ! Initialize
245  ALLOCATE (polar_tmp(3, 3))
246  polar_tmp(:, :) = 0.0_dp
247 
248  DO i = 1, 3 ! directions of electric field
249  DO z = 1, 3 !dipole directions
250  DO ispin = 1, dft_control%nspins
251  !SL compute trace
252  ptmp = 0.0_dp
253  CALL cp_fm_trace(psi1_dberry(i, ispin), dberry_psi0(z, ispin), ptmp)
254  polar_tmp(i, z) = polar_tmp(i, z) - 2.0_dp*ptmp
255  END DO
256  END DO
257  END DO !spin
258 
259  IF (do_periodic) THEN
260  polar(:, :) = matmul(matmul(cell%hmat, polar_tmp), transpose(cell%hmat))/(twopi*twopi)
261  ELSE
262  polar(:, :) = polar_tmp(:, :)
263  END IF
264  !SL evtl maxocc instead?
265  IF (dft_control%nspins == 1) THEN
266  polar(:, :) = 2.0_dp*polar(:, :)
267  END IF
268 
269  IF (ASSOCIATED(polar_tmp)) THEN
270  DEALLOCATE (polar_tmp)
271  END IF
272 
273  END IF ! do_raman
274 
275  CALL timestop(handle)
276 
277  END SUBROUTINE polar_polar
278 
279 ! **************************************************************************************************
280 !> \brief Print information related to the polarisability tensor
281 !> \param qs_env ...
282 !> \par History
283 !> 06.2018 polar_env integrated into qs_env (MK)
284 ! **************************************************************************************************
285  SUBROUTINE polar_print(qs_env)
286 
287  TYPE(qs_environment_type), POINTER :: qs_env
288 
289  CHARACTER(LEN=default_string_length) :: description
290  INTEGER :: iounit, unit_p
291  LOGICAL :: do_raman, run_stopped
292  REAL(dp), DIMENSION(:, :), POINTER :: polar
293  TYPE(cp_logger_type), POINTER :: logger
294  TYPE(cp_result_type), POINTER :: results
295  TYPE(dft_control_type), POINTER :: dft_control
296  TYPE(mp_para_env_type), POINTER :: para_env
297  TYPE(polar_env_type), POINTER :: polar_env
298  TYPE(section_vals_type), POINTER :: polar_section
299 
300  NULLIFY (logger, dft_control, para_env, results)
301 
302  CALL get_qs_env(qs_env=qs_env, &
303  dft_control=dft_control, &
304  polar_env=polar_env, &
305  results=results, &
306  para_env=para_env)
307 
308  logger => cp_get_default_logger()
309  iounit = cp_logger_get_default_io_unit(logger)
310 
311  polar_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES%POLAR")
312 
313  CALL get_polar_env(polar_env=polar_env, &
314  polar=polar, &
315  do_raman=do_raman, &
316  run_stopped=run_stopped)
317 
318  IF (.NOT. run_stopped .AND. do_raman) THEN
319 
320  description = "[POLAR]"
321  CALL cp_results_erase(results, description=description)
322  CALL put_results(results, description=description, values=polar(:, :))
323 
324  IF (btest(cp_print_key_should_output(logger%iter_info, polar_section, &
325  "PRINT%POLAR_MATRIX"), cp_p_file)) THEN
326 
327  unit_p = cp_print_key_unit_nr(logger, polar_section, "PRINT%POLAR_MATRIX", &
328  extension=".data", middle_name="raman", log_filename=.false.)
329  IF (unit_p > 0) THEN
330  IF (unit_p /= iounit) THEN
331  WRITE (unit_p, *)
332  WRITE (unit_p, '(T10,A)') 'POLARIZABILITY TENSOR (atomic units):'
333  WRITE (unit_p, '(T10,A,3F15.5)') "xx,yy,zz", polar(1, 1), polar(2, 2), polar(3, 3)
334  WRITE (unit_p, '(T10,A,3F15.5)') "xy,xz,yz", polar(1, 2), polar(1, 3), polar(2, 3)
335  WRITE (unit_p, '(T10,A,3F15.5)') "yx,zx,zy", polar(2, 1), polar(3, 1), polar(3, 2)
336  WRITE (unit_p, '(T10,A)') 'POLARIZABILITY TENSOR (Angstrom^3):'
337  WRITE (unit_p, '(T10,A,3F15.5)') "xx,yy,zz", polar(1, 1)*angstrom**3, &
338  polar(2, 2)*angstrom**3, polar(3, 3)*angstrom**3
339  WRITE (unit_p, '(T10,A,3F15.5)') "xy,xz,yz", polar(1, 2)*angstrom**3, &
340  polar(1, 3)*angstrom**3, polar(2, 3)*angstrom**3
341  WRITE (unit_p, '(T10,A,3F15.5)') "yx,zx,zy", polar(2, 1)*angstrom**3, &
342  polar(3, 1)*angstrom**3, polar(3, 2)*angstrom**3
343  CALL cp_print_key_finished_output(unit_p, logger, polar_section, &
344  "PRINT%POLAR_MATRIX")
345  END IF
346  END IF
347  END IF
348  IF (iounit > 0) THEN
349  WRITE (iounit, '(/,T2,A)') &
350  'POLAR| Polarizability tensor [a.u.]'
351  WRITE (iounit, '(T2,A,T24,3(1X,F18.12))') &
352  'POLAR| xx,yy,zz', polar(1, 1), polar(2, 2), polar(3, 3)
353  WRITE (iounit, '(T2,A,T24,3(1X,F18.12))') &
354  'POLAR| xy,xz,yz', polar(1, 2), polar(1, 3), polar(2, 3)
355  WRITE (iounit, '(T2,A,T24,3(1X,F18.12))') &
356  'POLAR| yx,zx,zy', polar(2, 1), polar(3, 1), polar(3, 2)
357  WRITE (iounit, '(/,T2,A)') &
358  'POLAR| Polarizability tensor [ang^3]'
359  WRITE (iounit, '(T2,A,T24,3(1X,F18.12))') &
360  'POLAR| xx,yy,zz', polar(1, 1)*angstrom**3, polar(2, 2)*angstrom**3, polar(3, 3)*angstrom**3
361  WRITE (iounit, '(T2,A,T24,3(1X,F18.12))') &
362  'POLAR| xy,xz,yz', polar(1, 2)*angstrom**3, polar(1, 3)*angstrom**3, polar(2, 3)*angstrom**3
363  WRITE (iounit, '(T2,A,T24,3(1X,F18.12))') &
364  'POLAR| yx,zx,zy', polar(2, 1)*angstrom**3, polar(3, 1)*angstrom**3, polar(3, 2)*angstrom**3
365  END IF
366  END IF
367 
368  END SUBROUTINE polar_print
369 
370 ! **************************************************************************************************
371 !> \brief Calculate the polarisability tensor using response theory
372 !> \param p_env ...
373 !> \param qs_env ...
374 !> \par History
375 !> 06.2018 polar_env integrated into qs_env (MK)
376 ! **************************************************************************************************
377  SUBROUTINE polar_response(p_env, qs_env)
378 
379  TYPE(qs_p_env_type) :: p_env
380  TYPE(qs_environment_type), POINTER :: qs_env
381 
382  CHARACTER(LEN=*), PARAMETER :: routinen = 'polar_response'
383 
384  INTEGER :: handle, idir, iounit, ispin, nao, nmo, &
385  nspins
386  LOGICAL :: do_periodic, do_raman, should_stop
387  TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
388  TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:) :: h1_psi0, psi0_order, psi1
389  TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: dberry_psi0, psi1_dberry
390  TYPE(cp_fm_type), POINTER :: mo_coeff
391  TYPE(cp_logger_type), POINTER :: logger
392  TYPE(dft_control_type), POINTER :: dft_control
393  TYPE(linres_control_type), POINTER :: linres_control
394  TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
395  TYPE(mp_para_env_type), POINTER :: para_env
396  TYPE(polar_env_type), POINTER :: polar_env
397  TYPE(qs_matrix_pools_type), POINTER :: mpools
398  TYPE(section_vals_type), POINTER :: lr_section, polar_section
399 
400  CALL timeset(routinen, handle)
401 
402  NULLIFY (dft_control, linres_control, lr_section, polar_section)
403  NULLIFY (logger, mpools, mo_coeff, para_env)
404  NULLIFY (tmp_fm_struct, psi1_dberry, dberry_psi0)
405 
406  logger => cp_get_default_logger()
407  lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
408  polar_section => section_vals_get_subs_vals(qs_env%input, &
409  "PROPERTIES%LINRES%POLAR")
410 
411  iounit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
412  extension=".linresLog")
413  IF (iounit > 0) THEN
414  WRITE (unit=iounit, fmt="(T2,A,/)") &
415  "POLAR| Self consistent optimization of the response wavefunctions"
416  END IF
417 
418  CALL get_qs_env(qs_env=qs_env, &
419  dft_control=dft_control, &
420  mpools=mpools, &
421  linres_control=linres_control, &
422  mos=mos, &
423  polar_env=polar_env, &
424  para_env=para_env)
425 
426  nspins = dft_control%nspins
427 
428  CALL get_polar_env(polar_env=polar_env, do_raman=do_raman, do_periodic=do_periodic)
429 
430  ! Allocate the vectors
431  ALLOCATE (psi0_order(nspins))
432  ALLOCATE (psi1(nspins), h1_psi0(nspins))
433  DO ispin = 1, nspins
434  CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
435  psi0_order(ispin) = mo_coeff
436  CALL cp_fm_get_info(mo_coeff, ncol_global=nmo, nrow_global=nao)
437  NULLIFY (tmp_fm_struct)
438  CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
439  ncol_global=nmo, &
440  context=mo_coeff%matrix_struct%context)
441  CALL cp_fm_create(psi1(ispin), tmp_fm_struct)
442  CALL cp_fm_create(h1_psi0(ispin), tmp_fm_struct)
443  CALL cp_fm_struct_release(tmp_fm_struct)
444  END DO
445 
446  IF (do_raman) THEN
447  CALL get_polar_env(polar_env=polar_env, &
448  psi1_dberry=psi1_dberry, &
449  dberry_psi0=dberry_psi0)
450  ! Restart
451  IF (linres_control%linres_restart) THEN
452  DO idir = 1, 3
453  CALL linres_read_restart(qs_env, lr_section, psi1_dberry(idir, :), idir, "psi1_dBerry")
454  END DO
455  ELSE
456  DO idir = 1, 3
457  DO ispin = 1, nspins
458  CALL cp_fm_set_all(psi1_dberry(idir, ispin), 0.0_dp)
459  END DO
460  END DO
461  END IF
462  loop_idir: DO idir = 1, 3
463  IF (iounit > 0) THEN
464  IF (do_periodic) THEN
465  WRITE (iounit, "(/,T2,A)") &
466  "POLAR| Response to the perturbation operator Berry phase_"//achar(idir + 119)
467  ELSE
468  WRITE (iounit, "(/,T2,A)") &
469  "POLAR| Response to the perturbation operator R_"//achar(idir + 119)
470  END IF
471  END IF
472  ! Do scf cycle to optimize psi1
473  DO ispin = 1, nspins
474  CALL cp_fm_to_fm(psi1_dberry(idir, ispin), psi1(ispin))
475  CALL cp_fm_to_fm(dberry_psi0(idir, ispin), h1_psi0(ispin))
476  END DO
477  !
478  linres_control%lr_triplet = .false. ! we do singlet response
479  linres_control%do_kernel = .true. ! we do coupled response
480  linres_control%converged = .false.
481  CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, iounit, should_stop)
482 
483  ! Copy the response
484  DO ispin = 1, nspins
485  CALL cp_fm_to_fm(psi1(ispin), psi1_dberry(idir, ispin))
486  END DO
487  !
488  ! Write the new result to the restart file
489  IF (linres_control%linres_restart) THEN
490  CALL linres_write_restart(qs_env, lr_section, psi1_dberry(idir, :), idir, "psi1_dBerry")
491  END IF
492  END DO loop_idir
493  END IF ! do_raman
494 
495  CALL set_polar_env(polar_env, run_stopped=should_stop)
496 
497  ! Clean up
498  CALL cp_fm_release(psi1)
499  CALL cp_fm_release(h1_psi0)
500 
501  DEALLOCATE (psi0_order)
502 
503  CALL cp_print_key_finished_output(iounit, logger, lr_section, &
504  "PRINT%PROGRAM_RUN_INFO")
505 
506  CALL timestop(handle)
507 
508  END SUBROUTINE polar_response
509 
510 ! **************************************************************************************************
511 !> \brief Prints the polarisability tensor to a file during MD runs
512 !> \param force_env ...
513 !> \param motion_section ...
514 !> \param itimes ...
515 !> \param time ...
516 !> \param pos ...
517 !> \param act ...
518 !> \par History
519 !> 06.2018 Creation (MK)
520 !> \author Matthias Krack (MK)
521 ! **************************************************************************************************
522  SUBROUTINE write_polarisability_tensor(force_env, motion_section, itimes, time, pos, act)
523 
524  TYPE(force_env_type), POINTER :: force_env
525  TYPE(section_vals_type), POINTER :: motion_section
526  INTEGER, INTENT(IN) :: itimes
527  REAL(kind=dp), INTENT(IN) :: time
528  CHARACTER(LEN=default_string_length), INTENT(IN), &
529  OPTIONAL :: pos, act
530 
531  CHARACTER(LEN=default_string_length) :: my_act, my_pos
532  INTEGER :: iounit
533  LOGICAL :: do_raman, new_file, run_stopped
534  REAL(kind=dp), DIMENSION(:, :), POINTER :: polar
535  TYPE(cp_logger_type), POINTER :: logger
536  TYPE(polar_env_type), POINTER :: polar_env
537  TYPE(qs_environment_type), POINTER :: qs_env
538 
539  NULLIFY (qs_env)
540 
541  CALL force_env_get(force_env, qs_env=qs_env)
542  IF (ASSOCIATED(qs_env)) THEN
543  NULLIFY (polar_env)
544  CALL get_qs_env(qs_env=qs_env, polar_env=polar_env)
545  IF (ASSOCIATED(polar_env)) THEN
546  CALL get_polar_env(polar_env=polar_env, &
547  polar=polar, &
548  do_raman=do_raman, &
549  run_stopped=run_stopped)
550  IF (.NOT. run_stopped .AND. do_raman) THEN
551  NULLIFY (logger)
552  logger => cp_get_default_logger()
553  my_pos = "APPEND"
554  my_act = "WRITE"
555  IF (PRESENT(pos)) my_pos = pos
556  IF (PRESENT(act)) my_act = act
557  iounit = cp_print_key_unit_nr(logger, motion_section, "PRINT%POLAR_MATRIX", &
558  extension=".polar", file_position=my_pos, &
559  file_action=my_act, file_form="FORMATTED", &
560  is_new_file=new_file)
561  ELSE
562  iounit = 0
563  END IF
564  IF (iounit > 0) THEN
565  IF (new_file) THEN
566  WRITE (unit=iounit, fmt='(A,9(11X,A2," [a.u.]"),6X,A)') &
567  "# Step Time [fs]", "xx", "xy", "xz", "yx", "yy", "yz", "zx", "zy", "zz"
568  END IF
569  WRITE (unit=iounit, fmt='(I8,F12.3,9(1X,F19.8))') itimes, time, &
570  polar(1, 1), polar(1, 2), polar(1, 3), &
571  polar(2, 1), polar(2, 2), polar(2, 3), &
572  polar(3, 1), polar(3, 2), polar(3, 3)
573  CALL m_flush(iounit)
574  CALL cp_print_key_finished_output(iounit, logger, motion_section, "PRINT%POLAR_MATRIX")
575  END IF
576  END IF ! polar_env
577  END IF ! qs_env
578 
579  END SUBROUTINE write_polarisability_tensor
580 
581 END MODULE qs_linres_polar_utils
bibliography
collects all references to literature in CP2K as new algorithms / method are included from literature...
Definition: bibliography.F:28
bibliography::luber2014
integer, save, public luber2014
Definition: bibliography.F:43
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_fm_basic_linalg
basic linear algebra operations for full matrices
Definition: cp_fm_basic_linalg.F:14
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
cp_result_methods
set of type/routines to handle the storage of results in force_envs
Definition: cp_result_methods.F:17
cp_result_methods::cp_results_erase
subroutine, public cp_results_erase(results, description, nval)
erase a part of result_list
Definition: cp_result_methods.F:362
cp_result_types
set of type/routines to handle the storage of results in force_envs
Definition: cp_result_types.F:18
force_env_types
Interface for the force calculations.
Definition: force_env_types.F:18
force_env_types::force_env_get
recursive subroutine, public force_env_get(force_env, in_use, fist_env, qs_env, meta_env, fp_env, subsys, para_env, potential_energy, additional_potential, kinetic_energy, harmonic_shell, kinetic_shell, cell, sub_force_env, qmmm_env, qmmmx_env, eip_env, pwdft_env, globenv, input, force_env_section, method_name_id, root_section, mixed_env, nnp_env, embed_env)
returns various attributes about the force environment
Definition: force_env_types.F:329
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
machine
Machine interface based on Fortran 2003 and POSIX.
Definition: machine.F:17
machine::m_flush
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition: machine.F:106
mathconstants
Definition of mathematical constants and functions.
Definition: mathconstants.F:16
mathconstants::twopi
real(kind=dp), parameter, public twopi
Definition: mathconstants.F:112
message_passing
Interface to the message passing library MPI.
Definition: message_passing.F:23
physcon
Definition of physical constants:
Definition: physcon.F:68
physcon::angstrom
real(kind=dp), parameter, public angstrom
Definition: physcon.F:144
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_environment_types::set_qs_env
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, WannierCentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Set the QUICKSTEP environment.
Definition: qs_environment_types.F:1035
qs_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:808
qs_linres_methods::linres_read_restart
subroutine, public linres_read_restart(qs_env, linres_section, vec, ivec, tag, ind)
...
Definition: qs_linres_methods.F:922
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_polar_utils
Polarizability calculation by dfpt Initialization of the polar_env, Perturbation Hamiltonian by appli...
Definition: qs_linres_polar_utils.F:17
qs_linres_polar_utils::write_polarisability_tensor
subroutine, public write_polarisability_tensor(force_env, motion_section, itimes, time, pos, act)
Prints the polarisability tensor to a file during MD runs.
Definition: qs_linres_polar_utils.F:523
qs_linres_polar_utils::polar_polar
subroutine, public polar_polar(qs_env)
...
Definition: qs_linres_polar_utils.F:199
qs_linres_polar_utils::polar_print
subroutine, public polar_print(qs_env)
Print information related to the polarisability tensor.
Definition: qs_linres_polar_utils.F:286
qs_linres_polar_utils::polar_env_init
subroutine, public polar_env_init(qs_env)
Initialize the polar environment.
Definition: qs_linres_polar_utils.F:87
qs_linres_polar_utils::polar_response
subroutine, public polar_response(p_env, qs_env)
Calculate the polarisability tensor using response theory.
Definition: qs_linres_polar_utils.F:378
qs_linres_types
Type definitiona for linear response calculations.
Definition: qs_linres_types.F:12
qs_linres_types::get_polar_env
subroutine, public get_polar_env(polar_env, do_raman, do_periodic, dBerry_psi0, polar, psi1_dBerry, run_stopped)
...
Definition: qs_linres_types.F:1239
qs_linres_types::set_polar_env
subroutine, public set_polar_env(polar_env, do_raman, do_periodic, dBerry_psi0, polar, psi1_dBerry, run_stopped)
...
Definition: qs_linres_types.F:1270
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_p_env_types
basis types for the calculation of the perturbation of density theory.
Definition: qs_p_env_types.F:14