(git:51bfb72)
Loading...
Searching...
No Matches
mp2_ri_2c.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Framework for 2c-integrals for RI
10!> \par History
11!> 06.2012 created [Mauro Del Ben]
12!> 03.2019 separated from mp2_ri_gpw [Frederick Stein]
13! **************************************************************************************************
14MODULE mp2_ri_2c
15 USE atomic_kind_types, ONLY: atomic_kind_type,&
16 get_atomic_kind_set
17 USE basis_set_types, ONLY: gto_basis_set_p_type,&
18 gto_basis_set_type
19 USE cell_types, ONLY: cell_type,&
20 get_cell,&
21 plane_distance
22 USE constants_operator, ONLY: operator_coulomb
23 USE cp_blacs_env, ONLY: cp_blacs_env_create,&
24 cp_blacs_env_release,&
25 cp_blacs_env_type
26 USE cp_cfm_basic_linalg, ONLY: cp_cfm_scale_and_add_fm
27 USE cp_cfm_types, ONLY: cp_cfm_create,&
28 cp_cfm_release,&
29 cp_cfm_to_fm,&
30 cp_cfm_type
31 USE cp_control_types, ONLY: dft_control_type
32 USE cp_dbcsr_api, ONLY: &
33 dbcsr_copy, dbcsr_create, dbcsr_deallocate_matrix, dbcsr_desymmetrize, &
34 dbcsr_distribution_release, dbcsr_distribution_type, dbcsr_p_type, dbcsr_release, &
35 dbcsr_set, dbcsr_type, dbcsr_type_antisymmetric, dbcsr_type_no_symmetry, &
36 dbcsr_type_symmetric
37 USE cp_dbcsr_contrib, ONLY: dbcsr_reserve_all_blocks
38 USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
39 USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
40 cp_dbcsr_dist2d_to_dist,&
41 dbcsr_allocate_matrix_set,&
42 dbcsr_deallocate_matrix_set
43 USE cp_eri_mme_interface, ONLY: cp_eri_mme_param
44 USE cp_fm_basic_linalg, ONLY: cp_fm_column_scale,&
45 cp_fm_triangular_invert
46 USE cp_fm_cholesky, ONLY: cp_fm_cholesky_decompose
47 USE cp_fm_diag, ONLY: choose_eigv_solver,&
48 cp_fm_power,&
49 cp_fm_syevx
50 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
51 cp_fm_struct_release,&
52 cp_fm_struct_type
53 USE cp_fm_types, ONLY: cp_fm_copy_general,&
54 cp_fm_create,&
55 cp_fm_get_info,&
56 cp_fm_release,&
57 cp_fm_set_all,&
58 cp_fm_to_fm,&
59 cp_fm_type
60 USE distribution_2d_types, ONLY: distribution_2d_type
61 USE group_dist_types, ONLY: create_group_dist,&
62 get_group_dist,&
63 group_dist_d1_type,&
64 release_group_dist
65 USE input_constants, ONLY: do_eri_gpw,&
66 do_eri_mme,&
67 do_eri_os,&
68 do_potential_coulomb,&
69 do_potential_id,&
70 do_potential_long,&
71 do_potential_truncated
72 USE kinds, ONLY: dp
73 USE kpoint_coulomb_2c, ONLY: build_2c_coulomb_matrix_kp
74 USE kpoint_methods, ONLY: kpoint_init_cell_index,&
75 rskp_transform
76 USE kpoint_types, ONLY: get_kpoint_info,&
77 kpoint_type
78 USE libint_2c_3c, ONLY: compare_potential_types,&
79 libint_potential_type
80 USE machine, ONLY: m_flush
81 USE mathconstants, ONLY: gaussi,&
82 z_one,&
83 z_zero
84 USE message_passing, ONLY: mp_comm_type,&
85 mp_para_env_release,&
86 mp_para_env_type,&
87 mp_request_type
88 USE mp2_eri, ONLY: mp2_eri_2c_integrate
89 USE mp2_eri_gpw, ONLY: mp2_eri_2c_integrate_gpw
90 USE mp2_types, ONLY: integ_mat_buffer_type
91 USE parallel_gemm_api, ONLY: parallel_gemm
92 USE particle_methods, ONLY: get_particle_set
93 USE particle_types, ONLY: particle_type
94 USE qs_environment_types, ONLY: get_qs_env,&
95 qs_environment_type
96 USE qs_integral_utils, ONLY: basis_set_list_setup
97 USE qs_interactions, ONLY: init_interaction_radii
98 USE qs_kind_types, ONLY: get_qs_kind,&
99 qs_kind_type
100 USE qs_ks_types, ONLY: get_ks_env,&
101 set_ks_env
102 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type,&
103 release_neighbor_list_sets
104 USE qs_tensors, ONLY: build_2c_integrals,&
105 build_2c_neighbor_lists
106 USE rpa_communication, ONLY: communicate_buffer
107 USE rpa_gw_kpoints_util, ONLY: cp_cfm_power
108
109!$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num
110#include "./base/base_uses.f90"
111
112 IMPLICIT NONE
113
114 PRIVATE
115
116 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'mp2_ri_2c'
117
118 PUBLIC :: get_2c_integrals, trunc_coulomb_for_exchange, ri_2c_integral_mat, &
119 inversion_of_m_and_mult_with_chol_dec_of_v
120
121CONTAINS
122
123! **************************************************************************************************
124!> \brief ...
125!> \param qs_env ...
126!> \param eri_method ...
127!> \param eri_param ...
128!> \param para_env ...
129!> \param para_env_sub ...
130!> \param mp2_memory ...
131!> \param my_Lrows ...
132!> \param my_Vrows ...
133!> \param fm_matrix_PQ ...
134!> \param ngroup ...
135!> \param color_sub ...
136!> \param dimen_RI ...
137!> \param dimen_RI_red reduced RI dimension, needed if we perform SVD instead of Cholesky
138!> \param kpoints ...
139!> \param my_group_L_size ...
140!> \param my_group_L_start ...
141!> \param my_group_L_end ...
142!> \param gd_array ...
143!> \param calc_PQ_cond_num ...
144!> \param do_svd ...
145!> \param eps_svd ...
146!> \param potential ...
147!> \param ri_metric ...
148!> \param fm_matrix_L_kpoints ...
149!> \param fm_matrix_Minv_L_kpoints ...
150!> \param fm_matrix_Minv ...
151!> \param fm_matrix_Minv_Vtrunc_Minv ...
152!> \param do_im_time ...
153!> \param do_kpoints ...
154!> \param mp2_eps_pgf_orb_S ...
155!> \param qs_kind_set ...
156!> \param sab_orb_sub ...
157!> \param calc_forces ...
158!> \param unit_nr ...
159! **************************************************************************************************
160 SUBROUTINE get_2c_integrals(qs_env, eri_method, eri_param, para_env, para_env_sub, mp2_memory, &
161 my_Lrows, my_Vrows, fm_matrix_PQ, ngroup, color_sub, dimen_RI, dimen_RI_red, &
162 kpoints, my_group_L_size, my_group_L_start, my_group_L_end, &
163 gd_array, calc_PQ_cond_num, do_svd, eps_svd, potential, ri_metric, &
164 fm_matrix_L_kpoints, fm_matrix_Minv_L_kpoints, &
165 fm_matrix_Minv, fm_matrix_Minv_Vtrunc_Minv, &
166 do_im_time, do_kpoints, mp2_eps_pgf_orb_S, qs_kind_set, sab_orb_sub, calc_forces, unit_nr)
167
168 TYPE(qs_environment_type), POINTER :: qs_env
169 INTEGER, INTENT(IN) :: eri_method
170 TYPE(cp_eri_mme_param), POINTER :: eri_param
171 TYPE(mp_para_env_type), POINTER :: para_env, para_env_sub
172 REAL(kind=dp), INTENT(IN) :: mp2_memory
173 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :), &
174 INTENT(OUT) :: my_lrows, my_vrows
175 TYPE(cp_fm_type), INTENT(OUT) :: fm_matrix_pq
176 INTEGER, INTENT(IN) :: ngroup, color_sub
177 INTEGER, INTENT(OUT) :: dimen_ri, dimen_ri_red
178 TYPE(kpoint_type), POINTER :: kpoints
179 INTEGER, INTENT(OUT) :: my_group_l_size, my_group_l_start, &
180 my_group_l_end
181 TYPE(group_dist_d1_type), INTENT(OUT) :: gd_array
182 LOGICAL, INTENT(IN) :: calc_pq_cond_num, do_svd
183 REAL(kind=dp), INTENT(IN) :: eps_svd
184 TYPE(libint_potential_type) :: potential, ri_metric
185 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: fm_matrix_l_kpoints, &
186 fm_matrix_minv_l_kpoints, &
187 fm_matrix_minv, &
188 fm_matrix_minv_vtrunc_minv
189 LOGICAL, INTENT(IN) :: do_im_time, do_kpoints
190 REAL(kind=dp), INTENT(IN) :: mp2_eps_pgf_orb_s
191 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
192 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
193 POINTER :: sab_orb_sub
194 LOGICAL, INTENT(IN) :: calc_forces
195 INTEGER, INTENT(IN) :: unit_nr
196
197 CHARACTER(LEN=*), PARAMETER :: routinen = 'get_2c_integrals'
198
199 INTEGER :: handle, num_small_eigen
200 REAL(kind=dp) :: cond_num, eps_pgf_orb_old
201 TYPE(cp_fm_type) :: fm_matrix_l_work, fm_matrix_m_inv_work, &
202 fm_matrix_v
203 TYPE(dft_control_type), POINTER :: dft_control
204 TYPE(libint_potential_type) :: trunc_coulomb
205 TYPE(mp_para_env_type), POINTER :: para_env_l
206
207 CALL timeset(routinen, handle)
208
209 ! calculate V and store it in fm_matrix_L_work
210 CALL compute_2c_integrals(qs_env, eri_method, eri_param, para_env, para_env_sub, para_env_l, mp2_memory, &
211 fm_matrix_l_work, ngroup, color_sub, dimen_ri, &
212 my_group_l_size, my_group_l_start, my_group_l_end, &
213 gd_array, calc_pq_cond_num, cond_num, &
214 num_small_eigen, potential, sab_orb_sub, do_im_time=do_im_time)
215
216 IF (do_im_time .AND. calc_forces) THEN
217 !need a copy of the (P|Q) integrals
218 CALL cp_fm_create(fm_matrix_pq, fm_matrix_l_work%matrix_struct)
219 CALL cp_fm_to_fm(fm_matrix_l_work, fm_matrix_pq)
220 END IF
221
222 dimen_ri_red = dimen_ri
223
224 IF (compare_potential_types(ri_metric, potential) .AND. .NOT. do_im_time) THEN
225 CALL decomp_mat_l(fm_matrix_l_work, do_svd, eps_svd, num_small_eigen, cond_num, .true., gd_array, ngroup, &
226 dimen_ri, dimen_ri_red, para_env)
227 ELSE
228
229 ! RI-metric matrix (depending on RI metric: Coulomb, overlap or truncated Coulomb)
230 IF (do_im_time) THEN
231 CALL get_qs_env(qs_env, dft_control=dft_control)
232
233 ! re-init the radii to be able to generate pair lists with appropriate screening for overlap matrix
234 eps_pgf_orb_old = dft_control%qs_control%eps_pgf_orb
235 dft_control%qs_control%eps_pgf_orb = mp2_eps_pgf_orb_s
236 CALL init_interaction_radii(dft_control%qs_control, qs_kind_set)
237
238 CALL ri_2c_integral_mat(qs_env, fm_matrix_minv_l_kpoints, fm_matrix_l_work, dimen_ri, ri_metric, &
239 do_kpoints, kpoints, put_mat_ks_env=.true., &
240 regularization_ri=qs_env%mp2_env%ri_rpa_im_time%regularization_RI)
241
242 ! re-init the radii to previous values
243 dft_control%qs_control%eps_pgf_orb = eps_pgf_orb_old
244 CALL init_interaction_radii(dft_control%qs_control, qs_kind_set)
245
246 ! GPW overlap matrix
247 ELSE
248
249 CALL mp_para_env_release(para_env_l)
250 CALL release_group_dist(gd_array)
251
252 ALLOCATE (fm_matrix_minv_l_kpoints(1, 1))
253
254 ! Calculate matrix of RI operator (for overlap metric: S), store it in fm_matrix_Minv_L_kpoints
255 CALL compute_2c_integrals(qs_env, eri_method, eri_param, para_env, para_env_sub, para_env_l, mp2_memory, &
256 fm_matrix_minv_l_kpoints(1, 1), ngroup, color_sub, dimen_ri, &
257 my_group_l_size, my_group_l_start, my_group_l_end, &
258 gd_array, calc_pq_cond_num, cond_num, &
259 num_small_eigen, ri_metric, sab_orb_sub, &
260 fm_matrix_l_extern=fm_matrix_l_work)
261
262 END IF
263
264 IF (do_kpoints) THEN
265
266 ! k-dependent 1/r Coulomb matrix V_PQ(k)
267 CALL compute_v_by_lattice_sum(qs_env, fm_matrix_l_kpoints, fm_matrix_minv_l_kpoints, kpoints)
268
269 CALL inversion_of_m_and_mult_with_chol_dec_of_v(fm_matrix_minv_l_kpoints, fm_matrix_l_kpoints, dimen_ri, &
270 kpoints, qs_env%mp2_env%ri_rpa_im_time%eps_eigval_S)
271
272 CALL trunc_coulomb_for_exchange(qs_env, trunc_coulomb)
273
274 ! Gamma-only truncated Coulomb matrix V^tr with cutoff radius = half the unit cell size; for exchange self-energy
275 CALL ri_2c_integral_mat(qs_env, fm_matrix_minv_vtrunc_minv, fm_matrix_l_work, dimen_ri, trunc_coulomb, &
276 do_kpoints=.false., kpoints=kpoints, put_mat_ks_env=.false., regularization_ri=0.0_dp)
277
278 ! Gamma-only RI-metric matrix; for computing Gamma-only/MIC self-energy
279 CALL ri_2c_integral_mat(qs_env, fm_matrix_minv, fm_matrix_l_work, dimen_ri, ri_metric, &
280 do_kpoints=.false., kpoints=kpoints, put_mat_ks_env=.false., regularization_ri=0.0_dp)
281
282 CALL gamma_only_inversion_of_m_and_mult_with_chol_dec_of_vtrunc(fm_matrix_minv_vtrunc_minv, &
283 fm_matrix_minv, qs_env)
284 ELSE
285 IF (calc_forces .AND. (.NOT. do_im_time)) THEN
286 ! For the gradients, we make a copy of the inverse root of L
287 CALL cp_fm_create(fm_matrix_v, fm_matrix_l_work%matrix_struct)
288 CALL cp_fm_to_fm(fm_matrix_l_work, fm_matrix_v)
289
290 CALL decomp_mat_l(fm_matrix_v, do_svd, eps_svd, num_small_eigen, cond_num, .true., gd_array, ngroup, &
291 dimen_ri, dimen_ri_red, para_env)
292 END IF
293
294 CALL decomp_mat_l(fm_matrix_l_work, do_svd, eps_svd, num_small_eigen, cond_num, .false., gd_array, ngroup, &
295 dimen_ri, dimen_ri_red, para_env)
296
297 CALL decomp_mat_l(fm_matrix_minv_l_kpoints(1, 1), .false., 0.0_dp, num_small_eigen, cond_num, .true., &
298 gd_array, ngroup, dimen_ri, dimen_ri_red, para_env)
299
300 CALL cp_fm_create(fm_matrix_m_inv_work, fm_matrix_minv_l_kpoints(1, 1)%matrix_struct)
301 CALL cp_fm_set_all(fm_matrix_m_inv_work, 0.0_dp)
302
303 CALL parallel_gemm('N', 'T', dimen_ri, dimen_ri, dimen_ri, 1.0_dp, fm_matrix_minv_l_kpoints(1, 1), &
304 fm_matrix_minv_l_kpoints(1, 1), 0.0_dp, fm_matrix_m_inv_work)
305
306 IF (do_svd) THEN
307 ! We have to reset the size of fm_matrix_Minv_L_kpoints
308 CALL reset_size_matrix(fm_matrix_minv_l_kpoints(1, 1), dimen_ri_red, fm_matrix_l_work%matrix_struct)
309
310 ! L (=fm_matrix_Minv_L_kpoints) = S^(-1)*chol_dec(V)
311 CALL parallel_gemm('T', 'N', dimen_ri, dimen_ri_red, dimen_ri, 1.0_dp, fm_matrix_m_inv_work, &
312 fm_matrix_l_work, 0.0_dp, fm_matrix_minv_l_kpoints(1, 1))
313 ELSE
314 ! L (=fm_matrix_Minv_L_kpoints) = S^(-1)*chol_dec(V)
315 CALL parallel_gemm('T', 'T', dimen_ri, dimen_ri, dimen_ri, 1.0_dp, fm_matrix_m_inv_work, &
316 fm_matrix_l_work, 0.0_dp, fm_matrix_minv_l_kpoints(1, 1))
317 END IF
318
319 ! clean the S_inv matrix
320 CALL cp_fm_release(fm_matrix_m_inv_work)
321 END IF
322
323 IF (.NOT. do_im_time) THEN
324
325 CALL cp_fm_to_fm(fm_matrix_minv_l_kpoints(1, 1), fm_matrix_l_work)
326 CALL cp_fm_release(fm_matrix_minv_l_kpoints)
327
328 END IF
329
330 END IF
331
332 ! If the group distribution changed because of an SVD, we get the new values here
333 CALL get_group_dist(gd_array, color_sub, my_group_l_start, my_group_l_end, my_group_l_size)
334
335 IF (.NOT. do_im_time) THEN
336 IF (unit_nr > 0) THEN
337 WRITE (unit=unit_nr, fmt="(T3,A,T75,i6)") "RI_INFO| Cholesky decomposition group size:", para_env_l%num_pe
338 WRITE (unit=unit_nr, fmt="(T3,A,T75,i6)") "RI_INFO| Number of groups for auxiliary basis functions", ngroup
339 IF (calc_pq_cond_num .OR. do_svd) THEN
340 WRITE (unit=unit_nr, fmt="(T3,A,T67,ES14.5)") &
341 "RI_INFO| Condition number of the (P|Q):", cond_num
342 WRITE (unit=unit_nr, fmt="(T3,A,T75,i6)") &
343 "RI_INFO| Number of non-positive Eigenvalues of (P|Q):", num_small_eigen
344 END IF
345 CALL m_flush(unit_nr)
346 END IF
347
348 ! replicate the necessary row of the L^{-1} matrix on each proc
349 CALL grep_lcols(fm_matrix_l_work, my_group_l_start, my_group_l_end, my_group_l_size, my_lrows)
350 IF (calc_forces .AND. .NOT. compare_potential_types(qs_env%mp2_env%ri_metric, &
351 qs_env%mp2_env%potential_parameter)) THEN
352 CALL grep_lcols(fm_matrix_v, my_group_l_start, my_group_l_end, my_group_l_size, my_vrows)
353 END IF
354 END IF
355 CALL cp_fm_release(fm_matrix_l_work)
356 IF (calc_forces .AND. .NOT. (do_im_time .OR. compare_potential_types(qs_env%mp2_env%ri_metric, &
357 qs_env%mp2_env%potential_parameter))) CALL cp_fm_release(fm_matrix_v)
358 CALL mp_para_env_release(para_env_l)
359
360 CALL timestop(handle)
361
362 END SUBROUTINE get_2c_integrals
363
364! **************************************************************************************************
365!> \brief ...
366!> \param fm_matrix_L ...
367!> \param do_svd ...
368!> \param eps_svd ...
369!> \param num_small_eigen ...
370!> \param cond_num ...
371!> \param do_inversion ...
372!> \param gd_array ...
373!> \param ngroup ...
374!> \param dimen_RI ...
375!> \param dimen_RI_red ...
376!> \param para_env ...
377! **************************************************************************************************
378 SUBROUTINE decomp_mat_l(fm_matrix_L, do_svd, eps_svd, num_small_eigen, cond_num, do_inversion, gd_array, ngroup, &
379 dimen_RI, dimen_RI_red, para_env)
380
381 TYPE(cp_fm_type), INTENT(INOUT) :: fm_matrix_l
382 LOGICAL, INTENT(IN) :: do_svd
383 REAL(kind=dp), INTENT(IN) :: eps_svd
384 INTEGER, INTENT(INOUT) :: num_small_eigen
385 REAL(kind=dp), INTENT(INOUT) :: cond_num
386 LOGICAL, INTENT(IN) :: do_inversion
387 TYPE(group_dist_d1_type), INTENT(INOUT) :: gd_array
388 INTEGER, INTENT(IN) :: ngroup, dimen_ri
389 INTEGER, INTENT(INOUT) :: dimen_ri_red
390 TYPE(mp_para_env_type), INTENT(IN) :: para_env
391
392 IF (do_svd) THEN
393 CALL matrix_root_with_svd(fm_matrix_l, num_small_eigen, cond_num, eps_svd, do_inversion, para_env)
394
395 dimen_ri_red = dimen_ri - num_small_eigen
396
397 ! We changed the size of fm_matrix_L in matrix_root_with_svd.
398 ! So, we have to get new group sizes
399 CALL release_group_dist(gd_array)
400
401 CALL create_group_dist(gd_array, ngroup, dimen_ri_red)
402 ELSE
403
404 ! calculate Cholesky decomposition L (V=LL^T)
405 CALL cholesky_decomp(fm_matrix_l, dimen_ri, do_inversion=do_inversion)
406
407 IF (do_inversion) CALL invert_mat(fm_matrix_l)
408 END IF
409
410 END SUBROUTINE decomp_mat_l
411
412! **************************************************************************************************
413!> \brief ...
414!> \param qs_env ...
415!> \param fm_matrix_L_kpoints ...
416!> \param fm_matrix_Minv_L_kpoints ...
417!> \param kpoints ...
418! **************************************************************************************************
419 SUBROUTINE compute_v_by_lattice_sum(qs_env, fm_matrix_L_kpoints, fm_matrix_Minv_L_kpoints, kpoints)
420 TYPE(qs_environment_type), POINTER :: qs_env
421 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: fm_matrix_l_kpoints, &
422 fm_matrix_minv_l_kpoints
423 TYPE(kpoint_type), POINTER :: kpoints
424
425 CHARACTER(LEN=*), PARAMETER :: routinen = 'compute_V_by_lattice_sum'
426
427 INTEGER :: handle, i_dim, i_real_imag, ikp, nkp, &
428 nkp_extra, nkp_orig
429 INTEGER, DIMENSION(3) :: nkp_grid_orig, periodic
430 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
431 TYPE(cell_type), POINTER :: cell
432 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_ri_aux_transl, matrix_v_ri_kp
433 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
434 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
435
436 CALL timeset(routinen, handle)
437
438 NULLIFY (matrix_s_ri_aux_transl, particle_set, cell, qs_kind_set)
439
440 CALL get_qs_env(qs_env=qs_env, &
441 matrix_s_ri_aux_kp=matrix_s_ri_aux_transl, &
442 particle_set=particle_set, &
443 cell=cell, &
444 qs_kind_set=qs_kind_set, &
445 atomic_kind_set=atomic_kind_set)
446
447 ! check that we have a 2n x 2m x 2k mesh to guarantee good convergence
448 CALL get_cell(cell=cell, periodic=periodic)
449 DO i_dim = 1, 3
450 IF (periodic(i_dim) == 1) THEN
451 cpassert(modulo(kpoints%nkp_grid(i_dim), 2) == 0)
452 END IF
453 END DO
454
455 nkp = kpoints%nkp
456
457 ALLOCATE (fm_matrix_l_kpoints(nkp, 2))
458 DO i_real_imag = 1, 2
459 DO ikp = 1, nkp
460 CALL cp_fm_create(fm_matrix_l_kpoints(ikp, i_real_imag), &
461 fm_matrix_minv_l_kpoints(1, i_real_imag)%matrix_struct)
462 CALL cp_fm_set_all(fm_matrix_l_kpoints(ikp, i_real_imag), 0.0_dp)
463 END DO
464 END DO
465
466 CALL allocate_matrix_v_ri_kp(matrix_v_ri_kp, matrix_s_ri_aux_transl, nkp)
467
468 IF (qs_env%mp2_env%ri_rpa_im_time%do_extrapolate_kpoints) THEN
469
470 nkp_orig = qs_env%mp2_env%ri_rpa_im_time%nkp_orig
471 nkp_extra = qs_env%mp2_env%ri_rpa_im_time%nkp_extra
472
473 CALL build_2c_coulomb_matrix_kp(matrix_v_ri_kp, kpoints, basis_type="RI_AUX", &
474 cell=cell, particle_set=particle_set, qs_kind_set=qs_kind_set, &
475 atomic_kind_set=atomic_kind_set, &
476 size_lattice_sum=qs_env%mp2_env%mp2_gpw%size_lattice_sum, &
477 operator_type=operator_coulomb, ikp_start=1, ikp_end=nkp_orig)
478
479 nkp_grid_orig = kpoints%nkp_grid
480 kpoints%nkp_grid(1:3) = qs_env%mp2_env%ri_rpa_im_time%kp_grid_extra(1:3)
481
482 CALL build_2c_coulomb_matrix_kp(matrix_v_ri_kp, kpoints, basis_type="RI_AUX", &
483 cell=cell, particle_set=particle_set, qs_kind_set=qs_kind_set, &
484 atomic_kind_set=atomic_kind_set, &
485 size_lattice_sum=qs_env%mp2_env%mp2_gpw%size_lattice_sum, &
486 operator_type=operator_coulomb, ikp_start=nkp_orig + 1, ikp_end=nkp)
487
488 kpoints%nkp_grid = nkp_grid_orig
489
490 ELSE
491
492 CALL build_2c_coulomb_matrix_kp(matrix_v_ri_kp, kpoints, basis_type="RI_AUX", &
493 cell=cell, particle_set=particle_set, qs_kind_set=qs_kind_set, &
494 atomic_kind_set=atomic_kind_set, &
495 size_lattice_sum=qs_env%mp2_env%mp2_gpw%size_lattice_sum, &
496 operator_type=operator_coulomb, ikp_start=1, ikp_end=nkp)
497
498 END IF
499
500 DO ikp = 1, nkp
501
502 CALL copy_dbcsr_to_fm(matrix_v_ri_kp(ikp, 1)%matrix, fm_matrix_l_kpoints(ikp, 1))
503 CALL copy_dbcsr_to_fm(matrix_v_ri_kp(ikp, 2)%matrix, fm_matrix_l_kpoints(ikp, 2))
504
505 END DO
506
507 CALL dbcsr_deallocate_matrix_set(matrix_v_ri_kp)
508
509 CALL timestop(handle)
510
511 END SUBROUTINE compute_v_by_lattice_sum
512
513! **************************************************************************************************
514!> \brief ...
515!> \param matrix_v_RI_kp ...
516!> \param matrix_s_RI_aux_transl ...
517!> \param nkp ...
518! **************************************************************************************************
519 SUBROUTINE allocate_matrix_v_ri_kp(matrix_v_RI_kp, matrix_s_RI_aux_transl, nkp)
520
521 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_v_ri_kp, matrix_s_ri_aux_transl
522 INTEGER :: nkp
523
524 INTEGER :: ikp
525
526 NULLIFY (matrix_v_ri_kp)
527 CALL dbcsr_allocate_matrix_set(matrix_v_ri_kp, nkp, 2)
528
529 DO ikp = 1, nkp
530
531 ALLOCATE (matrix_v_ri_kp(ikp, 1)%matrix)
532 CALL dbcsr_create(matrix_v_ri_kp(ikp, 1)%matrix, template=matrix_s_ri_aux_transl(1, 1)%matrix, &
533 matrix_type=dbcsr_type_no_symmetry)
534 CALL dbcsr_reserve_all_blocks(matrix_v_ri_kp(ikp, 1)%matrix)
535 CALL dbcsr_set(matrix_v_ri_kp(ikp, 1)%matrix, 0.0_dp)
536
537 ALLOCATE (matrix_v_ri_kp(ikp, 2)%matrix)
538 CALL dbcsr_create(matrix_v_ri_kp(ikp, 2)%matrix, template=matrix_s_ri_aux_transl(1, 1)%matrix, &
539 matrix_type=dbcsr_type_no_symmetry)
540 CALL dbcsr_reserve_all_blocks(matrix_v_ri_kp(ikp, 2)%matrix)
541 CALL dbcsr_set(matrix_v_ri_kp(ikp, 2)%matrix, 0.0_dp)
542
543 END DO
544
545 END SUBROUTINE allocate_matrix_v_ri_kp
546
547! **************************************************************************************************
548!> \brief ...
549!> \param qs_env ...
550!> \param fm_matrix_Minv_L_kpoints ...
551!> \param fm_matrix_L ...
552!> \param dimen_RI ...
553!> \param ri_metric ...
554!> \param do_kpoints ...
555!> \param kpoints ...
556!> \param put_mat_KS_env ...
557!> \param regularization_RI ...
558!> \param ikp_ext ...
559!> \param do_build_cell_index ...
560! **************************************************************************************************
561 SUBROUTINE ri_2c_integral_mat(qs_env, fm_matrix_Minv_L_kpoints, fm_matrix_L, dimen_RI, ri_metric, &
562 do_kpoints, kpoints, put_mat_KS_env, regularization_RI, ikp_ext, &
563 do_build_cell_index)
564
565 TYPE(qs_environment_type), POINTER :: qs_env
566 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: fm_matrix_minv_l_kpoints
567 TYPE(cp_fm_type), INTENT(IN) :: fm_matrix_l
568 INTEGER, INTENT(IN) :: dimen_ri
569 TYPE(libint_potential_type) :: ri_metric
570 LOGICAL, INTENT(IN) :: do_kpoints
571 TYPE(kpoint_type), OPTIONAL, POINTER :: kpoints
572 LOGICAL, OPTIONAL :: put_mat_ks_env
573 REAL(kind=dp), OPTIONAL :: regularization_ri
574 INTEGER, OPTIONAL :: ikp_ext
575 LOGICAL, OPTIONAL :: do_build_cell_index
576
577 CHARACTER(LEN=*), PARAMETER :: routinen = 'RI_2c_integral_mat'
578
579 INTEGER :: handle, i_real_imag, i_size, ikp, &
580 ikp_for_xkp, img, n_real_imag, natom, &
581 nimg, nkind, nkp
582 INTEGER, ALLOCATABLE, DIMENSION(:) :: sizes_ri
583 INTEGER, DIMENSION(:), POINTER :: col_bsize, row_bsize
584 INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
585 LOGICAL :: my_do_build_cell_index, my_put_mat_ks_env
586 REAL(kind=dp) :: my_regularization_ri
587 REAL(kind=dp), DIMENSION(:, :), POINTER :: xkp
588 TYPE(cp_blacs_env_type), POINTER :: blacs_env
589 TYPE(cp_fm_struct_type), POINTER :: fm_struct
590 TYPE(cp_fm_type) :: fm_matrix_s_global
591 TYPE(dbcsr_distribution_type) :: dbcsr_dist
592 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_ri_aux_transl
593 TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:) :: mat_2c
594 TYPE(dbcsr_type), POINTER :: cmatrix, matrix_s_ri_aux_desymm, &
595 rmatrix, tmpmat
596 TYPE(dft_control_type), POINTER :: dft_control
597 TYPE(distribution_2d_type), POINTER :: dist_2d
598 TYPE(gto_basis_set_p_type), DIMENSION(:), POINTER :: basis_set_ri
599 TYPE(mp_para_env_type), POINTER :: para_env
600 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
601 POINTER :: sab_ri
602 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
603 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
604
605 CALL timeset(routinen, handle)
606
607 NULLIFY (sab_ri, matrix_s_ri_aux_transl, dist_2d)
608
609 IF (PRESENT(regularization_ri)) THEN
610 my_regularization_ri = regularization_ri
611 ELSE
612 my_regularization_ri = 0.0_dp
613 END IF
614
615 IF (PRESENT(put_mat_ks_env)) THEN
616 my_put_mat_ks_env = put_mat_ks_env
617 ELSE
618 my_put_mat_ks_env = .false.
619 END IF
620
621 IF (PRESENT(do_build_cell_index)) THEN
622 my_do_build_cell_index = do_build_cell_index
623 ELSE
624 my_do_build_cell_index = .false.
625 END IF
626
627 CALL get_qs_env(qs_env=qs_env, &
628 para_env=para_env, &
629 blacs_env=blacs_env, &
630 nkind=nkind, &
631 distribution_2d=dist_2d, &
632 qs_kind_set=qs_kind_set, &
633 particle_set=particle_set, &
634 dft_control=dft_control, &
635 natom=natom)
636
637 ALLOCATE (sizes_ri(natom))
638 ALLOCATE (basis_set_ri(nkind))
639
640 CALL basis_set_list_setup(basis_set_ri, 'RI_AUX', qs_kind_set)
641 CALL get_particle_set(particle_set, qs_kind_set, nsgf=sizes_ri, basis=basis_set_ri)
642
643 CALL build_2c_neighbor_lists(sab_ri, basis_set_ri, basis_set_ri, ri_metric, "RPA_2c_nl_RI", qs_env, &
644 sym_ij=.true., dist_2d=dist_2d)
645
646 CALL cp_dbcsr_dist2d_to_dist(dist_2d, dbcsr_dist)
647 ALLOCATE (row_bsize(SIZE(sizes_ri)))
648 ALLOCATE (col_bsize(SIZE(sizes_ri)))
649 row_bsize(:) = sizes_ri
650 col_bsize(:) = sizes_ri
651
652 IF (do_kpoints) THEN
653 cpassert(PRESENT(kpoints))
654 IF (my_do_build_cell_index) THEN
655 CALL kpoint_init_cell_index(kpoints, sab_ri, para_env, dft_control)
656 END IF
657 CALL get_kpoint_info(kpoints, nkp=nkp, xkp=xkp, &
658 cell_to_index=cell_to_index)
659 n_real_imag = 2
660 nimg = dft_control%nimages
661 ELSE
662 nkp = 1
663 n_real_imag = 1
664 nimg = 1
665 END IF
666
667 ALLOCATE (mat_2c(nimg))
668 CALL dbcsr_create(mat_2c(1), "(RI|RI)", dbcsr_dist, dbcsr_type_symmetric, &
669 row_bsize, col_bsize)
670 DEALLOCATE (row_bsize, col_bsize)
671
672 DO img = 2, nimg
673 CALL dbcsr_create(mat_2c(img), template=mat_2c(1))
674 END DO
675
676 CALL build_2c_integrals(mat_2c, 0.0_dp, qs_env, sab_ri, basis_set_ri, basis_set_ri, &
677 ri_metric, do_kpoints=do_kpoints, ext_kpoints=kpoints, &
678 regularization_ri=my_regularization_ri)
679
680 CALL dbcsr_distribution_release(dbcsr_dist)
681 DEALLOCATE (basis_set_ri)
682
683 IF (my_put_mat_ks_env) THEN
684 CALL get_ks_env(qs_env%ks_env, matrix_s_ri_aux_kp=matrix_s_ri_aux_transl)
685 IF (ASSOCIATED(matrix_s_ri_aux_transl)) CALL dbcsr_deallocate_matrix_set(matrix_s_ri_aux_transl)
686 END IF
687 CALL dbcsr_allocate_matrix_set(matrix_s_ri_aux_transl, 1, nimg)
688 DO img = 1, nimg
689 ALLOCATE (matrix_s_ri_aux_transl(1, img)%matrix)
690 CALL dbcsr_copy(matrix_s_ri_aux_transl(1, img)%matrix, mat_2c(img))
691 CALL dbcsr_release(mat_2c(img))
692 END DO
693
694 IF (my_put_mat_ks_env) THEN
695 CALL set_ks_env(qs_env%ks_env, matrix_s_ri_aux_kp=matrix_s_ri_aux_transl)
696 END IF
697
698 IF (PRESENT(ikp_ext)) nkp = 1
699
700 ALLOCATE (fm_matrix_minv_l_kpoints(nkp, n_real_imag))
701 DO i_real_imag = 1, n_real_imag
702 DO i_size = 1, nkp
703 CALL cp_fm_create(fm_matrix_minv_l_kpoints(i_size, i_real_imag), fm_matrix_l%matrix_struct)
704 CALL cp_fm_set_all(fm_matrix_minv_l_kpoints(i_size, i_real_imag), 0.0_dp)
705 END DO
706 END DO
707
708 NULLIFY (fm_struct)
709 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=dimen_ri, &
710 ncol_global=dimen_ri, para_env=para_env)
711
712 CALL cp_fm_create(fm_matrix_s_global, fm_struct)
713 CALL cp_fm_set_all(fm_matrix_s_global, 0.0_dp)
714
715 IF (do_kpoints) THEN
716
717 ALLOCATE (rmatrix, cmatrix, tmpmat)
718 CALL dbcsr_create(rmatrix, template=matrix_s_ri_aux_transl(1, 1)%matrix, &
719 matrix_type=dbcsr_type_symmetric)
720 CALL dbcsr_create(cmatrix, template=matrix_s_ri_aux_transl(1, 1)%matrix, &
721 matrix_type=dbcsr_type_antisymmetric)
722 CALL dbcsr_create(tmpmat, template=matrix_s_ri_aux_transl(1, 1)%matrix, &
723 matrix_type=dbcsr_type_no_symmetry)
724 CALL cp_dbcsr_alloc_block_from_nbl(rmatrix, sab_ri)
725 CALL cp_dbcsr_alloc_block_from_nbl(cmatrix, sab_ri)
726
727 DO ikp = 1, nkp
728
729 ! s matrix is not spin dependent, double the work
730 CALL dbcsr_set(rmatrix, 0.0_dp)
731 CALL dbcsr_set(cmatrix, 0.0_dp)
732
733 IF (PRESENT(ikp_ext)) THEN
734 ikp_for_xkp = ikp_ext
735 ELSE
736 ikp_for_xkp = ikp
737 END IF
738
739 CALL rskp_transform(rmatrix=rmatrix, cmatrix=cmatrix, rsmat=matrix_s_ri_aux_transl, ispin=1, &
740 xkp=xkp(1:3, ikp_for_xkp), cell_to_index=cell_to_index, sab_nl=sab_ri)
741
742 CALL dbcsr_set(tmpmat, 0.0_dp)
743 CALL dbcsr_desymmetrize(rmatrix, tmpmat)
744
745 CALL cp_fm_set_all(fm_matrix_s_global, 0.0_dp)
746 CALL copy_dbcsr_to_fm(tmpmat, fm_matrix_s_global)
747 CALL cp_fm_copy_general(fm_matrix_s_global, fm_matrix_minv_l_kpoints(ikp, 1), para_env)
748
749 CALL dbcsr_set(tmpmat, 0.0_dp)
750 CALL dbcsr_desymmetrize(cmatrix, tmpmat)
751
752 CALL cp_fm_set_all(fm_matrix_s_global, 0.0_dp)
753 CALL copy_dbcsr_to_fm(tmpmat, fm_matrix_s_global)
754 CALL cp_fm_copy_general(fm_matrix_s_global, fm_matrix_minv_l_kpoints(ikp, 2), para_env)
755
756 END DO
757
758 CALL dbcsr_deallocate_matrix(rmatrix)
759 CALL dbcsr_deallocate_matrix(cmatrix)
760 CALL dbcsr_deallocate_matrix(tmpmat)
761
762 ELSE
763
764 NULLIFY (matrix_s_ri_aux_desymm)
765 ALLOCATE (matrix_s_ri_aux_desymm)
766 CALL dbcsr_create(matrix_s_ri_aux_desymm, template=matrix_s_ri_aux_transl(1, 1)%matrix, &
767 name='S_RI non_symm', matrix_type=dbcsr_type_no_symmetry)
768
769 CALL dbcsr_desymmetrize(matrix_s_ri_aux_transl(1, 1)%matrix, matrix_s_ri_aux_desymm)
770
771 CALL copy_dbcsr_to_fm(matrix_s_ri_aux_desymm, fm_matrix_s_global)
772
773 CALL cp_fm_copy_general(fm_matrix_s_global, fm_matrix_minv_l_kpoints(1, 1), para_env)
774
775 CALL dbcsr_deallocate_matrix(matrix_s_ri_aux_desymm)
776
777 END IF
778
779 CALL release_neighbor_list_sets(sab_ri)
780
781 CALL cp_fm_struct_release(fm_struct)
782
783 CALL cp_fm_release(fm_matrix_s_global)
784
785 IF (.NOT. my_put_mat_ks_env) THEN
786 CALL dbcsr_deallocate_matrix_set(matrix_s_ri_aux_transl)
787 END IF
788
789 CALL timestop(handle)
790
791 END SUBROUTINE ri_2c_integral_mat
792
793! **************************************************************************************************
794!> \brief ...
795!> \param qs_env ...
796!> \param eri_method ...
797!> \param eri_param ...
798!> \param para_env ...
799!> \param para_env_sub ...
800!> \param para_env_L ...
801!> \param mp2_memory ...
802!> \param fm_matrix_L ...
803!> \param ngroup ...
804!> \param color_sub ...
805!> \param dimen_RI ...
806!> \param my_group_L_size ...
807!> \param my_group_L_start ...
808!> \param my_group_L_end ...
809!> \param gd_array ...
810!> \param calc_PQ_cond_num ...
811!> \param cond_num ...
812!> \param num_small_eigen ...
813!> \param potential ...
814!> \param sab_orb_sub ...
815!> \param do_im_time ...
816!> \param fm_matrix_L_extern ...
817! **************************************************************************************************
818 SUBROUTINE compute_2c_integrals(qs_env, eri_method, eri_param, para_env, para_env_sub, para_env_L, mp2_memory, &
819 fm_matrix_L, ngroup, color_sub, dimen_RI, &
820 my_group_L_size, my_group_L_start, my_group_L_end, &
821 gd_array, calc_PQ_cond_num, cond_num, num_small_eigen, potential, &
822 sab_orb_sub, do_im_time, fm_matrix_L_extern)
823
824 TYPE(qs_environment_type), POINTER :: qs_env
825 INTEGER, INTENT(IN) :: eri_method
826 TYPE(cp_eri_mme_param), POINTER :: eri_param
827 TYPE(mp_para_env_type), INTENT(IN) :: para_env
828 TYPE(mp_para_env_type), INTENT(IN), POINTER :: para_env_sub
829 TYPE(mp_para_env_type), INTENT(OUT), POINTER :: para_env_l
830 REAL(kind=dp), INTENT(IN) :: mp2_memory
831 TYPE(cp_fm_type), INTENT(OUT) :: fm_matrix_l
832 INTEGER, INTENT(IN) :: ngroup, color_sub
833 INTEGER, INTENT(OUT) :: dimen_ri, my_group_l_size, &
834 my_group_l_start, my_group_l_end
835 TYPE(group_dist_d1_type), INTENT(OUT) :: gd_array
836 LOGICAL, INTENT(IN) :: calc_pq_cond_num
837 REAL(kind=dp), INTENT(OUT) :: cond_num
838 INTEGER, INTENT(OUT) :: num_small_eigen
839 TYPE(libint_potential_type), INTENT(IN) :: potential
840 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
841 POINTER :: sab_orb_sub
842 LOGICAL, INTENT(IN), OPTIONAL :: do_im_time
843 TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: fm_matrix_l_extern
844
845 CHARACTER(LEN=*), PARAMETER :: routinen = 'compute_2c_integrals'
846
847 INTEGER :: best_group_size, color_l, group_size, handle, handle2, i_global, iatom, iib, &
848 ikind, iproc, j_global, jjb, natom, ncol_local, nkind, nrow_local, nsgf, potential_type, &
849 proc_receive, proc_receive_static, proc_send_static, proc_shift, rec_l_end, rec_l_size, &
850 rec_l_start, strat_group_size
851 INTEGER, ALLOCATABLE, DIMENSION(:) :: kind_of
852 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
853 LOGICAL :: my_do_im_time
854 REAL(kind=dp) :: min_mem_for_qk
855 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: egen_l
856 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: l_external_col, l_local_col
857 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
858 TYPE(cp_blacs_env_type), POINTER :: blacs_env_l
859 TYPE(cp_fm_type) :: fm_matrix_l_diag
860 TYPE(group_dist_d1_type) :: gd_sub_array
861 TYPE(gto_basis_set_type), POINTER :: basis_set_a
862 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
863 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
864
865 CALL timeset(routinen, handle)
866
867 my_do_im_time = .false.
868 IF (PRESENT(do_im_time)) THEN
869 my_do_im_time = do_im_time
870 END IF
871
872 ! get stuff
873 CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set, qs_kind_set=qs_kind_set, &
874 particle_set=particle_set)
875
876 nkind = SIZE(qs_kind_set)
877 natom = SIZE(particle_set)
878
879 CALL get_atomic_kind_set(atomic_kind_set, kind_of=kind_of)
880
881 DO ikind = 1, nkind
882 CALL get_qs_kind(qs_kind=qs_kind_set(ikind), basis_set=basis_set_a, basis_type="RI_AUX")
883 cpassert(ASSOCIATED(basis_set_a))
884 END DO
885
886 dimen_ri = 0
887 DO iatom = 1, natom
888 ikind = kind_of(iatom)
889 CALL get_qs_kind(qs_kind=qs_kind_set(ikind), nsgf=nsgf, basis_type="RI_AUX")
890 dimen_ri = dimen_ri + nsgf
891 END DO
892
893 ! check that very small systems are not running on too many processes
894 IF (dimen_ri < ngroup) THEN
895 CALL cp_abort(__location__, "Product of block size and number "// &
896 "of RI functions should not exceed total number of processes")
897 END IF
898
899 CALL create_group_dist(gd_array, ngroup, dimen_ri)
900
901 CALL get_group_dist(gd_array, color_sub, my_group_l_start, my_group_l_end, my_group_l_size)
902
903 CALL timeset(routinen//"_loop_lm", handle2)
904
905 ALLOCATE (l_local_col(dimen_ri, my_group_l_size))
906 l_local_col = 0.0_dp
907
908 potential_type = potential%potential_type
909
910 IF ((eri_method .EQ. do_eri_mme .OR. eri_method .EQ. do_eri_os) &
911 .AND. potential_type .EQ. do_potential_coulomb .OR. &
912 (eri_method .EQ. do_eri_mme .AND. potential_type .EQ. do_potential_long)) THEN
913
914 CALL mp2_eri_2c_integrate(eri_param, potential, para_env_sub, qs_env, &
915 basis_type_a="RI_AUX", basis_type_b="RI_AUX", &
916 hab=l_local_col, first_b=my_group_l_start, last_b=my_group_l_end, &
917 eri_method=eri_method)
918
919 ELSEIF (eri_method .EQ. do_eri_gpw .OR. &
920 (potential_type == do_potential_long .AND. qs_env%mp2_env%eri_method == do_eri_os) &
921 .OR. (potential_type == do_potential_id .AND. qs_env%mp2_env%eri_method == do_eri_mme)) THEN
922
923 CALL mp2_eri_2c_integrate_gpw(qs_env, para_env_sub, my_group_l_start, my_group_l_end, &
924 natom, potential, sab_orb_sub, l_local_col, kind_of)
925
926 ELSE
927 cpabort("unknown ERI method")
928 END IF
929
930 CALL timestop(handle2)
931
932 ! split the total number of proc in a subgroup of the size of ~1/10 of the
933 ! total num of proc
934 best_group_size = para_env%num_pe
935
936 strat_group_size = max(1, para_env%num_pe/10)
937
938 min_mem_for_qk = real(dimen_ri, kind=dp)*dimen_ri*3.0_dp*8.0_dp/1024_dp/1024_dp
939
940 group_size = strat_group_size - 1
941 DO iproc = strat_group_size, para_env%num_pe
942 group_size = group_size + 1
943 ! check that group_size is a multiple of sub_group_size and a divisor of
944 ! the total num of proc
945 IF (mod(para_env%num_pe, group_size) /= 0 .OR. mod(group_size, para_env_sub%num_pe) /= 0) cycle
946
947 ! check for memory
948 IF (real(group_size, kind=dp)*mp2_memory < min_mem_for_qk) cycle
949
950 best_group_size = group_size
951 EXIT
952 END DO
953
954 IF (my_do_im_time) THEN
955 ! para_env_L is the para_env for im. time to avoid bug
956 best_group_size = para_env%num_pe
957 END IF
958
959 ! create the L group
960 color_l = para_env%mepos/best_group_size
961 ALLOCATE (para_env_l)
962 CALL para_env_l%from_split(para_env, color_l)
963
964 ! create the blacs_L
965 NULLIFY (blacs_env_l)
966 CALL cp_blacs_env_create(blacs_env=blacs_env_l, para_env=para_env_l)
967
968 ! create the full matrix L defined in the L group
969 CALL create_matrix_l(fm_matrix_l, blacs_env_l, dimen_ri, para_env_l, "fm_matrix_L", fm_matrix_l_extern)
970
971 IF (my_do_im_time .AND. para_env%num_pe > 1) THEN
972
973 CALL fill_fm_l_from_l_loc_non_blocking(fm_matrix_l, l_local_col, para_env, &
974 my_group_l_start, my_group_l_end, &
975 dimen_ri)
976
977 ELSE
978 block
979 TYPE(mp_comm_type) :: comm_exchange
980
981 CALL comm_exchange%from_split(para_env_l, para_env_sub%mepos)
982
983 CALL create_group_dist(gd_sub_array, my_group_l_start, &
984 my_group_l_end, my_group_l_size, comm_exchange)
985
986 CALL cp_fm_get_info(matrix=fm_matrix_l, &
987 nrow_local=nrow_local, &
988 ncol_local=ncol_local, &
989 row_indices=row_indices, &
990 col_indices=col_indices)
991
992 DO jjb = 1, ncol_local
993 j_global = col_indices(jjb)
994 IF (j_global >= my_group_l_start .AND. j_global <= my_group_l_end) THEN
995 DO iib = 1, nrow_local
996 i_global = row_indices(iib)
997 fm_matrix_l%local_data(iib, jjb) = l_local_col(i_global, j_global - my_group_l_start + 1)
998 END DO
999 END IF
1000 END DO
1001
1002 proc_send_static = modulo(comm_exchange%mepos + 1, comm_exchange%num_pe)
1003 proc_receive_static = modulo(comm_exchange%mepos - 1, comm_exchange%num_pe)
1004
1005 DO proc_shift = 1, comm_exchange%num_pe - 1
1006 proc_receive = modulo(comm_exchange%mepos - proc_shift, comm_exchange%num_pe)
1007
1008 CALL get_group_dist(gd_sub_array, proc_receive, rec_l_start, rec_l_end, rec_l_size)
1009
1010 ALLOCATE (l_external_col(dimen_ri, rec_l_size))
1011 l_external_col = 0.0_dp
1012
1013 CALL comm_exchange%sendrecv(l_local_col, proc_send_static, l_external_col, proc_receive_static)
1014
1015 DO jjb = 1, ncol_local
1016 j_global = col_indices(jjb)
1017 IF (j_global >= rec_l_start .AND. j_global <= rec_l_end) THEN
1018 DO iib = 1, nrow_local
1019 i_global = row_indices(iib)
1020 fm_matrix_l%local_data(iib, jjb) = l_external_col(i_global, j_global - rec_l_start + 1)
1021 END DO
1022 END IF
1023 END DO
1024
1025 CALL move_alloc(l_external_col, l_local_col)
1026 END DO
1027
1028 CALL release_group_dist(gd_sub_array)
1029 CALL comm_exchange%free()
1030 END block
1031 END IF
1032
1033 DEALLOCATE (l_local_col)
1034
1035 ! create the new group for the sum of the local data over all processes
1036 block
1037 TYPE(mp_comm_type) :: comm_exchange
1038 comm_exchange = fm_matrix_l%matrix_struct%context%interconnect(para_env)
1039 CALL comm_exchange%sum(fm_matrix_l%local_data)
1040 CALL comm_exchange%free()
1041 END block
1042
1043 cond_num = 1.0_dp
1044 num_small_eigen = 0
1045 IF (calc_pq_cond_num) THEN
1046 ! calculate the condition number of the (P|Q) matrix
1047 ! create a copy of the matrix
1048
1049 CALL create_matrix_l(fm_matrix_l_diag, blacs_env_l, dimen_ri, para_env_l, "fm_matrix_L_diag", fm_matrix_l_extern)
1050
1051 CALL cp_fm_to_fm(source=fm_matrix_l, destination=fm_matrix_l_diag)
1052
1053 ALLOCATE (egen_l(dimen_ri))
1054
1055 egen_l = 0.0_dp
1056 CALL cp_fm_syevx(matrix=fm_matrix_l_diag, eigenvalues=egen_l)
1057
1058 num_small_eigen = 0
1059 DO iib = 1, dimen_ri
1060 IF (abs(egen_l(iib)) < 0.001_dp) num_small_eigen = num_small_eigen + 1
1061 END DO
1062
1063 cond_num = maxval(abs(egen_l))/minval(abs(egen_l))
1064
1065 CALL cp_fm_release(fm_matrix_l_diag)
1066
1067 DEALLOCATE (egen_l)
1068 END IF
1069
1070 ! release blacs_env
1071 CALL cp_blacs_env_release(blacs_env_l)
1072
1073 CALL timestop(handle)
1074
1075 END SUBROUTINE compute_2c_integrals
1076
1077! **************************************************************************************************
1078!> \brief ...
1079!> \param matrix ...
1080!> \param num_small_evals ...
1081!> \param cond_num ...
1082!> \param eps_svd ...
1083!> \param do_inversion ...
1084!> \param para_env ...
1085! **************************************************************************************************
1086 SUBROUTINE matrix_root_with_svd(matrix, num_small_evals, cond_num, eps_svd, do_inversion, para_env)
1087 TYPE(cp_fm_type), INTENT(INOUT) :: matrix
1088 INTEGER, INTENT(OUT) :: num_small_evals
1089 REAL(kind=dp), INTENT(OUT) :: cond_num
1090 REAL(kind=dp), INTENT(IN) :: eps_svd
1091 LOGICAL, INTENT(IN) :: do_inversion
1092 TYPE(mp_para_env_type), INTENT(IN) :: para_env
1093
1094 CHARACTER(LEN=*), PARAMETER :: routinen = 'matrix_root_with_svd'
1095
1096 INTEGER :: group_size_l, handle, ii, needed_evals, &
1097 nrow, pos_max(1)
1098 INTEGER, ALLOCATABLE, DIMENSION(:) :: num_eval
1099 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: evals
1100 TYPE(cp_fm_type) :: evecs
1101 TYPE(mp_comm_type) :: comm_exchange
1102
1103 CALL timeset(routinen, handle)
1104
1105 ! Create arrays carrying eigenvalues and eigenvectors later
1106 CALL cp_fm_get_info(matrix=matrix, nrow_global=nrow)
1107
1108 ALLOCATE (evals(nrow))
1109 evals = 0
1110
1111 CALL cp_fm_create(evecs, matrix%matrix_struct)
1112
1113 ! Perform the EVD (=SVD of a positive semidefinite matrix)
1114 CALL choose_eigv_solver(matrix, evecs, evals)
1115
1116 ! Determine the number of neglectable eigenvalues assuming that the eigenvalues are in ascending order
1117 num_small_evals = 0
1118 DO ii = 1, nrow
1119 IF (evals(ii) > eps_svd) THEN
1120 num_small_evals = ii - 1
1121 EXIT
1122 END IF
1123 END DO
1124 needed_evals = nrow - num_small_evals
1125
1126 ! Get the condition number w.r.t. considered eigenvalues
1127 cond_num = evals(nrow)/evals(num_small_evals + 1)
1128
1129 ! Determine the eigenvalues of the request matrix root or its inverse
1130 evals(1:num_small_evals) = 0.0_dp
1131 IF (do_inversion) THEN
1132 evals(num_small_evals + 1:nrow) = 1.0_dp/sqrt(evals(num_small_evals + 1:nrow))
1133 ELSE
1134 evals(num_small_evals + 1:nrow) = sqrt(evals(num_small_evals + 1:nrow))
1135 END IF
1136
1137 CALL cp_fm_column_scale(evecs, evals)
1138
1139 ! As it turns out, the results in the subgroups differ.
1140 ! Thus, we have to equilibrate the results.
1141 ! Step 1: Get a communicator connecting processes with same id within subgroups
1142 ! use that group_size_L is divisible by the total number of processes (see above)
1143 group_size_l = para_env%num_pe/matrix%matrix_struct%para_env%num_pe
1144 comm_exchange = matrix%matrix_struct%context%interconnect(para_env)
1145
1146 ALLOCATE (num_eval(0:group_size_l - 1))
1147 CALL comm_exchange%allgather(num_small_evals, num_eval)
1148
1149 num_small_evals = minval(num_eval)
1150
1151 IF (num_small_evals /= maxval(num_eval)) THEN
1152 ! Step 2: Get position of maximum value
1153 pos_max = maxloc(num_eval)
1154 num_small_evals = num_eval(pos_max(1))
1155 needed_evals = nrow - num_small_evals
1156
1157 ! Step 3: Broadcast your local data to all other processes
1158 CALL comm_exchange%bcast(evecs%local_data, pos_max(1))
1159 CALL comm_exchange%bcast(cond_num, pos_max(1))
1160 END IF
1161
1162 DEALLOCATE (num_eval)
1163
1164 CALL comm_exchange%free()
1165
1166 CALL reset_size_matrix(matrix, needed_evals, matrix%matrix_struct)
1167
1168 ! Copy the needed eigenvectors
1169 CALL cp_fm_to_fm(evecs, matrix, needed_evals, num_small_evals + 1)
1170
1171 CALL cp_fm_release(evecs)
1172
1173 CALL timestop(handle)
1174
1175 END SUBROUTINE matrix_root_with_svd
1176
1177! **************************************************************************************************
1178!> \brief ...
1179!> \param matrix ...
1180!> \param new_size ...
1181!> \param fm_struct_template ...
1182! **************************************************************************************************
1183 SUBROUTINE reset_size_matrix(matrix, new_size, fm_struct_template)
1184 TYPE(cp_fm_type), INTENT(INOUT) :: matrix
1185 INTEGER, INTENT(IN) :: new_size
1186 TYPE(cp_fm_struct_type), POINTER :: fm_struct_template
1187
1188 CHARACTER(LEN=*), PARAMETER :: routinen = 'reset_size_matrix'
1189
1190 INTEGER :: handle
1191 TYPE(cp_fm_struct_type), POINTER :: fm_struct
1192
1193 CALL timeset(routinen, handle)
1194
1195 ! Choose the block sizes as large as in the template for the later steps
1196 NULLIFY (fm_struct)
1197 CALL cp_fm_struct_create(fm_struct, ncol_global=new_size, template_fmstruct=fm_struct_template, force_block=.true.)
1198
1199 CALL cp_fm_release(matrix)
1200
1201 CALL cp_fm_create(matrix, fm_struct)
1202 CALL cp_fm_set_all(matrix, 0.0_dp)
1203
1204 CALL cp_fm_struct_release(fm_struct)
1205
1206 CALL timestop(handle)
1207
1208 END SUBROUTINE reset_size_matrix
1209
1210! **************************************************************************************************
1211!> \brief ...
1212!> \param fm_matrix_L ...
1213!> \param L_local_col ...
1214!> \param para_env ...
1215!> \param my_group_L_start ...
1216!> \param my_group_L_end ...
1217!> \param dimen_RI ...
1218! **************************************************************************************************
1219 SUBROUTINE fill_fm_l_from_l_loc_non_blocking(fm_matrix_L, L_local_col, para_env, my_group_L_start, &
1220 my_group_L_end, dimen_RI)
1221 TYPE(cp_fm_type), INTENT(IN) :: fm_matrix_l
1222 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :), &
1223 INTENT(IN) :: l_local_col
1224 TYPE(mp_para_env_type), INTENT(IN) :: para_env
1225 INTEGER, INTENT(IN) :: my_group_l_start, my_group_l_end, &
1226 dimen_ri
1227
1228 CHARACTER(LEN=*), PARAMETER :: routinen = 'fill_fm_L_from_L_loc_non_blocking'
1229
1230 INTEGER :: dummy_proc, handle, handle2, i_entry_rec, i_row, iproc, j_col, lll, mmm, &
1231 ncol_local, nrow_local, proc_send, send_pcol, send_prow
1232 INTEGER, ALLOCATABLE, DIMENSION(:) :: entry_counter, num_entries_rec, &
1233 num_entries_send
1234 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
1235 TYPE(integ_mat_buffer_type), ALLOCATABLE, &
1236 DIMENSION(:) :: buffer_rec, buffer_send
1237 TYPE(mp_request_type), DIMENSION(:, :), POINTER :: req_array
1238
1239 CALL timeset(routinen, handle)
1240
1241 CALL timeset(routinen//"_1", handle2)
1242
1243 ! get info for the dest
1244 CALL cp_fm_get_info(matrix=fm_matrix_l, &
1245 nrow_local=nrow_local, &
1246 ncol_local=ncol_local, &
1247 row_indices=row_indices, &
1248 col_indices=col_indices)
1249
1250 ALLOCATE (num_entries_rec(0:para_env%num_pe - 1))
1251 ALLOCATE (num_entries_send(0:para_env%num_pe - 1))
1252
1253 num_entries_rec(:) = 0
1254 num_entries_send(:) = 0
1255
1256 dummy_proc = 0
1257
1258 ! get the process, where the elements from L_local_col have to be sent
1259 DO lll = 1, dimen_ri
1260
1261 send_prow = fm_matrix_l%matrix_struct%g2p_row(lll)
1262
1263 DO mmm = my_group_l_start, my_group_l_end
1264
1265 send_pcol = fm_matrix_l%matrix_struct%g2p_col(mmm)
1266
1267 proc_send = fm_matrix_l%matrix_struct%context%blacs2mpi(send_prow, send_pcol)
1268
1269 num_entries_send(proc_send) = num_entries_send(proc_send) + 1
1270
1271 END DO
1272
1273 END DO
1274
1275 CALL timestop(handle2)
1276
1277 CALL timeset(routinen//"_2", handle2)
1278
1279 CALL para_env%alltoall(num_entries_send, num_entries_rec, 1)
1280
1281 CALL timestop(handle2)
1282
1283 CALL timeset(routinen//"_3", handle2)
1284
1285 ! allocate buffers to send the entries and the information of the entries
1286 ALLOCATE (buffer_rec(0:para_env%num_pe - 1))
1287 ALLOCATE (buffer_send(0:para_env%num_pe - 1))
1288
1289 ! allocate data message and corresponding indices
1290 DO iproc = 0, para_env%num_pe - 1
1291
1292 ALLOCATE (buffer_rec(iproc)%msg(num_entries_rec(iproc)))
1293 buffer_rec(iproc)%msg = 0.0_dp
1294
1295 END DO
1296
1297 CALL timestop(handle2)
1298
1299 CALL timeset(routinen//"_4", handle2)
1300
1301 DO iproc = 0, para_env%num_pe - 1
1302
1303 ALLOCATE (buffer_send(iproc)%msg(num_entries_send(iproc)))
1304 buffer_send(iproc)%msg = 0.0_dp
1305
1306 END DO
1307
1308 CALL timestop(handle2)
1309
1310 CALL timeset(routinen//"_5", handle2)
1311
1312 DO iproc = 0, para_env%num_pe - 1
1313
1314 ALLOCATE (buffer_rec(iproc)%indx(num_entries_rec(iproc), 2))
1315 buffer_rec(iproc)%indx = 0
1316
1317 END DO
1318
1319 CALL timestop(handle2)
1320
1321 CALL timeset(routinen//"_6", handle2)
1322
1323 DO iproc = 0, para_env%num_pe - 1
1324
1325 ALLOCATE (buffer_send(iproc)%indx(num_entries_send(iproc), 2))
1326 buffer_send(iproc)%indx = 0
1327
1328 END DO
1329
1330 CALL timestop(handle2)
1331
1332 CALL timeset(routinen//"_7", handle2)
1333
1334 ALLOCATE (entry_counter(0:para_env%num_pe - 1))
1335 entry_counter(:) = 0
1336
1337 ! get the process, where the elements from L_local_col have to be sent and
1338 ! write the data and the info to buffer_send
1339 DO mmm = my_group_l_start, my_group_l_end
1340
1341 send_pcol = fm_matrix_l%matrix_struct%g2p_col(mmm)
1342
1343 DO lll = 1, dimen_ri
1344
1345 send_prow = fm_matrix_l%matrix_struct%g2p_row(lll)
1346
1347 proc_send = fm_matrix_l%matrix_struct%context%blacs2mpi(send_prow, send_pcol)
1348
1349 entry_counter(proc_send) = entry_counter(proc_send) + 1
1350
1351 buffer_send(proc_send)%msg(entry_counter(proc_send)) = &
1352 l_local_col(lll, mmm - my_group_l_start + 1)
1353
1354 buffer_send(proc_send)%indx(entry_counter(proc_send), 1) = lll
1355 buffer_send(proc_send)%indx(entry_counter(proc_send), 2) = mmm
1356
1357 END DO
1358
1359 END DO
1360
1361 ALLOCATE (req_array(1:para_env%num_pe, 4))
1362
1363 CALL timestop(handle2)
1364
1365 CALL timeset(routinen//"_8", handle2)
1366
1367 ! communicate the buffer
1368 CALL communicate_buffer(para_env, num_entries_rec, num_entries_send, buffer_rec, &
1369 buffer_send, req_array)
1370
1371 fm_matrix_l%local_data = 0.0_dp
1372
1373 CALL timestop(handle2)
1374
1375 CALL timeset(routinen//"_9", handle2)
1376
1377 ! fill fm_matrix_L with the entries from buffer_rec
1378 DO iproc = 0, para_env%num_pe - 1
1379 DO i_entry_rec = 1, num_entries_rec(iproc)
1380 DO j_col = 1, ncol_local
1381 IF (col_indices(j_col) == buffer_rec(iproc)%indx(i_entry_rec, 2)) THEN
1382 DO i_row = 1, nrow_local
1383 IF (row_indices(i_row) == buffer_rec(iproc)%indx(i_entry_rec, 1)) THEN
1384 fm_matrix_l%local_data(i_row, j_col) = buffer_rec(iproc)%msg(i_entry_rec)
1385 END IF
1386 END DO
1387 END IF
1388 END DO
1389 END DO
1390 END DO
1391
1392 CALL timestop(handle2)
1393
1394 CALL timeset(routinen//"_10", handle2)
1395
1396 DO iproc = 0, para_env%num_pe - 1
1397 DEALLOCATE (buffer_rec(iproc)%msg)
1398 DEALLOCATE (buffer_rec(iproc)%indx)
1399 DEALLOCATE (buffer_send(iproc)%msg)
1400 DEALLOCATE (buffer_send(iproc)%indx)
1401 END DO
1402
1403 DEALLOCATE (buffer_rec, buffer_send)
1404 DEALLOCATE (req_array)
1405 DEALLOCATE (entry_counter)
1406 DEALLOCATE (num_entries_rec, num_entries_send)
1407
1408 CALL timestop(handle2)
1409
1410 CALL timestop(handle)
1411
1412 END SUBROUTINE fill_fm_l_from_l_loc_non_blocking
1413
1414! **************************************************************************************************
1415!> \brief ...
1416!> \param fm_matrix_L ...
1417!> \param dimen_RI ...
1418!> \param do_inversion ...
1419! **************************************************************************************************
1420 SUBROUTINE cholesky_decomp(fm_matrix_L, dimen_RI, do_inversion)
1421
1422 TYPE(cp_fm_type), INTENT(IN) :: fm_matrix_l
1423 INTEGER, INTENT(IN) :: dimen_ri
1424 LOGICAL, INTENT(IN) :: do_inversion
1425
1426 CHARACTER(LEN=*), PARAMETER :: routinen = 'cholesky_decomp'
1427
1428 INTEGER :: handle, i_global, iib, info_chol, &
1429 j_global, jjb, ncol_local, nrow_local
1430 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
1431
1432 CALL timeset(routinen, handle)
1433
1434 ! do cholesky decomposition
1435 CALL cp_fm_cholesky_decompose(matrix=fm_matrix_l, n=dimen_ri, info_out=info_chol)
1436 cpassert(info_chol == 0)
1437
1438 IF (.NOT. do_inversion) THEN
1439 ! clean the lower part of the L^{-1} matrix (just to not have surprises afterwards)
1440 CALL cp_fm_get_info(matrix=fm_matrix_l, &
1441 nrow_local=nrow_local, &
1442 ncol_local=ncol_local, &
1443 row_indices=row_indices, &
1444 col_indices=col_indices)
1445 DO jjb = 1, ncol_local
1446 j_global = col_indices(jjb)
1447 DO iib = 1, nrow_local
1448 i_global = row_indices(iib)
1449 IF (j_global < i_global) fm_matrix_l%local_data(iib, jjb) = 0.0_dp
1450 END DO
1451 END DO
1452
1453 END IF
1454
1455 CALL timestop(handle)
1456
1457 END SUBROUTINE cholesky_decomp
1458
1459! **************************************************************************************************
1460!> \brief ...
1461!> \param fm_matrix_Minv_L_kpoints ...
1462!> \param fm_matrix_L_kpoints ...
1463!> \param dimen_RI ...
1464!> \param kpoints ...
1465!> \param eps_eigval_S ...
1466! **************************************************************************************************
1467 SUBROUTINE inversion_of_m_and_mult_with_chol_dec_of_v(fm_matrix_Minv_L_kpoints, fm_matrix_L_kpoints, &
1468 dimen_RI, kpoints, eps_eigval_S)
1469
1470 TYPE(cp_fm_type), DIMENSION(:, :), INTENT(IN) :: fm_matrix_minv_l_kpoints, &
1471 fm_matrix_l_kpoints
1472 INTEGER, INTENT(IN) :: dimen_ri
1473 TYPE(kpoint_type), POINTER :: kpoints
1474 REAL(kind=dp), INTENT(IN) :: eps_eigval_s
1475
1476 CHARACTER(LEN=*), PARAMETER :: routinen = 'inversion_of_M_and_mult_with_chol_dec_of_V'
1477
1478 INTEGER :: handle, ikp, nkp
1479 TYPE(cp_cfm_type) :: cfm_matrix_k_tmp, cfm_matrix_m_tmp, &
1480 cfm_matrix_v_tmp, cfm_matrix_vtrunc_tmp
1481 TYPE(cp_fm_struct_type), POINTER :: matrix_struct
1482
1483 CALL timeset(routinen, handle)
1484
1485 CALL cp_fm_get_info(fm_matrix_minv_l_kpoints(1, 1), matrix_struct=matrix_struct)
1486
1487 CALL cp_cfm_create(cfm_matrix_m_tmp, matrix_struct)
1488 CALL cp_cfm_create(cfm_matrix_v_tmp, matrix_struct)
1489 CALL cp_cfm_create(cfm_matrix_k_tmp, matrix_struct)
1490 CALL cp_cfm_create(cfm_matrix_vtrunc_tmp, matrix_struct)
1491
1492 CALL get_kpoint_info(kpoints, nkp=nkp)
1493
1494 DO ikp = 1, nkp
1495
1496 CALL cp_cfm_scale_and_add_fm(z_zero, cfm_matrix_m_tmp, z_one, fm_matrix_minv_l_kpoints(ikp, 1))
1497 CALL cp_cfm_scale_and_add_fm(z_one, cfm_matrix_m_tmp, gaussi, fm_matrix_minv_l_kpoints(ikp, 2))
1498
1499 CALL cp_cfm_scale_and_add_fm(z_zero, cfm_matrix_v_tmp, z_one, fm_matrix_l_kpoints(ikp, 1))
1500 CALL cp_cfm_scale_and_add_fm(z_one, cfm_matrix_v_tmp, gaussi, fm_matrix_l_kpoints(ikp, 2))
1501
1502 CALL cp_cfm_power(cfm_matrix_m_tmp, threshold=eps_eigval_s, exponent=-1.0_dp)
1503
1504 CALL cp_cfm_power(cfm_matrix_v_tmp, threshold=0.0_dp, exponent=0.5_dp)
1505
1506 ! save L(k) = SQRT(V(k))
1507 CALL cp_cfm_to_fm(cfm_matrix_v_tmp, fm_matrix_l_kpoints(ikp, 1), fm_matrix_l_kpoints(ikp, 2))
1508
1509 ! get K = M^(-1)*L, where L is the Cholesky decomposition of V = L^T*L
1510 CALL parallel_gemm("N", "C", dimen_ri, dimen_ri, dimen_ri, z_one, cfm_matrix_m_tmp, cfm_matrix_v_tmp, &
1511 z_zero, cfm_matrix_k_tmp)
1512
1513 ! move
1514 CALL cp_cfm_to_fm(cfm_matrix_k_tmp, fm_matrix_minv_l_kpoints(ikp, 1), fm_matrix_minv_l_kpoints(ikp, 2))
1515
1516 END DO
1517
1518 CALL cp_cfm_release(cfm_matrix_m_tmp)
1519 CALL cp_cfm_release(cfm_matrix_v_tmp)
1520 CALL cp_cfm_release(cfm_matrix_k_tmp)
1521 CALL cp_cfm_release(cfm_matrix_vtrunc_tmp)
1522
1523 CALL timestop(handle)
1524
1525 END SUBROUTINE inversion_of_m_and_mult_with_chol_dec_of_v
1526
1527! **************************************************************************************************
1528!> \brief ...
1529!> \param fm_matrix_Minv_Vtrunc_Minv ...
1530!> \param fm_matrix_Minv ...
1531!> \param qs_env ...
1532! **************************************************************************************************
1533 SUBROUTINE gamma_only_inversion_of_m_and_mult_with_chol_dec_of_vtrunc(fm_matrix_Minv_Vtrunc_Minv, &
1534 fm_matrix_Minv, qs_env)
1535
1536 TYPE(cp_fm_type), DIMENSION(:, :) :: fm_matrix_minv_vtrunc_minv, &
1537 fm_matrix_minv
1538 TYPE(qs_environment_type), POINTER :: qs_env
1539
1540 CHARACTER(LEN=*), PARAMETER :: &
1541 routinen = 'Gamma_only_inversion_of_M_and_mult_with_chol_dec_of_Vtrunc'
1542
1543 INTEGER :: dimen_ri, handle, ndep
1544 REAL(kind=dp) :: eps_eigval_s_gamma
1545 TYPE(cp_fm_type) :: fm_matrix_ri_metric_inv_work, fm_work
1546
1547 CALL timeset(routinen, handle)
1548
1549 CALL cp_fm_create(fm_work, fm_matrix_minv(1, 1)%matrix_struct)
1550 CALL cp_fm_set_all(fm_work, 0.0_dp)
1551
1552 CALL cp_fm_create(fm_matrix_ri_metric_inv_work, fm_matrix_minv(1, 1)%matrix_struct)
1553 CALL cp_fm_set_all(fm_matrix_ri_metric_inv_work, 0.0_dp)
1554
1555 CALL cp_fm_get_info(matrix=fm_matrix_minv(1, 1), nrow_global=dimen_ri)
1556
1557 eps_eigval_s_gamma = qs_env%mp2_env%ri_rpa_im_time%eps_eigval_S_Gamma
1558
1559 IF (eps_eigval_s_gamma > 1.0e-18) THEN
1560
1561 ! remove small eigenvalues
1562 CALL cp_fm_power(fm_matrix_minv(1, 1), fm_matrix_ri_metric_inv_work, -0.5_dp, &
1563 eps_eigval_s_gamma, ndep)
1564
1565 ELSE
1566
1567 CALL cholesky_decomp(fm_matrix_minv(1, 1), dimen_ri, do_inversion=.true.)
1568
1569 CALL invert_mat(fm_matrix_minv(1, 1))
1570
1571 END IF
1572
1573 CALL parallel_gemm('N', 'T', dimen_ri, dimen_ri, dimen_ri, 1.0_dp, fm_matrix_minv(1, 1), &
1574 fm_matrix_minv(1, 1), 0.0_dp, fm_matrix_ri_metric_inv_work)
1575
1576 CALL cp_fm_to_fm(fm_matrix_ri_metric_inv_work, fm_matrix_minv(1, 1))
1577
1578 CALL parallel_gemm('N', 'N', dimen_ri, dimen_ri, dimen_ri, 1.0_dp, fm_matrix_ri_metric_inv_work, &
1579 fm_matrix_minv_vtrunc_minv(1, 1), 0.0_dp, fm_work)
1580
1581 CALL parallel_gemm('N', 'N', dimen_ri, dimen_ri, dimen_ri, 1.0_dp, fm_work, &
1582 fm_matrix_ri_metric_inv_work, 0.0_dp, fm_matrix_minv_vtrunc_minv(1, 1))
1583
1584 CALL cp_fm_release(fm_work)
1585 CALL cp_fm_release(fm_matrix_ri_metric_inv_work)
1586
1587 CALL timestop(handle)
1588
1589 END SUBROUTINE gamma_only_inversion_of_m_and_mult_with_chol_dec_of_vtrunc
1590
1591! **************************************************************************************************
1592!> \brief ...
1593!> \param qs_env ...
1594!> \param trunc_coulomb ...
1595!> \param rel_cutoff_trunc_coulomb_ri_x ...
1596!> \param cell_grid ...
1597!> \param do_BvK_cell ...
1598! **************************************************************************************************
1599 SUBROUTINE trunc_coulomb_for_exchange(qs_env, trunc_coulomb, rel_cutoff_trunc_coulomb_ri_x, &
1600 cell_grid, do_BvK_cell)
1601 TYPE(qs_environment_type), POINTER :: qs_env
1602 TYPE(libint_potential_type), OPTIONAL :: trunc_coulomb
1603 REAL(kind=dp), OPTIONAL :: rel_cutoff_trunc_coulomb_ri_x
1604 INTEGER, DIMENSION(3), OPTIONAL :: cell_grid
1605 LOGICAL, OPTIONAL :: do_bvk_cell
1606
1607 CHARACTER(LEN=*), PARAMETER :: routinen = 'trunc_coulomb_for_exchange'
1608
1609 INTEGER :: handle, i_dim
1610 INTEGER, DIMENSION(3) :: periodic
1611 LOGICAL :: my_do_bvk_cell
1612 REAL(kind=dp) :: kp_fac, kp_fac_idim, my_rel_cutoff_trunc_coulomb_ri_x, &
1613 shortest_dist_cell_planes
1614 TYPE(cell_type), POINTER :: cell
1615 TYPE(kpoint_type), POINTER :: kpoints_scf
1616
1617 CALL timeset(routinen, handle)
1618
1619 NULLIFY (cell)
1620 CALL get_qs_env(qs_env, cell=cell, kpoints=kpoints_scf)
1621 CALL get_cell(cell=cell, periodic=periodic)
1622
1623 my_do_bvk_cell = .false.
1624 IF (PRESENT(do_bvk_cell)) my_do_bvk_cell = do_bvk_cell
1625 IF (my_do_bvk_cell) THEN
1626 kp_fac = 1.0e10_dp
1627 DO i_dim = 1, 3
1628 ! look for smallest k-point mesh in periodic direction
1629 IF (periodic(i_dim) == 1) THEN
1630 IF (PRESENT(cell_grid)) THEN
1631 kp_fac_idim = real(cell_grid(i_dim), kind=dp)
1632 ELSE
1633 kp_fac_idim = real(kpoints_scf%nkp_grid(i_dim), kind=dp)
1634 END IF
1635 IF (kp_fac > kp_fac_idim) kp_fac = kp_fac_idim
1636 END IF
1637 END DO
1638 ELSE
1639 kp_fac = 1.0_dp
1640 END IF
1641
1642 shortest_dist_cell_planes = 1.0e4_dp
1643 IF (periodic(1) == 1) THEN
1644 IF (shortest_dist_cell_planes > plane_distance(1, 0, 0, cell)) THEN
1645 shortest_dist_cell_planes = plane_distance(1, 0, 0, cell)
1646 END IF
1647 END IF
1648 IF (periodic(2) == 1) THEN
1649 IF (shortest_dist_cell_planes > plane_distance(0, 1, 0, cell)) THEN
1650 shortest_dist_cell_planes = plane_distance(0, 1, 0, cell)
1651 END IF
1652 END IF
1653 IF (periodic(3) == 1) THEN
1654 IF (shortest_dist_cell_planes > plane_distance(0, 0, 1, cell)) THEN
1655 shortest_dist_cell_planes = plane_distance(0, 0, 1, cell)
1656 END IF
1657 END IF
1658
1659 IF (PRESENT(rel_cutoff_trunc_coulomb_ri_x)) THEN
1660 my_rel_cutoff_trunc_coulomb_ri_x = rel_cutoff_trunc_coulomb_ri_x
1661 ELSE
1662 my_rel_cutoff_trunc_coulomb_ri_x = qs_env%mp2_env%ri_rpa_im_time%rel_cutoff_trunc_coulomb_ri_x
1663 END IF
1664
1665 IF (PRESENT(trunc_coulomb)) THEN
1666 trunc_coulomb%potential_type = do_potential_truncated
1667 trunc_coulomb%cutoff_radius = shortest_dist_cell_planes* &
1668 my_rel_cutoff_trunc_coulomb_ri_x* &
1669 kp_fac
1670 trunc_coulomb%filename = "t_c_g.dat"
1671 ! dummy
1672 trunc_coulomb%omega = 0.0_dp
1673 END IF
1674
1675 CALL timestop(handle)
1676
1677 END SUBROUTINE trunc_coulomb_for_exchange
1678
1679! **************************************************************************************************
1680!> \brief ...
1681!> \param fm_matrix_L ...
1682! **************************************************************************************************
1683 SUBROUTINE invert_mat(fm_matrix_L)
1684
1685 TYPE(cp_fm_type), INTENT(IN) :: fm_matrix_l
1686
1687 CHARACTER(LEN=*), PARAMETER :: routinen = 'invert_mat'
1688
1689 INTEGER :: handle, i_global, iib, j_global, jjb, &
1690 ncol_local, nrow_local
1691 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
1692
1693 CALL timeset(routinen, handle)
1694
1695 CALL cp_fm_triangular_invert(matrix_a=fm_matrix_l, uplo_tr='U')
1696
1697 ! clear the lower part of the L^{-1} matrix (just to have no surprises afterwards)
1698 CALL cp_fm_get_info(matrix=fm_matrix_l, &
1699 nrow_local=nrow_local, &
1700 ncol_local=ncol_local, &
1701 row_indices=row_indices, &
1702 col_indices=col_indices)
1703 ! assume upper triangular format
1704 DO jjb = 1, ncol_local
1705 j_global = col_indices(jjb)
1706 DO iib = 1, nrow_local
1707 i_global = row_indices(iib)
1708 IF (j_global < i_global) fm_matrix_l%local_data(iib, jjb) = 0.0_dp
1709 END DO
1710 END DO
1711
1712 CALL timestop(handle)
1713
1714 END SUBROUTINE invert_mat
1715
1716! **************************************************************************************************
1717!> \brief ...
1718!> \param fm_matrix_L ...
1719!> \param blacs_env_L ...
1720!> \param dimen_RI ...
1721!> \param para_env_L ...
1722!> \param name ...
1723!> \param fm_matrix_L_extern ...
1724! **************************************************************************************************
1725 SUBROUTINE create_matrix_l(fm_matrix_L, blacs_env_L, dimen_RI, para_env_L, name, fm_matrix_L_extern)
1726 TYPE(cp_fm_type), INTENT(OUT) :: fm_matrix_l
1727 TYPE(cp_blacs_env_type), POINTER :: blacs_env_l
1728 INTEGER, INTENT(IN) :: dimen_ri
1729 TYPE(mp_para_env_type), POINTER :: para_env_l
1730 CHARACTER(LEN=*), INTENT(IN) :: name
1731 TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: fm_matrix_l_extern
1732
1733 CHARACTER(LEN=*), PARAMETER :: routinen = 'create_matrix_L'
1734
1735 INTEGER :: handle
1736 TYPE(cp_fm_struct_type), POINTER :: fm_struct
1737
1738 CALL timeset(routinen, handle)
1739
1740 ! create the full matrix L defined in the L group
1741 IF (PRESENT(fm_matrix_l_extern)) THEN
1742 CALL cp_fm_create(fm_matrix_l, fm_matrix_l_extern%matrix_struct, name=name)
1743 ELSE
1744 NULLIFY (fm_struct)
1745 CALL cp_fm_struct_create(fm_struct, context=blacs_env_l, nrow_global=dimen_ri, &
1746 ncol_global=dimen_ri, para_env=para_env_l)
1747
1748 CALL cp_fm_create(fm_matrix_l, fm_struct, name=name)
1749
1750 CALL cp_fm_struct_release(fm_struct)
1751 END IF
1752
1753 CALL cp_fm_set_all(matrix=fm_matrix_l, alpha=0.0_dp)
1754
1755 CALL timestop(handle)
1756
1757 END SUBROUTINE create_matrix_l
1758
1759! **************************************************************************************************
1760!> \brief ...
1761!> \param fm_matrix_L ...
1762!> \param my_group_L_start ...
1763!> \param my_group_L_end ...
1764!> \param my_group_L_size ...
1765!> \param my_Lrows ...
1766! **************************************************************************************************
1767 SUBROUTINE grep_lcols(fm_matrix_L, &
1768 my_group_L_start, my_group_L_end, my_group_L_size, my_Lrows)
1769 TYPE(cp_fm_type), INTENT(IN) :: fm_matrix_l
1770 INTEGER, INTENT(IN) :: my_group_l_start, my_group_l_end, &
1771 my_group_l_size
1772 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :), &
1773 INTENT(OUT) :: my_lrows
1774
1775 CHARACTER(LEN=*), PARAMETER :: routinen = 'grep_Lcols'
1776
1777 INTEGER :: dimen_ri, handle, handle2, i_global, iib, j_global, jjb, max_row_col_local, &
1778 ncol_local, ncol_rec, nrow_local, nrow_rec, proc_receive_static, proc_send_static, &
1779 proc_shift
1780 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: local_col_row_info, rec_col_row_info
1781 INTEGER, DIMENSION(:), POINTER :: col_indices, col_indices_rec, &
1782 row_indices, row_indices_rec
1783 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: local_l, rec_l
1784 REAL(kind=dp), CONTIGUOUS, DIMENSION(:, :), &
1785 POINTER :: local_l_internal
1786 TYPE(mp_para_env_type), POINTER :: para_env
1787
1788 CALL timeset(routinen, handle)
1789
1790 CALL cp_fm_get_info(matrix=fm_matrix_l, &
1791 nrow_local=nrow_local, &
1792 ncol_local=ncol_local, &
1793 row_indices=row_indices, &
1794 col_indices=col_indices, &
1795 nrow_global=dimen_ri, &
1796 local_data=local_l_internal, &
1797 para_env=para_env)
1798
1799 ALLOCATE (my_lrows(dimen_ri, my_group_l_size))
1800 my_lrows = 0.0_dp
1801
1802 ALLOCATE (local_l(nrow_local, ncol_local))
1803 local_l(:, :) = local_l_internal(1:nrow_local, 1:ncol_local)
1804
1805 max_row_col_local = max(nrow_local, ncol_local)
1806 CALL para_env%max(max_row_col_local)
1807
1808 ALLOCATE (local_col_row_info(0:max_row_col_local, 2))
1809 local_col_row_info = 0
1810 ! 0,1 nrows
1811 local_col_row_info(0, 1) = nrow_local
1812 local_col_row_info(1:nrow_local, 1) = row_indices(1:nrow_local)
1813 ! 0,2 ncols
1814 local_col_row_info(0, 2) = ncol_local
1815 local_col_row_info(1:ncol_local, 2) = col_indices(1:ncol_local)
1816
1817 ALLOCATE (rec_col_row_info(0:max_row_col_local, 2))
1818
1819 ! accumulate data on my_Lrows starting from myself
1820 DO jjb = 1, ncol_local
1821 j_global = col_indices(jjb)
1822 IF (j_global >= my_group_l_start .AND. j_global <= my_group_l_end) THEN
1823 DO iib = 1, nrow_local
1824 i_global = row_indices(iib)
1825 my_lrows(i_global, j_global - my_group_l_start + 1) = local_l(iib, jjb)
1826 END DO
1827 END IF
1828 END DO
1829
1830 proc_send_static = modulo(para_env%mepos + 1, para_env%num_pe)
1831 proc_receive_static = modulo(para_env%mepos - 1, para_env%num_pe)
1832
1833 CALL timeset(routinen//"_comm", handle2)
1834
1835 DO proc_shift = 1, para_env%num_pe - 1
1836 ! first exchange information on the local data
1837 rec_col_row_info = 0
1838 CALL para_env%sendrecv(local_col_row_info, proc_send_static, rec_col_row_info, proc_receive_static)
1839 nrow_rec = rec_col_row_info(0, 1)
1840 ncol_rec = rec_col_row_info(0, 2)
1841
1842 ALLOCATE (row_indices_rec(nrow_rec))
1843 row_indices_rec = rec_col_row_info(1:nrow_rec, 1)
1844
1845 ALLOCATE (col_indices_rec(ncol_rec))
1846 col_indices_rec = rec_col_row_info(1:ncol_rec, 2)
1847
1848 ALLOCATE (rec_l(nrow_rec, ncol_rec))
1849 rec_l = 0.0_dp
1850
1851 ! then send and receive the real data
1852 CALL para_env%sendrecv(local_l, proc_send_static, rec_l, proc_receive_static)
1853
1854 ! accumulate the received data on my_Lrows
1855 DO jjb = 1, ncol_rec
1856 j_global = col_indices_rec(jjb)
1857 IF (j_global >= my_group_l_start .AND. j_global <= my_group_l_end) THEN
1858 DO iib = 1, nrow_rec
1859 i_global = row_indices_rec(iib)
1860 my_lrows(i_global, j_global - my_group_l_start + 1) = rec_l(iib, jjb)
1861 END DO
1862 END IF
1863 END DO
1864
1865 local_col_row_info(:, :) = rec_col_row_info
1866 CALL move_alloc(rec_l, local_l)
1867
1868 DEALLOCATE (col_indices_rec)
1869 DEALLOCATE (row_indices_rec)
1870 END DO
1871 CALL timestop(handle2)
1872
1873 DEALLOCATE (local_col_row_info)
1874 DEALLOCATE (rec_col_row_info)
1875 DEALLOCATE (local_l)
1876
1877 CALL timestop(handle)
1878
1879 END SUBROUTINE grep_lcols
1880
1881END MODULE mp2_ri_2c
modulo
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Definition grid_common.h:120
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:194
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition cp_dbcsr_operations.F:77
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:85
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:87
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:82
group_dist_types::create_group_dist
Definition group_dist_types.F:39
group_dist_types::get_group_dist
Definition group_dist_types.F:47
group_dist_types::release_group_dist
Definition group_dist_types.F:53
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:38
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind_set
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
Definition atomic_kind_types.F:223
basis_set_types
Definition basis_set_types.F:15
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cell_types::get_cell
subroutine, public get_cell(cell, alpha, beta, gamma, deth, orthorhombic, abc, periodic, h, h_inv, symmetry_id, tag)
Get informations about a simulation cell.
Definition cell_types.F:195
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
constants_operator
constants for the different operators of the 2c-integrals
Definition constants_operator.F:14
constants_operator::operator_coulomb
integer, parameter, public operator_coulomb
Definition constants_operator.F:22
cp_blacs_env
methods related to the blacs parallel environment
Definition cp_blacs_env.F:15
cp_blacs_env::cp_blacs_env_release
subroutine, public cp_blacs_env_release(blacs_env)
releases the given blacs_env
Definition cp_blacs_env.F:274
cp_blacs_env::cp_blacs_env_create
subroutine, public cp_blacs_env_create(blacs_env, para_env, blacs_grid_layout, blacs_repeatable, row_major, grid_2d)
allocates and initializes a type that represent a blacs context
Definition cp_blacs_env.F:123
cp_cfm_basic_linalg
Basic linear algebra operations for complex full matrices.
Definition cp_cfm_basic_linalg.F:16
cp_cfm_basic_linalg::cp_cfm_scale_and_add_fm
subroutine, public cp_cfm_scale_and_add_fm(alpha, matrix_a, beta, matrix_b)
Scale and add two BLACS matrices (a = alpha*a + beta*b). where b is a real matrix (adapted from cp_cf...
Definition cp_cfm_basic_linalg.F:351
cp_cfm_types
Represents a complex full matrix distributed on many processors.
Definition cp_cfm_types.F:12
cp_cfm_types::cp_cfm_create
subroutine, public cp_cfm_create(matrix, matrix_struct, name)
Creates a new full matrix with the given structure.
Definition cp_cfm_types.F:121
cp_cfm_types::cp_cfm_release
subroutine, public cp_cfm_release(matrix)
Releases a full matrix.
Definition cp_cfm_types.F:155
cp_cfm_types::cp_cfm_to_fm
subroutine, public cp_cfm_to_fm(msource, mtargetr, mtargeti)
Copy real and imaginary parts of a complex full matrix into separate real-value full matrices.
Definition cp_cfm_types.F:761
cp_control_types
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Definition cp_control_types.F:12
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_distribution_release
subroutine, public dbcsr_distribution_release(dist)
...
Definition cp_dbcsr_api.F:642
cp_dbcsr_api::dbcsr_deallocate_matrix
subroutine, public dbcsr_deallocate_matrix(matrix)
...
Definition cp_dbcsr_api.F:235
cp_dbcsr_api::dbcsr_desymmetrize
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
Definition cp_dbcsr_api.F:524
cp_dbcsr_api::dbcsr_copy
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
Definition cp_dbcsr_api.F:370
cp_dbcsr_api::dbcsr_set
subroutine, public dbcsr_set(matrix, alpha)
...
Definition cp_dbcsr_api.F:1194
cp_dbcsr_api::dbcsr_release
subroutine, public dbcsr_release(matrix)
...
Definition cp_dbcsr_api.F:1132
cp_dbcsr_contrib
Definition cp_dbcsr_contrib.F:8
cp_dbcsr_contrib::dbcsr_reserve_all_blocks
subroutine, public dbcsr_reserve_all_blocks(matrix)
Reserves all blocks.
Definition cp_dbcsr_contrib.F:295
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_dist2d_to_dist
subroutine, public cp_dbcsr_dist2d_to_dist(dist2d, dist)
Creates a DBCSR distribution from a distribution_2d.
Definition cp_dbcsr_operations.F:419
cp_dbcsr_operations::copy_dbcsr_to_fm
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:214
cp_eri_mme_interface
Interface to Minimax-Ewald method for periodic ERI's to be used in CP2K.
Definition cp_eri_mme_interface.F:12
cp_fm_basic_linalg
Basic linear algebra operations for full matrices.
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_column_scale
subroutine, public cp_fm_column_scale(matrixa, scaling)
scales column i of matrix a with scaling(i)
Definition cp_fm_basic_linalg.F:1879
cp_fm_basic_linalg::cp_fm_triangular_invert
subroutine, public cp_fm_triangular_invert(matrix_a, uplo_tr)
inverts a triangular matrix
Definition cp_fm_basic_linalg.F:2191
cp_fm_cholesky
various cholesky decomposition related routines
Definition cp_fm_cholesky.F:14
cp_fm_cholesky::cp_fm_cholesky_decompose
subroutine, public cp_fm_cholesky_decompose(matrix, n, info_out)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U,...
Definition cp_fm_cholesky.F:63
cp_fm_diag
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
cp_fm_diag::cp_fm_power
subroutine, public cp_fm_power(matrix, work, exponent, threshold, n_dependent, verbose, eigvals)
...
Definition cp_fm_diag.F:997
cp_fm_diag::choose_eigv_solver
subroutine, public choose_eigv_solver(matrix, eigenvectors, eigenvalues, info)
Choose the Eigensolver depending on which library is available ELPA seems to be unstable for small sy...
Definition cp_fm_diag.F:234
cp_fm_diag::cp_fm_syevx
subroutine, public cp_fm_syevx(matrix, eigenvectors, eigenvalues, neig, work_syevx)
compute eigenvalues and optionally eigenvectors of a real symmetric matrix using scalapack....
Definition cp_fm_diag.F:667
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:351
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_copy_general
subroutine, public cp_fm_copy_general(source, destination, para_env)
General copy of a fm matrix to another fm matrix. Uses non-blocking MPI rather than ScaLAPACK.
Definition cp_fm_types.F:1537
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:1009
cp_fm_types::cp_fm_set_all
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
Definition cp_fm_types.F:528
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Definition cp_fm_types.F:164
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
group_dist_types
Types to describe group distributions.
Definition group_dist_types.F:14
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_eri_os
integer, parameter, public do_eri_os
Definition input_constants.F:1258
input_constants::do_eri_mme
integer, parameter, public do_eri_mme
Definition input_constants.F:1258
input_constants::do_potential_truncated
integer, parameter, public do_potential_truncated
Definition input_constants.F:811
input_constants::do_potential_id
integer, parameter, public do_potential_id
Definition input_constants.F:811
input_constants::do_potential_coulomb
integer, parameter, public do_potential_coulomb
Definition input_constants.F:811
input_constants::do_potential_long
integer, parameter, public do_potential_long
Definition input_constants.F:811
input_constants::do_eri_gpw
integer, parameter, public do_eri_gpw
Definition input_constants.F:1258
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kpoint_coulomb_2c
Routines to compute the Coulomb integral V_(alpha beta)(k) for a k-point k using lattice summation in...
Definition kpoint_coulomb_2c.F:15
kpoint_coulomb_2c::build_2c_coulomb_matrix_kp
subroutine, public build_2c_coulomb_matrix_kp(matrix_v_kp, kpoints, basis_type, cell, particle_set, qs_kind_set, atomic_kind_set, size_lattice_sum, operator_type, ikp_start, ikp_end)
...
Definition kpoint_coulomb_2c.F:76
kpoint_methods
Routines needed for kpoint calculation.
Definition kpoint_methods.F:15
kpoint_methods::rskp_transform
subroutine, public rskp_transform(rmatrix, cmatrix, rsmat, ispin, xkp, cell_to_index, sab_nl, is_complex, rs_sign)
Transformation of real space matrices to a kpoint.
Definition kpoint_methods.F:812
kpoint_methods::kpoint_init_cell_index
subroutine, public kpoint_init_cell_index(kpoint, sab_nl, para_env, dft_control)
Generates the mapping of cell indices and linear RS index CELL (0,0,0) is always mapped to index 1.
Definition kpoint_methods.F:679
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:365
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::compare_potential_types
pure logical function, public compare_potential_types(potential1, potential2)
Helper function to compare libint_potential_types.
Definition libint_2c_3c.F:1062
machine
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
machine::m_flush
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition machine.F:136
mathconstants
Definition of mathematical constants and functions.
Definition mathconstants.F:16
mathconstants::z_one
complex(kind=dp), parameter, public z_one
Definition mathconstants.F:143
mathconstants::gaussi
complex(kind=dp), parameter, public gaussi
Definition mathconstants.F:142
mathconstants::z_zero
complex(kind=dp), parameter, public z_zero
Definition mathconstants.F:143
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
message_passing::mp_para_env_release
subroutine, public mp_para_env_release(para_env)
releases the para object (to be called when you don't want anymore the shared copy of this object)
Definition message_passing.F:2026
mp2_eri_gpw
Routines to calculate 2c- and 3c-integrals for RI with GPW.
Definition mp2_eri_gpw.F:13
mp2_eri_gpw::mp2_eri_2c_integrate_gpw
subroutine, public mp2_eri_2c_integrate_gpw(qs_env, para_env_sub, my_group_l_start, my_group_l_end, natom, potential_parameter, sab_orb_sub, l_local_col, kind_of)
Integrates the potential of an RI function.
Definition mp2_eri_gpw.F:167
mp2_eri
Interface to direct methods for electron repulsion integrals for MP2.
Definition mp2_eri.F:12
mp2_eri::mp2_eri_2c_integrate
subroutine, public mp2_eri_2c_integrate(param, potential_parameter, para_env, qs_env, basis_type_a, basis_type_b, hab, first_b, last_b, eri_method, pab, force_a, force_b, hdab, hadb, reflection_z_a, reflection_z_b, do_reflection_a, do_reflection_b)
high-level integration routine for 2c integrals over CP2K basis sets. Contiguous column-wise distribu...
Definition mp2_eri.F:115
mp2_ri_2c
Framework for 2c-integrals for RI.
Definition mp2_ri_2c.F:14
mp2_ri_2c::get_2c_integrals
subroutine, public get_2c_integrals(qs_env, eri_method, eri_param, para_env, para_env_sub, mp2_memory, my_lrows, my_vrows, fm_matrix_pq, ngroup, color_sub, dimen_ri, dimen_ri_red, kpoints, my_group_l_size, my_group_l_start, my_group_l_end, gd_array, calc_pq_cond_num, do_svd, eps_svd, potential, ri_metric, fm_matrix_l_kpoints, fm_matrix_minv_l_kpoints, fm_matrix_minv, fm_matrix_minv_vtrunc_minv, do_im_time, do_kpoints, mp2_eps_pgf_orb_s, qs_kind_set, sab_orb_sub, calc_forces, unit_nr)
...
Definition mp2_ri_2c.F:167
mp2_ri_2c::trunc_coulomb_for_exchange
subroutine, public trunc_coulomb_for_exchange(qs_env, trunc_coulomb, rel_cutoff_trunc_coulomb_ri_x, cell_grid, do_bvk_cell)
...
Definition mp2_ri_2c.F:1601
mp2_ri_2c::inversion_of_m_and_mult_with_chol_dec_of_v
subroutine, public inversion_of_m_and_mult_with_chol_dec_of_v(fm_matrix_minv_l_kpoints, fm_matrix_l_kpoints, dimen_ri, kpoints, eps_eigval_s)
...
Definition mp2_ri_2c.F:1469
mp2_ri_2c::ri_2c_integral_mat
subroutine, public ri_2c_integral_mat(qs_env, fm_matrix_minv_l_kpoints, fm_matrix_l, dimen_ri, ri_metric, do_kpoints, kpoints, put_mat_ks_env, regularization_ri, ikp_ext, do_build_cell_index)
...
Definition mp2_ri_2c.F:564
mp2_types
Types needed for MP2 calculations.
Definition mp2_types.F:14
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
particle_methods
Define methods related to particle_type.
Definition particle_methods.F:14
particle_methods::get_particle_set
subroutine, public get_particle_set(particle_set, qs_kind_set, first_sgf, last_sgf, nsgf, nmao, basis)
Get the components of a particle set.
Definition particle_methods.F:98
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
qs_environment_types
Definition qs_environment_types.F:14
qs_environment_types::get_qs_env
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, 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, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
Definition qs_environment_types.F:526
qs_integral_utils
Some utility functions for the calculation of integrals.
Definition qs_integral_utils.F:14
qs_integral_utils::basis_set_list_setup
subroutine, public basis_set_list_setup(basis_set_list, basis_type, qs_kind_set)
Set up an easy accessible list of the basis sets for all kinds.
Definition qs_integral_utils.F:149
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:76
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, 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_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition qs_kind_types.F:466
qs_ks_types
Definition qs_ks_types.F:15
qs_ks_types::set_ks_env
subroutine, public set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, kpoints, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
Definition qs_ks_types.F:578
qs_ks_types::get_ks_env
subroutine, public get_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, rho, rho_xc, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, kpoints, do_kpoints, atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, particle_set, energy, force, local_particles, local_molecules, molecule_kind_set, molecule_set, subsys, cp_subsys, virial, results, atprop, nkind, natom, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env, nelectron_total, nelectron_spin)
...
Definition qs_ks_types.F:336
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:1095
qs_tensors
Utility methods to build 3-center integral tensors of various types.
Definition qs_tensors.F:11
qs_tensors::build_2c_integrals
subroutine, public build_2c_integrals(t2c, filter_eps, qs_env, nl_2c, basis_i, basis_j, potential_parameter, do_kpoints, do_hfx_kpoints, ext_kpoints, regularization_ri)
...
Definition qs_tensors.F:3156
qs_tensors::build_2c_neighbor_lists
subroutine, public build_2c_neighbor_lists(ij_list, basis_i, basis_j, potential_parameter, name, qs_env, sym_ij, molecular, dist_2d, pot_to_rad)
Build 2-center neighborlists adapted to different operators This mainly wraps build_neighbor_lists fo...
Definition qs_tensors.F:144
rpa_communication
Auxiliary routines necessary to redistribute an fm_matrix from a given blacs_env to another.
Definition rpa_communication.F:14
rpa_communication::communicate_buffer
subroutine, public communicate_buffer(para_env, num_entries_rec, num_entries_send, buffer_rec, buffer_send, req_array, do_indx, do_msg)
...
Definition rpa_communication.F:434
rpa_gw_kpoints_util
Routines treating GW and RPA calculations with kpoints.
Definition rpa_gw_kpoints_util.F:13
rpa_gw_kpoints_util::cp_cfm_power
subroutine, public cp_cfm_power(matrix, threshold, exponent, min_eigval)
...
Definition rpa_gw_kpoints_util.F:274
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
basis_set_types::gto_basis_set_p_type
Definition basis_set_types.F:94
basis_set_types::gto_basis_set_type
Definition basis_set_types.F:72
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
cp_blacs_env::cp_blacs_env_type
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Definition cp_blacs_env.F:53
cp_cfm_types::cp_cfm_type
Represent a complex full matrix.
Definition cp_cfm_types.F:68
cp_control_types::dft_control_type
Definition cp_control_types.F:602
cp_dbcsr_api::dbcsr_distribution_type
Definition cp_dbcsr_api.F:182
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
cp_dbcsr_api::dbcsr_type
Definition cp_dbcsr_api.F:176
cp_eri_mme_interface::cp_eri_mme_param
Definition cp_eri_mme_interface.F:79
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:113
distribution_2d_types::distribution_2d_type
distributes pairs on a 2d grid of processors
Definition distribution_2d_types.F:62
group_dist_types::group_dist_d1_type
Definition group_dist_types.F:35
kpoint_types::kpoint_type
Contains information about kpoints.
Definition kpoint_types.F:150
libint_2c_3c::libint_potential_type
Definition libint_2c_3c.F:67
message_passing::mp_comm_type
Definition message_passing.F:189
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:763
message_passing::mp_request_type
Definition message_passing.F:623
mp2_types::integ_mat_buffer_type
Definition mp2_types.F:409
particle_types::particle_type
Definition particle_types.F:35
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:220
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:177
qs_neighbor_list_types::neighbor_list_set_p_type
Definition qs_neighbor_list_types.F:67