(git:f2099e5)
Loading...
Searching...
No Matches
qs_dcdr.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Calculate the derivatives of the MO coefficients wrt nuclear coordinates
10!> \author Sandra Luber, Edward Ditler
11! **************************************************************************************************
12
13MODULE qs_dcdr
14
15 USE atomic_kind_types, ONLY: get_atomic_kind
16 USE cell_types, ONLY: cell_type,&
17 pbc
18 USE cp_array_utils, ONLY: cp_2d_r_p_type
19 USE cp_dbcsr_api, ONLY: dbcsr_add,&
20 dbcsr_copy,&
21 dbcsr_desymmetrize,&
22 dbcsr_p_type,&
23 dbcsr_set
24 USE cp_dbcsr_operations, ONLY: cp_dbcsr_sm_fm_multiply,&
25 dbcsr_allocate_matrix_set,&
26 dbcsr_deallocate_matrix_set
27 USE cp_fm_basic_linalg, ONLY: cp_fm_scale,&
28 cp_fm_scale_and_add,&
29 cp_fm_trace
30 USE cp_fm_types, ONLY: cp_fm_create,&
31 cp_fm_get_diag,&
32 cp_fm_release,&
33 cp_fm_set_all,&
34 cp_fm_to_fm,&
35 cp_fm_type
36 USE cp_log_handling, ONLY: cp_get_default_logger,&
37 cp_logger_type
38 USE cp_output_handling, ONLY: cp_print_key_finished_output,&
39 cp_print_key_unit_nr
40 USE input_section_types, ONLY: section_vals_get_subs_vals,&
41 section_vals_type
42 USE kinds, ONLY: dp
43 USE molecule_types, ONLY: molecule_of_atom,&
44 molecule_type
45 USE parallel_gemm_api, ONLY: parallel_gemm
46 USE particle_types, ONLY: particle_type
47 USE qs_dcdr_ao, ONLY: apply_op_constant_term,&
48 core_dr,&
49 d_core_charge_density_dr,&
50 d_vhxc_dr,&
51 hr_mult_by_delta_1d,&
52 vhxc_r_perturbed_basis_functions
53 USE qs_dcdr_utils, ONLY: dcdr_read_restart,&
54 dcdr_write_restart,&
55 multiply_localization,&
56 shift_wannier_into_cell
57 USE qs_environment_types, ONLY: get_qs_env,&
58 qs_environment_type
59 USE qs_kind_types, ONLY: get_qs_kind,&
60 qs_kind_type
61 USE qs_linres_methods, ONLY: linres_solver
62 USE qs_linres_types, ONLY: dcdr_env_type,&
63 get_polar_env,&
64 linres_control_type,&
65 polar_env_type
66 USE qs_mo_types, ONLY: get_mo_set,&
67 mo_set_type
68 USE qs_moments, ONLY: build_local_moment_matrix,&
69 build_local_moments_der_matrix
70 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
71 USE qs_p_env_types, ONLY: qs_p_env_type
72#include "./base/base_uses.f90"
73
74 IMPLICIT NONE
75
76 PRIVATE
77 PUBLIC :: prepare_per_atom, dcdr_response_dr, dcdr_build_op_dr, apt_dr, apt_dr_localization
78
79 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_dcdr'
80
81CONTAINS
82
83! **************************************************************************************************
84!> \brief Prepare the environment for a choice of lambda
85!> \param dcdr_env ...
86!> \param qs_env ...
87!> \author Edward Ditler
88! **************************************************************************************************
89 SUBROUTINE prepare_per_atom(dcdr_env, qs_env)
90 TYPE(dcdr_env_type) :: dcdr_env
91 TYPE(qs_environment_type), POINTER :: qs_env
92
93 CHARACTER(LEN=*), PARAMETER :: routinen = 'prepare_per_atom'
94
95 INTEGER :: handle, i, ispin, j, natom
96 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
97 POINTER :: sab_all
98 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
99 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
100
101 CALL timeset(routinen, handle)
102
103 NULLIFY (sab_all, qs_kind_set, particle_set)
104 CALL get_qs_env(qs_env=qs_env, &
105 sab_all=sab_all, &
106 qs_kind_set=qs_kind_set, &
107 particle_set=particle_set)
108
109 natom = SIZE(particle_set)
110 IF (dcdr_env%distributed_origin) dcdr_env%ref_point(:) = particle_set(dcdr_env%lambda)%r(:)
111
112 dcdr_env%delta_basis_function = 0._dp
113 dcdr_env%delta_basis_function(:, dcdr_env%lambda) = 1._dp
114
115 ! S matrix
116 ! S1 = - < da/dr | b > * delta_a - < a | db/dr > * delta_b
117
118 ! matrix_s(2:4) are anti-symmetric matrices and contain derivatives wrt. to < a |
119 ! = < da/dR | b > = - < da/dr | b > = < a | db/dr >
120 ! matrix_s1(2:4) = d/dR < a | b >
121 ! and it's built as
122 ! = - matrix_s * delta_b + matrix_s * delta_a
123 ! = - < da/dR | b > * delta_b + < da/dR | b > * delta_a
124 ! = + < da/dr | b > * delta_b - < da/dr | b > * delta_a
125 ! = - < a | db/dr > * delta_b - < da/dr | b > * delta_a
126
127 DO i = 1, 3
128 ! S matrix
129 CALL dbcsr_set(dcdr_env%matrix_nosym_temp(i)%matrix, 0._dp)
130 CALL dbcsr_desymmetrize(dcdr_env%matrix_s(1 + i)%matrix, dcdr_env%matrix_s1(1 + i)%matrix)
131 CALL dbcsr_desymmetrize(dcdr_env%matrix_s(1 + i)%matrix, dcdr_env%matrix_nosym_temp(i)%matrix)
132
133 CALL hr_mult_by_delta_1d(dcdr_env%matrix_s1(1 + i)%matrix, qs_kind_set, "ORB", &
134 sab_all, dcdr_env%lambda, direction_or=.true.)
135 CALL hr_mult_by_delta_1d(dcdr_env%matrix_nosym_temp(i)%matrix, qs_kind_set, "ORB", &
136 sab_all, dcdr_env%lambda, direction_or=.false.)
137
138 CALL dbcsr_add(dcdr_env%matrix_s1(1 + i)%matrix, dcdr_env%matrix_nosym_temp(i)%matrix, -1._dp, +1._dp)
139 CALL dbcsr_set(dcdr_env%matrix_nosym_temp(i)%matrix, 0._dp)
140
141 ! T matrix
142 CALL dbcsr_set(dcdr_env%matrix_nosym_temp(i)%matrix, 0._dp)
143 CALL dbcsr_desymmetrize(dcdr_env%matrix_t(1 + i)%matrix, dcdr_env%matrix_t1(1 + i)%matrix)
144 CALL dbcsr_desymmetrize(dcdr_env%matrix_t(1 + i)%matrix, dcdr_env%matrix_nosym_temp(i)%matrix)
145
146 CALL hr_mult_by_delta_1d(dcdr_env%matrix_t1(1 + i)%matrix, qs_kind_set, "ORB", &
147 sab_all, dcdr_env%lambda, direction_or=.true.)
148 CALL hr_mult_by_delta_1d(dcdr_env%matrix_nosym_temp(i)%matrix, qs_kind_set, "ORB", &
149 sab_all, dcdr_env%lambda, direction_or=.false.)
150
151 CALL dbcsr_add(dcdr_env%matrix_t1(1 + i)%matrix, dcdr_env%matrix_nosym_temp(i)%matrix, -1._dp, +1._dp)
152 CALL dbcsr_set(dcdr_env%matrix_nosym_temp(i)%matrix, 0._dp)
153 END DO
154
155 ! Operator:
156 DO ispin = 1, dcdr_env%nspins
157 DO i = 1, 3
158 CALL dbcsr_set(dcdr_env%matrix_ppnl_1(i)%matrix, 0.0_dp)
159 CALL dbcsr_set(dcdr_env%matrix_hc(i)%matrix, 0.0_dp)
160 CALL dbcsr_set(dcdr_env%matrix_vhxc_perturbed_basis(ispin, i)%matrix, 0.0_dp)
161 CALL dbcsr_set(dcdr_env%matrix_vhxc_perturbed_basis(ispin, i + 3)%matrix, 0.0_dp)
162 CALL dbcsr_set(dcdr_env%matrix_d_vhxc_dR(i, ispin)%matrix, 0.0_dp)
163 CALL dbcsr_set(dcdr_env%matrix_core_charge_1(i)%matrix, 0.0_dp)
164 END DO
165 END DO
166
167 CALL core_dr(qs_env, dcdr_env) ! dcdr_env%matrix_ppnl_1, hc
168 CALL d_vhxc_dr(qs_env, dcdr_env) ! dcdr_env%matrix_d_vhxc_dR
169 CALL d_core_charge_density_dr(qs_env, dcdr_env) ! dcdr_env%matrix_core_charge_1
170 CALL vhxc_r_perturbed_basis_functions(qs_env, dcdr_env) ! dcdr_env%matrix_vhxc_perturbed_basis
171
172 ! APT:
173 DO i = 1, 3
174 DO j = 1, 3
175 CALL dbcsr_set(dcdr_env%matrix_difdip(i, j)%matrix, 0._dp)
176 END DO
177 END DO
178
179 CALL build_local_moments_der_matrix(qs_env, dcdr_env%matrix_difdip, 1, 0, &
180 ref_point=dcdr_env%ref_point, ordered=.true., &
181 deltar=dcdr_env%delta_basis_function, neighbor_image=.false.)
182
183 CALL timestop(handle)
184 END SUBROUTINE prepare_per_atom
185
186! **************************************************************************************************
187!> \brief Build the operator for the position perturbation
188!> \param dcdr_env ...
189!> \param qs_env ...
190!> \authors Sandra Luber
191!> Edward Ditler
192!> Ravi Kumar
193!> Rangsiman Ketkaew
194! **************************************************************************************************
195 SUBROUTINE dcdr_build_op_dr(dcdr_env, qs_env)
196
197 TYPE(dcdr_env_type) :: dcdr_env
198 TYPE(qs_environment_type), POINTER :: qs_env
199
200 CHARACTER(LEN=*), PARAMETER :: routinen = 'dcdr_build_op_dR'
201 REAL(kind=dp), PARAMETER :: one = 1.0_dp, zero = 0.0_dp
202
203 INTEGER :: handle, ispin, nao, nmo
204 TYPE(cp_fm_type) :: buf
205 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: opdr_sym
206
207 CALL timeset(routinen, handle)
208
209 nao = dcdr_env%nao
210
211 ! allocate matrix for the sum of the perturbation terms of the operator (dbcsr matrix)
212 NULLIFY (opdr_sym)
213 CALL dbcsr_allocate_matrix_set(opdr_sym, 1)
214 ALLOCATE (opdr_sym(1)%matrix)
215 CALL dbcsr_copy(opdr_sym(1)%matrix, dcdr_env%matrix_s1(1)%matrix) ! symmetric
216 CALL dbcsr_set(opdr_sym(1)%matrix, 0.0_dp)
217
218 DO ispin = 1, dcdr_env%nspins
219 nmo = dcdr_env%nmo(ispin)
220
221 CALL apply_op_constant_term(qs_env, dcdr_env) ! dcdr_env%matrix_apply_op_constant
222 ! Hartree and Exchange-Correlation contributions
223 CALL dbcsr_add(opdr_sym(1)%matrix, dcdr_env%matrix_core_charge_1(dcdr_env%beta)%matrix, zero, one)
224 CALL dbcsr_add(opdr_sym(1)%matrix, dcdr_env%matrix_d_vhxc_dR(dcdr_env%beta, ispin)%matrix, one, one)
225 CALL dbcsr_add(opdr_sym(1)%matrix, dcdr_env%matrix_vhxc_perturbed_basis(ispin, dcdr_env%beta)%matrix, one, one)
226
227 ! Core Hamiltonian contributions
228 CALL dbcsr_add(opdr_sym(1)%matrix, dcdr_env%matrix_hc(dcdr_env%beta)%matrix, one, one)
229 CALL dbcsr_add(opdr_sym(1)%matrix, dcdr_env%matrix_ppnl_1(dcdr_env%beta)%matrix, one, one)
230 CALL dbcsr_add(opdr_sym(1)%matrix, dcdr_env%matrix_apply_op_constant(ispin)%matrix, one, one)
231
232 CALL dbcsr_desymmetrize(opdr_sym(1)%matrix, dcdr_env%hamiltonian1(1)%matrix)
233 CALL dbcsr_add(dcdr_env%hamiltonian1(1)%matrix, dcdr_env%matrix_t1(dcdr_env%beta + 1)%matrix, one, one)
234
235 CALL cp_dbcsr_sm_fm_multiply(dcdr_env%hamiltonian1(1)%matrix, dcdr_env%mo_coeff(ispin), &
236 dcdr_env%op_dR(ispin), ncol=nmo)
237
238 ! The overlap derivative terms for the Sternheimer equation
239 ! buf = mo * (-mo * matrix_ks * mo)
240 CALL cp_fm_create(buf, dcdr_env%likemos_fm_struct(ispin)%struct)
241 CALL parallel_gemm('N', 'N', nao, nmo, nmo, &
242 -1.0_dp, dcdr_env%mo_coeff(ispin), dcdr_env%chc(ispin), &
243 0.0_dp, buf)
244
245 CALL cp_dbcsr_sm_fm_multiply(dcdr_env%matrix_s1(dcdr_env%beta + 1)%matrix, buf, dcdr_env%op_dR(ispin), &
246 nmo, alpha=1.0_dp, beta=1.0_dp)
247 CALL cp_fm_release(buf)
248
249 ! SL multiply by -1 for response solver (H-S<H> C + dR_coupled= - (op_dR)
250 CALL cp_fm_scale(-1.0_dp, dcdr_env%op_dR(ispin))
251
252 IF (dcdr_env%z_matrix_method) THEN
253 CALL cp_fm_to_fm(dcdr_env%op_dR(ispin), dcdr_env%matrix_m_alpha(dcdr_env%beta, ispin))
254 END IF
255
256 END DO
257
258 CALL dbcsr_deallocate_matrix_set(opdr_sym)
259
260 CALL timestop(handle)
261 END SUBROUTINE dcdr_build_op_dr
262
263! **************************************************************************************************
264!> \brief Get the dC/dR by solving the Sternheimer equation, using the op_dR matrix
265!> \param dcdr_env ...
266!> \param p_env ...
267!> \param qs_env ...
268!> \authors SL, ED
269! **************************************************************************************************
270 SUBROUTINE dcdr_response_dr(dcdr_env, p_env, qs_env)
271
272 TYPE(dcdr_env_type) :: dcdr_env
273 TYPE(qs_p_env_type) :: p_env
274 TYPE(qs_environment_type), POINTER :: qs_env
275
276 CHARACTER(LEN=*), PARAMETER :: routinen = 'dcdr_response_dR'
277
278 INTEGER :: handle, ispin, output_unit
279 LOGICAL :: should_stop
280 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:) :: h1_psi0, psi0_order, psi1
281 TYPE(cp_fm_type), POINTER :: mo_coeff
282 TYPE(cp_logger_type), POINTER :: logger
283 TYPE(linres_control_type), POINTER :: linres_control
284 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
285 TYPE(section_vals_type), POINTER :: lr_section
286
287 CALL timeset(routinen, handle)
288 NULLIFY (linres_control, lr_section, logger)
289
290 CALL get_qs_env(qs_env=qs_env, &
291 linres_control=linres_control, &
292 mos=mos)
293
294 logger => cp_get_default_logger()
295 lr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%LINRES")
296
297 output_unit = cp_print_key_unit_nr(logger, lr_section, "PRINT%PROGRAM_RUN_INFO", &
298 extension=".linresLog")
299 IF (output_unit > 0) THEN
300 WRITE (unit=output_unit, fmt="(T10,A,/)") &
301 "*** Self consistent optimization of the response wavefunction ***"
302 END IF
303
304 ! allocate the vectors
305 ALLOCATE (psi0_order(dcdr_env%nspins))
306 ALLOCATE (psi1(dcdr_env%nspins))
307 ALLOCATE (h1_psi0(dcdr_env%nspins))
308
309 DO ispin = 1, dcdr_env%nspins
310 CALL cp_fm_create(psi1(ispin), dcdr_env%likemos_fm_struct(ispin)%struct)
311 CALL cp_fm_create(h1_psi0(ispin), dcdr_env%likemos_fm_struct(ispin)%struct)
312 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
313 psi0_order(ispin) = mo_coeff
314 END DO
315
316 DO ispin = 1, dcdr_env%nspins
317 CALL cp_fm_set_all(psi1(ispin), 0.0_dp)
318 CALL cp_fm_set_all(h1_psi0(ispin), 0.0_dp)
319 END DO
320 ! Restart
321 IF (linres_control%linres_restart) THEN
322 CALL dcdr_read_restart(qs_env, lr_section, psi1, dcdr_env%lambda, dcdr_env%beta, "dCdR")
323 END IF
324
325 IF (output_unit > 0) THEN
326 WRITE (output_unit, "(T10,A,I4,A)") &
327 "Response to the perturbation operator referring to atom ", dcdr_env%lambda, &
328 " displaced in "//achar(dcdr_env%beta + 119)
329 END IF
330 DO ispin = 1, dcdr_env%nspins
331 CALL cp_fm_set_all(dcdr_env%dCR(ispin), 0.0_dp)
332 CALL cp_fm_to_fm(dcdr_env%op_dR(ispin), h1_psi0(ispin))
333 END DO
334
335 linres_control%lr_triplet = .false. ! we do singlet response
336 linres_control%do_kernel = .true.
337 linres_control%converged = .false.
338
339 ! Position perturbation to get dCR
340 ! (H0-E0) psi1 = (H1-E1) psi0
341 ! psi1 = the perturbed wavefunction
342 ! h1_psi0 = (H1-E1-S1*\varepsilon)
343 ! psi0_order = the unperturbed wavefunction
344 CALL linres_solver(p_env, qs_env, psi1, h1_psi0, psi0_order, &
345 output_unit, should_stop)
346 DO ispin = 1, dcdr_env%nspins
347 CALL cp_fm_to_fm(psi1(ispin), dcdr_env%dCR(ispin))
348 END DO
349
350 ! Write the new result to the restart file
351 IF (linres_control%linres_restart) THEN
352 CALL dcdr_write_restart(qs_env, lr_section, psi1, dcdr_env%lambda, dcdr_env%beta, "dCdR")
353 END IF
354
355 ! clean up
356 DO ispin = 1, dcdr_env%nspins
357 CALL cp_fm_release(psi1(ispin))
358 CALL cp_fm_release(h1_psi0(ispin))
359 END DO
360 DEALLOCATE (psi1, h1_psi0, psi0_order)
361 CALL cp_print_key_finished_output(output_unit, logger, lr_section, &
362 "PRINT%PROGRAM_RUN_INFO")
363
364 CALL timestop(handle)
365
366 END SUBROUTINE dcdr_response_dr
367
368! **************************************************************************************************
369!> \brief Calculate atomic polar tensor
370!> \param qs_env ...
371!> \param dcdr_env ...
372!> \authors Sandra Luber
373!> Edward Ditler
374!> Ravi Kumar
375!> Rangsiman Ketkaew
376! **************************************************************************************************
377 SUBROUTINE apt_dr(qs_env, dcdr_env)
378 TYPE(qs_environment_type), POINTER :: qs_env
379 TYPE(dcdr_env_type) :: dcdr_env
380
381 CHARACTER(LEN=*), PARAMETER :: routinen = 'apt_dR'
382
383 INTEGER :: alpha, handle, ikind, ispin, nao, nmo
384 LOGICAL :: ghost
385 REAL(dp) :: apt_basis_derivative, &
386 apt_coeff_derivative, charge, f_spin, &
387 temp1, temp2
388 REAL(dp), DIMENSION(:, :, :), POINTER :: apt_el, apt_nuc
389 TYPE(cp_fm_type) :: overlap1_mo, tmp_fm_like_mos
390 TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: dberry_psi0, psi1_dberry
391 TYPE(cp_fm_type), POINTER :: mo_coeff
392 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
393 TYPE(polar_env_type), POINTER :: polar_env
394 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
395
396 apt_basis_derivative = 0._dp
397 apt_coeff_derivative = 0._dp
398
399 CALL timeset(routinen, handle)
400
401 NULLIFY (qs_kind_set, particle_set)
402 CALL get_qs_env(qs_env=qs_env, &
403 qs_kind_set=qs_kind_set, &
404 particle_set=particle_set)
405
406 nao = dcdr_env%nao
407 apt_el => dcdr_env%apt_el_dcdr
408 apt_nuc => dcdr_env%apt_nuc_dcdr
409
410 f_spin = 2._dp/dcdr_env%nspins
411
412 DO ispin = 1, dcdr_env%nspins
413 ! Compute S^(1,R)_(ij)
414 CALL cp_fm_create(tmp_fm_like_mos, dcdr_env%likemos_fm_struct(ispin)%struct)
415 CALL cp_fm_create(overlap1_mo, dcdr_env%momo_fm_struct(ispin)%struct)
416 nmo = dcdr_env%nmo(ispin)
417 mo_coeff => dcdr_env%mo_coeff(ispin)
418 CALL cp_fm_set_all(tmp_fm_like_mos, 0.0_dp)
419 CALL cp_fm_scale_and_add(0._dp, dcdr_env%dCR_prime(ispin), 1._dp, dcdr_env%dCR(ispin))
420 CALL cp_dbcsr_sm_fm_multiply(dcdr_env%matrix_s1(dcdr_env%beta + 1)%matrix, mo_coeff, &
421 tmp_fm_like_mos, ncol=nmo)
422 CALL parallel_gemm("T", "N", nmo, nmo, nao, &
423 1.0_dp, mo_coeff, tmp_fm_like_mos, &
424 0.0_dp, overlap1_mo)
425
426 ! C^1 <- -dCR - 0.5 * mo_coeff @ S1_ij
427 ! We get the negative of the coefficients out of the linres solver
428 ! And apply the constant correction due to the overlap derivative.
429 CALL parallel_gemm("N", "N", nao, nmo, nmo, &
430 -0.5_dp, mo_coeff, overlap1_mo, &
431 -1.0_dp, dcdr_env%dCR_prime(ispin))
432 CALL cp_fm_release(overlap1_mo)
433
434 DO alpha = 1, 3
435 IF (.NOT. dcdr_env%z_matrix_method) THEN
436
437 ! FIRST CONTRIBUTION: dCR * moments * mo
438 CALL cp_fm_set_all(tmp_fm_like_mos, 0._dp)
439 CALL dbcsr_desymmetrize(dcdr_env%matrix_s1(1)%matrix, dcdr_env%matrix_nosym_temp(1)%matrix)
440 CALL dbcsr_desymmetrize(dcdr_env%moments(alpha)%matrix, dcdr_env%matrix_nosym_temp(2)%matrix)
441 CALL dbcsr_add(dcdr_env%matrix_nosym_temp(1)%matrix, dcdr_env%matrix_nosym_temp(2)%matrix, &
442 -dcdr_env%ref_point(alpha), 1._dp)
443
444 CALL cp_dbcsr_sm_fm_multiply(dcdr_env%matrix_nosym_temp(1)%matrix, dcdr_env%dCR_prime(ispin), &
445 tmp_fm_like_mos, ncol=nmo)
446 CALL cp_fm_trace(mo_coeff, tmp_fm_like_mos, apt_coeff_derivative)
447
448 apt_coeff_derivative = (-2._dp)*f_spin*apt_coeff_derivative
449 apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) &
450 = apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) + apt_coeff_derivative
451 ELSE
452 CALL get_qs_env(qs_env=qs_env, polar_env=polar_env)
453 CALL get_polar_env(polar_env=polar_env, psi1_dberry=psi1_dberry, &
454 dberry_psi0=dberry_psi0)
455
456 ! Note that here dcdr_env%dCR_prime contains only occ-occ block contribution,
457 ! dcdr_env%dCR(ispin) is zero because we didn't run response calculation for dcdR.
458
459 CALL cp_fm_trace(dberry_psi0(alpha, ispin), &
460 dcdr_env%dCR_prime(ispin), &
461 temp1)
462
463 CALL cp_fm_trace(dcdr_env%matrix_m_alpha(dcdr_env%beta, ispin), &
464 psi1_dberry(alpha, ispin), &
465 temp2)
466
467 apt_coeff_derivative = temp1 - temp2
468
469 ! !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
470 ! - apt_coeff_derivative , here the trace is negative to compensate the
471 ! -ve sign in APTs= - 2 Z. M_alpha
472
473 apt_coeff_derivative = (-2._dp)*f_spin*apt_coeff_derivative
474 apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) &
475 = apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) + apt_coeff_derivative
476 END IF
477
478 ! SECOND CONTRIBUTION: We assemble all combinations of r_i, d(chi)/d(idir)
479 ! difdip contains derivatives with respect to atom dcdr_env%lambda
480 ! difdip(alpha, beta): < a | r_alpha | db/dR_beta >
481 ! Multiply by the MO coefficients
482 CALL cp_fm_set_all(tmp_fm_like_mos, 0.0_dp)
483 CALL cp_dbcsr_sm_fm_multiply(dcdr_env%matrix_difdip(alpha, dcdr_env%beta)%matrix, mo_coeff, &
484 tmp_fm_like_mos, ncol=nmo)
485 CALL cp_fm_trace(mo_coeff, tmp_fm_like_mos, apt_basis_derivative)
486
487 ! The negative sign compensates for the nuclear-coordinate derivative convention.
488 apt_basis_derivative = -f_spin*apt_basis_derivative
489 apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) = &
490 apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) + apt_basis_derivative
491
492 END DO ! alpha
493
494 CALL cp_fm_release(tmp_fm_like_mos)
495 END DO !ispin
496
497 ! Finally the nuclear contribution: nuclear charge * Kronecker_delta_{dcdr_env%beta,i}
498 CALL get_atomic_kind(particle_set(dcdr_env%lambda)%atomic_kind, kind_number=ikind)
499 CALL get_qs_kind(qs_kind_set(ikind), core_charge=charge, ghost=ghost)
500 IF (.NOT. ghost) THEN
501 apt_nuc(dcdr_env%beta, dcdr_env%beta, dcdr_env%lambda) = &
502 apt_nuc(dcdr_env%beta, dcdr_env%beta, dcdr_env%lambda) + charge
503 END IF
504
505 ! And deallocate all the things!
506 CALL cp_fm_release(tmp_fm_like_mos)
507 CALL cp_fm_release(overlap1_mo)
508
509 CALL timestop(handle)
510 END SUBROUTINE apt_dr
511
512! **************************************************************************************************
513!> \brief Calculate atomic polar tensor using the localized dipole operator
514!> \param qs_env ...
515!> \param dcdr_env ...
516!> \authors Edward Ditler
517!> Ravi Kumar
518!> Rangsiman Ketkaew
519! **************************************************************************************************
520 SUBROUTINE apt_dr_localization(qs_env, dcdr_env)
521 TYPE(qs_environment_type), POINTER :: qs_env
522 TYPE(dcdr_env_type) :: dcdr_env
523
524 CHARACTER(LEN=*), PARAMETER :: routinen = 'apt_dR_localization'
525
526 INTEGER :: alpha, handle, i, icenter, ikind, ispin, &
527 map_atom, map_molecule, &
528 max_nbr_center, nao, natom, nmo, &
529 nsubset
530 INTEGER, ALLOCATABLE, DIMENSION(:) :: mapping_atom_molecule
531 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: mapping_wannier_atom
532 LOGICAL :: ghost
533 REAL(dp) :: apt_basis_derivative, &
534 apt_coeff_derivative, charge, f_spin, &
535 smallest_r, this_factor, tmp_aptcontr, &
536 tmp_r
537 REAL(dp), ALLOCATABLE, DIMENSION(:) :: diagonal_elements, diagonal_elements2
538 REAL(dp), DIMENSION(3) :: distance, r_shifted
539 REAL(dp), DIMENSION(:, :, :), POINTER :: apt_el, apt_nuc
540 REAL(dp), DIMENSION(:, :, :, :), POINTER :: apt_center, apt_subset
541 TYPE(cell_type), POINTER :: cell
542 TYPE(cp_2d_r_p_type), DIMENSION(:), POINTER :: centers_set
543 TYPE(cp_fm_type), DIMENSION(:, :), POINTER :: dberry_psi0, psi1_dberry
544 TYPE(cp_fm_type), POINTER :: mo_coeff, overlap1_mo, tmp_fm, &
545 tmp_fm_like_mos, tmp_fm_momo, &
546 tmp_fm_momo2
547 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
548 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
549 TYPE(polar_env_type), POINTER :: polar_env
550 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
551
552 CALL timeset(routinen, handle)
553
554 NULLIFY (qs_kind_set, particle_set, molecule_set, cell)
555
556 CALL get_qs_env(qs_env=qs_env, &
557 qs_kind_set=qs_kind_set, &
558 particle_set=particle_set, &
559 molecule_set=molecule_set, &
560 cell=cell)
561
562 nsubset = SIZE(molecule_set)
563 natom = SIZE(particle_set)
564 apt_el => dcdr_env%apt_el_dcdr
565 apt_nuc => dcdr_env%apt_nuc_dcdr
566 apt_subset => dcdr_env%apt_el_dcdr_per_subset
567 apt_center => dcdr_env%apt_el_dcdr_per_center
568
569 ! Map wannier functions to atoms
570 IF (dcdr_env%nspins == 1) THEN
571 max_nbr_center = dcdr_env%nbr_center(1)
572 ELSE
573 max_nbr_center = max(dcdr_env%nbr_center(1), dcdr_env%nbr_center(2))
574 END IF
575 ALLOCATE (mapping_wannier_atom(max_nbr_center, dcdr_env%nspins))
576 ALLOCATE (mapping_atom_molecule(natom))
577 centers_set => dcdr_env%centers_set
578
579 DO ispin = 1, dcdr_env%nspins
580 DO icenter = 1, dcdr_env%nbr_center(ispin)
581 ! For every center we check which atom is closest
582 CALL shift_wannier_into_cell(r=centers_set(ispin)%array(1:3, icenter), &
583 cell=cell, &
584 r_shifted=r_shifted)
585
586 smallest_r = huge(0._dp)
587 DO i = 1, natom
588 distance = pbc(r_shifted, particle_set(i)%r(1:3), cell)
589 tmp_r = sum(distance**2)
590 IF (tmp_r < smallest_r) THEN
591 mapping_wannier_atom(icenter, ispin) = i
592 smallest_r = tmp_r
593 END IF
594 END DO
595 END DO
596
597 ! Map atoms to molecules
598 CALL molecule_of_atom(molecule_set, atom_to_mol=mapping_atom_molecule)
599 IF (dcdr_env%lambda == 1 .AND. dcdr_env%beta == 1) THEN
600 DO icenter = 1, dcdr_env%nbr_center(ispin)
601 map_atom = mapping_wannier_atom(icenter, ispin)
602 map_molecule = mapping_atom_molecule(map_atom)
603 END DO
604 END IF
605 END DO !ispin
606
607 nao = dcdr_env%nao
608 f_spin = 2._dp/dcdr_env%nspins
609
610 DO ispin = 1, dcdr_env%nspins
611 ! Compute S^(1,R)_(ij)
612
613 ALLOCATE (tmp_fm_like_mos)
614 ALLOCATE (overlap1_mo)
615 CALL cp_fm_create(tmp_fm_like_mos, dcdr_env%likemos_fm_struct(ispin)%struct)
616 CALL cp_fm_create(overlap1_mo, dcdr_env%momo_fm_struct(ispin)%struct)
617 nmo = dcdr_env%nmo(ispin)
618 mo_coeff => dcdr_env%mo_coeff(ispin)
619 CALL cp_fm_set_all(tmp_fm_like_mos, 0.0_dp)
620 CALL cp_fm_scale_and_add(0._dp, dcdr_env%dCR_prime(ispin), 1._dp, dcdr_env%dCR(ispin))
621 CALL cp_dbcsr_sm_fm_multiply(dcdr_env%matrix_s1(dcdr_env%beta + 1)%matrix, mo_coeff, &
622 tmp_fm_like_mos, ncol=nmo)
623 CALL parallel_gemm("T", "N", nmo, nmo, nao, &
624 1.0_dp, mo_coeff, tmp_fm_like_mos, &
625 0.0_dp, overlap1_mo)
626
627 ! C^1 <- -dCR - 0.5 * mo_coeff @ S1_ij
628 ! We get the negative of the coefficients out of the linres solver
629 ! And apply the constant correction due to the overlap derivative.
630 CALL parallel_gemm("N", "N", nao, nmo, nmo, &
631 -0.5_dp, mo_coeff, overlap1_mo, &
632 -1.0_dp, dcdr_env%dCR_prime(ispin))
633 CALL cp_fm_release(overlap1_mo)
634
635 ALLOCATE (diagonal_elements(nmo))
636 ALLOCATE (diagonal_elements2(nmo))
637
638 ! Allocate temporary matrices
639 ALLOCATE (tmp_fm)
640 ALLOCATE (tmp_fm_momo)
641 ALLOCATE (tmp_fm_momo2)
642 CALL cp_fm_create(tmp_fm, dcdr_env%likemos_fm_struct(ispin)%struct)
643 CALL cp_fm_create(tmp_fm_momo, dcdr_env%momo_fm_struct(ispin)%struct)
644 CALL cp_fm_create(tmp_fm_momo2, dcdr_env%momo_fm_struct(ispin)%struct)
645
646 ! FIRST CONTRIBUTION: dCR * moments * mo
647 this_factor = -2._dp*f_spin
648 DO alpha = 1, 3
649 IF (.NOT. dcdr_env%z_matrix_method) THEN
650
651 DO icenter = 1, dcdr_env%nbr_center(ispin)
652 CALL dbcsr_set(dcdr_env%moments(alpha)%matrix, 0.0_dp)
653 CALL build_local_moment_matrix(qs_env, dcdr_env%moments, 1, &
654 ref_point=centers_set(ispin)%array(1:3, icenter))
655 CALL multiply_localization(ao_matrix=dcdr_env%moments(alpha)%matrix, &
656 mo_coeff=dcdr_env%dCR_prime(ispin), work=tmp_fm, nmo=nmo, &
657 icenter=icenter, &
658 res=tmp_fm_like_mos)
659 END DO
660
661 CALL parallel_gemm("T", "N", nmo, nmo, nao, &
662 1.0_dp, mo_coeff, tmp_fm_like_mos, &
663 0.0_dp, tmp_fm_momo)
664 CALL cp_fm_get_diag(tmp_fm_momo, diagonal_elements)
665
666 ELSE
667 CALL get_qs_env(qs_env=qs_env, polar_env=polar_env)
668 CALL get_polar_env(polar_env=polar_env, psi1_dberry=psi1_dberry, &
669 dberry_psi0=dberry_psi0)
670
671 CALL parallel_gemm("T", "N", nmo, nmo, nao, &
672 1.0_dp, dcdr_env%dCR_prime(ispin), dberry_psi0(alpha, ispin), &
673 0.0_dp, tmp_fm_momo)
674 CALL cp_fm_get_diag(tmp_fm_momo, diagonal_elements)
675
676 CALL parallel_gemm("T", "N", nmo, nmo, nao, &
677 1.0_dp, dcdr_env%matrix_m_alpha(dcdr_env%beta, ispin), &
678 psi1_dberry(alpha, ispin), 0.0_dp, tmp_fm_momo2)
679 CALL cp_fm_get_diag(tmp_fm_momo2, diagonal_elements2)
680
681 diagonal_elements(:) = diagonal_elements(:) - diagonal_elements2(:)
682 END IF
683
684 DO icenter = 1, dcdr_env%nbr_center(ispin)
685 map_atom = mapping_wannier_atom(icenter, ispin)
686 map_molecule = mapping_atom_molecule(map_atom)
687 tmp_aptcontr = this_factor*diagonal_elements(icenter)
688
689 apt_subset(dcdr_env%beta, alpha, dcdr_env%lambda, map_molecule) &
690 = apt_subset(dcdr_env%beta, alpha, dcdr_env%lambda, map_molecule) + tmp_aptcontr
691
692 apt_center(dcdr_env%beta, alpha, dcdr_env%lambda, icenter) &
693 = apt_center(dcdr_env%beta, alpha, dcdr_env%lambda, icenter) + tmp_aptcontr
694 END DO
695
696 apt_coeff_derivative = this_factor*sum(diagonal_elements)
697 apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) &
698 = apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) + apt_coeff_derivative
699 END DO
700
701 ! SECOND CONTRIBUTION: We assemble all combinations of r_i, dphi/d(idir)
702 ! build part with AOs differentiated with respect to nuclear coordinates
703 ! difdip contains derivatives with respect to atom dcdr_env%lambda
704 ! difdip(alpha, beta): < a | r_alpha | d b/dR_beta >
705 this_factor = -f_spin
706 DO alpha = 1, 3
707 DO icenter = 1, dcdr_env%nbr_center(ispin)
708 ! Build the AO matrix with the right wannier center as reference point
709 CALL dbcsr_set(dcdr_env%matrix_difdip(1, dcdr_env%beta)%matrix, 0._dp)
710 CALL dbcsr_set(dcdr_env%matrix_difdip(2, dcdr_env%beta)%matrix, 0._dp)
711 CALL dbcsr_set(dcdr_env%matrix_difdip(3, dcdr_env%beta)%matrix, 0._dp)
712 CALL build_local_moments_der_matrix(qs_env, dcdr_env%matrix_difdip, 1, 0, &
713 ref_point=centers_set(ispin)%array(1:3, icenter), &
714 ordered=.true., deltar=dcdr_env%delta_basis_function, &
715 neighbor_image=.false.)
716 CALL multiply_localization(ao_matrix=dcdr_env%matrix_difdip(alpha, dcdr_env%beta)%matrix, &
717 mo_coeff=mo_coeff, work=tmp_fm, nmo=nmo, &
718 icenter=icenter, &
719 res=tmp_fm_like_mos)
720 END DO ! icenter
721
722 CALL parallel_gemm("T", "N", nmo, nmo, nao, &
723 1.0_dp, mo_coeff, tmp_fm_like_mos, &
724 0.0_dp, tmp_fm_momo)
725 CALL cp_fm_get_diag(tmp_fm_momo, diagonal_elements)
726
727 DO icenter = 1, dcdr_env%nbr_center(ispin)
728 map_atom = mapping_wannier_atom(icenter, ispin)
729 map_molecule = mapping_atom_molecule(map_atom)
730 tmp_aptcontr = this_factor*diagonal_elements(icenter)
731
732 apt_subset(dcdr_env%beta, alpha, dcdr_env%lambda, map_molecule) &
733 = apt_subset(dcdr_env%beta, alpha, dcdr_env%lambda, map_molecule) + tmp_aptcontr
734
735 apt_center(dcdr_env%beta, alpha, dcdr_env%lambda, icenter) &
736 = apt_center(dcdr_env%beta, alpha, dcdr_env%lambda, icenter) + tmp_aptcontr
737 END DO
738
739 ! The negative sign compensates for the nuclear-coordinate derivative convention.
740 apt_basis_derivative = this_factor*sum(diagonal_elements)
741
742 apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) &
743 = apt_el(dcdr_env%beta, alpha, dcdr_env%lambda) + apt_basis_derivative
744
745 END DO ! alpha
746 DEALLOCATE (diagonal_elements)
747 DEALLOCATE (diagonal_elements2)
748
749 CALL cp_fm_release(tmp_fm)
750 CALL cp_fm_release(tmp_fm_like_mos)
751 CALL cp_fm_release(tmp_fm_momo)
752 CALL cp_fm_release(tmp_fm_momo2)
753 DEALLOCATE (overlap1_mo)
754 DEALLOCATE (tmp_fm)
755 DEALLOCATE (tmp_fm_like_mos)
756 DEALLOCATE (tmp_fm_momo)
757 DEALLOCATE (tmp_fm_momo2)
758 END DO !ispin
759
760 ! Finally the nuclear contribution: nuclear charge * Kronecker_delta_{dcdr_env%beta,i}
761 CALL get_atomic_kind(particle_set(dcdr_env%lambda)%atomic_kind, kind_number=ikind)
762 CALL get_qs_kind(qs_kind_set(ikind), core_charge=charge, ghost=ghost)
763 IF (.NOT. ghost) THEN
764 apt_nuc(dcdr_env%beta, dcdr_env%beta, dcdr_env%lambda) = &
765 apt_nuc(dcdr_env%beta, dcdr_env%beta, dcdr_env%lambda) + charge
766
767 map_molecule = mapping_atom_molecule(dcdr_env%lambda)
768 apt_subset(dcdr_env%beta, dcdr_env%beta, dcdr_env%lambda, map_molecule) &
769 = apt_subset(dcdr_env%beta, dcdr_env%beta, dcdr_env%lambda, map_molecule) + charge
770 END IF
771
772 ! And deallocate all the things!
773
774 CALL timestop(handle)
775 END SUBROUTINE apt_dr_localization
776
777END MODULE qs_dcdr
cell_types::pbc
Definition cell_types.F:104
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition cp_dbcsr_operations.F:95
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:103
cp_fm_basic_linalg::cp_fm_trace
Definition cp_fm_basic_linalg.F:78
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:88
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:83
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:39
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_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_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_desymmetrize
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
Definition cp_dbcsr_api.F:525
cp_dbcsr_api::dbcsr_copy
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
Definition cp_dbcsr_api.F:371
cp_dbcsr_api::dbcsr_set
subroutine, public dbcsr_set(matrix, alpha)
...
Definition cp_dbcsr_api.F:1204
cp_dbcsr_api::dbcsr_add
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
Definition cp_dbcsr_api.F:254
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:799
cp_fm_basic_linalg
Basic linear algebra operations for full matrices.
Definition cp_fm_basic_linalg.F:14
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:170
cp_fm_basic_linalg::cp_fm_scale
subroutine, public cp_fm_scale(alpha, matrix_a)
scales a matrix matrix_a = alpha * matrix_b
Definition cp_fm_basic_linalg.F:1502
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_diag
subroutine, public cp_fm_get_diag(matrix, diag)
returns the diagonal elements of a fm
Definition cp_fm_types.F:564
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:534
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
creates a new full matrix with the given structure
Definition cp_fm_types.F:165
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:1073
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:735
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
molecule_types
Define the data structure for the molecule information.
Definition molecule_types.F:16
molecule_types::molecule_of_atom
subroutine, public molecule_of_atom(molecule_set, atom_to_mol)
finds for each atom the molecule it belongs to
Definition molecule_types.F:478
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
qs_dcdr_ao
Calculate the derivatives of the MO coefficients wrt nuclear coordinates.
Definition qs_dcdr_ao.F:13
qs_dcdr_ao::core_dr
subroutine, public core_dr(qs_env, dcdr_env)
Core Hamiltonian contributions to the operator (the pseudopotentials).
Definition qs_dcdr_ao.F:381
qs_dcdr_ao::apply_op_constant_term
subroutine, public apply_op_constant_term(qs_env, dcdr_env, overlap1)
Build the perturbed density matrix correction depending on the overlap derivative.
Definition qs_dcdr_ao.F:108
qs_dcdr_ao::hr_mult_by_delta_1d
subroutine, public hr_mult_by_delta_1d(matrix, qs_kind_set, basis_type, sab_nl, lambda, direction_or)
Enforce that one of the basis functions in < a | O | b > is centered on atom lambda.
Definition qs_dcdr_ao.F:631
qs_dcdr_ao::vhxc_r_perturbed_basis_functions
subroutine, public vhxc_r_perturbed_basis_functions(qs_env, dcdr_env)
The derivatives of the basis functions over which the HXC potential is integrated,...
Definition qs_dcdr_ao.F:558
qs_dcdr_ao::d_vhxc_dr
subroutine, public d_vhxc_dr(qs_env, dcdr_env)
The derivatives of the basis functions going into the HXC potential wrt nuclear positions.
Definition qs_dcdr_ao.F:421
qs_dcdr_ao::d_core_charge_density_dr
subroutine, public d_core_charge_density_dr(qs_env, dcdr_env)
Calculate the derivative of the Hartree term due to the core charge density.
Definition qs_dcdr_ao.F:309
qs_dcdr_utils
Calculate the derivatives of the MO coefficients wrt nuclear coordinates.
Definition qs_dcdr_utils.F:13
qs_dcdr_utils::multiply_localization
subroutine, public multiply_localization(ao_matrix, mo_coeff, work, nmo, icenter, res)
Multiply (ao_matrix @ mo_coeff) and store the column icenter in res.
Definition qs_dcdr_utils.F:113
qs_dcdr_utils::dcdr_read_restart
subroutine, public dcdr_read_restart(qs_env, linres_section, vec, lambda, beta, tag)
Copied from linres_read_restart.
Definition qs_dcdr_utils.F:149
qs_dcdr_utils::shift_wannier_into_cell
subroutine, public shift_wannier_into_cell(r, cell, r_shifted)
...
Definition qs_dcdr_utils.F:528
qs_dcdr_utils::dcdr_write_restart
subroutine, public dcdr_write_restart(qs_env, linres_section, vec, lambda, beta, tag)
Copied from linres_write_restart.
Definition qs_dcdr_utils.F:300
qs_dcdr
Calculate the derivatives of the MO coefficients wrt nuclear coordinates.
Definition qs_dcdr.F:13
qs_dcdr::prepare_per_atom
subroutine, public prepare_per_atom(dcdr_env, qs_env)
Prepare the environment for a choice of lambda.
Definition qs_dcdr.F:90
qs_dcdr::apt_dr_localization
subroutine, public apt_dr_localization(qs_env, dcdr_env)
Calculate atomic polar tensor using the localized dipole operator.
Definition qs_dcdr.F:521
qs_dcdr::apt_dr
subroutine, public apt_dr(qs_env, dcdr_env)
Calculate atomic polar tensor.
Definition qs_dcdr.F:378
qs_dcdr::dcdr_response_dr
subroutine, public dcdr_response_dr(dcdr_env, p_env, qs_env)
Get the dC/dR by solving the Sternheimer equation, using the op_dR matrix.
Definition qs_dcdr.F:271
qs_dcdr::dcdr_build_op_dr
subroutine, public dcdr_build_op_dr(dcdr_env, qs_env)
Build the operator for the position perturbation.
Definition qs_dcdr.F:196
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, mimic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_vhxc, 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, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
Definition qs_environment_types.F:545
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, hund_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, proj_shell_charge, lr_atom, do_mtlr, u_j_loop, ao_coef, init_u_ramping_each_scf, reltmat, ghost, monovalent, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition qs_kind_types.F:491
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, silent)
scf loop to optimize the first order wavefunctions (psi1) given a perturbation as an operator applied...
Definition qs_linres_methods.F:218
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:1242
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, cmo_coeff)
Get the components of a MO set data structure.
Definition qs_mo_types.F:445
qs_moments
Calculates the moment integrals <a|r^m|b> and <a|r x d/dr|b>.
Definition qs_moments.F:14
qs_moments::build_local_moment_matrix
subroutine, public build_local_moment_matrix(qs_env, moments, nmoments, ref_point, ref_points, basis_type, all_images, minimum_image, neighbor_image, first_component)
...
Definition qs_moments.F:166
qs_moments::build_local_moments_der_matrix
subroutine, public build_local_moments_der_matrix(qs_env, moments_der, nmoments_der, nmoments, ref_point, moments, basis_type, minimum_image, ordered, lambda, deltar, neighbor_image)
Calculate right-hand sided derivatives of multipole moments, e. g. < a | xy d/dz | b > Optionally sto...
Definition qs_moments.F:487
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
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:60
cp_array_utils::cp_2d_r_p_type
represent a pointer to a 2d array
Definition cp_array_utils.F:109
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:114
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
molecule_types::molecule_type
Definition molecule_types.F:159
particle_types::particle_type
Definition particle_types.F:35
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:224
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:190
qs_linres_types::dcdr_env_type
Definition qs_linres_types.F:247
qs_linres_types::linres_control_type
General settings for linear response calculations.
Definition qs_linres_types.F:53
qs_linres_types::polar_env_type
Definition qs_linres_types.F:136
qs_mo_types::mo_set_type
Definition qs_mo_types.F:53
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