35 USE dbcsr_api,
ONLY: dbcsr_add,&
109#include "./base/base_uses.f90"
115 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'rpa_main'
170 SUBROUTINE rpa_ri_compute_en(qs_env, Erpa, mp2_env, BIb_C, BIb_C_gw, BIb_C_bse_ij, BIb_C_bse_ab, &
171 para_env, para_env_sub, color_sub, &
172 gd_array, gd_B_virtual, gd_B_all, gd_B_occ_bse, gd_B_virt_bse, &
173 mo_coeff, fm_matrix_PQ, fm_matrix_L_kpoints, fm_matrix_Minv_L_kpoints, &
174 fm_matrix_Minv, fm_matrix_Minv_Vtrunc_Minv, kpoints, &
175 Eigenval, nmo, homo, dimen_RI, dimen_RI_red, gw_corr_lev_occ, gw_corr_lev_virt, &
176 unit_nr, do_ri_sos_laplace_mp2, my_do_gw, do_im_time, do_bse, matrix_s, &
177 mat_munu, mat_P_global, &
178 t_3c_M, t_3c_O, t_3c_O_compressed, t_3c_O_ind, &
179 starts_array_mc, ends_array_mc, &
180 starts_array_mc_block, ends_array_mc_block, calc_forces)
183 REAL(kind=
dp),
INTENT(OUT) :: erpa
184 TYPE(
mp2_type),
INTENT(INOUT) :: mp2_env
186 INTENT(INOUT) :: bib_c, bib_c_gw
187 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :), &
188 INTENT(INOUT) :: bib_c_bse_ij, bib_c_bse_ab
190 INTEGER,
INTENT(INOUT) :: color_sub
193 INTENT(INOUT) :: gd_b_virtual
195 TYPE(
cp_fm_type),
DIMENSION(:),
INTENT(IN) :: mo_coeff
197 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_matrix_l_kpoints, &
198 fm_matrix_minv_l_kpoints, &
200 fm_matrix_minv_vtrunc_minv
202 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :), &
203 INTENT(INOUT) :: eigenval
204 INTEGER,
INTENT(IN) :: nmo
205 INTEGER,
DIMENSION(:),
INTENT(IN) :: homo
206 INTEGER,
INTENT(IN) :: dimen_ri, dimen_ri_red
207 INTEGER,
DIMENSION(:),
INTENT(IN) :: gw_corr_lev_occ, gw_corr_lev_virt
208 INTEGER,
INTENT(IN) :: unit_nr
209 LOGICAL,
INTENT(IN) :: do_ri_sos_laplace_mp2, my_do_gw, &
211 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s
212 TYPE(dbcsr_p_type),
INTENT(IN) :: mat_munu
213 TYPE(dbcsr_p_type),
INTENT(INOUT) :: mat_p_global
214 TYPE(dbt_type) :: t_3c_m
215 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:, :) :: t_3c_o
217 DIMENSION(:, :, :),
INTENT(INOUT) :: t_3c_o_compressed
219 DIMENSION(:, :, :),
INTENT(INOUT) :: t_3c_o_ind
220 INTEGER,
ALLOCATABLE,
DIMENSION(:),
INTENT(IN) :: starts_array_mc, ends_array_mc, &
221 starts_array_mc_block, &
223 LOGICAL,
INTENT(IN) :: calc_forces
225 CHARACTER(LEN=*),
PARAMETER :: routinen =
'rpa_ri_compute_en'
227 INTEGER :: best_integ_group_size, best_num_integ_point, color_rpa_group, dimen_homo_square, &
228 dimen_nm_gw, dimen_virt_square, handle, handle2, handle3, ierr, iib, &
229 input_num_integ_groups, integ_group_size, ispin, jjb, min_integ_group_size, &
230 my_ab_comb_bse_end, my_ab_comb_bse_size, my_ab_comb_bse_start, my_group_l_end, &
231 my_group_l_size, my_group_l_start, my_ij_comb_bse_end, my_ij_comb_bse_size, &
232 my_ij_comb_bse_start, my_nm_gw_end, my_nm_gw_size, my_nm_gw_start, ncol_block_mat, &
233 ngroup, nrow_block_mat, nspins, num_integ_group, num_integ_points, pos_integ_group
234 INTEGER(KIND=int_8) :: mem
235 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: dimen_ia, my_ia_end, my_ia_size, &
237 LOGICAL :: do_kpoints_from_gamma, do_minimax_quad, &
238 my_open_shell, skip_integ_group_opt
239 REAL(kind=
dp) :: allowed_memory, avail_mem, e_range, emax, emin, mem_for_iak, mem_for_qk, &
240 mem_min, mem_per_group, mem_per_rank, mem_per_repl, mem_real
241 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: eigenval_kp
242 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:) :: fm_mat_q, fm_mat_q_gemm, fm_mat_s, &
244 TYPE(
cp_fm_type),
DIMENSION(1) :: fm_mat_r_gw, fm_mat_s_ab_bse, &
248 DIMENSION(:) :: bib_c_2d, bib_c_2d_gw
251 CALL timeset(routinen, handle)
256 IF (mp2_env%ri_rpa%do_ri_axk)
THEN
259 IF (mp2_env%ri_rpa%do_rse)
THEN
275 my_open_shell = (nspins == 2)
276 ALLOCATE (virtual(nspins), dimen_ia(nspins), my_ia_end(nspins), my_ia_start(nspins), my_ia_size(nspins))
277 virtual(:) = nmo - homo(:)
278 dimen_ia(:) = virtual(:)*homo(:)
280 ALLOCATE (eigenval_kp(nmo, 1, nspins))
281 eigenval_kp(:, 1, :) = eigenval(:, :)
283 IF (do_im_time) mp2_env%ri_rpa%minimax_quad = .true.
284 do_minimax_quad = mp2_env%ri_rpa%minimax_quad
286 IF (do_ri_sos_laplace_mp2)
THEN
287 num_integ_points = mp2_env%ri_laplace%n_quadrature
288 input_num_integ_groups = mp2_env%ri_laplace%num_integ_groups
291 e_range = mp2_env%e_range
292 IF (mp2_env%e_range <= 1.0_dp .OR. mp2_env%e_gap <= 0.0_dp)
THEN
296 IF (homo(ispin) > 0)
THEN
297 emin = min(emin, 2.0_dp*(eigenval(homo(ispin) + 1, ispin) - eigenval(homo(ispin), ispin)))
298 emax = max(emax, 2.0_dp*(maxval(eigenval(:, ispin)) - minval(eigenval(:, ispin))))
303 IF (e_range < 2.0_dp) e_range = 2.0_dp
307 jjb = num_integ_points - 1
309 num_integ_points = num_integ_points - 1
315 cpassert(num_integ_points >= 1)
317 num_integ_points = mp2_env%ri_rpa%rpa_num_quad_points
318 input_num_integ_groups = mp2_env%ri_rpa%rpa_num_integ_groups
320 IF (my_do_gw .AND. do_minimax_quad)
THEN
321 IF (num_integ_points > 34)
THEN
323 CALL cp_warn(__location__, &
324 "The required number of quadrature point exceeds the maximum possible in the "// &
325 "Minimax quadrature scheme. The number of quadrature point has been reset to 30.")
326 num_integ_points = 30
329 IF (do_minimax_quad .AND. num_integ_points > 20)
THEN
331 CALL cp_warn(__location__, &
332 "The required number of quadrature point exceeds the maximum possible in the "// &
333 "Minimax quadrature scheme. The number of quadrature point has been reset to 20.")
334 num_integ_points = 20
338 allowed_memory = mp2_env%mp2_memory
340 CALL get_group_dist(gd_array, color_sub, my_group_l_start, my_group_l_end, my_group_l_size)
342 ngroup = para_env%num_pe/para_env_sub%num_pe
345 IF (do_im_time .OR. mp2_env%ri_g0w0%do_periodic .OR. do_bse)
THEN
347 integ_group_size = ngroup
348 best_num_integ_point = num_integ_points
354 mem_for_iak = real(sum(dimen_ia), kind=
dp)*dimen_ri_red*8.0_dp/(1024_dp**2)
355 mem_for_qk = real(dimen_ri_red, kind=
dp)*nspins*dimen_ri_red*8.0_dp/(1024_dp**2)
358 mem_real = (mem + 1024*1024 - 1)/(1024*1024)
359 CALL para_env%min(mem_real)
361 mem_per_rank = 0.0_dp
364 mem_per_repl = mem_for_iak
366 mem_per_repl = mem_per_repl + 2.0_dp*mem_for_qk
368 IF (calc_forces)
CALL rpa_grad_needed_mem(homo, virtual, dimen_ri_red, mem_per_rank, mem_per_repl, do_ri_sos_laplace_mp2)
370 mem_min = mem_per_repl/para_env%num_pe + mem_per_rank
372 IF (unit_nr > 0)
THEN
373 WRITE (unit_nr,
'(T3,A,T68,F9.2,A4)')
'RI_INFO| Minimum required memory per MPI process:', mem_min,
' MiB'
374 WRITE (unit_nr,
'(T3,A,T68,F9.2,A4)')
'RI_INFO| Available memory per MPI process:', mem_real,
' MiB'
378 mem_real = min(mem_real, allowed_memory)
380 mem_real = mem_real - mem_per_rank
382 mem_per_group = mem_real*para_env_sub%num_pe
385 skip_integ_group_opt = .false.
388 IF (input_num_integ_groups > 0)
THEN
389 IF (mod(num_integ_points, input_num_integ_groups) == 0)
THEN
390 IF (mod(ngroup, input_num_integ_groups) == 0)
THEN
391 best_integ_group_size = ngroup/input_num_integ_groups
392 best_num_integ_point = num_integ_points/input_num_integ_groups
393 skip_integ_group_opt = .true.
395 IF (unit_nr > 0)
WRITE (unit_nr,
'(T3,A)')
'Total number of groups not multiple of NUM_INTEG_GROUPS'
398 IF (unit_nr > 0)
WRITE (unit_nr,
'(T3,A)')
'NUM_QUAD_POINTS not multiple of NUM_INTEG_GROUPS'
402 IF (.NOT. skip_integ_group_opt)
THEN
403 best_integ_group_size = ngroup
404 best_num_integ_point = num_integ_points
406 min_integ_group_size = max(1, ngroup/num_integ_points)
408 integ_group_size = min_integ_group_size - 1
409 DO iib = min_integ_group_size + 1, ngroup
410 integ_group_size = integ_group_size + 1
413 IF (mod(ngroup, integ_group_size) /= 0) cycle
416 avail_mem = integ_group_size*mem_per_group
417 IF (avail_mem < mem_per_repl) cycle
420 num_integ_group = ngroup/integ_group_size
421 IF (mod(num_integ_points, num_integ_group) /= 0) cycle
423 best_num_integ_point = num_integ_points/num_integ_group
424 best_integ_group_size = integ_group_size
431 integ_group_size = best_integ_group_size
435 IF (unit_nr > 0 .AND. .NOT. do_im_time)
THEN
436 IF (do_ri_sos_laplace_mp2)
THEN
437 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
438 "RI_INFO| Group size for laplace numerical integration:", integ_group_size*para_env_sub%num_pe
439 WRITE (unit=unit_nr, fmt=
"(T3,A)") &
440 "INTEG_INFO| MINIMAX approximation"
441 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
442 "INTEG_INFO| Number of integration points:", num_integ_points
443 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
444 "INTEG_INFO| Number of integration points per Laplace group:", best_num_integ_point
446 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
447 "RI_INFO| Group size for frequency integration:", integ_group_size*para_env_sub%num_pe
448 IF (do_minimax_quad)
THEN
449 WRITE (unit=unit_nr, fmt=
"(T3,A)") &
450 "INTEG_INFO| MINIMAX quadrature"
452 WRITE (unit=unit_nr, fmt=
"(T3,A)") &
453 "INTEG_INFO| Clenshaw-Curtius quadrature"
455 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
456 "INTEG_INFO| Number of integration points:", num_integ_points
457 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
458 "INTEG_INFO| Number of integration points per RPA group:", best_num_integ_point
463 num_integ_group = ngroup/integ_group_size
465 pos_integ_group = mod(color_sub, integ_group_size)
466 color_rpa_group = color_sub/integ_group_size
468 CALL timeset(routinen//
"_reorder", handle2)
472 ALLOCATE (bib_c_2d(nspins))
473 IF (.NOT. do_im_time)
THEN
478 CALL calculate_bib_c_2d(bib_c_2d(ispin)%array, bib_c(ispin)%array, para_env_sub, dimen_ia(ispin), &
479 homo(ispin), virtual(ispin), gd_b_virtual(ispin), &
480 my_ia_size(ispin), my_ia_start(ispin), my_ia_end(ispin), my_group_l_size)
482 DEALLOCATE (bib_c(ispin)%array)
489 ALLOCATE (bib_c_2d_gw(nspins))
491 CALL timeset(routinen//
"_reorder_gw", handle3)
493 dimen_nm_gw = nmo*(gw_corr_lev_occ(1) + gw_corr_lev_virt(1))
497 CALL calculate_bib_c_2d(bib_c_2d_gw(ispin)%array, bib_c_gw(ispin)%array, para_env_sub, dimen_nm_gw, &
498 gw_corr_lev_occ(ispin) + gw_corr_lev_virt(ispin), nmo, gd_b_all, &
499 my_nm_gw_size, my_nm_gw_start, my_nm_gw_end, my_group_l_size)
500 DEALLOCATE (bib_c_gw(ispin)%array)
505 CALL timestop(handle3)
512 CALL timeset(routinen//
"_reorder_bse1", handle3)
514 dimen_homo_square = homo(1)**2
516 CALL calculate_bib_c_2d(bib_c_2d_bse_ij(1)%array, bib_c_bse_ij, para_env_sub, dimen_homo_square, &
517 homo(1), homo(1), gd_b_occ_bse, &
518 my_ij_comb_bse_size, my_ij_comb_bse_start, my_ij_comb_bse_end, my_group_l_size)
520 DEALLOCATE (bib_c_bse_ij)
523 CALL timestop(handle3)
525 CALL timeset(routinen//
"_reorder_bse2", handle3)
527 dimen_virt_square = virtual(1)**2
529 CALL calculate_bib_c_2d(bib_c_2d_bse_ab(1)%array, bib_c_bse_ab, para_env_sub, dimen_virt_square, &
530 virtual(1), virtual(1), gd_b_virt_bse, &
531 my_ab_comb_bse_size, my_ab_comb_bse_start, my_ab_comb_bse_end, my_group_l_size)
533 DEALLOCATE (bib_c_bse_ab)
536 CALL timestop(handle3)
540 CALL timestop(handle2)
542 IF (num_integ_group > 1)
THEN
543 ALLOCATE (para_env_rpa)
544 CALL para_env_rpa%from_split(para_env, color_rpa_group)
546 para_env_rpa => para_env
553 ALLOCATE (fm_mat_q(nspins), fm_mat_q_gemm(1), fm_mat_s(1))
555 ALLOCATE (fm_mat_q(nspins), fm_mat_q_gemm(nspins), fm_mat_s(nspins))
558 CALL create_integ_mat(bib_c_2d, para_env, para_env_sub, color_sub, ngroup, integ_group_size, &
559 dimen_ri_red, dimen_ia, color_rpa_group, &
560 mp2_env%block_size_row, mp2_env%block_size_col, unit_nr, &
561 my_ia_size, my_ia_start, my_ia_end, &
562 my_group_l_size, my_group_l_start, my_group_l_end, &
563 para_env_rpa, fm_mat_s, nrow_block_mat, ncol_block_mat, &
564 dimen_ia_for_block_size=dimen_ia(1), &
565 do_im_time=do_im_time, fm_mat_q_gemm=fm_mat_q_gemm, fm_mat_q=fm_mat_q, qs_env=qs_env)
567 DEALLOCATE (bib_c_2d, my_ia_end, my_ia_size, my_ia_start)
570 ALLOCATE (fm_mat_s_gw(nspins))
571 IF (my_do_gw .AND. .NOT. do_im_time)
THEN
573 CALL create_integ_mat(bib_c_2d_gw, para_env, para_env_sub, color_sub, ngroup, integ_group_size, &
574 dimen_ri_red, [dimen_nm_gw, dimen_nm_gw], color_rpa_group, &
575 mp2_env%block_size_row, mp2_env%block_size_col, unit_nr, &
576 [my_nm_gw_size, my_nm_gw_size], [my_nm_gw_start, my_nm_gw_start], [my_nm_gw_end, my_nm_gw_end], &
577 my_group_l_size, my_group_l_start, my_group_l_end, &
578 para_env_rpa, fm_mat_s_gw, nrow_block_mat, ncol_block_mat, &
579 fm_mat_q(1)%matrix_struct%context, fm_mat_q(1)%matrix_struct%context, &
580 fm_mat_q=fm_mat_r_gw)
581 DEALLOCATE (bib_c_2d_gw)
587 CALL create_integ_mat(bib_c_2d_bse_ij, para_env, para_env_sub, color_sub, ngroup, integ_group_size, &
588 dimen_ri_red, [dimen_homo_square], color_rpa_group, &
589 mp2_env%block_size_row, mp2_env%block_size_col, unit_nr, &
590 [my_ij_comb_bse_size], [my_ij_comb_bse_start], [my_ij_comb_bse_end], &
591 my_group_l_size, my_group_l_start, my_group_l_end, &
592 para_env_rpa, fm_mat_s_ij_bse, nrow_block_mat, ncol_block_mat, &
593 fm_mat_q(1)%matrix_struct%context, fm_mat_q(1)%matrix_struct%context)
595 CALL create_integ_mat(bib_c_2d_bse_ab, para_env, para_env_sub, color_sub, ngroup, integ_group_size, &
596 dimen_ri_red, [dimen_virt_square], color_rpa_group, &
597 mp2_env%block_size_row, mp2_env%block_size_col, unit_nr, &
598 [my_ab_comb_bse_size], [my_ab_comb_bse_start], [my_ab_comb_bse_end], &
599 my_group_l_size, my_group_l_start, my_group_l_end, &
600 para_env_rpa, fm_mat_s_ab_bse, nrow_block_mat, ncol_block_mat, &
601 fm_mat_q(1)%matrix_struct%context, fm_mat_q(1)%matrix_struct%context)
605 do_kpoints_from_gamma = qs_env%mp2_env%ri_rpa_im_time%do_kpoints_from_Gamma
606 IF (do_kpoints_from_gamma)
THEN
612 CALL rpa_num_int(qs_env, erpa, mp2_env, para_env, para_env_rpa, para_env_sub, unit_nr, &
613 homo, virtual, dimen_ri, dimen_ri_red, dimen_ia, dimen_nm_gw, &
614 eigenval_kp, num_integ_points, num_integ_group, color_rpa_group, &
615 fm_matrix_pq, fm_mat_s, fm_mat_q_gemm, fm_mat_q, fm_mat_s_gw, fm_mat_r_gw(1), &
616 fm_mat_s_ij_bse(1), fm_mat_s_ab_bse(1), &
617 my_do_gw, do_bse, gw_corr_lev_occ, gw_corr_lev_virt, &
619 do_im_time, mo_coeff, &
620 fm_matrix_l_kpoints, fm_matrix_minv_l_kpoints, &
621 fm_matrix_minv, fm_matrix_minv_vtrunc_minv, &
622 mat_munu, mat_p_global, &
623 t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, &
624 starts_array_mc, ends_array_mc, &
625 starts_array_mc_block, ends_array_mc_block, &
626 matrix_s, do_kpoints_from_gamma, kpoints, gd_array, color_sub, &
627 do_ri_sos_laplace_mp2=do_ri_sos_laplace_mp2, calc_forces=calc_forces)
633 IF (.NOT. do_im_time)
THEN
639 IF (my_do_gw .AND. .NOT. do_im_time)
THEN
649 IF (mp2_env%ri_rpa%do_ri_axk)
THEN
650 CALL dbcsr_release(mp2_env%ri_rpa%mo_coeff_o)
651 DEALLOCATE (mp2_env%ri_rpa%mo_coeff_o)
652 CALL dbcsr_release(mp2_env%ri_rpa%mo_coeff_v)
653 DEALLOCATE (mp2_env%ri_rpa%mo_coeff_v)
656 CALL timestop(handle)
677 SUBROUTINE calculate_bib_c_2d(BIb_C_2D, BIb_C, para_env_sub, dimen_ia, homo, virtual, &
679 my_ia_size, my_ia_start, my_ia_end, my_group_L_size)
681 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :), &
682 INTENT(OUT) :: bib_c_2d
683 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :), &
686 INTEGER,
INTENT(IN) :: dimen_ia, homo, virtual
688 INTEGER :: my_ia_size, my_ia_start, my_ia_end, &
691 INTEGER,
PARAMETER :: occ_chunk = 128
693 INTEGER :: ia_global, iib, itmp(2), jjb, my_b_size, my_b_virtual_start, occ_high, occ_low, &
694 proc_receive, proc_send, proc_shift, rec_b_size, rec_b_virtual_end, rec_b_virtual_start
695 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:),
TARGET :: bib_c_rec_1d
696 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: bib_c_rec
698 itmp =
get_limit(dimen_ia, para_env_sub%num_pe, para_env_sub%mepos)
699 my_ia_start = itmp(1)
701 my_ia_size = my_ia_end - my_ia_start + 1
703 CALL get_group_dist(gd_b_virtual, para_env_sub%mepos, sizes=my_b_size, starts=my_b_virtual_start)
706 ALLOCATE (bib_c_2d(my_group_l_size, my_ia_size))
712 DO jjb = 1, my_b_size
713 ia_global = (iib - 1)*virtual + my_b_virtual_start + jjb - 1
714 IF (ia_global >= my_ia_start .AND. ia_global <= my_ia_end)
THEN
715 bib_c_2d(1:my_group_l_size, ia_global - my_ia_start + 1) = bib_c(1:my_group_l_size, jjb, iib)
720 IF (para_env_sub%num_pe > 1)
THEN
721 ALLOCATE (bib_c_rec_1d(int(my_group_l_size,
int_8)*
maxsize(gd_b_virtual)*min(homo, occ_chunk)))
722 DO proc_shift = 1, para_env_sub%num_pe - 1
723 proc_send =
modulo(para_env_sub%mepos + proc_shift, para_env_sub%num_pe)
724 proc_receive =
modulo(para_env_sub%mepos - proc_shift, para_env_sub%num_pe)
726 CALL get_group_dist(gd_b_virtual, proc_receive, rec_b_virtual_start, rec_b_virtual_end, rec_b_size)
729 DO occ_low = 1, homo, occ_chunk
730 occ_high = min(homo, occ_low + occ_chunk - 1)
731 bib_c_rec(1:my_group_l_size, 1:rec_b_size, 1:occ_high - occ_low + 1) => &
732 bib_c_rec_1d(1:int(my_group_l_size,
int_8)*rec_b_size*(occ_high - occ_low + 1))
733 CALL para_env_sub%sendrecv(bib_c(:, :, occ_low:occ_high), proc_send, &
734 bib_c_rec(:, :, 1:occ_high - occ_low + 1), proc_receive)
738 DO iib = occ_low, occ_high
739 DO jjb = 1, rec_b_size
740 ia_global = (iib - 1)*virtual + rec_b_virtual_start + jjb - 1
741 IF (ia_global >= my_ia_start .AND. ia_global <= my_ia_end)
THEN
742 bib_c_2d(1:my_group_l_size, ia_global - my_ia_start + 1) = bib_c_rec(1:my_group_l_size, jjb, iib - occ_low + 1)
749 DEALLOCATE (bib_c_rec_1d)
752 END SUBROUTINE calculate_bib_c_2d
786 SUBROUTINE create_integ_mat(BIb_C_2D, para_env, para_env_sub, color_sub, ngroup, integ_group_size, &
787 dimen_RI, dimen_ia, color_rpa_group, &
788 ext_row_block_size, ext_col_block_size, unit_nr, &
789 my_ia_size, my_ia_start, my_ia_end, &
790 my_group_L_size, my_group_L_start, my_group_L_end, &
791 para_env_RPA, fm_mat_S, nrow_block_mat, ncol_block_mat, &
792 blacs_env_ext, blacs_env_ext_S, dimen_ia_for_block_size, &
793 do_im_time, fm_mat_Q_gemm, fm_mat_Q, qs_env)
796 INTENT(INOUT) :: bib_c_2d
798 INTEGER,
INTENT(IN) :: color_sub, ngroup, integ_group_size, &
800 INTEGER,
DIMENSION(:),
INTENT(IN) :: dimen_ia
801 INTEGER,
INTENT(IN) :: color_rpa_group, ext_row_block_size, &
802 ext_col_block_size, unit_nr
803 INTEGER,
DIMENSION(:),
INTENT(IN) :: my_ia_size, my_ia_start, my_ia_end
804 INTEGER,
INTENT(IN) :: my_group_l_size, my_group_l_start, &
807 TYPE(
cp_fm_type),
DIMENSION(:),
INTENT(INOUT) :: fm_mat_s
808 INTEGER,
INTENT(INOUT) :: nrow_block_mat, ncol_block_mat
810 INTEGER,
INTENT(IN),
OPTIONAL :: dimen_ia_for_block_size
811 LOGICAL,
INTENT(IN),
OPTIONAL :: do_im_time
812 TYPE(
cp_fm_type),
DIMENSION(:),
OPTIONAL :: fm_mat_q_gemm, fm_mat_q
816 CHARACTER(LEN=*),
PARAMETER :: routinen =
'create_integ_mat'
818 INTEGER :: col_row_proc_ratio, grid_2d(2), handle, &
819 iproc, iproc_col, iproc_row, ispin, &
821 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: group_grid_2_mepos
822 LOGICAL :: my_blacs_ext, my_blacs_s_ext, &
828 CALL timeset(routinen, handle)
830 cpassert(
PRESENT(blacs_env_ext) .OR.
PRESENT(dimen_ia_for_block_size))
832 my_blacs_ext = .false.
833 IF (
PRESENT(blacs_env_ext)) my_blacs_ext = .true.
835 my_blacs_s_ext = .false.
836 IF (
PRESENT(blacs_env_ext_s)) my_blacs_s_ext = .true.
838 my_do_im_time = .false.
839 IF (
PRESENT(do_im_time)) my_do_im_time = do_im_time
843 IF (my_blacs_s_ext)
THEN
844 blacs_env => blacs_env_ext_s
846 IF (para_env_rpa%num_pe > 1)
THEN
847 col_row_proc_ratio = max(1, dimen_ia_for_block_size/dimen_ri)
849 iproc_col = min(max(int(sqrt(real(para_env_rpa%num_pe*col_row_proc_ratio, kind=
dp))), 1), para_env_rpa%num_pe) + 1
850 DO iproc = 1, para_env_rpa%num_pe
851 iproc_col = iproc_col - 1
852 IF (mod(para_env_rpa%num_pe, iproc_col) == 0)
EXIT
855 iproc_row = para_env_rpa%num_pe/iproc_col
856 grid_2d(1) = iproc_row
857 grid_2d(2) = iproc_col
863 IF (unit_nr > 0 .AND. .NOT. my_do_im_time)
THEN
864 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
865 "MATRIX_INFO| Number row processes:", grid_2d(1)
866 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
867 "MATRIX_INFO| Number column processes:", grid_2d(2)
871 IF (ext_row_block_size > 0)
THEN
872 nrow_block_mat = ext_row_block_size
874 nrow_block_mat = max(1, dimen_ri/grid_2d(1)/2)
878 IF (ext_col_block_size > 0)
THEN
879 ncol_block_mat = ext_col_block_size
881 ncol_block_mat = max(1, dimen_ia_for_block_size/grid_2d(2)/2)
884 IF (unit_nr > 0 .AND. .NOT. my_do_im_time)
THEN
885 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
886 "MATRIX_INFO| Row block size:", nrow_block_mat
887 WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
888 "MATRIX_INFO| Column block size:", ncol_block_mat
892 IF (.NOT. my_do_im_time)
THEN
893 DO ispin = 1,
SIZE(bib_c_2d)
895 IF (my_blacs_ext)
THEN
897 ncol_global=dimen_ia(ispin), para_env=para_env_rpa)
900 ncol_global=dimen_ia(ispin), para_env=para_env_rpa, &
901 nrow_block=nrow_block_mat, ncol_block=ncol_block_mat, force_block=.true.)
905 CALL create_group_dist(gd_ia, my_ia_start(ispin), my_ia_end(ispin), my_ia_size(ispin), para_env_rpa)
906 CALL create_group_dist(gd_l, my_group_l_start, my_group_l_end, my_group_l_size, para_env_rpa)
910 mepos_in_rpa_group = mod(color_sub, integ_group_size)
911 ALLOCATE (group_grid_2_mepos(0:integ_group_size - 1, 0:para_env_sub%num_pe - 1))
912 group_grid_2_mepos = 0
913 group_grid_2_mepos(mepos_in_rpa_group, para_env_sub%mepos) = para_env_rpa%mepos
914 CALL para_env_rpa%sum(group_grid_2_mepos)
916 CALL array2fm(bib_c_2d(ispin)%array, fm_struct, my_group_l_start, my_group_l_end, &
917 my_ia_start(ispin), my_ia_end(ispin), gd_l, gd_ia, &
918 group_grid_2_mepos, ngroup, para_env_sub%num_pe, fm_mat_s(ispin), &
919 integ_group_size, color_rpa_group)
921 DEALLOCATE (group_grid_2_mepos)
929 IF (para_env_rpa%num_pe /= para_env%num_pe)
THEN
932 comm_exchange = fm_mat_s(ispin)%matrix_struct%context%interconnect(para_env)
933 CALL comm_exchange%sum(fm_mat_s(ispin)%local_data)
934 CALL comm_exchange%free()
940 IF (
PRESENT(fm_mat_q_gemm) .AND. .NOT. my_do_im_time)
THEN
944 ncol_global=dimen_ri, para_env=para_env_rpa, &
945 nrow_block=nrow_block_mat, ncol_block=ncol_block_mat, force_block=.true.)
946 DO ispin = 1,
SIZE(fm_mat_q_gemm)
947 CALL cp_fm_create(fm_mat_q_gemm(ispin), fm_struct, name=
"fm_mat_Q_gemm")
952 IF (
PRESENT(fm_mat_q))
THEN
953 NULLIFY (blacs_env_q)
954 IF (my_blacs_ext)
THEN
955 blacs_env_q => blacs_env_ext
956 ELSE IF (para_env_rpa%num_pe == para_env%num_pe .AND.
PRESENT(qs_env))
THEN
957 CALL get_qs_env(qs_env, blacs_env=blacs_env_q)
963 ncol_global=dimen_ri, para_env=para_env_rpa)
964 DO ispin = 1,
SIZE(fm_mat_q)
965 CALL cp_fm_create(fm_mat_q(ispin), fm_struct, name=
"fm_mat_Q")
970 IF (.NOT. (my_blacs_ext .OR. (para_env_rpa%num_pe == para_env%num_pe .AND.
PRESENT(qs_env)))) &
977 CALL timestop(handle)
979 END SUBROUTINE create_integ_mat
1037 SUBROUTINE rpa_num_int(qs_env, Erpa, mp2_env, para_env, para_env_RPA, para_env_sub, unit_nr, &
1038 homo, virtual, dimen_RI, dimen_RI_red, dimen_ia, dimen_nm_gw, &
1039 Eigenval, num_integ_points, num_integ_group, color_rpa_group, &
1040 fm_matrix_PQ, fm_mat_S, fm_mat_Q_gemm, fm_mat_Q, fm_mat_S_gw, fm_mat_R_gw, &
1041 fm_mat_S_ij_bse, fm_mat_S_ab_bse, &
1042 my_do_gw, do_bse, gw_corr_lev_occ, gw_corr_lev_virt, &
1043 do_minimax_quad, do_im_time, mo_coeff, &
1044 fm_matrix_L_kpoints, fm_matrix_Minv_L_kpoints, &
1045 fm_matrix_Minv, fm_matrix_Minv_Vtrunc_Minv, &
1046 mat_munu, mat_P_global, &
1047 t_3c_M, t_3c_O, t_3c_O_compressed, t_3c_O_ind, &
1048 starts_array_mc, ends_array_mc, &
1049 starts_array_mc_block, ends_array_mc_block, &
1050 matrix_s, do_kpoints_from_Gamma, kpoints, gd_array, color_sub, &
1051 do_ri_sos_laplace_mp2, calc_forces)
1053 TYPE(qs_environment_type),
POINTER :: qs_env
1054 REAL(kind=dp),
INTENT(OUT) :: erpa
1055 TYPE(mp2_type) :: mp2_env
1056 TYPE(mp_para_env_type),
POINTER :: para_env, para_env_rpa, para_env_sub
1057 INTEGER,
INTENT(IN) :: unit_nr
1058 INTEGER,
DIMENSION(:),
INTENT(IN) :: homo, virtual
1059 INTEGER,
INTENT(IN) :: dimen_ri, dimen_ri_red
1060 INTEGER,
DIMENSION(:),
INTENT(IN) :: dimen_ia
1061 INTEGER,
INTENT(IN) :: dimen_nm_gw
1062 REAL(kind=dp),
ALLOCATABLE,
DIMENSION(:, :, :), &
1063 INTENT(INOUT) :: eigenval
1064 INTEGER,
INTENT(IN) :: num_integ_points, num_integ_group, &
1066 TYPE(cp_fm_type),
INTENT(IN) :: fm_matrix_pq
1067 TYPE(cp_fm_type),
DIMENSION(:),
INTENT(IN) :: fm_mat_s, fm_mat_q_gemm, fm_mat_q, &
1069 TYPE(cp_fm_type),
INTENT(IN) :: fm_mat_r_gw, fm_mat_s_ij_bse, &
1071 LOGICAL,
INTENT(IN) :: my_do_gw, do_bse
1072 INTEGER,
DIMENSION(:),
INTENT(IN) :: gw_corr_lev_occ, gw_corr_lev_virt
1073 LOGICAL,
INTENT(IN) :: do_minimax_quad, do_im_time
1074 TYPE(cp_fm_type),
DIMENSION(:),
INTENT(IN) :: mo_coeff
1075 TYPE(cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_matrix_l_kpoints, &
1076 fm_matrix_minv_l_kpoints, &
1078 fm_matrix_minv_vtrunc_minv
1079 TYPE(dbcsr_p_type),
INTENT(IN) :: mat_munu
1080 TYPE(dbcsr_p_type),
INTENT(INOUT) :: mat_p_global
1081 TYPE(dbt_type),
INTENT(INOUT) :: t_3c_m
1082 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:, :), &
1083 INTENT(INOUT) :: t_3c_o
1084 TYPE(hfx_compression_type),
ALLOCATABLE, &
1085 DIMENSION(:, :, :),
INTENT(INOUT) :: t_3c_o_compressed
1086 TYPE(block_ind_type),
ALLOCATABLE, &
1087 DIMENSION(:, :, :),
INTENT(INOUT) :: t_3c_o_ind
1088 INTEGER,
ALLOCATABLE,
DIMENSION(:),
INTENT(IN) :: starts_array_mc, ends_array_mc, &
1089 starts_array_mc_block, &
1091 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s
1092 LOGICAL :: do_kpoints_from_gamma
1093 TYPE(kpoint_type),
POINTER :: kpoints
1094 TYPE(group_dist_d1_type),
INTENT(IN) :: gd_array
1095 INTEGER,
INTENT(IN) :: color_sub
1096 LOGICAL,
INTENT(IN) :: do_ri_sos_laplace_mp2, calc_forces
1098 CHARACTER(LEN=*),
PARAMETER :: routinen =
'rpa_num_int'
1100 COMPLEX(KIND=dp),
ALLOCATABLE, &
1101 DIMENSION(:, :, :, :) :: vec_sigma_c_gw
1102 INTEGER :: count_ev_sc_gw, cut_memory, group_size_p, gw_corr_lev_tot, handle, handle3, i, &
1103 ikp_local, ispin, iter_evgw, iter_sc_gw0, j, jquad, min_bsize, mm_style, nkp, &
1104 nkp_self_energy, nmo, nspins, num_3c_repl, num_cells_dm, num_fit_points, pspin, qspin, &
1106 INTEGER(int_8) :: dbcsr_nflop
1107 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: index_to_cell_3c
1108 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :) :: cell_to_index_3c
1109 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, prim_blk_sizes, &
1111 LOGICAL :: do_apply_ic_corr_to_gw, do_gw_im_time, do_ic_model, do_kpoints_cubic_rpa, &
1112 do_periodic, do_print, do_ri_sigma_x, exit_ev_gw, first_cycle, &
1113 first_cycle_periodic_correction, my_open_shell, print_ic_values
1114 LOGICAL,
ALLOCATABLE,
DIMENSION(:, :, :, :, :) :: has_mat_p_blocks
1115 REAL(kind=dp) :: a_scaling, alpha, dbcsr_time, e_axk, e_axk_corr, eps_filter, &
1116 eps_filter_im_time, ext_scaling, fermi_level_offset, fermi_level_offset_input, &
1117 my_flop_rate, omega, omega_max_fit, omega_old, tau, tau_old
1118 REAL(kind=dp),
ALLOCATABLE,
DIMENSION(:) :: delta_corr, e_fermi, tau_tj, tau_wj, tj, &
1119 trace_qomega, vec_omega_fit_gw, wj, &
1121 REAL(kind=dp),
ALLOCATABLE,
DIMENSION(:, :) :: vec_w_gw, weights_cos_tf_t_to_w, &
1122 weights_cos_tf_w_to_t, &
1123 weights_sin_tf_t_to_w
1124 REAL(kind=dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: eigenval_last, eigenval_scf, &
1126 TYPE(cp_cfm_type) :: cfm_mat_q
1127 TYPE(cp_fm_type) :: fm_mat_q_static_bse, fm_mat_q_static_bse_gemm, fm_mat_ri_global_work, &
1128 fm_mat_s_ia_bse, fm_mat_work, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1129 fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau
1130 TYPE(cp_fm_type),
ALLOCATABLE,
DIMENSION(:) :: fm_mat_s_gw_work, fm_mat_w, &
1131 fm_mo_coeff_occ, fm_mo_coeff_virt
1132 TYPE(cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_mat_l_kpoints, fm_mat_minv_l_kpoints
1133 TYPE(dbcsr_p_type) :: mat_dm, mat_l, mat_m_p_munu_occ, &
1134 mat_m_p_munu_virt, mat_minvvminv
1135 TYPE(dbcsr_p_type),
ALLOCATABLE, &
1136 DIMENSION(:, :, :) :: mat_p_omega
1137 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_berry_im_mo_mo, &
1138 matrix_berry_re_mo_mo
1139 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: mat_p_omega_kp
1140 TYPE(dbcsr_type),
POINTER :: mat_w, mat_work
1141 TYPE(dbt_type) :: t_3c_overl_int_ao_mo
1142 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:) :: t_3c_overl_int_gw_ao, &
1143 t_3c_overl_int_gw_ri, &
1144 t_3c_overl_nnp_ic, &
1145 t_3c_overl_nnp_ic_reflected
1146 TYPE(dgemm_counter_type) :: dgemm_counter
1147 TYPE(hfx_compression_type),
ALLOCATABLE, &
1148 DIMENSION(:) :: t_3c_o_mo_compressed
1149 TYPE(im_time_force_type) :: force_data
1151 TYPE(two_dim_int_array),
ALLOCATABLE,
DIMENSION(:) :: t_3c_o_mo_ind
1153 CALL timeset(routinen, handle)
1156 nmo = homo(1) + virtual(1)
1157 my_open_shell = (nspins == 2)
1159 do_gw_im_time = my_do_gw .AND. do_im_time
1160 do_ri_sigma_x = mp2_env%ri_g0w0%do_ri_Sigma_x
1161 do_ic_model = mp2_env%ri_g0w0%do_ic_model
1162 print_ic_values = mp2_env%ri_g0w0%print_ic_values
1163 do_periodic = mp2_env%ri_g0w0%do_periodic
1164 do_kpoints_cubic_rpa = mp2_env%ri_rpa_im_time%do_im_time_kpoints
1167 mm_style = wfc_mm_style_gemm
1168 IF (.NOT. do_ri_sos_laplace_mp2) mm_style = mp2_env%ri_rpa%mm_style
1171 ext_scaling = 0.2_dp
1172 omega_max_fit = mp2_env%ri_g0w0%omega_max_fit
1173 fermi_level_offset_input = mp2_env%ri_g0w0%fermi_level_offset
1174 iter_evgw = mp2_env%ri_g0w0%iter_evGW
1175 iter_sc_gw0 = mp2_env%ri_g0w0%iter_sc_GW0
1176 IF ((.NOT. do_im_time))
THEN
1177 IF (iter_sc_gw0 .NE. 1 .AND. iter_evgw .NE. 1) cpabort(
"Mixed scGW0/evGW not implemented.")
1179 IF (iter_sc_gw0 .NE. 1) iter_evgw = iter_sc_gw0
1182 ext_scaling = 0.0_dp
1187 IF (do_kpoints_cubic_rpa .AND. do_ri_sos_laplace_mp2)
THEN
1188 cpabort(
"RI-SOS-Laplace-MP2 with k-point-sampling is not implemented.")
1191 do_apply_ic_corr_to_gw = .false.
1192 IF (mp2_env%ri_g0w0%ic_corr_list(1)%array(1) > 0.0_dp) do_apply_ic_corr_to_gw = .true.
1194 IF (do_im_time)
THEN
1195 cpassert(do_minimax_quad .OR. do_ri_sos_laplace_mp2)
1197 group_size_p = mp2_env%ri_rpa_im_time%group_size_P
1198 cut_memory = mp2_env%ri_rpa_im_time%cut_memory
1199 eps_filter = mp2_env%ri_rpa_im_time%eps_filter
1200 eps_filter_im_time = mp2_env%ri_rpa_im_time%eps_filter* &
1201 mp2_env%ri_rpa_im_time%eps_filter_factor
1203 min_bsize = mp2_env%ri_rpa_im_time%min_bsize
1205 CALL alloc_im_time(qs_env, para_env, dimen_ri, dimen_ri_red, &
1206 num_integ_points, nspins, fm_mat_q(1), fm_mo_coeff_occ, fm_mo_coeff_virt, &
1207 fm_matrix_minv_l_kpoints, fm_matrix_l_kpoints, mat_p_global, &
1208 t_3c_o, matrix_s, kpoints, eps_filter_im_time, &
1209 cut_memory, nkp, num_cells_dm, num_3c_repl, &
1210 size_p, ikp_local, &
1214 do_ic_model, do_kpoints_cubic_rpa, &
1215 do_kpoints_from_gamma, do_ri_sigma_x, my_open_shell, &
1216 has_mat_p_blocks, wkp_w, &
1217 cfm_mat_q, fm_mat_minv_l_kpoints, fm_mat_l_kpoints, &
1218 fm_mat_ri_global_work, fm_mat_work, fm_mo_coeff_occ_scaled, &
1219 fm_mo_coeff_virt_scaled, mat_dm, mat_l, mat_m_p_munu_occ, mat_m_p_munu_virt, &
1220 mat_minvvminv, mat_p_omega, mat_p_omega_kp, mat_work, mo_coeff, &
1221 fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, homo, nmo)
1223 IF (calc_forces)
CALL init_im_time_forces(force_data, fm_matrix_pq, t_3c_m, unit_nr, mp2_env, qs_env)
1227 CALL dbcsr_get_info(mat_p_global%matrix, &
1228 row_blk_size=ri_blk_sizes)
1230 CALL dbcsr_get_info(matrix_s(1)%matrix, &
1231 row_blk_size=prim_blk_sizes)
1233 gw_corr_lev_tot = gw_corr_lev_occ(1) + gw_corr_lev_virt(1)
1235 IF (.NOT. do_kpoints_cubic_rpa)
THEN
1236 CALL allocate_matrices_gw_im_time(gw_corr_lev_occ, gw_corr_lev_virt, homo, nmo, &
1237 num_integ_points, unit_nr, &
1238 ri_blk_sizes, do_ic_model, &
1239 para_env, fm_mat_w, fm_mat_q(1), &
1241 t_3c_overl_int_ao_mo, t_3c_o_mo_compressed, t_3c_o_mo_ind, &
1242 t_3c_overl_int_gw_ri, t_3c_overl_int_gw_ao, &
1243 starts_array_mc, ends_array_mc, &
1244 t_3c_overl_nnp_ic, t_3c_overl_nnp_ic_reflected, &
1245 matrix_s, mat_w, t_3c_o, &
1246 t_3c_o_compressed, t_3c_o_ind, &
1254 IF (do_ic_model)
THEN
1256 cpassert(do_gw_im_time)
1257 cpassert(.NOT. do_periodic)
1258 IF (cut_memory .NE. 1) cpabort(
"For IC, use MEMORY_CUT 1 in the LOW_SCALING section.")
1261 ALLOCATE (e_fermi(nspins))
1262 IF (do_minimax_quad .OR. do_ri_sos_laplace_mp2)
THEN
1263 do_print = .NOT. do_ic_model
1264 CALL get_minimax_grid(para_env, unit_nr, homo, eigenval, num_integ_points, do_im_time, &
1265 do_ri_sos_laplace_mp2, do_print, &
1266 tau_tj, tau_wj, qs_env, do_gw_im_time, &
1267 do_kpoints_cubic_rpa, e_fermi(1), tj, wj, &
1268 weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, weights_sin_tf_t_to_w, &
1269 qs_env%mp2_env%ri_g0w0%regularization_minimax)
1272 IF (qs_env%mp2_env%ri_rpa_im_time%keep_quad)
THEN
1273 CALL keep_initial_quad(tj, wj, tau_tj, tau_wj, weights_cos_tf_t_to_w, &
1274 weights_cos_tf_w_to_t, do_ri_sos_laplace_mp2, do_im_time, &
1275 num_integ_points, unit_nr, qs_env)
1278 IF (calc_forces) cpabort(
"Forces with Clenshaw-Curtis grid not implemented.")
1279 CALL get_clenshaw_grid(para_env, para_env_rpa, unit_nr, homo, virtual, eigenval, num_integ_points, &
1280 num_integ_group, color_rpa_group, fm_mat_s, my_do_gw, &
1281 ext_scaling, a_scaling, tj, wj)
1285 IF (.NOT. do_ri_sos_laplace_mp2)
THEN
1286 ALLOCATE (trace_qomega(dimen_ri_red))
1289 IF (do_ri_sos_laplace_mp2 .AND. .NOT. do_im_time)
THEN
1291 ELSE IF (my_open_shell .OR. do_ri_sos_laplace_mp2)
THEN
1298 CALL allocate_matrices_gw(vec_sigma_c_gw, color_rpa_group, dimen_nm_gw, &
1299 gw_corr_lev_occ, gw_corr_lev_virt, homo, &
1300 nmo, num_integ_group, num_integ_points, unit_nr, &
1301 gw_corr_lev_tot, num_fit_points, omega_max_fit, &
1302 do_minimax_quad, do_periodic, do_ri_sigma_x,.NOT. do_im_time, &
1303 first_cycle_periodic_correction, &
1304 a_scaling, eigenval, tj, vec_omega_fit_gw, vec_sigma_x_gw, &
1305 delta_corr, eigenval_last, eigenval_scf, vec_w_gw, &
1306 fm_mat_s_gw, fm_mat_s_gw_work, &
1307 para_env, mp2_env, kpoints, nkp, nkp_self_energy, &
1308 do_kpoints_cubic_rpa, do_kpoints_from_gamma)
1312 CALL cp_fm_create(fm_mat_q_static_bse_gemm, fm_mat_q_gemm(1)%matrix_struct)
1313 CALL cp_fm_to_fm(fm_mat_q_gemm(1), fm_mat_q_static_bse_gemm)
1314 CALL cp_fm_set_all(fm_mat_q_static_bse_gemm, 0.0_dp)
1316 CALL cp_fm_create(fm_mat_q_static_bse, fm_mat_q(1)%matrix_struct)
1317 CALL cp_fm_to_fm(fm_mat_q(1), fm_mat_q_static_bse)
1318 CALL cp_fm_set_all(fm_mat_q_static_bse, 0.0_dp)
1324 IF (calc_forces .AND. .NOT. do_im_time)
CALL rpa_grad_create(
rpa_grad, fm_mat_q(1), &
1325 fm_mat_s, homo, virtual, mp2_env, eigenval(:, 1, :), &
1326 unit_nr, do_ri_sos_laplace_mp2)
1329 IF (mp2_env%ri_rpa%do_ri_axk) e_axk = 0.0_dp
1330 first_cycle = .true.
1332 CALL dgemm_counter_init(dgemm_counter, unit_nr, mp2_env%ri_rpa%print_dgemm_info)
1334 DO count_ev_sc_gw = 1, iter_evgw
1339 IF (do_ic_model) cycle
1344 vec_sigma_c_gw = z_zero
1345 first_cycle = .true.
1349 IF (do_im_time)
THEN
1351 IF (.NOT. do_kpoints_cubic_rpa)
THEN
1352 DO ispin = 1, nspins
1353 e_fermi(ispin) = (eigenval(homo(ispin), 1, ispin) + eigenval(homo(ispin) + 1, 1, ispin))*0.5_dp
1360 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(/T3,A,T66,i15)") &
1361 "MEMORY_INFO| Memory cut:", cut_memory
1362 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3,A,T66,ES15.2)") &
1363 "SPARSITY_INFO| Eps filter for M virt/occ tensors:", eps_filter
1364 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3,A,T66,ES15.2)") &
1365 "SPARSITY_INFO| Eps filter for P matrix:", eps_filter_im_time
1366 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3,A,T66,i15)") &
1367 "SPARSITY_INFO| Minimum tensor block size:", min_bsize
1370 CALL zero_mat_p_omega(mat_p_omega(:, :, 1))
1373 CALL compute_mat_p_omega(mat_p_omega(:, :, 1), fm_scaled_dm_occ_tau, &
1374 fm_scaled_dm_virt_tau, fm_mo_coeff_occ(1), fm_mo_coeff_virt(1), &
1375 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1376 mat_p_global, matrix_s, 1, &
1377 t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, &
1378 starts_array_mc, ends_array_mc, &
1379 starts_array_mc_block, ends_array_mc_block, &
1380 weights_cos_tf_t_to_w, tj, tau_tj, e_fermi(1), eps_filter, alpha, &
1381 eps_filter_im_time, eigenval(:, 1, 1), nmo, &
1382 num_integ_points, cut_memory, &
1383 unit_nr, mp2_env, para_env, &
1384 qs_env, do_kpoints_from_gamma, &
1385 index_to_cell_3c, cell_to_index_3c, &
1386 has_mat_p_blocks, do_ri_sos_laplace_mp2, &
1387 dbcsr_time, dbcsr_nflop)
1390 IF (my_open_shell)
THEN
1391 CALL zero_mat_p_omega(mat_p_omega(:, :, 2))
1392 CALL compute_mat_p_omega(mat_p_omega(:, :, 2), fm_scaled_dm_occ_tau, &
1393 fm_scaled_dm_virt_tau, fm_mo_coeff_occ(2), &
1394 fm_mo_coeff_virt(2), &
1395 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1396 mat_p_global, matrix_s, 2, &
1397 t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, &
1398 starts_array_mc, ends_array_mc, &
1399 starts_array_mc_block, ends_array_mc_block, &
1400 weights_cos_tf_t_to_w, tj, tau_tj, e_fermi(2), eps_filter, alpha, &
1401 eps_filter_im_time, eigenval(:, 1, 2), nmo, &
1402 num_integ_points, cut_memory, &
1403 unit_nr, mp2_env, para_env, &
1404 qs_env, do_kpoints_from_gamma, &
1405 index_to_cell_3c, cell_to_index_3c, &
1406 has_mat_p_blocks, do_ri_sos_laplace_mp2, &
1407 dbcsr_time, dbcsr_nflop)
1410 IF (.NOT. do_ri_sos_laplace_mp2)
THEN
1411 DO j = 1,
SIZE(mat_p_omega, 2)
1412 DO i = 1,
SIZE(mat_p_omega, 1)
1413 CALL dbcsr_add(mat_p_omega(i, j, 1)%matrix, mat_p_omega(i, j, 2)%matrix, 1.0_dp, 1.0_dp)
1414 IF (.NOT. calc_forces)
CALL dbcsr_clear(mat_p_omega(i, j, 2)%matrix)
1422 DO jquad = 1, num_integ_points
1424 IF (
modulo(jquad, num_integ_group) /= color_rpa_group) cycle
1426 CALL timeset(routinen//
"_RPA_matrix_operations", handle3)
1428 IF (do_ri_sos_laplace_mp2)
THEN
1429 omega = tau_tj(jquad)
1431 IF (do_minimax_quad)
THEN
1434 omega = a_scaling/tan(tj(jquad))
1438 IF (do_im_time)
THEN
1441 IF (.NOT. (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma))
THEN
1443 DO ispin = 1,
SIZE(mat_p_omega, 3)
1444 CALL contract_p_omega_with_mat_l(mat_p_omega(jquad, 1, ispin)%matrix, mat_l%matrix, mat_work, &
1445 eps_filter_im_time, fm_mat_work, dimen_ri, dimen_ri_red, &
1446 fm_mat_minv_l_kpoints(1, 1), fm_mat_q(ispin))
1452 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3, A, 1X, I3, 1X, A, 1X, I3)") &
1453 "INTEG_INFO| Started with Integration point", jquad,
"of", num_integ_points
1455 IF (first_cycle .AND. count_ev_sc_gw > 1)
THEN
1456 IF (iter_sc_gw0 == 1)
THEN
1457 DO ispin = 1, nspins
1458 CALL remove_scaling_factor_rpa(fm_mat_s(ispin), virtual(ispin), &
1459 eigenval_last(:, 1, ispin), homo(ispin), omega_old)
1462 DO ispin = 1, nspins
1463 CALL remove_scaling_factor_rpa(fm_mat_s(ispin), virtual(ispin), &
1464 eigenval_scf(:, 1, ispin), homo(ispin), omega_old)
1469 CALL calc_mat_q(fm_mat_s(1), do_ri_sos_laplace_mp2, first_cycle, iter_sc_gw0, virtual(1), &
1470 eigenval(:, 1, 1), eigenval_scf(:, 1, 1), &
1471 homo(1), omega, omega_old, jquad, mm_style, dimen_ri_red, dimen_ia(1), alpha, fm_mat_q(1), &
1472 fm_mat_q_gemm(1), do_bse, fm_mat_q_static_bse_gemm, dgemm_counter, &
1475 IF (my_open_shell)
THEN
1476 CALL calc_mat_q(fm_mat_s(2), do_ri_sos_laplace_mp2, first_cycle, iter_sc_gw0, virtual(2), &
1477 eigenval(:, 1, 2), eigenval_scf(:, 1, 2), &
1478 homo(2), omega, omega_old, jquad, mm_style, &
1479 dimen_ri_red, dimen_ia(2), alpha, fm_mat_q(2), fm_mat_q_gemm(2), do_bse, &
1480 fm_mat_q_static_bse_gemm, dgemm_counter, &
1484 IF (.NOT. do_ri_sos_laplace_mp2)
THEN
1485 CALL cp_fm_scale_and_add(alpha=1.0_dp, matrix_a=fm_mat_q(1), beta=1.0_dp, matrix_b=fm_mat_q(2))
1493 IF (mp2_env%ri_rpa%do_ri_axk)
THEN
1494 CALL compute_axk_ener(qs_env, fm_mat_q(1), fm_mat_q_gemm(1), dimen_ri_red, dimen_ia(1), para_env_sub, &
1495 para_env_rpa, eigenval(:, 1, 1), fm_mat_s(1), homo(1), virtual(1), omega, &
1496 mp2_env, mat_munu, unit_nr, e_axk_corr)
1499 e_axk = e_axk + e_axk_corr*wj(jquad)
1502 IF (do_ri_sos_laplace_mp2)
THEN
1503 CALL sos_mp2_postprocessing(fm_mat_q, erpa, tau_wj(jquad))
1505 IF (calc_forces .AND. .NOT. do_im_time)
CALL rpa_grad_matrix_operations(mp2_env,
rpa_grad, do_ri_sos_laplace_mp2, &
1506 fm_mat_q, fm_mat_q_gemm, dgemm_counter, fm_mat_s, omega, homo, virtual, &
1507 eigenval(:, 1, :), tau_wj(jquad), unit_nr)
1509 IF (calc_forces .AND. .NOT. do_im_time)
CALL rpa_grad_copy_q(fm_mat_q(1),
rpa_grad)
1511 CALL q_trace_and_add_unit_matrix(dimen_ri_red, trace_qomega, fm_mat_q(1), para_env_rpa)
1513 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
1514 CALL invert_eps_compute_w_and_erpa_kp(dimen_ri, num_integ_points, jquad, nkp, count_ev_sc_gw, para_env, &
1515 erpa, tau_tj, tj, wj, weights_cos_tf_w_to_t, &
1516 wkp_w, do_gw_im_time, do_ri_sigma_x, do_kpoints_from_gamma, &
1517 cfm_mat_q, ikp_local, &
1518 mat_p_omega(:, :, 1), mat_p_omega_kp, qs_env, eps_filter_im_time, unit_nr, &
1519 kpoints, fm_mat_minv_l_kpoints, fm_mat_l_kpoints, &
1520 fm_mat_w, fm_mat_ri_global_work, mat_minvvminv, &
1521 fm_matrix_minv, fm_matrix_minv_vtrunc_minv)
1523 CALL compute_erpa_by_freq_int(dimen_ri_red, trace_qomega, fm_mat_q(1), para_env_rpa, erpa, wj(jquad))
1526 IF (calc_forces .AND. .NOT. do_im_time)
CALL rpa_grad_matrix_operations(mp2_env,
rpa_grad, do_ri_sos_laplace_mp2, &
1527 fm_mat_q, fm_mat_q_gemm, dgemm_counter, fm_mat_s, omega, homo, virtual, &
1528 eigenval(:, 1, :), wj(jquad), unit_nr)
1532 first_cycle = .false.
1535 CALL timestop(handle3)
1539 CALL get_fermi_level_offset(fermi_level_offset, fermi_level_offset_input, eigenval(:, 1, :), homo)
1542 IF (do_im_time)
THEN
1544 IF (.NOT. (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma))
THEN
1545 CALL compute_w_cubic_gw(fm_mat_w, fm_mat_q(1), fm_mat_work, dimen_ri, fm_mat_minv_l_kpoints, num_integ_points, &
1546 tj, tau_tj, weights_cos_tf_w_to_t, jquad, omega)
1549 CALL compute_gw_self_energy(vec_sigma_c_gw, dimen_nm_gw, dimen_ri_red, gw_corr_lev_occ, &
1550 gw_corr_lev_virt, homo, jquad, nmo, num_fit_points, num_integ_points, &
1551 do_bse, do_im_time, do_periodic, first_cycle_periodic_correction, &
1552 fermi_level_offset, &
1553 omega, eigenval(:, 1, :), delta_corr, vec_omega_fit_gw, vec_w_gw, wj, &
1554 fm_mat_q(1), fm_mat_q_static_bse, fm_mat_r_gw, fm_mat_s_gw, &
1555 fm_mat_s_gw_work, mo_coeff(1), para_env, &
1556 para_env_rpa, matrix_berry_im_mo_mo, matrix_berry_re_mo_mo, &
1557 kpoints, qs_env, mp2_env)
1561 IF (unit_nr > 0)
CALL m_flush(unit_nr)
1562 CALL para_env%sync()
1566 CALL para_env%sum(erpa)
1568 IF (.NOT. do_ri_sos_laplace_mp2)
THEN
1569 erpa = erpa/(pi*2.0_dp)
1570 IF (do_minimax_quad) erpa = erpa/2.0_dp
1573 IF (mp2_env%ri_rpa%do_ri_axk)
THEN
1574 CALL para_env%sum(e_axk)
1575 e_axk = e_axk/(pi*2.0_dp)
1576 IF (do_minimax_quad) e_axk = e_axk/2.0_dp
1577 mp2_env%ri_rpa%ener_axk = e_axk
1580 IF (calc_forces .AND. do_ri_sos_laplace_mp2 .AND. do_im_time)
THEN
1581 IF (my_open_shell)
THEN
1584 CALL calc_laplace_loop_forces(force_data, mat_p_omega(:, 1, :), t_3c_m, t_3c_o(1, 1), &
1585 t_3c_o_compressed(1, 1, :), t_3c_o_ind(1, 1, :), fm_scaled_dm_occ_tau, &
1586 fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, &
1587 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1588 starts_array_mc, ends_array_mc, starts_array_mc_block, &
1589 ends_array_mc_block, num_integ_points, nmo, eigenval(:, 1, :), &
1590 tau_tj, tau_wj, cut_memory, pspin, qspin, my_open_shell, &
1591 unit_nr, dbcsr_time, dbcsr_nflop, mp2_env, qs_env)
1594 CALL calc_laplace_loop_forces(force_data, mat_p_omega(:, 1, :), t_3c_m, t_3c_o(1, 1), &
1595 t_3c_o_compressed(1, 1, :), t_3c_o_ind(1, 1, :), fm_scaled_dm_occ_tau, &
1596 fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, &
1597 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1598 starts_array_mc, ends_array_mc, starts_array_mc_block, &
1599 ends_array_mc_block, num_integ_points, nmo, eigenval(:, 1, :), &
1600 tau_tj, tau_wj, cut_memory, pspin, qspin, my_open_shell, &
1601 unit_nr, dbcsr_time, dbcsr_nflop, mp2_env, qs_env)
1606 CALL calc_laplace_loop_forces(force_data, mat_p_omega(:, 1, :), t_3c_m, t_3c_o(1, 1), &
1607 t_3c_o_compressed(1, 1, :), t_3c_o_ind(1, 1, :), fm_scaled_dm_occ_tau, &
1608 fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, &
1609 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1610 starts_array_mc, ends_array_mc, starts_array_mc_block, &
1611 ends_array_mc_block, num_integ_points, nmo, eigenval(:, 1, :), &
1612 tau_tj, tau_wj, cut_memory, pspin, qspin, my_open_shell, &
1613 unit_nr, dbcsr_time, dbcsr_nflop, mp2_env, qs_env)
1615 CALL calc_post_loop_forces(force_data, unit_nr, qs_env)
1618 IF (calc_forces .AND. do_im_time .AND. .NOT. do_ri_sos_laplace_mp2)
THEN
1619 DO ispin = 1, nspins
1620 CALL calc_rpa_loop_forces(force_data, mat_p_omega(:, 1, :), t_3c_m, t_3c_o(1, 1), &
1621 t_3c_o_compressed(1, 1, :), t_3c_o_ind(1, 1, :), fm_scaled_dm_occ_tau, &
1622 fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, &
1623 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, &
1624 starts_array_mc, ends_array_mc, starts_array_mc_block, &
1625 ends_array_mc_block, num_integ_points, nmo, eigenval(:, 1, :), &
1626 e_fermi(ispin), weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, tj, &
1627 wj, tau_tj, cut_memory, ispin, my_open_shell, unit_nr, dbcsr_time, &
1628 dbcsr_nflop, mp2_env, qs_env)
1630 CALL calc_post_loop_forces(force_data, unit_nr, qs_env)
1633 IF (do_im_time)
THEN
1635 my_flop_rate = real(dbcsr_nflop, dp)/(1.0e09_dp*dbcsr_time)
1636 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(/T3,A,T73,ES8.2)") &
1637 "PERFORMANCE| DBCSR total number of flops:", real(dbcsr_nflop*para_env%num_pe, dp)
1638 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3,A,T66,F15.2)") &
1639 "PERFORMANCE| DBCSR total execution time:", dbcsr_time
1640 IF (unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3,A,T66,F15.2)") &
1641 "PERFORMANCE| DBCSR flop rate (Gflops / MPI rank):", my_flop_rate
1645 CALL dgemm_counter_write(dgemm_counter, para_env)
1654 CALL compute_qp_energies(vec_sigma_c_gw, count_ev_sc_gw, gw_corr_lev_occ, &
1655 gw_corr_lev_tot, gw_corr_lev_virt, homo, &
1656 nmo, num_fit_points, num_integ_points, &
1657 unit_nr, do_apply_ic_corr_to_gw, do_im_time, &
1658 do_periodic, do_ri_sigma_x, first_cycle_periodic_correction, &
1659 e_fermi, eps_filter, fermi_level_offset, &
1660 delta_corr, eigenval, &
1661 eigenval_last, eigenval_scf, iter_sc_gw0, exit_ev_gw, tau_tj, tj, &
1662 vec_omega_fit_gw, vec_sigma_x_gw, mp2_env%ri_g0w0%ic_corr_list, &
1663 weights_cos_tf_t_to_w, weights_sin_tf_t_to_w, &
1664 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, fm_mo_coeff_occ, &
1665 fm_mo_coeff_virt, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, &
1666 mo_coeff(1), fm_mat_w, para_env, para_env_rpa, mat_dm, mat_minvvminv, &
1667 t_3c_o, t_3c_m, t_3c_overl_int_ao_mo, t_3c_o_compressed, t_3c_o_mo_compressed, &
1668 t_3c_o_ind, t_3c_o_mo_ind, &
1669 t_3c_overl_int_gw_ri, t_3c_overl_int_gw_ao, &
1670 matrix_berry_im_mo_mo, matrix_berry_re_mo_mo, mat_w, matrix_s, &
1671 kpoints, mp2_env, qs_env, nkp_self_energy, do_kpoints_cubic_rpa, &
1672 starts_array_mc, ends_array_mc)
1675 IF (exit_ev_gw)
EXIT
1681 IF (do_ic_model)
THEN
1683 IF (my_open_shell)
THEN
1685 CALL calculate_ic_correction(eigenval(:, 1, 1), mat_minvvminv%matrix, &
1686 t_3c_overl_nnp_ic(1), t_3c_overl_nnp_ic_reflected(1), &
1688 gw_corr_lev_occ(1), gw_corr_lev_virt(1), homo(1), unit_nr, &
1689 print_ic_values, para_env, do_alpha=.true.)
1691 CALL calculate_ic_correction(eigenval(:, 1, 2), mat_minvvminv%matrix, &
1692 t_3c_overl_nnp_ic(2), t_3c_overl_nnp_ic_reflected(2), &
1694 gw_corr_lev_occ(2), gw_corr_lev_virt(2), homo(2), unit_nr, &
1695 print_ic_values, para_env, do_beta=.true.)
1699 CALL calculate_ic_correction(eigenval(:, 1, 1), mat_minvvminv%matrix, &
1700 t_3c_overl_nnp_ic(1), t_3c_overl_nnp_ic_reflected(1), &
1702 gw_corr_lev_occ(1), gw_corr_lev_virt(1), homo(1), unit_nr, &
1703 print_ic_values, para_env)
1712 CALL cp_fm_create(fm_mat_s_ia_bse, fm_mat_s(1)%matrix_struct)
1713 CALL cp_fm_to_fm(fm_mat_s(1), fm_mat_s_ia_bse)
1715 IF (iter_sc_gw0 == 1)
THEN
1716 CALL remove_scaling_factor_rpa(fm_mat_s_ia_bse, virtual(1), &
1717 eigenval_last(:, 1, 1), homo(1), omega)
1719 CALL remove_scaling_factor_rpa(fm_mat_s_ia_bse, virtual(1), &
1720 eigenval_scf(:, 1, 1), homo(1), omega)
1723 CALL start_bse_calculation(fm_mat_s_ia_bse, fm_mat_s_ij_bse, fm_mat_s_ab_bse, &
1724 fm_mat_q_static_bse, fm_mat_q_static_bse_gemm, &
1725 eigenval, homo, virtual, dimen_ri, dimen_ri_red, &
1726 gd_array, color_sub, mp2_env, unit_nr)
1728 CALL cp_fm_release(fm_mat_s_ia_bse)
1732 CALL deallocate_matrices_gw(fm_mat_s_gw_work, vec_w_gw, vec_sigma_c_gw, vec_omega_fit_gw, &
1733 mp2_env%ri_g0w0%vec_Sigma_x_minus_vxc_gw, &
1734 eigenval_last, eigenval_scf, do_periodic, matrix_berry_re_mo_mo, matrix_berry_im_mo_mo, &
1735 kpoints, vec_sigma_x_gw,.NOT. do_im_time)
1738 IF (do_im_time)
THEN
1740 CALL dealloc_im_time(fm_mo_coeff_occ, fm_mo_coeff_virt, &
1741 fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, index_to_cell_3c, &
1742 cell_to_index_3c, do_ic_model, &
1743 do_kpoints_cubic_rpa, do_kpoints_from_gamma, do_ri_sigma_x, &
1745 wkp_w, cfm_mat_q, fm_mat_minv_l_kpoints, fm_mat_l_kpoints, &
1746 fm_matrix_minv, fm_matrix_minv_vtrunc_minv, fm_mat_ri_global_work, fm_mat_work, &
1747 fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, mat_dm, mat_l, &
1748 mat_minvvminv, mat_p_omega, mat_p_omega_kp, &
1749 t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, mat_work, qs_env)
1752 CALL deallocate_matrices_gw_im_time(weights_cos_tf_w_to_t, weights_sin_tf_t_to_w, do_ic_model, &
1753 do_kpoints_cubic_rpa, fm_mat_w, &
1754 t_3c_overl_int_ao_mo, t_3c_o_mo_compressed, t_3c_o_mo_ind, &
1755 t_3c_overl_int_gw_ri, t_3c_overl_int_gw_ao, &
1756 t_3c_overl_nnp_ic, t_3c_overl_nnp_ic_reflected, &
1762 IF (.NOT. do_ri_sos_laplace_mp2)
THEN
1765 DEALLOCATE (trace_qomega)
1768 IF (do_im_time .OR. do_ri_sos_laplace_mp2)
THEN
1773 IF (do_im_time .AND. calc_forces)
THEN
1774 CALL im_time_force_release(force_data)
1777 IF (calc_forces .AND. .NOT. do_im_time)
CALL rpa_grad_finalize(
rpa_grad, mp2_env, para_env_sub, para_env, &
1778 qs_env, gd_array, color_sub, do_ri_sos_laplace_mp2, &
1781 CALL timestop(handle)
1783 END SUBROUTINE rpa_num_int
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public wilhelm2016b
integer, save, public delben2015
integer, save, public wilhelm2016a
integer, save, public wilhelm2018
integer, save, public delben2013
integer, save, public ren2013
integer, save, public wilhelm2017
integer, save, public ren2011
integer, save, public bates2013
Main routines for GW + Bethe-Salpeter for computing electronic excitations.
subroutine, public start_bse_calculation(fm_mat_s_ia_bse, fm_mat_s_ij_bse, fm_mat_s_ab_bse, fm_mat_q_static_bse, fm_mat_q_static_bse_gemm, eigenval, homo, virtual, dimen_ri, dimen_ri_red, gd_array, color_sub, mp2_env, unit_nr)
Main subroutine managing BSE calculations.
methods related to the blacs parallel environment
subroutine, public cp_blacs_env_release(blacs_env)
releases the given blacs_env
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
Represents a complex full matrix distributed on many processors.
basic linear algebra operations for full matrices
subroutine, public cp_fm_scale_and_add(alpha, matrix_a, beta, matrix_b)
calc A <- alpha*A + beta*B optimized for alpha == 1.0 (just add beta*B) and beta == 0....
represent the structure of a full matrix
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
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
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
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
This is the start of a dbt_api, all publically needed functions are exported here....
Counters to determine the performance of parallel DGEMMs.
elemental subroutine, public dgemm_counter_init(dgemm_counter, unit_nr, print_info)
Initialize a dgemm_counter.
subroutine, public dgemm_counter_write(dgemm_counter, para_env)
calculate and print flop rates
Types to describe group distributions.
elemental integer function, public maxsize(this)
...
Types and set/get functions for HFX.
Defines the basic variable types.
integer, parameter, public int_8
integer, parameter, public dp
Types and basic routines needed for a kpoint calculation.
Machine interface based on Fortran 2003 and POSIX.
subroutine, public m_memory(mem)
Returns the total amount of memory [bytes] in use, if known, zero otherwise.
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
complex(kind=dp), parameter, public z_zero
Interface to the message passing library MPI.
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)
Routines to calculate the minimax coefficients in order to approximate 1/x as a sum over exponential ...
subroutine, public check_exp_minimax_range(k, rc, ierr)
Check that a minimax approximation is available for given input k, Rc. ierr == 0: everything ok ierr ...
Routines to calculate frequency and time grids (integration points and weights) for correlation metho...
subroutine, public get_minimax_grid(para_env, unit_nr, homo, eigenval, num_integ_points, do_im_time, do_ri_sos_laplace_mp2, do_print, tau_tj, tau_wj, qs_env, do_gw_im_time, do_kpoints_cubic_rpa, e_fermi, tj, wj, weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, weights_sin_tf_t_to_w, regularization)
...
subroutine, public get_clenshaw_grid(para_env, para_env_rpa, unit_nr, homo, virtual, eigenval, num_integ_points, num_integ_group, color_rpa_group, fm_mat_s, my_do_gw, ext_scaling, a_scaling, tj, wj)
...
Routines to calculate MP2 energy with laplace approach.
subroutine, public sos_mp2_postprocessing(fm_mat_q, erpa, tau_wjquad)
...
Routines for calculating RI-MP2 gradients.
subroutine, public array2fm(mat2d, fm_struct, my_start_row, my_end_row, my_start_col, my_end_col, gd_row, gd_col, group_grid_2_mepos, ngroup_row, ngroup_col, fm_mat, integ_group_size, color_group, do_release_mat)
redistribute local part of array to fm
Types needed for MP2 calculations.
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Auxiliary routines needed for RPA-AXK given blacs_env to another.
subroutine, public compute_axk_ener(qs_env, fm_mat_q, fm_mat_q_gemm, dimen_ri, dimen_ia, para_env_sub, para_env_rpa, eig, fm_mat_s, homo, virtual, omega, mp2_env, mat_munu, unit_nr, e_axk_corr)
Main driver for RPA-AXK energies.
Routines to calculate RI-RPA and SOS-MP2 gradients.
subroutine, public rpa_grad_finalize(rpa_grad, mp2_env, para_env_sub, para_env, qs_env, gd_array, color_sub, do_ri_sos_laplace_mp2, homo, virtual)
...
subroutine, public rpa_grad_copy_q(fm_mat_q, rpa_grad)
...
pure subroutine, public rpa_grad_needed_mem(homo, virtual, dimen_ri, mem_per_rank, mem_per_repl, do_ri_sos_laplace_mp2)
Calculates the necessary minimum memory for the Gradient code ion MiB.
subroutine, public rpa_grad_create(rpa_grad, fm_mat_q, fm_mat_s, homo, virtual, mp2_env, eigenval, unit_nr, do_ri_sos_laplace_mp2)
Creates the arrays of a rpa_grad_type.
subroutine, public rpa_grad_matrix_operations(mp2_env, rpa_grad, do_ri_sos_laplace_mp2, fm_mat_q, fm_mat_q_gemm, dgemm_counter, fm_mat_s, omega, homo, virtual, eigenval, weight, unit_nr)
...
Routines to calculate image charge corrections.
subroutine, public calculate_ic_correction(eigenval, mat_sinvvsinv, t_3c_overl_nnp_ic, t_3c_overl_nnp_ic_reflected, gw_corr_lev_tot, gw_corr_lev_occ, gw_corr_lev_virt, homo, unit_nr, print_ic_values, para_env, do_alpha, do_beta)
...
Routines treating GW and RPA calculations with kpoints.
subroutine, public invert_eps_compute_w_and_erpa_kp(dimen_ri, num_integ_points, jquad, nkp, count_ev_sc_gw, para_env, erpa, tau_tj, tj, wj, weights_cos_tf_w_to_t, wkp_w, do_gw_im_time, do_ri_sigma_x, do_kpoints_from_gamma, cfm_mat_q, ikp_local, mat_p_omega, mat_p_omega_kp, qs_env, eps_filter_im_time, unit_nr, kpoints, fm_mat_minv_l_kpoints, fm_matrix_l_kpoints, fm_mat_w, fm_mat_ri_global_work, mat_minvvminv, fm_matrix_minv, fm_matrix_minv_vtrunc_minv)
...
subroutine, public get_bandstruc_and_k_dependent_mos(qs_env, eigenval_kp)
...
Routines for GW, continuous development [Jan Wilhelm].
subroutine, public deallocate_matrices_gw(fm_mat_s_gw_work, vec_w_gw, vec_sigma_c_gw, vec_omega_fit_gw, vec_sigma_x_minus_vxc_gw, eigenval_last, eigenval_scf, do_periodic, matrix_berry_re_mo_mo, matrix_berry_im_mo_mo, kpoints, vec_sigma_x_gw, my_do_gw)
...
subroutine, public allocate_matrices_gw(vec_sigma_c_gw, color_rpa_group, dimen_nm_gw, gw_corr_lev_occ, gw_corr_lev_virt, homo, nmo, num_integ_group, num_integ_points, unit_nr, gw_corr_lev_tot, num_fit_points, omega_max_fit, do_minimax_quad, do_periodic, do_ri_sigma_x, my_do_gw, first_cycle_periodic_correction, a_scaling, eigenval, tj, vec_omega_fit_gw, vec_sigma_x_gw, delta_corr, eigenval_last, eigenval_scf, vec_w_gw, fm_mat_s_gw, fm_mat_s_gw_work, para_env, mp2_env, kpoints, nkp, nkp_self_energy, do_kpoints_cubic_rpa, do_kpoints_from_gamma)
...
subroutine, public get_fermi_level_offset(fermi_level_offset, fermi_level_offset_input, eigenval, homo)
...
subroutine, public allocate_matrices_gw_im_time(gw_corr_lev_occ, gw_corr_lev_virt, homo, nmo, num_integ_points, unit_nr, ri_blk_sizes, do_ic_model, para_env, fm_mat_w, fm_mat_q, mo_coeff, t_3c_overl_int_ao_mo, t_3c_o_mo_compressed, t_3c_o_mo_ind, t_3c_overl_int_gw_ri, t_3c_overl_int_gw_ao, starts_array_mc, ends_array_mc, t_3c_overl_nnp_ic, t_3c_overl_nnp_ic_reflected, matrix_s, mat_w, t_3c_overl_int, t_3c_o_compressed, t_3c_o_ind, qs_env)
...
subroutine, public deallocate_matrices_gw_im_time(weights_cos_tf_w_to_t, weights_sin_tf_t_to_w, do_ic_model, do_kpoints_cubic_rpa, fm_mat_w, t_3c_overl_int_ao_mo, t_3c_o_mo_compressed, t_3c_o_mo_ind, t_3c_overl_int_gw_ri, t_3c_overl_int_gw_ao, t_3c_overl_nnp_ic, t_3c_overl_nnp_ic_reflected, mat_w, qs_env)
...
subroutine, public compute_qp_energies(vec_sigma_c_gw, count_ev_sc_gw, gw_corr_lev_occ, gw_corr_lev_tot, gw_corr_lev_virt, homo, nmo, num_fit_points, num_integ_points, unit_nr, do_apply_ic_corr_to_gw, do_im_time, do_periodic, do_ri_sigma_x, first_cycle_periodic_correction, e_fermi, eps_filter, fermi_level_offset, delta_corr, eigenval, eigenval_last, eigenval_scf, iter_sc_gw0, exit_ev_gw, tau_tj, tj, vec_omega_fit_gw, vec_sigma_x_gw, ic_corr_list, weights_cos_tf_t_to_w, weights_sin_tf_t_to_w, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, fm_mo_coeff_occ, fm_mo_coeff_virt, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, mo_coeff, fm_mat_w, para_env, para_env_rpa, mat_dm, mat_minvvminv, t_3c_o, t_3c_m, t_3c_overl_int_ao_mo, t_3c_o_compressed, t_3c_o_mo_compressed, t_3c_o_ind, t_3c_o_mo_ind, t_3c_overl_int_gw_ri, t_3c_overl_int_gw_ao, matrix_berry_im_mo_mo, matrix_berry_re_mo_mo, mat_w, matrix_s, kpoints, mp2_env, qs_env, nkp_self_energy, do_kpoints_cubic_rpa, starts_array_mc, ends_array_mc)
...
subroutine, public compute_w_cubic_gw(fm_mat_w, fm_mat_q, fm_mat_work, dimen_ri, fm_mat_l, num_integ_points, tj, tau_tj, weights_cos_tf_w_to_t, jquad, omega)
...
subroutine, public compute_gw_self_energy(vec_sigma_c_gw, dimen_nm_gw, dimen_ri, gw_corr_lev_occ, gw_corr_lev_virt, homo, jquad, nmo, num_fit_points, num_integ_points, do_bse, do_im_time, do_periodic, first_cycle_periodic_correction, fermi_level_offset, omega, eigenval, delta_corr, vec_omega_fit_gw, vec_w_gw, wj, fm_mat_q, fm_mat_q_static_bse, fm_mat_r_gw, fm_mat_s_gw, fm_mat_s_gw_work, mo_coeff, para_env, para_env_rpa, matrix_berry_im_mo_mo, matrix_berry_re_mo_mo, kpoints, qs_env, mp2_env)
...
Routines needed for cubic-scaling RPA and SOS-Laplace-MP2 forces.
subroutine, public calc_laplace_loop_forces(force_data, mat_p_omega, t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, starts_array_mc, ends_array_mc, starts_array_mc_block, ends_array_mc_block, num_integ_points, nmo, eigenval, tau_tj, tau_wj, cut_memory, pspin, qspin, open_shell, unit_nr, dbcsr_time, dbcsr_nflop, mp2_env, qs_env)
Updates the cubic-scaling SOS-Laplace-MP2 contribution to the forces at each quadrature point.
subroutine, public calc_post_loop_forces(force_data, unit_nr, qs_env)
All the forces that can be calculated after the loop on the Laplace quaradture, using terms collected...
subroutine, public calc_rpa_loop_forces(force_data, mat_p_omega, t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, starts_array_mc, ends_array_mc, starts_array_mc_block, ends_array_mc_block, num_integ_points, nmo, eigenval, e_fermi, weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, tj, wj, tau_tj, cut_memory, ispin, open_shell, unit_nr, dbcsr_time, dbcsr_nflop, mp2_env, qs_env)
Updates the cubic-scaling RPA contribution to the forces at each quadrature point....
subroutine, public keep_initial_quad(tj, wj, tau_tj, tau_wj, weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, do_laplace, do_im_time, num_integ_points, unit_nr, qs_env)
Overwrites the "optimal" Laplace quadrature with that of the first step.
subroutine, public init_im_time_forces(force_data, fm_matrix_pq, t_3c_m, unit_nr, mp2_env, qs_env)
Initializes and pre-calculates all needed tensors for the forces.
Types needed for cubic-scaling RPA and SOS-Laplace-MP2 forces.
subroutine, public im_time_force_release(force_data)
Cleans everything up.
Routines for low-scaling RPA/GW with imaginary time.
subroutine, public zero_mat_p_omega(mat_p_omega)
...
subroutine, public compute_mat_p_omega(mat_p_omega, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, fm_mo_coeff_occ, fm_mo_coeff_virt, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, mat_p_global, matrix_s, ispin, t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, starts_array_mc, ends_array_mc, starts_array_mc_block, ends_array_mc_block, weights_cos_tf_t_to_w, tj, tau_tj, e_fermi, eps_filter, alpha, eps_filter_im_time, eigenval, nmo, num_integ_points, cut_memory, unit_nr, mp2_env, para_env, qs_env, do_kpoints_from_gamma, index_to_cell_3c, cell_to_index_3c, has_mat_p_blocks, do_ri_sos_laplace_mp2, dbcsr_time, dbcsr_nflop)
...
Routines to calculate RI-RPA energy.
subroutine, public rpa_ri_compute_en(qs_env, erpa, mp2_env, bib_c, bib_c_gw, bib_c_bse_ij, bib_c_bse_ab, para_env, para_env_sub, color_sub, gd_array, gd_b_virtual, gd_b_all, gd_b_occ_bse, gd_b_virt_bse, mo_coeff, fm_matrix_pq, fm_matrix_l_kpoints, fm_matrix_minv_l_kpoints, fm_matrix_minv, fm_matrix_minv_vtrunc_minv, kpoints, eigenval, nmo, homo, dimen_ri, dimen_ri_red, gw_corr_lev_occ, gw_corr_lev_virt, unit_nr, do_ri_sos_laplace_mp2, my_do_gw, do_im_time, do_bse, matrix_s, mat_munu, mat_p_global, t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, starts_array_mc, ends_array_mc, starts_array_mc_block, ends_array_mc_block, calc_forces)
...
Utility functions for RPA calculations.
subroutine, public compute_erpa_by_freq_int(dimen_ri, trace_qomega, fm_mat_q, para_env_rpa, erpa, wjquad)
...
subroutine, public alloc_im_time(qs_env, para_env, dimen_ri, dimen_ri_red, num_integ_points, nspins, fm_mat_q, fm_mo_coeff_occ, fm_mo_coeff_virt, fm_matrix_minv_l_kpoints, fm_matrix_l_kpoints, mat_p_global, t_3c_o, matrix_s, kpoints, eps_filter_im_time, cut_memory, nkp, num_cells_dm, num_3c_repl, size_p, ikp_local, index_to_cell_3c, cell_to_index_3c, col_blk_size, do_ic_model, do_kpoints_cubic_rpa, do_kpoints_from_gamma, do_ri_sigma_x, my_open_shell, has_mat_p_blocks, wkp_w, cfm_mat_q, fm_mat_minv_l_kpoints, fm_mat_l_kpoints, fm_mat_ri_global_work, fm_mat_work, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, mat_dm, mat_l, mat_m_p_munu_occ, mat_m_p_munu_virt, mat_minvvminv, mat_p_omega, mat_p_omega_kp, mat_work, mo_coeff, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, homo, nmo)
...
subroutine, public dealloc_im_time(fm_mo_coeff_occ, fm_mo_coeff_virt, fm_scaled_dm_occ_tau, fm_scaled_dm_virt_tau, index_to_cell_3c, cell_to_index_3c, do_ic_model, do_kpoints_cubic_rpa, do_kpoints_from_gamma, do_ri_sigma_x, has_mat_p_blocks, wkp_w, cfm_mat_q, fm_mat_minv_l_kpoints, fm_mat_l_kpoints, fm_matrix_minv, fm_matrix_minv_vtrunc_minv, fm_mat_ri_global_work, fm_mat_work, fm_mo_coeff_occ_scaled, fm_mo_coeff_virt_scaled, mat_dm, mat_l, mat_minvvminv, mat_p_omega, mat_p_omega_kp, t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, mat_work, qs_env)
...
subroutine, public contract_p_omega_with_mat_l(mat_p_omega, mat_l, mat_work, eps_filter_im_time, fm_mat_work, dimen_ri, dimen_ri_red, fm_mat_l, fm_mat_q)
...
subroutine, public calc_mat_q(fm_mat_s, do_ri_sos_laplace_mp2, first_cycle, iter_sc_gw0, virtual, eigenval, eigenval_scf, homo, omega, omega_old, jquad, mm_style, dimen_ri, dimen_ia, alpha, fm_mat_q, fm_mat_q_gemm, do_bse, fm_mat_q_static_bse_gemm, dgemm_counter, num_integ_points)
...
subroutine, public remove_scaling_factor_rpa(fm_mat_s, virtual, eigenval_last, homo, omega_old)
...
subroutine, public q_trace_and_add_unit_matrix(dimen_ri, trace_qomega, fm_mat_q, para_env_rpa)
...
All kind of helpful little routines.
pure integer function, dimension(2), public get_limit(m, n, me)
divide m entries into n parts, return size of part me
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Represent a complex full matrix.
keeps the information about the structure of a full matrix
Contains information about kpoints.
stores all the informations relevant to an mpi environment