(git:34ef472)
hfx_admm_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 Utilities for hfx and admm methods
10 !>
11 !>
12 !> \par History
13 !> refactoring 03-2011 [MI]
14 !> Made GAPW compatible 12.2019 (A. Bussy)
15 !> \author MI
16 ! **************************************************************************************************
17 MODULE hfx_admm_utils
18  USE admm_dm_types, ONLY: admm_dm_create
19  USE admm_methods, ONLY: kpoint_calc_admm_matrices,&
20  scale_dm
21  USE admm_types, ONLY: admm_env_create,&
22  admm_gapw_r3d_rs_type,&
23  admm_type,&
24  get_admm_env,&
25  set_admm_env
26  USE atomic_kind_types, ONLY: atomic_kind_type
27  USE basis_set_container_types, ONLY: add_basis_set_to_container
28  USE basis_set_types, ONLY: copy_gto_basis_set,&
29  get_gto_basis_set,&
30  gto_basis_set_type
31  USE cell_types, ONLY: cell_type,&
32  plane_distance
33  USE cp_blacs_env, ONLY: cp_blacs_env_type
34  USE cp_control_types, ONLY: admm_control_type,&
35  dft_control_type
36  USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
37  USE cp_dbcsr_operations, ONLY: cp_dbcsr_m_by_n_from_row_template,&
38  dbcsr_allocate_matrix_set
39  USE cp_fm_pool_types, ONLY: cp_fm_pool_p_type
40  USE cp_fm_struct, ONLY: cp_fm_struct_create,&
41  cp_fm_struct_release,&
42  cp_fm_struct_type
43  USE cp_fm_types, ONLY: cp_fm_create,&
44  cp_fm_get_info,&
45  cp_fm_type
46  USE cp_log_handling, ONLY: cp_get_default_logger,&
47  cp_logger_get_default_io_unit,&
48  cp_logger_type,&
49  cp_to_string
50  USE dbcsr_api, ONLY: dbcsr_add,&
51  dbcsr_copy,&
52  dbcsr_create,&
53  dbcsr_init_p,&
54  dbcsr_p_type,&
55  dbcsr_set,&
56  dbcsr_type,&
57  dbcsr_type_no_symmetry
58  USE distribution_1d_types, ONLY: distribution_1d_type
59  USE distribution_2d_types, ONLY: distribution_2d_type
60  USE external_potential_types, ONLY: copy_potential
61  USE hfx_derivatives, ONLY: derivatives_four_center
62  USE hfx_energy_potential, ONLY: integrate_four_center
63  USE hfx_pw_methods, ONLY: pw_hfx
64  USE hfx_ri, ONLY: hfx_ri_update_forces,&
65  hfx_ri_update_ks
66  USE hfx_ri_kp, ONLY: hfx_ri_update_forces_kp,&
67  hfx_ri_update_ks_kp
68  USE hfx_types, ONLY: hfx_type
69  USE input_constants, ONLY: &
70  do_admm_aux_exch_func_bee, do_admm_aux_exch_func_bee_libxc, do_admm_aux_exch_func_default, &
71  do_admm_aux_exch_func_default_libxc, do_admm_aux_exch_func_none, &
72  do_admm_aux_exch_func_opt, do_admm_aux_exch_func_opt_libxc, do_admm_aux_exch_func_pbex, &
73  do_admm_aux_exch_func_pbex_libxc, do_admm_aux_exch_func_sx_libxc, &
74  do_admm_basis_projection, do_admm_charge_constrained_projection, do_admm_purify_none, &
75  do_potential_coulomb, do_potential_id, do_potential_long, do_potential_mix_cl, &
76  do_potential_mix_cl_trunc, do_potential_short, do_potential_truncated, &
77  xc_funct_no_shortcut, xc_none
78  USE input_section_types, ONLY: section_vals_duplicate,&
79  section_vals_get,&
80  section_vals_get_subs_vals,&
81  section_vals_get_subs_vals2,&
82  section_vals_remove_values,&
83  section_vals_type,&
84  section_vals_val_get,&
85  section_vals_val_set
86  USE kinds, ONLY: dp
87  USE kpoint_methods, ONLY: kpoint_initialize_mos
88  USE kpoint_transitional, ONLY: kpoint_transitional_release,&
89  set_2d_pointer
90  USE kpoint_types, ONLY: get_kpoint_info,&
91  kpoint_type
92  USE libint_2c_3c, ONLY: cutoff_screen_factor
93  USE mathlib, ONLY: erfc_cutoff
94  USE message_passing, ONLY: mp_para_env_type
95  USE molecule_types, ONLY: molecule_type
96  USE particle_types, ONLY: particle_type
97  USE paw_proj_set_types, ONLY: get_paw_proj_set,&
98  paw_proj_set_type
99  USE pw_env_types, ONLY: pw_env_get,&
100  pw_env_type
101  USE pw_poisson_types, ONLY: pw_poisson_type
102  USE pw_pool_types, ONLY: pw_pool_type
103  USE pw_types, ONLY: pw_r3d_rs_type
104  USE qs_energy_types, ONLY: qs_energy_type
105  USE qs_environment_types, ONLY: get_qs_env,&
106  qs_environment_type,&
107  set_qs_env
108  USE qs_interactions, ONLY: init_interaction_radii
109  USE qs_kind_types, ONLY: get_qs_kind,&
110  get_qs_kind_set,&
111  init_gapw_basis_set,&
112  init_gapw_nlcc,&
113  qs_kind_type
114  USE qs_ks_types, ONLY: qs_ks_env_type
115  USE qs_local_rho_types, ONLY: local_rho_set_create
116  USE qs_matrix_pools, ONLY: mpools_get
117  USE qs_mo_types, ONLY: allocate_mo_set,&
118  get_mo_set,&
119  init_mo_set,&
120  mo_set_type
121  USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type,&
122  release_neighbor_list_sets
123  USE qs_neighbor_lists, ONLY: atom2d_build,&
124  atom2d_cleanup,&
125  build_neighbor_lists,&
126  local_atoms_type,&
127  pair_radius_setup,&
128  write_neighbor_lists
129  USE qs_oce_methods, ONLY: build_oce_matrices
130  USE qs_oce_types, ONLY: allocate_oce_set,&
131  create_oce_set
132  USE qs_overlap, ONLY: build_overlap_matrix
133  USE qs_rho_atom_methods, ONLY: init_rho_atom
134  USE qs_rho_methods, ONLY: qs_rho_rebuild
135  USE qs_rho_types, ONLY: qs_rho_create,&
136  qs_rho_get,&
137  qs_rho_type
138  USE rt_propagation_types, ONLY: rt_prop_type
139  USE task_list_methods, ONLY: generate_qs_task_list
140  USE task_list_types, ONLY: allocate_task_list,&
141  deallocate_task_list
142  USE virial_types, ONLY: virial_type
143  USE xc_adiabatic_utils, ONLY: rescale_xc_potential
144 #include "./base/base_uses.f90"
145 
146  IMPLICIT NONE
147 
148  PRIVATE
149 
150  ! *** Public subroutines ***
151  PUBLIC :: hfx_ks_matrix, hfx_admm_init, aux_admm_init, create_admm_xc_section, &
152  tddft_hfx_matrix, hfx_ks_matrix_kp
153 
154  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'hfx_admm_utils'
155 
156 CONTAINS
157 
158 ! **************************************************************************************************
159 !> \brief ...
160 !> \param qs_env ...
161 !> \param calculate_forces ...
162 ! **************************************************************************************************
163  SUBROUTINE hfx_admm_init(qs_env, calculate_forces)
164 
165  TYPE(qs_environment_type), POINTER :: qs_env
166  LOGICAL, INTENT(IN), OPTIONAL :: calculate_forces
167 
168  CHARACTER(LEN=*), PARAMETER :: routinen = 'hfx_admm_init'
169 
170  INTEGER :: handle, ispin, n_rep_hf, nao_aux_fit, &
171  natoms, nelectron, nmo
172  LOGICAL :: calc_forces, do_kpoints, &
173  s_mstruct_changed, use_virial
174  REAL(dp) :: maxocc
175  TYPE(admm_type), POINTER :: admm_env
176  TYPE(cp_blacs_env_type), POINTER :: blacs_env
177  TYPE(cp_fm_struct_type), POINTER :: aux_fit_fm_struct
178  TYPE(cp_fm_type), POINTER :: mo_coeff_aux_fit
179  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_aux_fit_kp
180  TYPE(dbcsr_type), POINTER :: mo_coeff_b
181  TYPE(dft_control_type), POINTER :: dft_control
182  TYPE(mo_set_type), DIMENSION(:), POINTER :: mos, mos_aux_fit
183  TYPE(mp_para_env_type), POINTER :: para_env
184  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
185  TYPE(qs_ks_env_type), POINTER :: ks_env
186  TYPE(section_vals_type), POINTER :: hfx_sections, input, xc_section
187  TYPE(virial_type), POINTER :: virial
188 
189  CALL timeset(routinen, handle)
190 
191  NULLIFY (admm_env, hfx_sections, mos, mos_aux_fit, para_env, virial, &
192  mo_coeff_aux_fit, xc_section, ks_env, dft_control, input, &
193  qs_kind_set, mo_coeff_b, aux_fit_fm_struct, blacs_env)
194 
195  CALL get_qs_env(qs_env, &
196  mos=mos, &
197  admm_env=admm_env, &
198  para_env=para_env, &
199  blacs_env=blacs_env, &
200  s_mstruct_changed=s_mstruct_changed, &
201  ks_env=ks_env, &
202  dft_control=dft_control, &
203  input=input, &
204  virial=virial, &
205  do_kpoints=do_kpoints)
206 
207  calc_forces = .false.
208  IF (PRESENT(calculate_forces)) calc_forces = .true.
209 
210  hfx_sections => section_vals_get_subs_vals(input, "DFT%XC%HF")
211 
212  CALL section_vals_get(hfx_sections, n_repetition=n_rep_hf)
213  IF (n_rep_hf > 1) &
214  cpabort("ADMM can handle only one HF section.")
215 
216  IF (.NOT. ASSOCIATED(admm_env)) THEN
217  ! setup admm environment
218  CALL get_qs_env(qs_env, input=input, natom=natoms, qs_kind_set=qs_kind_set)
219  CALL get_qs_kind_set(qs_kind_set, nsgf=nao_aux_fit, basis_type="AUX_FIT")
220  CALL admm_env_create(admm_env, dft_control%admm_control, mos, para_env, natoms, nao_aux_fit)
221  CALL set_qs_env(qs_env, admm_env=admm_env)
222  xc_section => section_vals_get_subs_vals(input, "DFT%XC")
223  CALL create_admm_xc_section(x_data=qs_env%x_data, xc_section=xc_section, &
224  admm_env=admm_env)
225 
226  ! Initialize the GAPW data types
227  IF (dft_control%qs_control%gapw .OR. dft_control%qs_control%gapw_xc) &
228  CALL init_admm_gapw(qs_env)
229 
230  ! ADMM neighbor lists and overlap matrices
231  CALL admm_init_hamiltonians(admm_env, qs_env, "AUX_FIT")
232 
233  !The aux_fit task list and densities
234  ALLOCATE (admm_env%rho_aux_fit)
235  CALL qs_rho_create(admm_env%rho_aux_fit)
236  ALLOCATE (admm_env%rho_aux_fit_buffer)
237  CALL qs_rho_create(admm_env%rho_aux_fit_buffer)
238  CALL admm_update_s_mstruct(admm_env, qs_env, "AUX_FIT")
239  IF (admm_env%do_gapw) CALL update_admm_gapw(qs_env)
240 
241  !The ADMM KS matrices
242  CALL admm_alloc_ks_matrices(admm_env, qs_env)
243 
244  !The aux_fit MOs and derivatives
245  ALLOCATE (mos_aux_fit(dft_control%nspins))
246  DO ispin = 1, dft_control%nspins
247  CALL get_mo_set(mo_set=mos(ispin), nmo=nmo, nelectron=nelectron, maxocc=maxocc)
248  CALL allocate_mo_set(mo_set=mos_aux_fit(ispin), &
249  nao=nao_aux_fit, &
250  nmo=nmo, &
251  nelectron=nelectron, &
252  n_el_f=real(nelectron, dp), &
253  maxocc=maxocc, &
254  flexible_electron_count=dft_control%relax_multiplicity)
255  END DO
256  admm_env%mos_aux_fit => mos_aux_fit
257 
258  DO ispin = 1, dft_control%nspins
259  CALL get_mo_set(mo_set=mos(ispin), nmo=nmo)
260  CALL cp_fm_struct_create(aux_fit_fm_struct, context=blacs_env, para_env=para_env, &
261  nrow_global=nao_aux_fit, ncol_global=nmo)
262  CALL get_mo_set(mos_aux_fit(ispin), mo_coeff=mo_coeff_aux_fit, mo_coeff_b=mo_coeff_b)
263  IF (.NOT. ASSOCIATED(mo_coeff_aux_fit)) THEN
264  CALL init_mo_set(mos_aux_fit(ispin), fm_struct=aux_fit_fm_struct, &
265  name="qs_env%mo_aux_fit"//trim(adjustl(cp_to_string(ispin))))
266  END IF
267  CALL cp_fm_struct_release(aux_fit_fm_struct)
268 
269  IF (.NOT. ASSOCIATED(mo_coeff_b)) THEN
270  CALL cp_fm_get_info(mos_aux_fit(ispin)%mo_coeff, ncol_global=nmo)
271  CALL dbcsr_init_p(mos_aux_fit(ispin)%mo_coeff_b)
272  CALL get_admm_env(admm_env, matrix_s_aux_fit_kp=matrix_s_aux_fit_kp)
273  CALL cp_dbcsr_m_by_n_from_row_template(mos_aux_fit(ispin)%mo_coeff_b, &
274  template=matrix_s_aux_fit_kp(1, 1)%matrix, &
275  n=nmo, sym=dbcsr_type_no_symmetry)
276  END IF
277  END DO
278 
279  IF (qs_env%requires_mo_derivs) THEN
280  ALLOCATE (admm_env%mo_derivs_aux_fit(dft_control%nspins))
281  DO ispin = 1, dft_control%nspins
282  CALL get_mo_set(admm_env%mos_aux_fit(ispin), mo_coeff=mo_coeff_aux_fit)
283  CALL cp_fm_create(admm_env%mo_derivs_aux_fit(ispin), mo_coeff_aux_fit%matrix_struct)
284  END DO
285  END IF
286 
287  IF (do_kpoints) THEN
288  block
289  TYPE(kpoint_type), POINTER :: kpoints
290  TYPE(mo_set_type), DIMENSION(:, :), POINTER :: mos_aux_fit_kp
291  TYPE(cp_fm_pool_p_type), DIMENSION(:), POINTER :: ao_mo_fm_pools_aux_fit
292  TYPE(cp_fm_struct_type), POINTER :: ao_ao_fm_struct
293  INTEGER :: ic, ik, ikk, is
294  INTEGER, PARAMETER :: nwork1 = 4
295  LOGICAL :: use_real_wfn
296 
297  NULLIFY (ao_mo_fm_pools_aux_fit, mos_aux_fit_kp)
298 
299  CALL get_qs_env(qs_env=qs_env, kpoints=kpoints)
300  CALL get_kpoint_info(kpoints, use_real_wfn=use_real_wfn)
301 
302  !Test combinations of input values. So far, only ADMM2 is availavle
303  IF (.NOT. admm_env%purification_method == do_admm_purify_none) &
304  cpabort("Only ADMM_PURIFICATION_METHOD NONE implemeted for ADMM K-points")
305  IF (.NOT. (dft_control%admm_control%method == do_admm_basis_projection &
306  .OR. dft_control%admm_control%method == do_admm_charge_constrained_projection)) &
307  cpabort("Only BASIS_PROJECTION and CHARGE_CONSTRAINED_PROJECTION implemented for KP")
308  IF (admm_env%do_admms .OR. admm_env%do_admmp .OR. admm_env%do_admmq) THEN
309  IF (use_real_wfn) cpabort("Only KP-HFX ADMM2 is implemented with REAL wavefunctions")
310  END IF
311 
312  CALL kpoint_initialize_mos(kpoints, admm_env%mos_aux_fit, for_aux_fit=.true.)
313 
314  CALL mpools_get(kpoints%mpools_aux_fit, ao_mo_fm_pools=ao_mo_fm_pools_aux_fit)
315  DO ik = 1, SIZE(kpoints%kp_aux_env)
316  mos_aux_fit_kp => kpoints%kp_aux_env(ik)%kpoint_env%mos
317  ikk = kpoints%kp_range(1) + ik - 1
318  DO ispin = 1, SIZE(mos_aux_fit_kp, 2)
319  DO ic = 1, SIZE(mos_aux_fit_kp, 1)
320  CALL get_mo_set(mos_aux_fit_kp(ic, ispin), mo_coeff=mo_coeff_aux_fit, mo_coeff_b=mo_coeff_b)
321 
322  ! no sparse matrix representation of kpoint MO vectors
323  cpassert(.NOT. ASSOCIATED(mo_coeff_b))
324 
325  IF (.NOT. ASSOCIATED(mo_coeff_aux_fit)) THEN
326  CALL init_mo_set(mos_aux_fit_kp(ic, ispin), &
327  fm_pool=ao_mo_fm_pools_aux_fit(ispin)%pool, &
328  name="kpoints_"//trim(adjustl(cp_to_string(ikk)))// &
329  "%mo_aux_fit"//trim(adjustl(cp_to_string(ispin))))
330  END IF
331  END DO
332  END DO
333  END DO
334 
335  ALLOCATE (admm_env%scf_work_aux_fit(nwork1))
336 
337  ! create an ao_ao parallel matrix structure
338  CALL cp_fm_struct_create(ao_ao_fm_struct, context=blacs_env, para_env=para_env, &
339  nrow_global=nao_aux_fit, &
340  ncol_global=nao_aux_fit)
341 
342  DO is = 1, nwork1
343  CALL cp_fm_create(admm_env%scf_work_aux_fit(is), &
344  matrix_struct=ao_ao_fm_struct, &
345  name="SCF-WORK_MATRIX-AUX-"//trim(adjustl(cp_to_string(is))))
346  END DO
347  CALL cp_fm_struct_release(ao_ao_fm_struct)
348 
349  ! Create and populate the internal ADMM overlap matrices at each KP
350  CALL kpoint_calc_admm_matrices(qs_env, calc_forces)
351 
352  END block
353  END IF
354 
355  ELSE IF (s_mstruct_changed) THEN
356  CALL admm_init_hamiltonians(admm_env, qs_env, "AUX_FIT")
357  CALL admm_update_s_mstruct(admm_env, qs_env, "AUX_FIT")
358  CALL admm_alloc_ks_matrices(admm_env, qs_env)
359  IF (admm_env%do_gapw) CALL update_admm_gapw(qs_env)
360  IF (do_kpoints) CALL kpoint_calc_admm_matrices(qs_env, calc_forces)
361  END IF
362 
363  IF (admm_env%do_gapw .AND. dft_control%do_admm_dm) THEN
364  cpabort("GAPW ADMM not implemented for MCWEENY or NONE_DM purification.")
365  END IF
366 
367  !ADMMS and ADMMP stress tensors only available for close-shell systesms, because virial cannot
368  !be scaled by gsi spin component wise
369  use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
370  IF (use_virial .AND. admm_env%do_admms .AND. dft_control%nspins == 2) THEN
371  cpabort("ADMMS stress tensor is only available for closed-shell systems")
372  END IF
373  IF (use_virial .AND. admm_env%do_admmp .AND. dft_control%nspins == 2) THEN
374  cpabort("ADMMP stress tensor is only available for closed-shell systems")
375  END IF
376 
377  IF (dft_control%do_admm_dm .AND. .NOT. ASSOCIATED(admm_env%admm_dm)) THEN
378  CALL admm_dm_create(admm_env%admm_dm, dft_control%admm_control, nspins=dft_control%nspins, natoms=natoms)
379  END IF
380 
381  CALL timestop(handle)
382 
383  END SUBROUTINE hfx_admm_init
384 
385 ! **************************************************************************************************
386 !> \brief Minimal setup routine for admm_env
387 !> No forces
388 !> No k-points
389 !> No DFT correction terms
390 !> \param qs_env ...
391 !> \param mos ...
392 !> \param admm_env ...
393 !> \param admm_control ...
394 !> \param basis_type ...
395 ! **************************************************************************************************
396  SUBROUTINE aux_admm_init(qs_env, mos, admm_env, admm_control, basis_type)
397 
398  TYPE(qs_environment_type), POINTER :: qs_env
399  TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
400  TYPE(admm_type), POINTER :: admm_env
401  TYPE(admm_control_type), POINTER :: admm_control
402  CHARACTER(LEN=*) :: basis_type
403 
404  CHARACTER(LEN=*), PARAMETER :: routinen = 'aux_admm_init'
405 
406  INTEGER :: handle, ispin, nao_aux_fit, natoms, &
407  nelectron, nmo
408  LOGICAL :: do_kpoints
409  REAL(dp) :: maxocc
410  TYPE(cp_blacs_env_type), POINTER :: blacs_env
411  TYPE(cp_fm_struct_type), POINTER :: aux_fit_fm_struct
412  TYPE(cp_fm_type), POINTER :: mo_coeff_aux_fit
413  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_aux_fit_kp
414  TYPE(dbcsr_type), POINTER :: mo_coeff_b
415  TYPE(dft_control_type), POINTER :: dft_control
416  TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_aux_fit
417  TYPE(mp_para_env_type), POINTER :: para_env
418  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
419  TYPE(qs_ks_env_type), POINTER :: ks_env
420 
421  CALL timeset(routinen, handle)
422 
423  cpassert(.NOT. ASSOCIATED(admm_env))
424 
425  CALL get_qs_env(qs_env, &
426  para_env=para_env, &
427  blacs_env=blacs_env, &
428  ks_env=ks_env, &
429  dft_control=dft_control, &
430  do_kpoints=do_kpoints)
431 
432  cpassert(.NOT. do_kpoints)
433  IF (dft_control%qs_control%gapw .OR. dft_control%qs_control%gapw_xc) THEN
434  cpabort("AUX ADMM not possible with GAPW")
435  END IF
436 
437  ! setup admm environment
438  CALL get_qs_env(qs_env, natom=natoms, qs_kind_set=qs_kind_set)
439  CALL get_qs_kind_set(qs_kind_set, nsgf=nao_aux_fit, basis_type=basis_type)
440  !
441  CALL admm_env_create(admm_env, admm_control, mos, para_env, natoms, nao_aux_fit)
442  ! no XC correction used
443  NULLIFY (admm_env%xc_section_aux, admm_env%xc_section_primary)
444  ! ADMM neighbor lists and overlap matrices
445  CALL admm_init_hamiltonians(admm_env, qs_env, basis_type)
446  NULLIFY (admm_env%rho_aux_fit, admm_env%rho_aux_fit_buffer)
447  !The ADMM KS matrices
448  CALL admm_alloc_ks_matrices(admm_env, qs_env)
449  !The aux_fit MOs and derivatives
450  ALLOCATE (mos_aux_fit(dft_control%nspins))
451  DO ispin = 1, dft_control%nspins
452  CALL get_mo_set(mo_set=mos(ispin), nmo=nmo, nelectron=nelectron, maxocc=maxocc)
453  CALL allocate_mo_set(mo_set=mos_aux_fit(ispin), nao=nao_aux_fit, nmo=nmo, &
454  nelectron=nelectron, n_el_f=real(nelectron, dp), &
455  maxocc=maxocc, flexible_electron_count=0.0_dp)
456  END DO
457  admm_env%mos_aux_fit => mos_aux_fit
458 
459  DO ispin = 1, dft_control%nspins
460  CALL get_mo_set(mo_set=mos(ispin), nmo=nmo)
461  CALL cp_fm_struct_create(aux_fit_fm_struct, context=blacs_env, para_env=para_env, &
462  nrow_global=nao_aux_fit, ncol_global=nmo)
463  CALL get_mo_set(mos_aux_fit(ispin), mo_coeff=mo_coeff_aux_fit, mo_coeff_b=mo_coeff_b)
464  IF (.NOT. ASSOCIATED(mo_coeff_aux_fit)) THEN
465  CALL init_mo_set(mos_aux_fit(ispin), fm_struct=aux_fit_fm_struct, &
466  name="mo_aux_fit"//trim(adjustl(cp_to_string(ispin))))
467  END IF
468  CALL cp_fm_struct_release(aux_fit_fm_struct)
469 
470  IF (.NOT. ASSOCIATED(mo_coeff_b)) THEN
471  CALL cp_fm_get_info(mos_aux_fit(ispin)%mo_coeff, ncol_global=nmo)
472  CALL dbcsr_init_p(mos_aux_fit(ispin)%mo_coeff_b)
473  CALL get_admm_env(admm_env, matrix_s_aux_fit_kp=matrix_s_aux_fit_kp)
474  CALL cp_dbcsr_m_by_n_from_row_template(mos_aux_fit(ispin)%mo_coeff_b, &
475  template=matrix_s_aux_fit_kp(1, 1)%matrix, &
476  n=nmo, sym=dbcsr_type_no_symmetry)
477  END IF
478  END DO
479 
480  CALL timestop(handle)
481 
482  END SUBROUTINE aux_admm_init
483 
484 ! **************************************************************************************************
485 !> \brief Sets up the admm_gapw env
486 !> \param qs_env ...
487 ! **************************************************************************************************
488  SUBROUTINE init_admm_gapw(qs_env)
489 
490  TYPE(qs_environment_type), POINTER :: qs_env
491 
492  INTEGER :: ikind, nkind
493  TYPE(admm_gapw_r3d_rs_type), POINTER :: admm_gapw_env
494  TYPE(admm_type), POINTER :: admm_env
495  TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
496  TYPE(dft_control_type), POINTER :: dft_control
497  TYPE(gto_basis_set_type), POINTER :: aux_fit_basis, aux_fit_soft_basis, &
498  orb_basis, soft_basis
499  TYPE(mp_para_env_type), POINTER :: para_env
500  TYPE(qs_kind_type), DIMENSION(:), POINTER :: admm_kind_set, qs_kind_set
501  TYPE(section_vals_type), POINTER :: input
502 
503  NULLIFY (admm_kind_set, aux_fit_basis, atomic_kind_set, aux_fit_soft_basis, &
504  dft_control, input, orb_basis, para_env, qs_kind_set, soft_basis)
505 
506  CALL get_qs_env(qs_env, admm_env=admm_env, &
507  atomic_kind_set=atomic_kind_set, &
508  dft_control=dft_control, &
509  input=input, &
510  para_env=para_env, &
511  qs_kind_set=qs_kind_set)
512 
513  admm_env%do_gapw = .true.
514  ALLOCATE (admm_env%admm_gapw_env)
515  admm_gapw_env => admm_env%admm_gapw_env
516  NULLIFY (admm_gapw_env%local_rho_set)
517  NULLIFY (admm_gapw_env%admm_kind_set)
518  NULLIFY (admm_gapw_env%task_list)
519 
520  !Create a new kind set for the ADMM stuff (paw_proj soft AUX_FIT basis, etc)
521  nkind = SIZE(qs_kind_set)
522  ALLOCATE (admm_gapw_env%admm_kind_set(nkind))
523  admm_kind_set => admm_gapw_env%admm_kind_set
524 
525  !In this new kind set, we want the AUX_FIT basis to be known as ORB, such that GAPW routines work
526  DO ikind = 1, nkind
527  !copying over simple data of interest from qs_kind_set
528  admm_kind_set(ikind)%name = qs_kind_set(ikind)%name
529  admm_kind_set(ikind)%element_symbol = qs_kind_set(ikind)%element_symbol
530  admm_kind_set(ikind)%natom = qs_kind_set(ikind)%natom
531  admm_kind_set(ikind)%hard_radius = qs_kind_set(ikind)%hard_radius
532  admm_kind_set(ikind)%max_rad_local = qs_kind_set(ikind)%max_rad_local
533  admm_kind_set(ikind)%gpw_r3d_rs_type_forced = qs_kind_set(ikind)%gpw_r3d_rs_type_forced
534  admm_kind_set(ikind)%ngrid_rad = qs_kind_set(ikind)%ngrid_rad
535  admm_kind_set(ikind)%ngrid_ang = qs_kind_set(ikind)%ngrid_ang
536 
537  !copying potentials of interest from qs_kind_set
538  IF (ASSOCIATED(qs_kind_set(ikind)%all_potential)) THEN
539  CALL copy_potential(qs_kind_set(ikind)%all_potential, admm_kind_set(ikind)%all_potential)
540  END IF
541  IF (ASSOCIATED(qs_kind_set(ikind)%gth_potential)) THEN
542  CALL copy_potential(qs_kind_set(ikind)%gth_potential, admm_kind_set(ikind)%gth_potential)
543  END IF
544  IF (ASSOCIATED(qs_kind_set(ikind)%sgp_potential)) THEN
545  CALL copy_potential(qs_kind_set(ikind)%sgp_potential, admm_kind_set(ikind)%sgp_potential)
546  END IF
547 
548  NULLIFY (orb_basis)
549  CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_fit_basis, basis_type="AUX_FIT")
550  CALL copy_gto_basis_set(aux_fit_basis, orb_basis)
551  CALL add_basis_set_to_container(admm_kind_set(ikind)%basis_sets, orb_basis, "ORB")
552  END DO
553 
554  !Create the corresponding soft basis set (and projectors)
555  CALL init_gapw_basis_set(admm_kind_set, dft_control%qs_control, input, &
556  modify_qs_control=.false.)
557 
558  !Make sure the basis and the projectors are well initialized
559  CALL init_interaction_radii(dft_control%qs_control, admm_kind_set)
560 
561  !We also init the atomic grids and harmonics
562  CALL local_rho_set_create(admm_gapw_env%local_rho_set)
563  CALL init_rho_atom(admm_gapw_env%local_rho_set%rho_atom_set, &
564  atomic_kind_set, admm_kind_set, dft_control, para_env)
565 
566  !Make sure that any NLCC potential is well initialized
567  CALL init_gapw_nlcc(admm_kind_set)
568 
569  !Need to have access to the soft AUX_FIT basis from the qs_env => add it to the qs_kinds
570  DO ikind = 1, nkind
571  NULLIFY (aux_fit_soft_basis)
572  CALL get_qs_kind(admm_kind_set(ikind), basis_set=soft_basis, basis_type="ORB_SOFT")
573  CALL copy_gto_basis_set(soft_basis, aux_fit_soft_basis)
574  CALL add_basis_set_to_container(qs_kind_set(ikind)%basis_sets, aux_fit_soft_basis, "AUX_FIT_SOFT")
575  END DO
576 
577  END SUBROUTINE init_admm_gapw
578 
579 ! **************************************************************************************************
580 !> \brief Builds the ADMM nmeighbor lists and overlap matrix on the model of qs_energies_init_hamiltonians()
581 !> \param admm_env ...
582 !> \param qs_env ...
583 !> \param aux_basis_type ...
584 ! **************************************************************************************************
585  SUBROUTINE admm_init_hamiltonians(admm_env, qs_env, aux_basis_type)
586 
587  TYPE(admm_type), POINTER :: admm_env
588  TYPE(qs_environment_type), POINTER :: qs_env
589  CHARACTER(len=*) :: aux_basis_type
590 
591  CHARACTER(len=*), PARAMETER :: routinen = 'admm_init_hamiltonians'
592 
593  INTEGER :: handle, hfx_pot, ikind, nkind
594  LOGICAL :: do_kpoints, mic, molecule_only
595  LOGICAL, ALLOCATABLE, DIMENSION(:) :: aux_fit_present, orb_present
596  REAL(dp) :: eps_schwarz, omega, pdist, roperator, &
597  subcells
598  REAL(dp), ALLOCATABLE, DIMENSION(:) :: aux_fit_radius, orb_radius
599  REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: pair_radius
600  TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
601  TYPE(cell_type), POINTER :: cell
602  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_aux_fit_kp, &
603  matrix_s_aux_fit_vs_orb_kp
604  TYPE(dft_control_type), POINTER :: dft_control
605  TYPE(distribution_1d_type), POINTER :: distribution_1d
606  TYPE(distribution_2d_type), POINTER :: distribution_2d
607  TYPE(gto_basis_set_type), POINTER :: aux_fit_basis_set, orb_basis_set
608  TYPE(kpoint_type), POINTER :: kpoints
609  TYPE(local_atoms_type), ALLOCATABLE, DIMENSION(:) :: atom2d
610  TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
611  TYPE(mp_para_env_type), POINTER :: para_env
612  TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
613  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
614  TYPE(qs_ks_env_type), POINTER :: ks_env
615  TYPE(section_vals_type), POINTER :: hfx_sections, neighbor_list_section
616 
617  NULLIFY (particle_set, cell, kpoints, distribution_1d, distribution_2d, molecule_set, &
618  atomic_kind_set, dft_control, neighbor_list_section, aux_fit_basis_set, orb_basis_set, &
619  ks_env, para_env, qs_kind_set, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb_kp)
620 
621  CALL timeset(routinen, handle)
622 
623  CALL get_qs_env(qs_env, nkind=nkind, particle_set=particle_set, cell=cell, kpoints=kpoints, &
624  local_particles=distribution_1d, distribution_2d=distribution_2d, &
625  molecule_set=molecule_set, atomic_kind_set=atomic_kind_set, do_kpoints=do_kpoints, &
626  dft_control=dft_control, para_env=para_env, qs_kind_set=qs_kind_set)
627  ALLOCATE (orb_present(nkind), aux_fit_present(nkind))
628  ALLOCATE (orb_radius(nkind), aux_fit_radius(nkind), pair_radius(nkind, nkind))
629  aux_fit_radius(:) = 0.0_dp
630 
631  molecule_only = .false.
632  IF (dft_control%qs_control%do_kg) molecule_only = .true.
633  mic = molecule_only
634  IF (kpoints%nkp > 0) THEN
635  mic = .false.
636  ELSEIF (dft_control%qs_control%semi_empirical) THEN
637  mic = .true.
638  END IF
639 
640  pdist = dft_control%qs_control%pairlist_radius
641 
642  CALL section_vals_val_get(qs_env%input, "DFT%SUBCELLS", r_val=subcells)
643  neighbor_list_section => section_vals_get_subs_vals(qs_env%input, "DFT%PRINT%NEIGHBOR_LISTS")
644 
645  ALLOCATE (atom2d(nkind))
646  CALL atom2d_build(atom2d, distribution_1d, distribution_2d, atomic_kind_set, &
647  molecule_set, molecule_only, particle_set=particle_set)
648 
649  DO ikind = 1, nkind
650  CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type="ORB")
651  IF (ASSOCIATED(orb_basis_set)) THEN
652  orb_present(ikind) = .true.
653  CALL get_gto_basis_set(gto_basis_set=orb_basis_set, kind_radius=orb_radius(ikind))
654  ELSE
655  orb_present(ikind) = .false.
656  END IF
657 
658  CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_fit_basis_set, basis_type=aux_basis_type)
659  IF (ASSOCIATED(aux_fit_basis_set)) THEN
660  aux_fit_present(ikind) = .true.
661  CALL get_gto_basis_set(gto_basis_set=aux_fit_basis_set, kind_radius=aux_fit_radius(ikind))
662  ELSE
663  aux_fit_present(ikind) = .false.
664  END IF
665  END DO
666 
667  IF (pdist < 0.0_dp) THEN
668  pdist = max(plane_distance(1, 0, 0, cell), &
669  plane_distance(0, 1, 0, cell), &
670  plane_distance(0, 0, 1, cell))
671  END IF
672 
673  !In case of K-points, we need to add the HFX potential range to sab_aux_fit, because it is used
674  !to populate AUX density and KS matrices
675  roperator = 0.0_dp
676  IF (do_kpoints) THEN
677  hfx_sections => section_vals_get_subs_vals(qs_env%input, "DFT%XC%HF")
678  CALL section_vals_val_get(hfx_sections, "INTERACTION_POTENTIAL%POTENTIAL_TYPE", i_val=hfx_pot)
679 
680  SELECT CASE (hfx_pot)
681  CASE (do_potential_id)
682  roperator = 0.0_dp
683  CASE (do_potential_truncated)
684  CALL section_vals_val_get(hfx_sections, "INTERACTION_POTENTIAL%CUTOFF_RADIUS", r_val=roperator)
685  CASE (do_potential_short)
686  CALL section_vals_val_get(hfx_sections, "INTERACTION_POTENTIAL%OMEGA", r_val=omega)
687  CALL section_vals_val_get(hfx_sections, "SCREENING%EPS_SCHWARZ", r_val=eps_schwarz)
688  CALL erfc_cutoff(eps_schwarz, omega, roperator)
689  CASE DEFAULT
690  cpabort("HFX potential not available for K-points (NYI)")
691  END SELECT
692  END IF
693 
694  CALL pair_radius_setup(aux_fit_present, aux_fit_present, aux_fit_radius, aux_fit_radius, pair_radius, pdist)
695  pair_radius = pair_radius + cutoff_screen_factor*roperator
696  CALL build_neighbor_lists(admm_env%sab_aux_fit, particle_set, atom2d, cell, pair_radius, &
697  mic=mic, molecular=molecule_only, subcells=subcells, nlname="sab_aux_fit")
698  CALL build_neighbor_lists(admm_env%sab_aux_fit_asymm, particle_set, atom2d, cell, pair_radius, &
699  mic=mic, symmetric=.false., molecular=molecule_only, subcells=subcells, &
700  nlname="sab_aux_fit_asymm")
701  CALL pair_radius_setup(aux_fit_present, orb_present, aux_fit_radius, orb_radius, pair_radius)
702  CALL build_neighbor_lists(admm_env%sab_aux_fit_vs_orb, particle_set, atom2d, cell, pair_radius, &
703  mic=mic, symmetric=.false., molecular=molecule_only, subcells=subcells, &
704  nlname="sab_aux_fit_vs_orb")
705 
706  CALL write_neighbor_lists(admm_env%sab_aux_fit, particle_set, cell, para_env, neighbor_list_section, &
707  "/SAB_AUX_FIT", "sab_aux_fit", "AUX_FIT_ORBITAL AUX_FIT_ORBITAL")
708  CALL write_neighbor_lists(admm_env%sab_aux_fit_vs_orb, particle_set, cell, para_env, neighbor_list_section, &
709  "/SAB_AUX_FIT_VS_ORB", "sab_aux_fit_vs_orb", "ORBITAL AUX_FIT_ORBITAL")
710 
711  CALL atom2d_cleanup(atom2d)
712 
713  !The ADMM overlap matrices (initially in qs_core_hamiltonian.F)
714  CALL get_qs_env(qs_env, ks_env=ks_env)
715 
716  CALL kpoint_transitional_release(admm_env%matrix_s_aux_fit)
717  CALL build_overlap_matrix(ks_env, matrixkp_s=matrix_s_aux_fit_kp, &
718  matrix_name="AUX_FIT_OVERLAP", &
719  basis_type_a=aux_basis_type, &
720  basis_type_b=aux_basis_type, &
721  sab_nl=admm_env%sab_aux_fit)
722  CALL set_2d_pointer(admm_env%matrix_s_aux_fit, matrix_s_aux_fit_kp)
723  CALL kpoint_transitional_release(admm_env%matrix_s_aux_fit_vs_orb)
724  CALL build_overlap_matrix(ks_env, matrixkp_s=matrix_s_aux_fit_vs_orb_kp, &
725  matrix_name="MIXED_OVERLAP", &
726  basis_type_a=aux_basis_type, &
727  basis_type_b="ORB", &
728  sab_nl=admm_env%sab_aux_fit_vs_orb)
729  CALL set_2d_pointer(admm_env%matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp)
730 
731  CALL timestop(handle)
732 
733  END SUBROUTINE admm_init_hamiltonians
734 
735 ! **************************************************************************************************
736 !> \brief Updates the ADMM task_list and density based on the model of qs_env_update_s_mstruct()
737 !> \param admm_env ...
738 !> \param qs_env ...
739 !> \param aux_basis_type ...
740 ! **************************************************************************************************
741  SUBROUTINE admm_update_s_mstruct(admm_env, qs_env, aux_basis_type)
742 
743  TYPE(admm_type), POINTER :: admm_env
744  TYPE(qs_environment_type), POINTER :: qs_env
745  CHARACTER(len=*) :: aux_basis_type
746 
747  CHARACTER(len=*), PARAMETER :: routinen = 'admm_update_s_mstruct'
748 
749  INTEGER :: handle
750  LOGICAL :: skip_load_balance_distributed
751  TYPE(dft_control_type), POINTER :: dft_control
752  TYPE(qs_ks_env_type), POINTER :: ks_env
753 
754  NULLIFY (ks_env, dft_control)
755 
756  CALL timeset(routinen, handle)
757 
758  CALL get_qs_env(qs_env, ks_env=ks_env, dft_control=dft_control)
759 
760  !The aux_fit task_list
761  skip_load_balance_distributed = dft_control%qs_control%skip_load_balance_distributed
762  IF (ASSOCIATED(admm_env%task_list_aux_fit)) CALL deallocate_task_list(admm_env%task_list_aux_fit)
763  CALL allocate_task_list(admm_env%task_list_aux_fit)
764  CALL generate_qs_task_list(ks_env, admm_env%task_list_aux_fit, &
765  reorder_rs_grid_ranks=.false., soft_valid=.false., &
766  basis_type=aux_basis_type, &
767  skip_load_balance_distributed=skip_load_balance_distributed, &
768  sab_orb_external=admm_env%sab_aux_fit)
769 
770  !The aux_fit densities
771  CALL qs_rho_rebuild(admm_env%rho_aux_fit, qs_env=qs_env, admm=.true.)
772  CALL qs_rho_rebuild(admm_env%rho_aux_fit_buffer, qs_env=qs_env, admm=.true.)
773 
774  CALL timestop(handle)
775 
776  END SUBROUTINE admm_update_s_mstruct
777 
778 ! **************************************************************************************************
779 !> \brief Update the admm_gapw_env internals to the current qs_env (i.e. atomic positions)
780 !> \param qs_env ...
781 ! **************************************************************************************************
782  SUBROUTINE update_admm_gapw(qs_env)
783 
784  TYPE(qs_environment_type), POINTER :: qs_env
785 
786  CHARACTER(len=*), PARAMETER :: routinen = 'update_admm_gapw'
787 
788  INTEGER :: handle, ikind, nkind
789  LOGICAL :: paw_atom
790  LOGICAL, ALLOCATABLE, DIMENSION(:) :: aux_present, oce_present
791  REAL(dp) :: subcells
792  REAL(dp), ALLOCATABLE, DIMENSION(:) :: aux_radius, oce_radius
793  REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: pair_radius
794  TYPE(admm_gapw_r3d_rs_type), POINTER :: admm_gapw_env
795  TYPE(admm_type), POINTER :: admm_env
796  TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
797  TYPE(cell_type), POINTER :: cell
798  TYPE(dft_control_type), POINTER :: dft_control
799  TYPE(distribution_1d_type), POINTER :: distribution_1d
800  TYPE(distribution_2d_type), POINTER :: distribution_2d
801  TYPE(gto_basis_set_type), POINTER :: aux_fit_basis
802  TYPE(local_atoms_type), ALLOCATABLE, DIMENSION(:) :: atom2d
803  TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
804  TYPE(neighbor_list_set_p_type), DIMENSION(:), &
805  POINTER :: sap_oce
806  TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
807  TYPE(paw_proj_set_type), POINTER :: paw_proj
808  TYPE(qs_kind_type), DIMENSION(:), POINTER :: admm_kind_set, qs_kind_set
809  TYPE(qs_ks_env_type), POINTER :: ks_env
810 
811  NULLIFY (ks_env, qs_kind_set, admm_kind_set, aux_fit_basis, cell, distribution_1d)
812  NULLIFY (distribution_2d, paw_proj, particle_set, molecule_set, admm_env, admm_gapw_env)
813  NULLIFY (dft_control, atomic_kind_set, sap_oce)
814 
815  CALL timeset(routinen, handle)
816 
817  CALL get_qs_env(qs_env, ks_env=ks_env, qs_kind_set=qs_kind_set, admm_env=admm_env, &
818  dft_control=dft_control)
819  admm_gapw_env => admm_env%admm_gapw_env
820  admm_kind_set => admm_gapw_env%admm_kind_set
821  nkind = SIZE(qs_kind_set)
822 
823  !Update the task lisft for the AUX_FIT_SOFT basis
824  IF (ASSOCIATED(admm_gapw_env%task_list)) CALL deallocate_task_list(admm_gapw_env%task_list)
825  CALL allocate_task_list(admm_gapw_env%task_list)
826 
827  !note: we set soft_valid to .FALSE. want to use AUX_FIT_SOFT and not the normal ORB SOFT basis
828  CALL generate_qs_task_list(ks_env, admm_gapw_env%task_list, reorder_rs_grid_ranks=.false., &
829  soft_valid=.false., basis_type="AUX_FIT_SOFT", &
830  skip_load_balance_distributed=dft_control%qs_control%skip_load_balance_distributed, &
831  sab_orb_external=admm_env%sab_aux_fit)
832 
833  !Update the precomputed oce integrals
834  !a sap_oce neighbor list is required => build it here
835  ALLOCATE (aux_present(nkind), oce_present(nkind))
836  aux_present = .false.; oce_present = .false.
837  ALLOCATE (aux_radius(nkind), oce_radius(nkind))
838  aux_radius = 0.0_dp; oce_radius = 0.0_dp
839 
840  DO ikind = 1, nkind
841  CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_fit_basis, basis_type="AUX_FIT")
842  IF (ASSOCIATED(aux_fit_basis)) THEN
843  aux_present(ikind) = .true.
844  CALL get_gto_basis_set(aux_fit_basis, kind_radius=aux_radius(ikind))
845  END IF
846 
847  !note: get oce info from admm_kind_set
848  CALL get_qs_kind(admm_kind_set(ikind), paw_atom=paw_atom, paw_proj_set=paw_proj)
849  IF (paw_atom) THEN
850  oce_present(ikind) = .true.
851  CALL get_paw_proj_set(paw_proj, rcprj=oce_radius(ikind))
852  END IF
853  END DO
854 
855  ALLOCATE (pair_radius(nkind, nkind))
856  pair_radius = 0.0_dp
857  CALL pair_radius_setup(aux_present, oce_present, aux_radius, oce_radius, pair_radius)
858 
859  CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set, cell=cell, &
860  distribution_2d=distribution_2d, local_particles=distribution_1d, &
861  particle_set=particle_set, molecule_set=molecule_set)
862  CALL section_vals_val_get(qs_env%input, "DFT%SUBCELLS", r_val=subcells)
863 
864  ALLOCATE (atom2d(nkind))
865  CALL atom2d_build(atom2d, distribution_1d, distribution_2d, atomic_kind_set, &
866  molecule_set, .false., particle_set)
867  CALL build_neighbor_lists(sap_oce, particle_set, atom2d, cell, pair_radius, &
868  subcells=subcells, operator_type="ABBA", nlname="AUX_PAW-PRJ")
869  CALL atom2d_cleanup(atom2d)
870 
871  !actually compute the oce matrices
872  CALL create_oce_set(admm_gapw_env%oce)
873  CALL allocate_oce_set(admm_gapw_env%oce, nkind)
874 
875  !always compute the derivative, cheap anyways
876  CALL build_oce_matrices(admm_gapw_env%oce%intac, calculate_forces=.true., nder=1, &
877  qs_kind_set=admm_kind_set, particle_set=particle_set, &
878  sap_oce=sap_oce, eps_fit=dft_control%qs_control%gapw_control%eps_fit)
879 
880  CALL release_neighbor_list_sets(sap_oce)
881 
882  CALL timestop(handle)
883 
884  END SUBROUTINE update_admm_gapw
885 
886 ! **************************************************************************************************
887 !> \brief Allocates the various ADMM KS matrices
888 !> \param admm_env ...
889 !> \param qs_env ...
890 ! **************************************************************************************************
891  SUBROUTINE admm_alloc_ks_matrices(admm_env, qs_env)
892 
893  TYPE(admm_type), POINTER :: admm_env
894  TYPE(qs_environment_type), POINTER :: qs_env
895 
896  CHARACTER(len=*), PARAMETER :: routinen = 'admm_alloc_ks_matrices'
897 
898  INTEGER :: handle, ic, ispin
899  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_ks_aux_fit_dft_kp, &
900  matrix_ks_aux_fit_hfx_kp, &
901  matrix_ks_aux_fit_kp, &
902  matrix_s_aux_fit_kp
903  TYPE(dft_control_type), POINTER :: dft_control
904 
905  NULLIFY (dft_control, matrix_s_aux_fit_kp, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp)
906 
907  CALL timeset(routinen, handle)
908 
909  CALL get_qs_env(qs_env, dft_control=dft_control)
910  CALL get_admm_env(admm_env, matrix_s_aux_fit_kp=matrix_s_aux_fit_kp)
911 
912  CALL kpoint_transitional_release(admm_env%matrix_ks_aux_fit)
913  CALL kpoint_transitional_release(admm_env%matrix_ks_aux_fit_dft)
914  CALL kpoint_transitional_release(admm_env%matrix_ks_aux_fit_hfx)
915 
916  CALL dbcsr_allocate_matrix_set(matrix_ks_aux_fit_kp, dft_control%nspins, dft_control%nimages)
917  CALL dbcsr_allocate_matrix_set(matrix_ks_aux_fit_dft_kp, dft_control%nspins, dft_control%nimages)
918  CALL dbcsr_allocate_matrix_set(matrix_ks_aux_fit_hfx_kp, dft_control%nspins, dft_control%nimages)
919 
920  DO ispin = 1, dft_control%nspins
921  DO ic = 1, dft_control%nimages
922  ALLOCATE (matrix_ks_aux_fit_kp(ispin, ic)%matrix)
923  CALL dbcsr_create(matrix_ks_aux_fit_kp(ispin, ic)%matrix, template=matrix_s_aux_fit_kp(1, ic)%matrix, &
924  name="KOHN-SHAM_MATRIX for ADMM")
925  CALL cp_dbcsr_alloc_block_from_nbl(matrix_ks_aux_fit_kp(ispin, ic)%matrix, admm_env%sab_aux_fit)
926  CALL dbcsr_set(matrix_ks_aux_fit_kp(ispin, ic)%matrix, 0.0_dp)
927 
928  ALLOCATE (matrix_ks_aux_fit_dft_kp(ispin, ic)%matrix)
929  CALL dbcsr_create(matrix_ks_aux_fit_dft_kp(ispin, ic)%matrix, template=matrix_s_aux_fit_kp(1, 1)%matrix, &
930  name="KOHN-SHAM_MATRIX for ADMM")
931  CALL cp_dbcsr_alloc_block_from_nbl(matrix_ks_aux_fit_dft_kp(ispin, ic)%matrix, admm_env%sab_aux_fit)
932  CALL dbcsr_set(matrix_ks_aux_fit_dft_kp(ispin, ic)%matrix, 0.0_dp)
933 
934  ALLOCATE (matrix_ks_aux_fit_hfx_kp(ispin, ic)%matrix)
935  CALL dbcsr_create(matrix_ks_aux_fit_hfx_kp(ispin, ic)%matrix, template=matrix_s_aux_fit_kp(1, 1)%matrix, &
936  name="KOHN-SHAM_MATRIX for ADMM")
937  CALL cp_dbcsr_alloc_block_from_nbl(matrix_ks_aux_fit_hfx_kp(ispin, ic)%matrix, admm_env%sab_aux_fit)
938  CALL dbcsr_set(matrix_ks_aux_fit_hfx_kp(ispin, ic)%matrix, 0.0_dp)
939  END DO
940  END DO
941 
942  CALL set_admm_env(admm_env, &
943  matrix_ks_aux_fit_kp=matrix_ks_aux_fit_kp, &
944  matrix_ks_aux_fit_dft_kp=matrix_ks_aux_fit_dft_kp, &
945  matrix_ks_aux_fit_hfx_kp=matrix_ks_aux_fit_hfx_kp)
946 
947  CALL timestop(handle)
948 
949  END SUBROUTINE admm_alloc_ks_matrices
950 
951 ! **************************************************************************************************
952 !> \brief Add the HFX K-point contribution to the real-space Hamiltonians
953 !> \param qs_env ...
954 !> \param matrix_ks ...
955 !> \param energy ...
956 !> \param calculate_forces ...
957 ! **************************************************************************************************
958  SUBROUTINE hfx_ks_matrix_kp(qs_env, matrix_ks, energy, calculate_forces)
959  TYPE(qs_environment_type), POINTER :: qs_env
960  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_ks
961  TYPE(qs_energy_type), POINTER :: energy
962  LOGICAL, INTENT(in) :: calculate_forces
963 
964  CHARACTER(LEN=*), PARAMETER :: routinen = 'hfx_ks_matrix_kp'
965 
966  INTEGER :: handle, img, irep, ispin, n_rep_hf, &
967  nimages, nspins
968  LOGICAL :: do_adiabatic_rescaling, &
969  s_mstruct_changed, use_virial
970  REAL(dp) :: eh1, ehfx, eold
971  REAL(dp), ALLOCATABLE, DIMENSION(:) :: hf_energy
972  TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks_aux_fit_im, matrix_ks_im
973  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h, matrix_ks_aux_fit_hfx_kp, &
974  matrix_ks_aux_fit_kp, matrix_ks_orb, &
975  rho_ao_orb
976  TYPE(dft_control_type), POINTER :: dft_control
977  TYPE(hfx_type), DIMENSION(:, :), POINTER :: x_data
978  TYPE(mp_para_env_type), POINTER :: para_env
979  TYPE(pw_env_type), POINTER :: pw_env
980  TYPE(pw_poisson_type), POINTER :: poisson_env
981  TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
982  TYPE(qs_rho_type), POINTER :: rho_orb
983  TYPE(section_vals_type), POINTER :: adiabatic_rescaling_section, &
984  hfx_sections, input
985  TYPE(virial_type), POINTER :: virial
986 
987  CALL timeset(routinen, handle)
988 
989  NULLIFY (auxbas_pw_pool, dft_control, hfx_sections, input, &
990  para_env, poisson_env, pw_env, virial, matrix_ks_im, &
991  matrix_ks_orb, rho_ao_orb, matrix_h, matrix_ks_aux_fit_kp, &
992  matrix_ks_aux_fit_im, matrix_ks_aux_fit_hfx_kp)
993 
994  CALL get_qs_env(qs_env=qs_env, &
995  dft_control=dft_control, &
996  input=input, &
997  matrix_h_kp=matrix_h, &
998  para_env=para_env, &
999  pw_env=pw_env, &
1000  virial=virial, &
1001  matrix_ks_im=matrix_ks_im, &
1002  s_mstruct_changed=s_mstruct_changed, &
1003  x_data=x_data)
1004 
1005  ! No RTP
1006  IF (qs_env%run_rtp) cpabort("No RTP implementation with K-points HFX")
1007 
1008  ! No adiabatic rescaling
1009  adiabatic_rescaling_section => section_vals_get_subs_vals(input, "DFT%XC%ADIABATIC_RESCALING")
1010  CALL section_vals_get(adiabatic_rescaling_section, explicit=do_adiabatic_rescaling)
1011  IF (do_adiabatic_rescaling) cpabort("No adiabatic rescaling implementation with K-points HFX")
1012 
1013  IF (dft_control%do_admm) THEN
1014  CALL get_admm_env(qs_env%admm_env, matrix_ks_aux_fit_kp=matrix_ks_aux_fit_kp, &
1015  matrix_ks_aux_fit_im=matrix_ks_aux_fit_im, &
1016  matrix_ks_aux_fit_hfx_kp=matrix_ks_aux_fit_hfx_kp)
1017  END IF
1018 
1019  nspins = dft_control%nspins
1020  nimages = dft_control%nimages
1021 
1022  use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1023  IF (use_virial .AND. calculate_forces) virial%pv_fock_4c = 0.0_dp
1024 
1025  hfx_sections => section_vals_get_subs_vals(input, "DFT%XC%HF")
1026  CALL section_vals_get(hfx_sections, n_repetition=n_rep_hf)
1027 
1028  ! *** Initialize the auxiliary ks matrix to zero if required
1029  IF (dft_control%do_admm) THEN
1030  DO ispin = 1, nspins
1031  DO img = 1, nimages
1032  CALL dbcsr_set(matrix_ks_aux_fit_kp(ispin, img)%matrix, 0.0_dp)
1033  END DO
1034  END DO
1035  END IF
1036  DO ispin = 1, nspins
1037  DO img = 1, nimages
1038  CALL dbcsr_set(matrix_ks(ispin, img)%matrix, 0.0_dp)
1039  END DO
1040  END DO
1041 
1042  ALLOCATE (hf_energy(n_rep_hf))
1043 
1044  eold = 0.0_dp
1045 
1046  DO irep = 1, n_rep_hf
1047 
1048  ! fetch the correct matrices for normal HFX or ADMM
1049  IF (dft_control%do_admm) THEN
1050  CALL get_admm_env(qs_env%admm_env, matrix_ks_aux_fit_kp=matrix_ks_orb, rho_aux_fit=rho_orb)
1051  ELSE
1052  CALL get_qs_env(qs_env=qs_env, matrix_ks_kp=matrix_ks_orb, rho=rho_orb)
1053  END IF
1054  CALL qs_rho_get(rho_struct=rho_orb, rho_ao_kp=rho_ao_orb)
1055 
1056  ! Finally the real hfx calulation
1057  ehfx = 0.0_dp
1058 
1059  IF (.NOT. x_data(irep, 1)%do_hfx_ri) THEN
1060  cpabort("Only RI-HFX is implemented for K-points")
1061  END IF
1062 
1063  CALL hfx_ri_update_ks_kp(qs_env, x_data(irep, 1)%ri_data, matrix_ks_orb, ehfx, &
1064  rho_ao_orb, s_mstruct_changed, nspins, &
1065  x_data(irep, 1)%general_parameter%fraction)
1066 
1067  IF (calculate_forces) THEN
1068  !Scale auxiliary density matrix for ADMMP (see Merlot2014) with gsi(ispin) to scale force
1069  IF (dft_control%do_admm) THEN
1070  CALL scale_dm(qs_env, rho_ao_orb, scale_back=.false.)
1071  END IF
1072 
1073  CALL hfx_ri_update_forces_kp(qs_env, x_data(irep, 1)%ri_data, nspins, &
1074  x_data(irep, 1)%general_parameter%fraction, &
1075  rho_ao_orb, use_virial=use_virial)
1076 
1077  IF (dft_control%do_admm) THEN
1078  CALL scale_dm(qs_env, rho_ao_orb, scale_back=.true.)
1079  END IF
1080  END IF
1081 
1082  CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, poisson_env=poisson_env)
1083  eh1 = ehfx - eold
1084  CALL pw_hfx(qs_env, eh1, hfx_sections, poisson_env, auxbas_pw_pool, irep)
1085  eold = ehfx
1086 
1087  END DO
1088 
1089  ! *** Set the total HFX energy
1090  energy%ex = ehfx
1091 
1092  ! *** Add Core-Hamiltonian-Matrix ***
1093  DO ispin = 1, nspins
1094  DO img = 1, nimages
1095  CALL dbcsr_add(matrix_ks(ispin, img)%matrix, matrix_h(1, img)%matrix, &
1096  1.0_dp, 1.0_dp)
1097  END DO
1098  END DO
1099  IF (use_virial .AND. calculate_forces) THEN
1100  virial%pv_exx = virial%pv_exx - virial%pv_fock_4c
1101  virial%pv_virial = virial%pv_virial - virial%pv_fock_4c
1102  virial%pv_calculate = .false.
1103  END IF
1104 
1105  !update the hfx aux_fit matrix
1106  IF (dft_control%do_admm) THEN
1107  DO ispin = 1, nspins
1108  DO img = 1, nimages
1109  CALL dbcsr_add(matrix_ks_aux_fit_hfx_kp(ispin, img)%matrix, matrix_ks_aux_fit_kp(ispin, img)%matrix, &
1110  0.0_dp, 1.0_dp)
1111  END DO
1112  END DO
1113  END IF
1114 
1115  CALL timestop(handle)
1116 
1117  END SUBROUTINE hfx_ks_matrix_kp
1118 
1119 ! **************************************************************************************************
1120 !> \brief Add the hfx contributions to the Hamiltonian
1121 !>
1122 !> \param qs_env ...
1123 !> \param matrix_ks ...
1124 !> \param rho ...
1125 !> \param energy ...
1126 !> \param calculate_forces ...
1127 !> \param just_energy ...
1128 !> \param v_rspace_new ...
1129 !> \param v_tau_rspace ...
1130 !> \par History
1131 !> refactoring 03-2011 [MI]
1132 ! **************************************************************************************************
1133 
1134  SUBROUTINE hfx_ks_matrix(qs_env, matrix_ks, rho, energy, calculate_forces, &
1135  just_energy, v_rspace_new, v_tau_rspace)
1136 
1137  TYPE(qs_environment_type), POINTER :: qs_env
1138  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_ks
1139  TYPE(qs_rho_type), POINTER :: rho
1140  TYPE(qs_energy_type), POINTER :: energy
1141  LOGICAL, INTENT(in) :: calculate_forces, just_energy
1142  TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: v_rspace_new, v_tau_rspace
1143 
1144  CHARACTER(LEN=*), PARAMETER :: routinen = 'hfx_ks_matrix'
1145 
1146  INTEGER :: handle, img, irep, ispin, mspin, &
1147  n_rep_hf, nimages, ns, nspins
1148  LOGICAL :: distribute_fock_matrix, &
1149  do_adiabatic_rescaling, &
1150  hfx_treat_lsd_in_core, &
1151  s_mstruct_changed, use_virial
1152  REAL(dp) :: eh1, ehfx, ehfxrt, eold
1153  REAL(dp), ALLOCATABLE, DIMENSION(:) :: hf_energy
1154  TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks_1d, matrix_ks_aux_fit, &
1155  matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_im, matrix_ks_im, rho_ao_1d, rho_ao_resp
1156  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h, matrix_h_im, matrix_ks_orb, &
1157  rho_ao_orb
1158  TYPE(dft_control_type), POINTER :: dft_control
1159  TYPE(hfx_type), DIMENSION(:, :), POINTER :: x_data
1160  TYPE(mo_set_type), DIMENSION(:), POINTER :: mo_array
1161  TYPE(mp_para_env_type), POINTER :: para_env
1162  TYPE(pw_env_type), POINTER :: pw_env
1163  TYPE(pw_poisson_type), POINTER :: poisson_env
1164  TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
1165  TYPE(qs_rho_type), POINTER :: rho_orb
1166  TYPE(rt_prop_type), POINTER :: rtp
1167  TYPE(section_vals_type), POINTER :: adiabatic_rescaling_section, &
1168  hfx_sections, input
1169  TYPE(virial_type), POINTER :: virial
1170 
1171  CALL timeset(routinen, handle)
1172 
1173  NULLIFY (auxbas_pw_pool, dft_control, hfx_sections, input, &
1174  para_env, poisson_env, pw_env, virial, matrix_ks_im, &
1175  matrix_ks_orb, rho_ao_orb, matrix_h, matrix_h_im, matrix_ks_aux_fit, &
1176  matrix_ks_aux_fit_im, matrix_ks_aux_fit_hfx)
1177 
1178  CALL get_qs_env(qs_env=qs_env, &
1179  dft_control=dft_control, &
1180  input=input, &
1181  matrix_h_kp=matrix_h, &
1182  matrix_h_im_kp=matrix_h_im, &
1183  para_env=para_env, &
1184  pw_env=pw_env, &
1185  virial=virial, &
1186  matrix_ks_im=matrix_ks_im, &
1187  s_mstruct_changed=s_mstruct_changed, &
1188  x_data=x_data)
1189 
1190  IF (dft_control%do_admm) THEN
1191  CALL get_admm_env(qs_env%admm_env, mos_aux_fit=mo_array, matrix_ks_aux_fit=matrix_ks_aux_fit, &
1192  matrix_ks_aux_fit_im=matrix_ks_aux_fit_im, matrix_ks_aux_fit_hfx=matrix_ks_aux_fit_hfx)
1193  ELSE
1194  CALL get_qs_env(qs_env=qs_env, mos=mo_array)
1195  END IF
1196 
1197  nspins = dft_control%nspins
1198  nimages = dft_control%nimages
1199 
1200  use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1201 
1202  IF (use_virial .AND. calculate_forces) virial%pv_fock_4c = 0.0_dp
1203 
1204  hfx_sections => section_vals_get_subs_vals(input, "DFT%XC%HF")
1205  CALL section_vals_get(hfx_sections, n_repetition=n_rep_hf)
1206  CALL section_vals_val_get(hfx_sections, "TREAT_LSD_IN_CORE", l_val=hfx_treat_lsd_in_core, &
1207  i_rep_section=1)
1208  adiabatic_rescaling_section => section_vals_get_subs_vals(input, "DFT%XC%ADIABATIC_RESCALING")
1209  CALL section_vals_get(adiabatic_rescaling_section, explicit=do_adiabatic_rescaling)
1210 
1211  ! *** Initialize the auxiliary ks matrix to zero if required
1212  IF (dft_control%do_admm) THEN
1213  DO ispin = 1, nspins
1214  CALL dbcsr_set(matrix_ks_aux_fit(ispin)%matrix, 0.0_dp)
1215  END DO
1216  END IF
1217  DO ispin = 1, nspins
1218  DO img = 1, nimages
1219  CALL dbcsr_set(matrix_ks(ispin, img)%matrix, 0.0_dp)
1220  END DO
1221  END DO
1222 
1223  CALL section_vals_get(hfx_sections, n_repetition=n_rep_hf)
1224 
1225  ALLOCATE (hf_energy(n_rep_hf))
1226 
1227  eold = 0.0_dp
1228 
1229  DO irep = 1, n_rep_hf
1230  ! Remember: Vhfx is added, energy is calclulated from total Vhfx,
1231  ! so energy of last iteration is correct
1232 
1233  IF (do_adiabatic_rescaling .AND. hfx_treat_lsd_in_core) &
1234  cpabort("HFX_TREAT_LSD_IN_CORE not implemented for adiabatically rescaled hybrids")
1235  ! everything is calculated with adiabatic rescaling but the potential is not added in a first step
1236  distribute_fock_matrix = .NOT. do_adiabatic_rescaling
1237 
1238  mspin = 1
1239  IF (hfx_treat_lsd_in_core) mspin = nspins
1240 
1241  ! fetch the correct matrices for normal HFX or ADMM
1242  IF (dft_control%do_admm) THEN
1243  CALL get_admm_env(qs_env%admm_env, matrix_ks_aux_fit=matrix_ks_1d, rho_aux_fit=rho_orb)
1244  ns = SIZE(matrix_ks_1d)
1245  matrix_ks_orb(1:ns, 1:1) => matrix_ks_1d(1:ns)
1246  ELSE
1247  CALL get_qs_env(qs_env=qs_env, matrix_ks_kp=matrix_ks_orb, rho=rho_orb)
1248  END IF
1249  CALL qs_rho_get(rho_struct=rho_orb, rho_ao_kp=rho_ao_orb)
1250  ! Finally the real hfx calulation
1251  ehfx = 0.0_dp
1252 
1253  IF (x_data(irep, 1)%do_hfx_ri) THEN
1254 
1255  CALL hfx_ri_update_ks(qs_env, x_data(irep, 1)%ri_data, matrix_ks_orb, ehfx, &
1256  mo_array, rho_ao_orb, &
1257  s_mstruct_changed, nspins, &
1258  x_data(irep, 1)%general_parameter%fraction)
1259  IF (dft_control%do_admm) THEN
1260  !for ADMMS, we need the exchange matrix k(d) for both spins
1261  DO ispin = 1, nspins
1262  CALL dbcsr_copy(matrix_ks_aux_fit_hfx(ispin)%matrix, matrix_ks_orb(ispin, 1)%matrix, &
1263  name="HF exch. part of matrix_ks_aux_fit for ADMMS")
1264  END DO
1265  END IF
1266 
1267  ELSE
1268 
1269  DO ispin = 1, mspin
1270  CALL integrate_four_center(qs_env, x_data, matrix_ks_orb, eh1, rho_ao_orb, hfx_sections, &
1271  para_env, s_mstruct_changed, irep, distribute_fock_matrix, &
1272  ispin=ispin)
1273  ehfx = ehfx + eh1
1274  END DO
1275  END IF
1276 
1277  IF (calculate_forces .AND. .NOT. do_adiabatic_rescaling) THEN
1278  !Scale auxiliary density matrix for ADMMP (see Merlot2014) with gsi(ispin) to scale force
1279  IF (dft_control%do_admm) THEN
1280  CALL scale_dm(qs_env, rho_ao_orb, scale_back=.false.)
1281  END IF
1282  NULLIFY (rho_ao_resp)
1283 
1284  IF (x_data(irep, 1)%do_hfx_ri) THEN
1285 
1286  CALL hfx_ri_update_forces(qs_env, x_data(irep, 1)%ri_data, nspins, &
1287  x_data(irep, 1)%general_parameter%fraction, &
1288  rho_ao=rho_ao_orb, mos=mo_array, &
1289  rho_ao_resp=rho_ao_resp, &
1290  use_virial=use_virial)
1291 
1292  ELSE
1293 
1294  CALL derivatives_four_center(qs_env, rho_ao_orb, rho_ao_resp, hfx_sections, &
1295  para_env, irep, use_virial)
1296 
1297  END IF
1298 
1299  !Scale auxiliary density matrix for ADMMP back with 1/gsi(ispin)
1300  IF (dft_control%do_admm) THEN
1301  CALL scale_dm(qs_env, rho_ao_orb, scale_back=.true.)
1302  END IF
1303  END IF
1304 
1305  !! If required, the calculation of the forces will be done later with adiabatic rescaling
1306  IF (do_adiabatic_rescaling) hf_energy(irep) = ehfx
1307 
1308  ! special case RTP/EMD we have a full complex density and HFX has a contribution from the imaginary part
1309  ehfxrt = 0.0_dp
1310  IF (qs_env%run_rtp) THEN
1311 
1312  CALL get_qs_env(qs_env=qs_env, rtp=rtp)
1313  DO ispin = 1, nspins
1314  CALL dbcsr_set(matrix_ks_im(ispin)%matrix, 0.0_dp)
1315  END DO
1316  IF (dft_control%do_admm) THEN
1317  ! matrix_ks_orb => matrix_ks_aux_fit_im
1318  ns = SIZE(matrix_ks_aux_fit_im)
1319  matrix_ks_orb(1:ns, 1:1) => matrix_ks_aux_fit_im(1:ns)
1320  DO ispin = 1, nspins
1321  CALL dbcsr_set(matrix_ks_aux_fit_im(ispin)%matrix, 0.0_dp)
1322  END DO
1323  ELSE
1324  ! matrix_ks_orb => matrix_ks_im
1325  ns = SIZE(matrix_ks_im)
1326  matrix_ks_orb(1:ns, 1:1) => matrix_ks_im(1:ns)
1327  END IF
1328 
1329  CALL qs_rho_get(rho_orb, rho_ao_im=rho_ao_1d)
1330  ns = SIZE(rho_ao_1d)
1331  rho_ao_orb(1:ns, 1:1) => rho_ao_1d(1:ns)
1332 
1333  ehfxrt = 0.0_dp
1334 
1335  IF (x_data(irep, 1)%do_hfx_ri) THEN
1336  CALL hfx_ri_update_ks(qs_env, x_data(irep, 1)%ri_data, matrix_ks_orb, ehfx, &
1337  mo_array, rho_ao_orb, &
1338  .false., nspins, &
1339  x_data(irep, 1)%general_parameter%fraction)
1340  IF (dft_control%do_admm) THEN
1341  !for ADMMS, we need the exchange matrix k(d) for both spins
1342  DO ispin = 1, nspins
1343  CALL dbcsr_copy(matrix_ks_aux_fit_hfx(ispin)%matrix, matrix_ks_orb(ispin, 1)%matrix, &
1344  name="HF exch. part of matrix_ks_aux_fit for ADMMS")
1345  END DO
1346  END IF
1347 
1348  ELSE
1349  DO ispin = 1, mspin
1350  CALL integrate_four_center(qs_env, x_data, matrix_ks_orb, eh1, rho_ao_orb, hfx_sections, &
1351  para_env, .false., irep, distribute_fock_matrix, &
1352  ispin=ispin)
1353  ehfxrt = ehfxrt + eh1
1354  END DO
1355  END IF
1356 
1357  IF (calculate_forces .AND. .NOT. do_adiabatic_rescaling) THEN
1358  NULLIFY (rho_ao_resp)
1359 
1360  IF (x_data(irep, 1)%do_hfx_ri) THEN
1361 
1362  CALL hfx_ri_update_forces(qs_env, x_data(irep, 1)%ri_data, nspins, &
1363  x_data(irep, 1)%general_parameter%fraction, &
1364  rho_ao=rho_ao_orb, mos=mo_array, &
1365  use_virial=use_virial)
1366 
1367  ELSE
1368  CALL derivatives_four_center(qs_env, rho_ao_orb, rho_ao_resp, hfx_sections, &
1369  para_env, irep, use_virial)
1370  END IF
1371  END IF
1372 
1373  !! If required, the calculation of the forces will be done later with adiabatic rescaling
1374  IF (do_adiabatic_rescaling) hf_energy(irep) = ehfx + ehfxrt
1375 
1376  IF (dft_control%rtp_control%velocity_gauge) THEN
1377  cpassert(ASSOCIATED(matrix_h_im))
1378  DO ispin = 1, nspins
1379  CALL dbcsr_add(matrix_ks_im(ispin)%matrix, matrix_h_im(1, 1)%matrix, &
1380  1.0_dp, 1.0_dp)
1381  END DO
1382  END IF
1383 
1384  END IF
1385 
1386  IF (.NOT. qs_env%run_rtp) THEN
1387  CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
1388  poisson_env=poisson_env)
1389  eh1 = ehfx - eold
1390  CALL pw_hfx(qs_env, eh1, hfx_sections, poisson_env, auxbas_pw_pool, irep)
1391  eold = ehfx
1392  END IF
1393 
1394  END DO
1395 
1396  ! *** Set the total HFX energy
1397  energy%ex = ehfx + ehfxrt
1398 
1399  ! *** Add Core-Hamiltonian-Matrix ***
1400  DO ispin = 1, nspins
1401  DO img = 1, nimages
1402  CALL dbcsr_add(matrix_ks(ispin, img)%matrix, matrix_h(1, img)%matrix, &
1403  1.0_dp, 1.0_dp)
1404  END DO
1405  END DO
1406  IF (use_virial .AND. calculate_forces) THEN
1407  virial%pv_exx = virial%pv_exx - virial%pv_fock_4c
1408  virial%pv_virial = virial%pv_virial - virial%pv_fock_4c
1409  virial%pv_calculate = .false.
1410  END IF
1411 
1412  !! If we perform adiabatic rescaling we are now able to rescale the xc-potential
1413  IF (do_adiabatic_rescaling) THEN
1414  CALL rescale_xc_potential(qs_env, matrix_ks, rho, energy, v_rspace_new, v_tau_rspace, &
1415  hf_energy, just_energy, calculate_forces, use_virial)
1416  END IF ! do_adiabatic_rescaling
1417 
1418  !update the hfx aux_fit matrixIF (dft_control%do_admm) THEN
1419  IF (dft_control%do_admm) THEN
1420  DO ispin = 1, nspins
1421  CALL dbcsr_add(matrix_ks_aux_fit_hfx(ispin)%matrix, matrix_ks_aux_fit(ispin)%matrix, &
1422  0.0_dp, 1.0_dp)
1423  END DO
1424  END IF
1425 
1426  CALL timestop(handle)
1427 
1428  END SUBROUTINE hfx_ks_matrix
1429 
1430 ! **************************************************************************************************
1431 !> \brief This routine modifies the xc section depending on the potential type
1432 !> used for the HF exchange and the resulting correction term. Currently
1433 !> three types of corrections are implemented:
1434 !>
1435 !> coulomb: Ex,hf = Ex,hf' + (PBEx-PBEx')
1436 !> shortrange: Ex,hf = Ex,hf' + (XWPBEX-XWPBEX')
1437 !> truncated: Ex,hf = Ex,hf' + ( (XWPBEX0-PBE_HOLE_TC_LR) -(XWPBEX0-PBE_HOLE_TC_LR)' )
1438 !>
1439 !> with ' denoting the auxiliary basis set and
1440 !>
1441 !> PBEx: PBE exchange functional
1442 !> XWPBEX: PBE exchange hole for short-range potential (erfc(omega*r)/r)
1443 !> XWPBEX0: PBE exchange hole for standard coulomb potential
1444 !> PBE_HOLE_TC_LR: PBE exchange hole for longrange truncated coulomb potential
1445 !>
1446 !> Above explanation is correct for the deafult case. If a specific functional is requested
1447 !> for the correction term (cfun), we get
1448 !> Ex,hf = Ex,hf' + (cfun-cfun')
1449 !> for all cases of operators.
1450 !>
1451 !> \param x_data ...
1452 !> \param xc_section the original xc_section
1453 !> \param admm_env the ADMM environment
1454 !> \par History
1455 !> 12.2009 created [Manuel Guidon]
1456 !> 05.2021 simplify for case of no correction [JGH]
1457 !> \author Manuel Guidon
1458 ! **************************************************************************************************
1459  SUBROUTINE create_admm_xc_section(x_data, xc_section, admm_env)
1460  TYPE(hfx_type), DIMENSION(:, :), POINTER :: x_data
1461  TYPE(section_vals_type), POINTER :: xc_section
1462  TYPE(admm_type), POINTER :: admm_env
1463 
1464  LOGICAL, PARAMETER :: debug_functional = .false.
1465 #if defined (__LIBXC)
1466  REAL(kind=dp), PARAMETER :: x_factor_c = 0.930525736349100025_dp
1467 #endif
1468 
1469  CHARACTER(LEN=20) :: name_x_func
1470  INTEGER :: hfx_potential_type, ifun, iounit, nfun
1471  LOGICAL :: funct_found
1472  REAL(dp) :: cutoff_radius, hfx_fraction, omega, &
1473  scale_coulomb, scale_longrange, scale_x
1474  TYPE(cp_logger_type), POINTER :: logger
1475  TYPE(section_vals_type), POINTER :: xc_fun, xc_fun_section
1476 
1477  logger => cp_get_default_logger()
1478  NULLIFY (admm_env%xc_section_aux, admm_env%xc_section_primary)
1479 
1480  !! ** Duplicate existing xc-section
1481  CALL section_vals_duplicate(xc_section, admm_env%xc_section_aux)
1482  CALL section_vals_duplicate(xc_section, admm_env%xc_section_primary)
1483  !** Now modify the auxiliary basis
1484  !** First remove all functionals
1485  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_aux, "XC_FUNCTIONAL")
1486 
1487  !* Overwrite possible shortcut
1488  CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", &
1489  i_val=xc_funct_no_shortcut)
1490 
1491  !** Get number of Functionals in the list
1492  ifun = 0
1493  nfun = 0
1494  DO
1495  ifun = ifun + 1
1496  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1497  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1498  nfun = nfun + 1
1499  END DO
1500 
1501  ifun = 0
1502  DO ifun = 1, nfun
1503  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=1)
1504  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1505  CALL section_vals_remove_values(xc_fun)
1506  END DO
1507 
1508  IF (ASSOCIATED(x_data)) THEN
1509  hfx_potential_type = x_data(1, 1)%potential_parameter%potential_type
1510  hfx_fraction = x_data(1, 1)%general_parameter%fraction
1511  ELSE
1512  cpwarn("ADMM requested without a DFT%XC%HF section. It will be ignored for the SCF.")
1513  admm_env%aux_exch_func = do_admm_aux_exch_func_none
1514  END IF
1515 
1516  !in case of no admm exchange corr., no auxiliary exchange functional needed
1517  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_none) THEN
1518  CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", &
1519  i_val=xc_none)
1520  hfx_fraction = 0.0_dp
1521  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_default) THEN
1522  ! default PBE Functional
1523  !! ** Add functionals evaluated with auxiliary basis
1524  SELECT CASE (hfx_potential_type)
1525  CASE (do_potential_coulomb)
1526  CALL section_vals_val_set(xc_fun_section, "PBE%_SECTION_PARAMETERS_", &
1527  l_val=.true.)
1528  CALL section_vals_val_set(xc_fun_section, "PBE%SCALE_X", &
1529  r_val=-hfx_fraction)
1530  CALL section_vals_val_set(xc_fun_section, "PBE%SCALE_C", &
1531  r_val=0.0_dp)
1532  CASE (do_potential_short)
1533  omega = x_data(1, 1)%potential_parameter%omega
1534  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1535  l_val=.true.)
1536  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1537  r_val=-hfx_fraction)
1538  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1539  r_val=0.0_dp)
1540  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1541  r_val=omega)
1542  CASE (do_potential_truncated)
1543  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
1544  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
1545  l_val=.true.)
1546  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1547  r_val=hfx_fraction)
1548  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1549  r_val=cutoff_radius)
1550  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1551  l_val=.true.)
1552  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1553  r_val=0.0_dp)
1554  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1555  r_val=-hfx_fraction)
1556  CASE (do_potential_long)
1557  omega = x_data(1, 1)%potential_parameter%omega
1558  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1559  l_val=.true.)
1560  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1561  r_val=hfx_fraction)
1562  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1563  r_val=-hfx_fraction)
1564  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1565  r_val=omega)
1566  CASE (do_potential_mix_cl)
1567  omega = x_data(1, 1)%potential_parameter%omega
1568  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
1569  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
1570  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1571  l_val=.true.)
1572  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1573  r_val=hfx_fraction*scale_longrange)
1574  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1575  r_val=-hfx_fraction*(scale_longrange + scale_coulomb))
1576  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1577  r_val=omega)
1578  CASE (do_potential_mix_cl_trunc)
1579  omega = x_data(1, 1)%potential_parameter%omega
1580  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
1581  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
1582  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
1583  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
1584  l_val=.true.)
1585  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1586  r_val=hfx_fraction*(scale_longrange + scale_coulomb))
1587  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1588  r_val=cutoff_radius)
1589  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1590  l_val=.true.)
1591  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1592  r_val=hfx_fraction*scale_longrange)
1593  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1594  r_val=-hfx_fraction*(scale_longrange + scale_coulomb))
1595  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1596  r_val=omega)
1597  CASE DEFAULT
1598  cpabort("Unknown potential operator!")
1599  END SELECT
1600 
1601  !** Now modify the functionals for the primary basis
1602  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_primary, "XC_FUNCTIONAL")
1603  !* Overwrite possible shortcut
1604  CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", &
1605  i_val=xc_funct_no_shortcut)
1606 
1607  SELECT CASE (hfx_potential_type)
1608  CASE (do_potential_coulomb)
1609  ifun = 0
1610  funct_found = .false.
1611  DO
1612  ifun = ifun + 1
1613  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1614  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1615  IF (xc_fun%section%name == "PBE") THEN
1616  funct_found = .true.
1617  END IF
1618  END DO
1619  IF (.NOT. funct_found) THEN
1620  CALL section_vals_val_set(xc_fun_section, "PBE%_SECTION_PARAMETERS_", &
1621  l_val=.true.)
1622  CALL section_vals_val_set(xc_fun_section, "PBE%SCALE_X", &
1623  r_val=hfx_fraction)
1624  CALL section_vals_val_set(xc_fun_section, "PBE%SCALE_C", &
1625  r_val=0.0_dp)
1626  ELSE
1627  CALL section_vals_val_get(xc_fun_section, "PBE%SCALE_X", &
1628  r_val=scale_x)
1629  scale_x = scale_x + hfx_fraction
1630  CALL section_vals_val_set(xc_fun_section, "PBE%SCALE_X", &
1631  r_val=scale_x)
1632  END IF
1633  CASE (do_potential_short)
1634  omega = x_data(1, 1)%potential_parameter%omega
1635  ifun = 0
1636  funct_found = .false.
1637  DO
1638  ifun = ifun + 1
1639  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1640  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1641  IF (xc_fun%section%name == "XWPBE") THEN
1642  funct_found = .true.
1643  END IF
1644  END DO
1645  IF (.NOT. funct_found) THEN
1646  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1647  l_val=.true.)
1648  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1649  r_val=hfx_fraction)
1650  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1651  r_val=0.0_dp)
1652  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1653  r_val=omega)
1654  ELSE
1655  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X", &
1656  r_val=scale_x)
1657  scale_x = scale_x + hfx_fraction
1658  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1659  r_val=scale_x)
1660  END IF
1661  CASE (do_potential_long)
1662  omega = x_data(1, 1)%potential_parameter%omega
1663  ifun = 0
1664  funct_found = .false.
1665  DO
1666  ifun = ifun + 1
1667  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1668  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1669  IF (xc_fun%section%name == "XWPBE") THEN
1670  funct_found = .true.
1671  END IF
1672  END DO
1673  IF (.NOT. funct_found) THEN
1674  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1675  l_val=.true.)
1676  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1677  r_val=-hfx_fraction)
1678  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1679  r_val=hfx_fraction)
1680  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1681  r_val=omega)
1682  ELSE
1683  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X", &
1684  r_val=scale_x)
1685  scale_x = scale_x - hfx_fraction
1686  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1687  r_val=scale_x)
1688  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X0", &
1689  r_val=scale_x)
1690  scale_x = scale_x + hfx_fraction
1691  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1692  r_val=scale_x)
1693 
1694  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1695  r_val=omega)
1696  END IF
1697  CASE (do_potential_truncated)
1698  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
1699  ifun = 0
1700  funct_found = .false.
1701  DO
1702  ifun = ifun + 1
1703  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1704  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1705  IF (xc_fun%section%name == "PBE_HOLE_T_C_LR") THEN
1706  funct_found = .true.
1707  END IF
1708  END DO
1709  IF (.NOT. funct_found) THEN
1710  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
1711  l_val=.true.)
1712  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1713  r_val=-hfx_fraction)
1714  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1715  r_val=cutoff_radius)
1716  ELSE
1717  CALL section_vals_val_get(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1718  r_val=scale_x)
1719  scale_x = scale_x - hfx_fraction
1720  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1721  r_val=scale_x)
1722  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1723  r_val=cutoff_radius)
1724  END IF
1725  ifun = 0
1726  funct_found = .false.
1727  DO
1728  ifun = ifun + 1
1729  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1730  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1731  IF (xc_fun%section%name == "XWPBE") THEN
1732  funct_found = .true.
1733  END IF
1734  END DO
1735  IF (.NOT. funct_found) THEN
1736  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1737  l_val=.true.)
1738  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1739  r_val=hfx_fraction)
1740  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1741  r_val=0.0_dp)
1742 
1743  ELSE
1744  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X0", &
1745  r_val=scale_x)
1746  scale_x = scale_x + hfx_fraction
1747  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1748  r_val=scale_x)
1749  END IF
1750  CASE (do_potential_mix_cl_trunc)
1751  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
1752  omega = x_data(1, 1)%potential_parameter%omega
1753  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
1754  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
1755  ifun = 0
1756  funct_found = .false.
1757  DO
1758  ifun = ifun + 1
1759  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1760  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1761  IF (xc_fun%section%name == "PBE_HOLE_T_C_LR") THEN
1762  funct_found = .true.
1763  END IF
1764  END DO
1765  IF (.NOT. funct_found) THEN
1766  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
1767  l_val=.true.)
1768  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1769  r_val=-hfx_fraction*(scale_coulomb + scale_longrange))
1770  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1771  r_val=cutoff_radius)
1772 
1773  ELSE
1774  CALL section_vals_val_get(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1775  r_val=scale_x)
1776  scale_x = scale_x - hfx_fraction*(scale_coulomb + scale_longrange)
1777  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1778  r_val=scale_x)
1779  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1780  r_val=cutoff_radius)
1781  END IF
1782  ifun = 0
1783  funct_found = .false.
1784  DO
1785  ifun = ifun + 1
1786  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1787  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1788  IF (xc_fun%section%name == "XWPBE") THEN
1789  funct_found = .true.
1790  END IF
1791  END DO
1792  IF (.NOT. funct_found) THEN
1793  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1794  l_val=.true.)
1795  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1796  r_val=hfx_fraction*(scale_coulomb + scale_longrange))
1797  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1798  r_val=-hfx_fraction*scale_longrange)
1799  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1800  r_val=omega)
1801 
1802  ELSE
1803  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X0", &
1804  r_val=scale_x)
1805  scale_x = scale_x + hfx_fraction*(scale_coulomb + scale_longrange)
1806  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1807  r_val=scale_x)
1808  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X", &
1809  r_val=scale_x)
1810  scale_x = scale_x - hfx_fraction*scale_longrange
1811  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1812  r_val=scale_x)
1813 
1814  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1815  r_val=omega)
1816  END IF
1817  CASE (do_potential_mix_cl)
1818  omega = x_data(1, 1)%potential_parameter%omega
1819  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
1820  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
1821  ifun = 0
1822  funct_found = .false.
1823  DO
1824  ifun = ifun + 1
1825  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1826  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1827  IF (xc_fun%section%name == "XWPBE") THEN
1828  funct_found = .true.
1829  END IF
1830  END DO
1831  IF (.NOT. funct_found) THEN
1832  CALL section_vals_val_set(xc_fun_section, "XWPBE%_SECTION_PARAMETERS_", &
1833  l_val=.true.)
1834  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1835  r_val=hfx_fraction*(scale_coulomb + scale_longrange))
1836  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1837  r_val=-hfx_fraction*scale_longrange)
1838  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1839  r_val=omega)
1840 
1841  ELSE
1842  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X0", &
1843  r_val=scale_x)
1844  scale_x = scale_x + hfx_fraction*(scale_coulomb + scale_longrange)
1845  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X0", &
1846  r_val=scale_x)
1847 
1848  CALL section_vals_val_get(xc_fun_section, "XWPBE%SCALE_X", &
1849  r_val=scale_x)
1850  scale_x = scale_x - hfx_fraction*scale_longrange
1851  CALL section_vals_val_set(xc_fun_section, "XWPBE%SCALE_X", &
1852  r_val=scale_x)
1853 
1854  CALL section_vals_val_set(xc_fun_section, "XWPBE%OMEGA", &
1855  r_val=omega)
1856  END IF
1857  END SELECT
1858  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_default_libxc) THEN
1859  ! default PBE Functional
1860  !! ** Add functionals evaluated with auxiliary basis
1861 #if defined (__LIBXC)
1862  SELECT CASE (hfx_potential_type)
1863  CASE (do_potential_coulomb)
1864  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
1865  l_val=.true.)
1866  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
1867  r_val=-hfx_fraction)
1868  CASE (do_potential_short)
1869  omega = x_data(1, 1)%potential_parameter%omega
1870  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
1871  l_val=.true.)
1872  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1873  r_val=-hfx_fraction)
1874  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
1875  r_val=omega)
1876  CASE (do_potential_truncated)
1877  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
1878  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
1879  l_val=.true.)
1880  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1881  r_val=hfx_fraction)
1882  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1883  r_val=cutoff_radius)
1884  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
1885  l_val=.true.)
1886  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1887  r_val=-hfx_fraction)
1888  CASE (do_potential_long)
1889  omega = x_data(1, 1)%potential_parameter%omega
1890  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
1891  l_val=.true.)
1892  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1893  r_val=hfx_fraction)
1894  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
1895  r_val=omega)
1896  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
1897  l_val=.true.)
1898  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
1899  r_val=-hfx_fraction)
1900  CASE (do_potential_mix_cl)
1901  omega = x_data(1, 1)%potential_parameter%omega
1902  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
1903  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
1904  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
1905  l_val=.true.)
1906  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1907  r_val=hfx_fraction*scale_longrange)
1908  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
1909  r_val=omega)
1910  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
1911  l_val=.true.)
1912  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
1913  r_val=-hfx_fraction*(scale_longrange + scale_coulomb))
1914  CASE (do_potential_mix_cl_trunc)
1915  omega = x_data(1, 1)%potential_parameter%omega
1916  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
1917  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
1918  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
1919  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
1920  l_val=.true.)
1921  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
1922  r_val=hfx_fraction*(scale_longrange + scale_coulomb))
1923  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
1924  r_val=cutoff_radius)
1925  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
1926  l_val=.true.)
1927  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1928  r_val=hfx_fraction*scale_longrange)
1929  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
1930  r_val=omega)
1931  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
1932  l_val=.true.)
1933  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
1934  r_val=-hfx_fraction*(scale_longrange + scale_coulomb))
1935  CASE DEFAULT
1936  cpabort("Unknown potential operator!")
1937  END SELECT
1938 
1939  !** Now modify the functionals for the primary basis
1940  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_primary, "XC_FUNCTIONAL")
1941  !* Overwrite possible shortcut
1942  CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", &
1943  i_val=xc_funct_no_shortcut)
1944 
1945  SELECT CASE (hfx_potential_type)
1946  CASE (do_potential_coulomb)
1947  ifun = 0
1948  funct_found = .false.
1949  DO
1950  ifun = ifun + 1
1951  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1952  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1953  IF (xc_fun%section%name == "GGA_X_PBE") THEN
1954  funct_found = .true.
1955  END IF
1956  END DO
1957  IF (.NOT. funct_found) THEN
1958  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
1959  l_val=.true.)
1960  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
1961  r_val=hfx_fraction)
1962  ELSE
1963  CALL section_vals_val_get(xc_fun_section, "GGA_X_PBE%SCALE", &
1964  r_val=scale_x)
1965  scale_x = scale_x + hfx_fraction
1966  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
1967  r_val=scale_x)
1968  END IF
1969  CASE (do_potential_short)
1970  omega = x_data(1, 1)%potential_parameter%omega
1971  ifun = 0
1972  funct_found = .false.
1973  DO
1974  ifun = ifun + 1
1975  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
1976  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
1977  IF (xc_fun%section%name == "GGA_X_WPBEH") THEN
1978  funct_found = .true.
1979  END IF
1980  END DO
1981  IF (.NOT. funct_found) THEN
1982  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
1983  l_val=.true.)
1984  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1985  r_val=hfx_fraction)
1986  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
1987  r_val=omega)
1988  ELSE
1989  CALL section_vals_val_get(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1990  r_val=scale_x)
1991  scale_x = scale_x + hfx_fraction
1992  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
1993  r_val=scale_x)
1994  END IF
1995  CASE (do_potential_long)
1996  omega = x_data(1, 1)%potential_parameter%omega
1997  ifun = 0
1998  funct_found = .false.
1999  DO
2000  ifun = ifun + 1
2001  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2002  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2003  IF (xc_fun%section%name == "GGA_X_WPBEH") THEN
2004  funct_found = .true.
2005  END IF
2006  END DO
2007  IF (.NOT. funct_found) THEN
2008  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
2009  l_val=.true.)
2010  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2011  r_val=-hfx_fraction)
2012  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
2013  r_val=omega)
2014  ELSE
2015  CALL section_vals_val_get(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2016  r_val=scale_x)
2017  scale_x = scale_x - hfx_fraction
2018  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2019  r_val=scale_x)
2020 
2021  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
2022  r_val=omega)
2023  END IF
2024  ifun = 0
2025  funct_found = .false.
2026  DO
2027  ifun = ifun + 1
2028  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2029  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2030  IF (xc_fun%section%name == "GGA_X_PBE") THEN
2031  funct_found = .true.
2032  END IF
2033  END DO
2034  IF (.NOT. funct_found) THEN
2035  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
2036  l_val=.true.)
2037  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2038  r_val=hfx_fraction)
2039  ELSE
2040  CALL section_vals_val_get(xc_fun_section, "GGA_X_PBE%SCALE", &
2041  r_val=scale_x)
2042  scale_x = scale_x + hfx_fraction
2043  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2044  r_val=scale_x)
2045  END IF
2046  CASE (do_potential_truncated)
2047  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
2048  ifun = 0
2049  funct_found = .false.
2050  DO
2051  ifun = ifun + 1
2052  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2053  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2054  IF (xc_fun%section%name == "PBE_HOLE_T_C_LR") THEN
2055  funct_found = .true.
2056  END IF
2057  END DO
2058  IF (.NOT. funct_found) THEN
2059  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
2060  l_val=.true.)
2061  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
2062  r_val=-hfx_fraction)
2063  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
2064  r_val=cutoff_radius)
2065 
2066  ELSE
2067  CALL section_vals_val_get(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
2068  r_val=scale_x)
2069  scale_x = scale_x - hfx_fraction
2070  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
2071  r_val=scale_x)
2072  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
2073  r_val=cutoff_radius)
2074  END IF
2075  ifun = 0
2076  funct_found = .false.
2077  DO
2078  ifun = ifun + 1
2079  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2080  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2081  IF (xc_fun%section%name == "GGA_X_PBE") THEN
2082  funct_found = .true.
2083  END IF
2084  END DO
2085  IF (.NOT. funct_found) THEN
2086  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
2087  l_val=.true.)
2088  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2089  r_val=hfx_fraction)
2090 
2091  ELSE
2092  CALL section_vals_val_get(xc_fun_section, "GGA_X_PBE%SCALE", &
2093  r_val=scale_x)
2094  scale_x = scale_x + hfx_fraction
2095  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2096  r_val=scale_x)
2097  END IF
2098  CASE (do_potential_mix_cl_trunc)
2099  cutoff_radius = x_data(1, 1)%potential_parameter%cutoff_radius
2100  omega = x_data(1, 1)%potential_parameter%omega
2101  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
2102  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
2103  ifun = 0
2104  funct_found = .false.
2105  DO
2106  ifun = ifun + 1
2107  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2108  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2109  IF (xc_fun%section%name == "PBE_HOLE_T_C_LR") THEN
2110  funct_found = .true.
2111  END IF
2112  END DO
2113  IF (.NOT. funct_found) THEN
2114  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%_SECTION_PARAMETERS_", &
2115  l_val=.true.)
2116  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
2117  r_val=-hfx_fraction*(scale_coulomb + scale_longrange))
2118  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
2119  r_val=cutoff_radius)
2120 
2121  ELSE
2122  CALL section_vals_val_get(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
2123  r_val=scale_x)
2124  scale_x = scale_x - hfx_fraction*(scale_coulomb + scale_longrange)
2125  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%SCALE_X", &
2126  r_val=scale_x)
2127  CALL section_vals_val_set(xc_fun_section, "PBE_HOLE_T_C_LR%CUTOFF_RADIUS", &
2128  r_val=cutoff_radius)
2129  END IF
2130  ifun = 0
2131  funct_found = .false.
2132  DO
2133  ifun = ifun + 1
2134  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2135  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2136  IF (xc_fun%section%name == "GGA_X_WPBEH") THEN
2137  funct_found = .true.
2138  END IF
2139  END DO
2140  IF (.NOT. funct_found) THEN
2141  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
2142  l_val=.true.)
2143  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2144  r_val=-hfx_fraction*scale_longrange)
2145  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
2146  r_val=omega)
2147 
2148  ELSE
2149  CALL section_vals_val_get(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2150  r_val=scale_x)
2151  scale_x = scale_x - hfx_fraction*scale_longrange
2152  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2153  r_val=scale_x)
2154 
2155  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
2156  r_val=omega)
2157  END IF
2158  ifun = 0
2159  funct_found = .false.
2160  DO
2161  ifun = ifun + 1
2162  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2163  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2164  IF (xc_fun%section%name == "GGA_X_PBE") THEN
2165  funct_found = .true.
2166  END IF
2167  END DO
2168  IF (.NOT. funct_found) THEN
2169  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
2170  l_val=.true.)
2171  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2172  r_val=hfx_fraction*(scale_coulomb + scale_longrange))
2173  ELSE
2174  CALL section_vals_val_get(xc_fun_section, "GGA_X_PBE%SCALE", &
2175  r_val=scale_x)
2176  scale_x = scale_x + hfx_fraction*(scale_coulomb + scale_longrange)
2177  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2178  r_val=scale_x)
2179  END IF
2180  CASE (do_potential_mix_cl)
2181  omega = x_data(1, 1)%potential_parameter%omega
2182  scale_coulomb = x_data(1, 1)%potential_parameter%scale_coulomb
2183  scale_longrange = x_data(1, 1)%potential_parameter%scale_longrange
2184  ifun = 0
2185  funct_found = .false.
2186  DO
2187  ifun = ifun + 1
2188  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2189  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2190  IF (xc_fun%section%name == "GGA_X_WPBEH") THEN
2191  funct_found = .true.
2192  END IF
2193  END DO
2194  IF (.NOT. funct_found) THEN
2195  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_SECTION_PARAMETERS_", &
2196  l_val=.true.)
2197  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2198  r_val=-hfx_fraction*scale_longrange)
2199  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
2200  r_val=omega)
2201 
2202  ELSE
2203  CALL section_vals_val_get(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2204  r_val=scale_x)
2205  scale_x = scale_x - hfx_fraction*scale_longrange
2206  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%SCALE", &
2207  r_val=scale_x)
2208 
2209  CALL section_vals_val_set(xc_fun_section, "GGA_X_WPBEH%_OMEGA", &
2210  r_val=omega)
2211  END IF
2212  ifun = 0
2213  funct_found = .false.
2214  DO
2215  ifun = ifun + 1
2216  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2217  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2218  IF (xc_fun%section%name == "GGA_X_PBE") THEN
2219  funct_found = .true.
2220  END IF
2221  END DO
2222  IF (.NOT. funct_found) THEN
2223  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%_SECTION_PARAMETERS_", &
2224  l_val=.true.)
2225  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2226  r_val=hfx_fraction*(scale_coulomb + scale_longrange))
2227  ELSE
2228  CALL section_vals_val_get(xc_fun_section, "GGA_X_PBE%SCALE", &
2229  r_val=scale_x)
2230  scale_x = scale_x + hfx_fraction*(scale_coulomb + scale_longrange)
2231  CALL section_vals_val_set(xc_fun_section, "GGA_X_PBE%SCALE", &
2232  r_val=scale_x)
2233  END IF
2234  END SELECT
2235 #else
2236  CALL cp_abort(__location__, "In order use a LibXC-based ADMM "// &
2237  "exchange correction functionals, you have to compile and link against LibXC!")
2238 #endif
2239 
2240  ! PBEX (always bare form), OPTX and Becke88 functional
2241  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_pbex .OR. &
2242  admm_env%aux_exch_func == do_admm_aux_exch_func_opt .OR. &
2243  admm_env%aux_exch_func == do_admm_aux_exch_func_bee) THEN
2244  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_pbex) THEN
2245  name_x_func = 'PBE'
2246  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt) THEN
2247  name_x_func = 'OPTX'
2248  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_bee) THEN
2249  name_x_func = 'BECKE88'
2250  END IF
2251  !primary basis
2252  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_SECTION_PARAMETERS_", &
2253  l_val=.true.)
2254  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE_X", &
2255  r_val=-hfx_fraction)
2256 
2257  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_pbex) THEN
2258  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE_C", r_val=0.0_dp)
2259  END IF
2260 
2261  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt) THEN
2262  IF (admm_env%aux_exch_func_param) THEN
2263  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%A1", &
2264  r_val=admm_env%aux_x_param(1))
2265  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%A2", &
2266  r_val=admm_env%aux_x_param(2))
2267  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%GAMMA", &
2268  r_val=admm_env%aux_x_param(3))
2269  END IF
2270  END IF
2271 
2272  !** Now modify the functionals for the primary basis
2273  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_primary, "XC_FUNCTIONAL")
2274  !* Overwrite possible L")
2275  !* Overwrite possible shortcut
2276  CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", &
2277  i_val=xc_funct_no_shortcut)
2278 
2279  ifun = 0
2280  funct_found = .false.
2281  DO
2282  ifun = ifun + 1
2283  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2284  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2285  IF (xc_fun%section%name == trim(name_x_func)) THEN
2286  funct_found = .true.
2287  END IF
2288  END DO
2289  IF (.NOT. funct_found) THEN
2290  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_SECTION_PARAMETERS_", &
2291  l_val=.true.)
2292  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE_X", &
2293  r_val=hfx_fraction)
2294  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_pbex) THEN
2295  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE_C", &
2296  r_val=0.0_dp)
2297  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt) THEN
2298  IF (admm_env%aux_exch_func_param) THEN
2299  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%A1", &
2300  r_val=admm_env%aux_x_param(1))
2301  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%A2", &
2302  r_val=admm_env%aux_x_param(2))
2303  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%GAMMA", &
2304  r_val=admm_env%aux_x_param(3))
2305  END IF
2306  END IF
2307 
2308  ELSE
2309  CALL section_vals_val_get(xc_fun_section, trim(name_x_func)//"%SCALE_X", &
2310  r_val=scale_x)
2311  scale_x = scale_x + hfx_fraction
2312  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE_X", &
2313  r_val=scale_x)
2314  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt) THEN
2315  cpassert(.NOT. admm_env%aux_exch_func_param)
2316  END IF
2317  END IF
2318 
2319  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_pbex_libxc .OR. &
2320  admm_env%aux_exch_func == do_admm_aux_exch_func_opt_libxc .OR. &
2321  admm_env%aux_exch_func == do_admm_aux_exch_func_sx_libxc .OR. &
2322  admm_env%aux_exch_func == do_admm_aux_exch_func_bee_libxc) THEN
2323 #if defined(__LIBXC)
2324  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_pbex_libxc) THEN
2325  name_x_func = 'GGA_X_PBE'
2326  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt_libxc) THEN
2327  name_x_func = 'GGA_X_OPTX'
2328  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_bee_libxc) THEN
2329  name_x_func = 'GGA_X_B88'
2330  ELSE IF (admm_env%aux_exch_func == do_admm_aux_exch_func_sx_libxc) THEN
2331  name_x_func = 'LDA_X'
2332  END IF
2333  !primary basis
2334  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_SECTION_PARAMETERS_", &
2335  l_val=.true.)
2336  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE", &
2337  r_val=-hfx_fraction)
2338 
2339  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt_libxc) THEN
2340  IF (admm_env%aux_exch_func_param) THEN
2341  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_A", &
2342  r_val=admm_env%aux_x_param(1))
2343  ! LibXC rescales the second parameter of the OPTX functional (see documentation there)
2344  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_B", &
2345  r_val=admm_env%aux_x_param(2)/x_factor_c)
2346  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_GAMMA", &
2347  r_val=admm_env%aux_x_param(3))
2348  END IF
2349  END IF
2350 
2351  !** Now modify the functionals for the primary basis
2352  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_primary, "XC_FUNCTIONAL")
2353  !* Overwrite possible L")
2354  !* Overwrite possible shortcut
2355  CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", &
2356  i_val=xc_funct_no_shortcut)
2357 
2358  ifun = 0
2359  funct_found = .false.
2360  DO
2361  ifun = ifun + 1
2362  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2363  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2364  IF (xc_fun%section%name == trim(name_x_func)) THEN
2365  funct_found = .true.
2366  END IF
2367  END DO
2368  IF (.NOT. funct_found) THEN
2369  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_SECTION_PARAMETERS_", &
2370  l_val=.true.)
2371  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE", &
2372  r_val=hfx_fraction)
2373  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt_libxc) THEN
2374  IF (admm_env%aux_exch_func_param) THEN
2375  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_A", &
2376  r_val=admm_env%aux_x_param(1))
2377  ! LibXC rescales the second parameter of the OPTX functional (see documentation there)
2378  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_B", &
2379  r_val=admm_env%aux_x_param(2)/x_factor_c)
2380  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%_GAMMA", &
2381  r_val=admm_env%aux_x_param(3))
2382  END IF
2383  END IF
2384 
2385  ELSE
2386  CALL section_vals_val_get(xc_fun_section, trim(name_x_func)//"%SCALE", &
2387  r_val=scale_x)
2388  scale_x = scale_x + hfx_fraction
2389  CALL section_vals_val_set(xc_fun_section, trim(name_x_func)//"%SCALE", &
2390  r_val=scale_x)
2391  IF (admm_env%aux_exch_func == do_admm_aux_exch_func_opt_libxc) THEN
2392  cpassert(.NOT. admm_env%aux_exch_func_param)
2393  END IF
2394  END IF
2395 #else
2396  CALL cp_abort(__location__, "In order use a LibXC-based ADMM "// &
2397  "exchange correction functionals, you have to compile and link against LibXC!")
2398 #endif
2399 
2400  ELSE
2401  cpabort("Unknown exchange correction functional!")
2402  END IF
2403 
2404  IF (debug_functional) THEN
2405  iounit = cp_logger_get_default_io_unit(logger)
2406  IF (iounit > 0) THEN
2407  WRITE (iounit, "(A)") " ADMM Primary Basis Set Functional"
2408  END IF
2409  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_primary, "XC_FUNCTIONAL")
2410  ifun = 0
2411  funct_found = .false.
2412  DO
2413  ifun = ifun + 1
2414  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2415  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2416 
2417  scale_x = -1000.0_dp
2418  IF (xc_fun%section%name /= "LYP" .AND. xc_fun%section%name /= "VWN") THEN
2419  CALL section_vals_val_get(xc_fun, "SCALE_X", r_val=scale_x)
2420  END IF
2421  IF (xc_fun%section%name == "XWPBE") THEN
2422  CALL section_vals_val_get(xc_fun, "SCALE_X0", r_val=hfx_fraction)
2423  IF (iounit > 0) THEN
2424  WRITE (iounit, "(T5,A,T25,2F10.3)") trim(xc_fun%section%name), scale_x, hfx_fraction
2425  END IF
2426  ELSE
2427  IF (iounit > 0) THEN
2428  WRITE (iounit, "(T5,A,T25,F10.3)") trim(xc_fun%section%name), scale_x
2429  END IF
2430  END IF
2431  END DO
2432 
2433  IF (iounit > 0) THEN
2434  WRITE (iounit, "(A)") " Auxiliary Basis Set Functional"
2435  END IF
2436  xc_fun_section => section_vals_get_subs_vals(admm_env%xc_section_aux, "XC_FUNCTIONAL")
2437  ifun = 0
2438  funct_found = .false.
2439  DO
2440  ifun = ifun + 1
2441  xc_fun => section_vals_get_subs_vals2(xc_fun_section, i_section=ifun)
2442  IF (.NOT. ASSOCIATED(xc_fun)) EXIT
2443  scale_x = -1000.0_dp
2444  IF (xc_fun%section%name /= "LYP" .AND. xc_fun%section%name /= "VWN") THEN
2445  CALL section_vals_val_get(xc_fun, "SCALE_X", r_val=scale_x)
2446  END IF
2447  IF (xc_fun%section%name == "XWPBE") THEN
2448  CALL section_vals_val_get(xc_fun, "SCALE_X0", r_val=hfx_fraction)
2449  IF (iounit > 0) THEN
2450  WRITE (iounit, "(T5,A,T25,2F10.3)") trim(xc_fun%section%name), scale_x, hfx_fraction
2451  END IF
2452  ELSE
2453  IF (iounit > 0) THEN
2454  WRITE (iounit, "(T5,A,T25,F10.3)") trim(xc_fun%section%name), scale_x
2455  END IF
2456  END IF
2457  END DO
2458  END IF
2459 
2460  END SUBROUTINE create_admm_xc_section
2461 
2462 ! **************************************************************************************************
2463 !> \brief Add the hfx contributions to the Hamiltonian
2464 !>
2465 !> \param matrix_ks Kohn-Sham matrix (updated on exit)
2466 !> \param rho_ao electron density expressed in terms of atomic orbitals
2467 !> \param qs_env Quickstep environment
2468 !> \param update_energy whether to update energy (default: yes)
2469 !> \param recalc_integrals whether to recalculate integrals (default: value of HF%TREAT_LSD_IN_CORE)
2470 !> \param external_hfx_sections ...
2471 !> \param external_x_data ...
2472 !> \note
2473 !> Simplified version of subroutine hfx_ks_matrix()
2474 ! **************************************************************************************************
2475  SUBROUTINE tddft_hfx_matrix(matrix_ks, rho_ao, qs_env, update_energy, recalc_integrals, &
2476  external_hfx_sections, external_x_data)
2477  TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
2478  TARGET :: matrix_ks, rho_ao
2479  TYPE(qs_environment_type), POINTER :: qs_env
2480  LOGICAL, INTENT(IN), OPTIONAL :: update_energy, recalc_integrals
2481  TYPE(section_vals_type), OPTIONAL, POINTER :: external_hfx_sections
2482  TYPE(hfx_type), DIMENSION(:, :), OPTIONAL, TARGET :: external_x_data
2483 
2484  CHARACTER(LEN=*), PARAMETER :: routinen = 'tddft_hfx_matrix'
2485 
2486  INTEGER :: handle, irep, ispin, mspin, n_rep_hf, &
2487  nspins
2488  LOGICAL :: distribute_fock_matrix, &
2489  hfx_treat_lsd_in_core, &
2490  my_update_energy, s_mstruct_changed
2491  REAL(kind=dp) :: eh1, ehfx
2492  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_ks_kp, rho_ao_kp
2493  TYPE(dft_control_type), POINTER :: dft_control
2494  TYPE(hfx_type), DIMENSION(:, :), POINTER :: x_data
2495  TYPE(mp_para_env_type), POINTER :: para_env
2496  TYPE(qs_energy_type), POINTER :: energy
2497  TYPE(section_vals_type), POINTER :: hfx_sections, input
2498 
2499  CALL timeset(routinen, handle)
2500 
2501  NULLIFY (dft_control, hfx_sections, input, para_env, matrix_ks_kp, rho_ao_kp)
2502 
2503  CALL get_qs_env(qs_env=qs_env, &
2504  dft_control=dft_control, &
2505  energy=energy, &
2506  input=input, &
2507  para_env=para_env, &
2508  s_mstruct_changed=s_mstruct_changed, &
2509  x_data=x_data)
2510 
2511  ! This should probably be the HF section from the TDDFPT XC section!
2512  hfx_sections => section_vals_get_subs_vals(input, "DFT%XC%HF")
2513 
2514  IF (PRESENT(external_hfx_sections)) hfx_sections => external_hfx_sections
2515  IF (PRESENT(external_x_data)) x_data => external_x_data
2516 
2517  my_update_energy = .true.
2518  IF (PRESENT(update_energy)) my_update_energy = update_energy
2519 
2520  IF (PRESENT(recalc_integrals)) s_mstruct_changed = recalc_integrals
2521 
2522  cpassert(dft_control%nimages == 1)
2523  nspins = dft_control%nspins
2524 
2525  CALL section_vals_get(hfx_sections, n_repetition=n_rep_hf)
2526  CALL section_vals_val_get(hfx_sections, "TREAT_LSD_IN_CORE", l_val=hfx_treat_lsd_in_core, &
2527  i_rep_section=1)
2528 
2529  CALL section_vals_get(hfx_sections, n_repetition=n_rep_hf)
2530  distribute_fock_matrix = .true.
2531 
2532  mspin = 1
2533  IF (hfx_treat_lsd_in_core) mspin = nspins
2534 
2535  matrix_ks_kp(1:nspins, 1:1) => matrix_ks(1:nspins)
2536  rho_ao_kp(1:nspins, 1:1) => rho_ao(1:nspins)
2537 
2538  DO irep = 1, n_rep_hf
2539  ! the real hfx calulation
2540  ehfx = 0.0_dp
2541 
2542  IF (x_data(irep, 1)%do_hfx_ri) THEN
2543  CALL hfx_ri_update_ks(qs_env, x_data(irep, 1)%ri_data, matrix_ks_kp, ehfx, &
2544  rho_ao=rho_ao_kp, geometry_did_change=s_mstruct_changed, &
2545  nspins=nspins, hf_fraction=x_data(irep, 1)%general_parameter%fraction)
2546 
2547  ELSE
2548  DO ispin = 1, mspin
2549  CALL integrate_four_center(qs_env, x_data, matrix_ks_kp, eh1, rho_ao_kp, hfx_sections, para_env, &
2550  s_mstruct_changed, irep, distribute_fock_matrix, ispin=ispin)
2551  ehfx = ehfx + eh1
2552  END DO
2553  END IF
2554  END DO
2555  IF (my_update_energy) energy%ex = ehfx
2556 
2557  CALL timestop(handle)
2558  END SUBROUTINE tddft_hfx_matrix
2559 
2560 END MODULE hfx_admm_utils
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition: cp_dbcsr_operations.F:81
admm_dm_types
Types and set/get functions for auxiliary density matrix methods.
Definition: admm_dm_types.F:14
admm_dm_types::admm_dm_create
subroutine, public admm_dm_create(admm_dm, admm_control, nspins, natoms)
Create a new admm_dm type.
Definition: admm_dm_types.F:63
admm_methods
Contains ADMM methods which require molecular orbitals.
Definition: admm_methods.F:15
admm_methods::scale_dm
subroutine, public scale_dm(qs_env, rho_ao_orb, scale_back)
Scale density matrix by gsi(ispin), is needed for force scaling in ADMMP.
Definition: admm_methods.F:2532
admm_methods::kpoint_calc_admm_matrices
subroutine, public kpoint_calc_admm_matrices(qs_env, calculate_forces)
Fill the ADMM overlp and basis change matrices in the KP env based on the real-space array.
Definition: admm_methods.F:3859
admm_types
Types and set/get functions for auxiliary density matrix methods.
Definition: admm_types.F:15
admm_types::get_admm_env
subroutine, public get_admm_env(admm_env, mo_derivs_aux_fit, mos_aux_fit, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, matrix_s_aux_fit, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp, task_list_aux_fit, matrix_ks_aux_fit, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp, rho_aux_fit, rho_aux_fit_buffer, admm_dm)
Get routine for the ADMM env.
Definition: admm_types.F:593
admm_types::set_admm_env
subroutine, public set_admm_env(admm_env, mo_derivs_aux_fit, mos_aux_fit, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, matrix_s_aux_fit, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp, task_list_aux_fit, matrix_ks_aux_fit, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp, rho_aux_fit, rho_aux_fit_buffer, admm_dm)
Set routine for the ADMM env.
Definition: admm_types.F:677
admm_types::admm_env_create
subroutine, public admm_env_create(admm_env, admm_control, mos, para_env, natoms, nao_aux_fit, blacs_env_ext)
creates ADMM environment, initializes the basic types
Definition: admm_types.F:220
atomic_kind_types
Define the atomic kind types and their sub types.
Definition: atomic_kind_types.F:23
basis_set_container_types
Definition: basis_set_container_types.F:13
basis_set_container_types::add_basis_set_to_container
subroutine, public add_basis_set_to_container(container, basis_set, basis_set_type)
...
Definition: basis_set_container_types.F:149
basis_set_types
Definition: basis_set_types.F:15
basis_set_types::copy_gto_basis_set
subroutine, public copy_gto_basis_set(basis_set_in, basis_set_out)
...
Definition: basis_set_types.F:212
basis_set_types::get_gto_basis_set
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius)
...
Definition: basis_set_types.F:615
cell_types
Handles all functions related to the CELL.
Definition: cell_types.F:15
cell_types::plane_distance
real(kind=dp) function, public plane_distance(h, k, l, cell)
Calculate the distance between two lattice planes as defined by a triple of Miller indices (hkl).
Definition: cell_types.F:252
cp_blacs_env
methods related to the blacs parallel environment
Definition: cp_blacs_env.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_m_by_n_from_row_template
subroutine, public cp_dbcsr_m_by_n_from_row_template(matrix, template, n, sym, data_type)
Utility function to create dbcsr matrix, m x n matrix (n arbitrary) with the same processor grid and ...
Definition: cp_dbcsr_operations.F:1156
cp_fm_pool_types
pool for for elements that are retained and released
Definition: cp_fm_pool_types.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_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
distribution_1d_types
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
Definition: distribution_1d_types.F:24
distribution_2d_types
stores a mapping of 2D info (e.g. matrix) on a 2D processor distribution (i.e. blacs grid) where cpus...
Definition: distribution_2d_types.F:15
external_potential_types
Definition of the atomic potential types.
Definition: external_potential_types.F:14
hfx_admm_utils
Utilities for hfx and admm methods.
Definition: hfx_admm_utils.F:17
hfx_admm_utils::hfx_ks_matrix
subroutine, public hfx_ks_matrix(qs_env, matrix_ks, rho, energy, calculate_forces, just_energy, v_rspace_new, v_tau_rspace)
Add the hfx contributions to the Hamiltonian.
Definition: hfx_admm_utils.F:1136
hfx_admm_utils::hfx_admm_init
subroutine, public hfx_admm_init(qs_env, calculate_forces)
...
Definition: hfx_admm_utils.F:164
hfx_admm_utils::aux_admm_init
subroutine, public aux_admm_init(qs_env, mos, admm_env, admm_control, basis_type)
Minimal setup routine for admm_env No forces No k-points No DFT correction terms.
Definition: hfx_admm_utils.F:397
hfx_admm_utils::create_admm_xc_section
subroutine, public create_admm_xc_section(x_data, xc_section, admm_env)
This routine modifies the xc section depending on the potential type used for the HF exchange and the...
Definition: hfx_admm_utils.F:1460
hfx_admm_utils::hfx_ks_matrix_kp
subroutine, public hfx_ks_matrix_kp(qs_env, matrix_ks, energy, calculate_forces)
Add the HFX K-point contribution to the real-space Hamiltonians.
Definition: hfx_admm_utils.F:959
hfx_admm_utils::tddft_hfx_matrix
subroutine, public tddft_hfx_matrix(matrix_ks, rho_ao, qs_env, update_energy, recalc_integrals, external_hfx_sections, external_x_data)
Add the hfx contributions to the Hamiltonian.
Definition: hfx_admm_utils.F:2477
hfx_derivatives
Routines to calculate derivatives with respect to basis function origin.
Definition: hfx_derivatives.F:14
hfx_derivatives::derivatives_four_center
subroutine, public derivatives_four_center(qs_env, rho_ao, rho_ao_resp, hfx_section, para_env, irep, use_virial, adiabatic_rescale_factor, resp_only, external_x_data)
computes four center derivatives for a full basis set and updates the forcesfock_4c arrays....
Definition: hfx_derivatives.F:117
hfx_energy_potential
Routines to calculate HFX energy and potential.
Definition: hfx_energy_potential.F:14
hfx_energy_potential::integrate_four_center
subroutine, public integrate_four_center(qs_env, x_data, ks_matrix, ehfx, rho_ao, hfx_section, para_env, geometry_did_change, irep, distribute_fock_matrix, ispin)
computes four center integrals for a full basis set and updates the Kohn-Sham-Matrix and energy....
Definition: hfx_energy_potential.F:137
hfx_pw_methods
Test routines for HFX caclulations using PW.
Definition: hfx_pw_methods.F:18
hfx_pw_methods::pw_hfx
subroutine, public pw_hfx(qs_env, ehfx, hfx_section, poisson_env, auxbas_pw_pool, irep)
computes the Hartree-Fock energy brute force in a pw basis
Definition: hfx_pw_methods.F:84
hfx_ri_kp
RI-methods for HFX and K-points. \auhtor Augustin Bussy (01.2023)
Definition: hfx_ri_kp.F:13
hfx_ri_kp::hfx_ri_update_forces_kp
subroutine, public hfx_ri_update_forces_kp(qs_env, ri_data, nspins, hf_fraction, rho_ao, use_virial)
Update the K-points RI-HFX forces.
Definition: hfx_ri_kp.F:771
hfx_ri_kp::hfx_ri_update_ks_kp
subroutine, public hfx_ri_update_ks_kp(qs_env, ri_data, ks_matrix, ehfx, rho_ao, geometry_did_change, nspins, hf_fraction)
Update the KS matrices for each real-space image.
Definition: hfx_ri_kp.F:420
hfx_ri
RI-methods for HFX.
Definition: hfx_ri.F:12
hfx_ri::hfx_ri_update_ks
subroutine, public hfx_ri_update_ks(qs_env, ri_data, ks_matrix, ehfx, mos, rho_ao, geometry_did_change, nspins, hf_fraction)
...
Definition: hfx_ri.F:1033
hfx_ri::hfx_ri_update_forces
subroutine, public hfx_ri_update_forces(qs_env, ri_data, nspins, hf_fraction, rho_ao, rho_ao_resp, mos, use_virial, resp_only, rescale_factor)
the general routine that calls the relevant force code
Definition: hfx_ri.F:3036
hfx_types
Types and set/get functions for HFX.
Definition: hfx_types.F:15
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition: input_constants.F:17
input_constants::do_admm_aux_exch_func_opt_libxc
integer, parameter, public do_admm_aux_exch_func_opt_libxc
Definition: input_constants.F:863
input_constants::do_admm_purify_none
integer, parameter, public do_admm_purify_none
Definition: input_constants.F:848
input_constants::xc_funct_no_shortcut
integer, parameter, public xc_funct_no_shortcut
Definition: input_constants.F:551
input_constants::do_admm_aux_exch_func_sx_libxc
integer, parameter, public do_admm_aux_exch_func_sx_libxc
Definition: input_constants.F:863
input_constants::do_admm_aux_exch_func_bee
integer, parameter, public do_admm_aux_exch_func_bee
Definition: input_constants.F:863
input_constants::do_potential_mix_cl
integer, parameter, public do_potential_mix_cl
Definition: input_constants.F:801
input_constants::do_admm_basis_projection
integer, parameter, public do_admm_basis_projection
Definition: input_constants.F:857
input_constants::do_admm_aux_exch_func_default_libxc
integer, parameter, public do_admm_aux_exch_func_default_libxc
Definition: input_constants.F:863
input_constants::do_admm_aux_exch_func_opt
integer, parameter, public do_admm_aux_exch_func_opt
Definition: input_constants.F:863
input_constants::do_admm_aux_exch_func_none
integer, parameter, public do_admm_aux_exch_func_none
Definition: input_constants.F:863
input_constants::do_admm_aux_exch_func_bee_libxc
integer, parameter, public do_admm_aux_exch_func_bee_libxc
Definition: input_constants.F:863
input_constants::do_admm_aux_exch_func_pbex_libxc
integer, parameter, public do_admm_aux_exch_func_pbex_libxc
Definition: input_constants.F:863
input_constants::do_admm_aux_exch_func_default
integer, parameter, public do_admm_aux_exch_func_default
Definition: input_constants.F:863
input_constants::do_potential_truncated
integer, parameter, public do_potential_truncated
Definition: input_constants.F:801
input_constants::do_admm_charge_constrained_projection
integer, parameter, public do_admm_charge_constrained_projection
Definition: input_constants.F:857
input_constants::do_potential_id
integer, parameter, public do_potential_id
Definition: input_constants.F:801
input_constants::do_potential_coulomb
integer, parameter, public do_potential_coulomb
Definition: input_constants.F:801
input_constants::do_potential_short
integer, parameter, public do_potential_short
Definition: input_constants.F:801
input_constants::do_potential_mix_cl_trunc
integer, parameter, public do_potential_mix_cl_trunc
Definition: input_constants.F:801
input_constants::xc_none
integer, parameter, public xc_none
Definition: input_constants.F:551
input_constants::do_potential_long
integer, parameter, public do_potential_long
Definition: input_constants.F:801
input_constants::do_admm_aux_exch_func_pbex
integer, parameter, public do_admm_aux_exch_func_pbex
Definition: input_constants.F:863
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_val_set
subroutine, public section_vals_val_set(section_vals, keyword_name, i_rep_section, i_rep_val, val, l_val, i_val, r_val, c_val, l_vals_ptr, i_vals_ptr, r_vals_ptr, c_vals_ptr)
sets the requested value
Definition: input_section_types.F:1216
input_section_types::section_vals_get_subs_vals2
type(section_vals_type) function, pointer, public section_vals_get_subs_vals2(section_vals, i_section, i_rep_section)
returns the values of the n-th non default subsection (null if no such section exists (not so many no...
Definition: input_section_types.F:779
input_section_types::section_vals_remove_values
subroutine, public section_vals_remove_values(section_vals)
removes the values of a repetition of the section
Definition: input_section_types.F:889
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_get
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
Definition: input_section_types.F:697
input_section_types::section_vals_duplicate
subroutine, public section_vals_duplicate(section_vals_in, section_vals_out, i_rep_start, i_rep_end)
creates a deep copy from section_vals_in to section_vals_out
Definition: input_section_types.F:1732
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
kpoint_methods
Routines needed for kpoint calculation.
Definition: kpoint_methods.F:15
kpoint_methods::kpoint_initialize_mos
subroutine, public kpoint_initialize_mos(kpoint, mos, added_mos, for_aux_fit)
Initialize a set of MOs and density matrix for each kpoint (kpoint group)
Definition: kpoint_methods.F:454
kpoint_transitional
Datatype to translate between k-points (2d) and gamma-point (1d) code.
Definition: kpoint_transitional.F:19
kpoint_transitional::kpoint_transitional_release
subroutine, public kpoint_transitional_release(this)
Release the matrix set, using the right pointer.
Definition: kpoint_transitional.F:119
kpoint_transitional::set_2d_pointer
subroutine, public set_2d_pointer(this, ptr_2d)
Assigns a 2D pointer.
Definition: kpoint_transitional.F:100
kpoint_types
Types and basic routines needed for a kpoint calculation.
Definition: kpoint_types.F:15
kpoint_types::get_kpoint_info
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym)
Retrieve information from a kpoint environment.
Definition: kpoint_types.F:333
libint_2c_3c
2- and 3-center electron repulsion integral routines based on libint2 Currently available operators: ...
Definition: libint_2c_3c.F:14
libint_2c_3c::cutoff_screen_factor
real(kind=dp), parameter, public cutoff_screen_factor
Definition: libint_2c_3c.F:48
mathlib
Collection of simple mathematical functions and subroutines.
Definition: mathlib.F:15
mathlib::erfc_cutoff
subroutine, public erfc_cutoff(eps, omg, r_cutoff)
compute a truncation radius for the shortrange operator
Definition: mathlib.F:1689
message_passing
Interface to the message passing library MPI.
Definition: message_passing.F:23
molecule_types
Define the data structure for the molecule information.
Definition: molecule_types.F:16
particle_types
Define the data structure for the particle information.
Definition: particle_types.F:19
paw_proj_set_types
Definition: paw_proj_set_types.F:13
paw_proj_set_types::get_paw_proj_set
subroutine, public get_paw_proj_set(paw_proj_set, csprj, chprj, first_prj, first_prjs, last_prj, local_oce_sphi_h, local_oce_sphi_s, maxl, ncgauprj, nsgauprj, nsatbas, nsotot, nprj, o2nindex, n2oindex, rcprj, rzetprj, zisomin, zetprj)
Get informations about a paw projectors set.
Definition: paw_proj_set_types.F:744
pw_env_types
container for various plainwaves related things
Definition: pw_env_types.F:14
pw_env_types::pw_env_get
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Definition: pw_env_types.F:113
pw_poisson_types
functions related to the poisson solver on regular grids
Definition: pw_poisson_types.F:22
pw_pool_types
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition: pw_pool_types.F:24
pw_types
Definition: pw_types.F:24
qs_energy_types
Definition: qs_energy_types.F:14
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_interactions
Calculate the interaction radii for the operator matrix calculation.
Definition: qs_interactions.F:17
qs_interactions::init_interaction_radii
subroutine, public init_interaction_radii(qs_control, qs_kind_set)
Initialize all the atomic kind radii for a given threshold value.
Definition: qs_interactions.F:75
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, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, 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_r3d_rs_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, 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, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition: qs_kind_types.F:451
qs_kind_types::get_qs_kind_set
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, basis_type, total_zeff_corr)
Get attributes of an atomic kind set.
Definition: qs_kind_types.F:864
qs_kind_types::init_gapw_nlcc
subroutine, public init_gapw_nlcc(qs_kind_set)
...
Definition: qs_kind_types.F:1393
qs_kind_types::init_gapw_basis_set
subroutine, public init_gapw_basis_set(qs_kind_set, qs_control, force_env_section, modify_qs_control)
...
Definition: qs_kind_types.F:1299
qs_ks_types
Definition: qs_ks_types.F:15
qs_local_rho_types
Definition: qs_local_rho_types.F:8
qs_local_rho_types::local_rho_set_create
subroutine, public local_rho_set_create(local_rho_set)
...
Definition: qs_local_rho_types.F:194
qs_matrix_pools
wrapper for the pools of matrixes
Definition: qs_matrix_pools.F:14
qs_matrix_pools::mpools_get
subroutine, public mpools_get(mpools, ao_mo_fm_pools, ao_ao_fm_pools, mo_mo_fm_pools, ao_mosub_fm_pools, mosub_mosub_fm_pools, maxao_maxmo_fm_pool, maxao_maxao_fm_pool, maxmo_maxmo_fm_pool)
returns various attributes of the mpools (notably the pools contained in it)
Definition: qs_matrix_pools.F:141
qs_mo_types
Definition and initialisation of the mo data type.
Definition: qs_mo_types.F:22
qs_mo_types::allocate_mo_set
subroutine, public allocate_mo_set(mo_set, nao, nmo, nelectron, n_el_f, maxocc, flexible_electron_count)
Allocates a mo set and partially initializes it (nao,nmo,nelectron, and flexible_electron_count are v...
Definition: qs_mo_types.F:206
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_mo_types::init_mo_set
subroutine, public init_mo_set(mo_set, fm_pool, fm_ref, fm_struct, name)
initializes an allocated mo_set. eigenvalues, mo_coeff, occupation_numbers are valid only after this ...
Definition: qs_mo_types.F:245
qs_neighbor_list_types
Define the neighbor list data types and the corresponding functionality.
Definition: qs_neighbor_list_types.F:17
qs_neighbor_list_types::release_neighbor_list_sets
subroutine, public release_neighbor_list_sets(nlists)
releases an array of neighbor_list_sets
Definition: qs_neighbor_list_types.F:1097
qs_neighbor_lists
Generate the atomic neighbor lists.
Definition: qs_neighbor_lists.F:19
qs_neighbor_lists::atom2d_cleanup
subroutine, public atom2d_cleanup(atom2d)
free the internals of atom2d
Definition: qs_neighbor_lists.F:131
qs_neighbor_lists::pair_radius_setup
subroutine, public pair_radius_setup(present_a, present_b, radius_a, radius_b, pair_radius, prmin)
...
Definition: qs_neighbor_lists.F:1618
qs_neighbor_lists::build_neighbor_lists
subroutine, public build_neighbor_lists(ab_list, particle_set, atom, cell, pair_radius, subcells, mic, symmetric, molecular, subset_of_mol, current_subset, operator_type, nlname, atomb_to_keep)
Build simple pair neighbor lists.
Definition: qs_neighbor_lists.F:1020
qs_neighbor_lists::write_neighbor_lists
subroutine, public write_neighbor_lists(ab, particle_set, cell, para_env, neighbor_list_section, nl_type, middle_name, nlname)
Write a set of neighbor lists to the output unit.
Definition: qs_neighbor_lists.F:1769
qs_neighbor_lists::atom2d_build
subroutine, public atom2d_build(atom2d, distribution_1d, distribution_2d, atomic_kind_set, molecule_set, molecule_only, particle_set)
Build some distribution structure of atoms, refactored from build_qs_neighbor_lists.
Definition: qs_neighbor_lists.F:173
qs_oce_methods
Routines for the construction of the coefficients for the expansion of the atomic densities rho1_hard...
Definition: qs_oce_methods.F:16
qs_oce_methods::build_oce_matrices
subroutine, public build_oce_matrices(intac, calculate_forces, nder, qs_kind_set, particle_set, sap_oce, eps_fit)
Set up the sparse matrix for the coefficients of one center expansions This routine uses the same log...
Definition: qs_oce_methods.F:506
qs_oce_types
Definition: qs_oce_types.F:8
qs_oce_types::allocate_oce_set
subroutine, public allocate_oce_set(oce_set, nkind)
Allocate and initialize the matrix set of oce coefficients.
Definition: qs_oce_types.F:60
qs_oce_types::create_oce_set
subroutine, public create_oce_set(oce_set)
...
Definition: qs_oce_types.F:79
qs_overlap
Calculation of overlap matrix, its derivatives and forces.
Definition: qs_overlap.F:19
qs_overlap::build_overlap_matrix
subroutine, public build_overlap_matrix(ks_env, matrix_s, matrixkp_s, matrix_name, nderivative, basis_type_a, basis_type_b, sab_nl, calculate_forces, matrix_p, matrixkp_p)
Calculation of the overlap matrix over Cartesian Gaussian functions.
Definition: qs_overlap.F:120
qs_rho_atom_methods
Definition: qs_rho_atom_methods.F:7
qs_rho_atom_methods::init_rho_atom
subroutine, public init_rho_atom(rho_atom_set, atomic_kind_set, qs_kind_set, dft_control, para_env)
...
Definition: qs_rho_atom_methods.F:754
qs_rho_methods
methods of the rho structure (defined in qs_rho_types)
Definition: qs_rho_methods.F:15
qs_rho_methods::qs_rho_rebuild
subroutine, public qs_rho_rebuild(rho, qs_env, rebuild_ao, rebuild_grids, admm, pw_env_external)
rebuilds rho (if necessary allocating and initializing it)
Definition: qs_rho_methods.F:96
qs_rho_types
superstucture that hold various representations of the density and keeps track of which ones are vali...
Definition: qs_rho_types.F:18
qs_rho_types::qs_rho_get
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
Definition: qs_rho_types.F:229
qs_rho_types::qs_rho_create
subroutine, public qs_rho_create(rho)
Allocates a new instance of rho.
Definition: qs_rho_types.F:99
rt_propagation_types
Types and set_get for real time propagation depending on runtype and diagonalization method different...
Definition: rt_propagation_types.F:21
task_list_methods
generate the tasks lists used by collocate and integrate routines
Definition: task_list_methods.F:14
task_list_methods::generate_qs_task_list
subroutine, public generate_qs_task_list(ks_env, task_list, reorder_rs_grid_ranks, skip_load_balance_distributed, soft_valid, basis_type, pw_env_external, sab_orb_external)
...
Definition: task_list_methods.F:119
task_list_types
types for task lists
Definition: task_list_types.F:14
task_list_types::deallocate_task_list
subroutine, public deallocate_task_list(task_list)
deallocates the components and the object itself
Definition: task_list_types.F:145
task_list_types::allocate_task_list
subroutine, public allocate_task_list(task_list)
allocates and initialised the components of the task_list_type
Definition: task_list_types.F:97
virial_types
Definition: virial_types.F:13
xc_adiabatic_utils
Definition: xc_adiabatic_utils.F:16
xc_adiabatic_utils::rescale_xc_potential
subroutine, public rescale_xc_potential(qs_env, ks_matrix, rho, energy, v_rspace_new, v_tau_rspace, hf_energy, just_energy, calculate_forces, use_virial)
Definition: xc_adiabatic_utils.F:70