46 USE dbt_api,
ONLY: dbt_destroy,&
68#include "./base/base_uses.f90"
74 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'rpa_util'
129 SUBROUTINE alloc_im_time(qs_env, para_env, dimen_RI, dimen_RI_red, num_integ_points, nspins, &
130 fm_mat_Q, cfm_mo_coeff, &
131 fm_matrix_Minv_L_kpoints, fm_matrix_L_kpoints, mat_P_global, &
132 t_3c_O, matrix_s, kpoints, eps_filter_im_time, &
133 cut_memory, nkp, num_cells_dm, num_3c_repl, &
138 do_ic_model, do_kpoints_cubic_RPA, &
139 do_kpoints_from_Gamma, do_ri_Sigma_x, my_open_shell, &
140 has_mat_P_blocks, wkp_W, &
141 cfm_mat_Q, fm_mat_Minv_L_kpoints, fm_mat_L_kpoints, &
142 fm_mat_RI_global_work, fm_mat_work, mat_dm, mat_L, mat_M_P_munu_occ, mat_M_P_munu_virt, &
143 mat_MinvVMinv, mat_P_omega, mat_P_omega_kp, &
148 INTEGER,
INTENT(IN) :: dimen_ri, dimen_ri_red, &
149 num_integ_points, nspins
151 TYPE(
cp_cfm_type),
ALLOCATABLE,
DIMENSION(:) :: cfm_mo_coeff
152 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_matrix_minv_l_kpoints, &
155 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:, :), &
156 INTENT(INOUT) :: t_3c_o
159 REAL(kind=
dp),
INTENT(IN) :: eps_filter_im_time
160 INTEGER,
INTENT(IN) :: cut_memory
161 INTEGER,
INTENT(OUT) :: nkp, num_cells_dm, num_3c_repl, size_p, &
163 INTEGER,
ALLOCATABLE,
DIMENSION(:, :),
INTENT(OUT) :: index_to_cell_3c
164 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :), &
165 INTENT(OUT) :: cell_to_index_3c
166 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size
167 LOGICAL,
INTENT(IN) :: do_ic_model, do_kpoints_cubic_rpa, &
168 do_kpoints_from_gamma, do_ri_sigma_x, &
170 LOGICAL,
ALLOCATABLE,
DIMENSION(:, :, :, :, :), &
171 INTENT(OUT) :: has_mat_p_blocks
172 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:), &
175 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_mat_minv_l_kpoints, fm_mat_l_kpoints
176 TYPE(
cp_fm_type),
INTENT(OUT) :: fm_mat_ri_global_work, fm_mat_work
177 TYPE(
dbcsr_p_type),
INTENT(OUT) :: mat_dm, mat_l, mat_m_p_munu_occ, &
178 mat_m_p_munu_virt, mat_minvvminv
180 DIMENSION(:, :, :) :: mat_p_omega
181 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: mat_p_omega_kp
183 TYPE(
cp_fm_type),
DIMENSION(:),
INTENT(IN) :: mo_coeff
185 CHARACTER(LEN=*),
PARAMETER :: routinen =
'alloc_im_time'
187 INTEGER :: cell_grid_dm(3), first_ikp_local, &
188 handle, i_dim, i_kp, ispin, jquad, &
189 nspins_p_omega, periodic(3)
190 INTEGER,
DIMENSION(:),
POINTER :: row_blk_size
191 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: wkp_v
195 CALL timeset(routinen, handle)
197 ALLOCATE (cfm_mo_coeff(nspins))
199 DO ispin = 1,
SIZE(mo_coeff)
200 CALL create_mo_coeff(cfm_mo_coeff(ispin), mo_coeff(ispin))
203 num_3c_repl =
SIZE(t_3c_o, 2)
205 IF (do_kpoints_cubic_rpa)
THEN
209 cell_grid_dm(i_dim) = (kpoints%nkp_grid(i_dim)/2)*2 - 1
211 num_cells_dm = cell_grid_dm(1)*cell_grid_dm(2)*cell_grid_dm(3)
212 ALLOCATE (index_to_cell_3c(3,
SIZE(kpoints%index_to_cell, 2)))
213 cpassert(
SIZE(kpoints%index_to_cell, 1) == 3)
214 index_to_cell_3c(:, :) = kpoints%index_to_cell(:, :)
215 ALLOCATE (cell_to_index_3c(lbound(kpoints%cell_to_index, 1):ubound(kpoints%cell_to_index, 1), &
216 lbound(kpoints%cell_to_index, 2):ubound(kpoints%cell_to_index, 2), &
217 lbound(kpoints%cell_to_index, 3):ubound(kpoints%cell_to_index, 3)))
218 cell_to_index_3c(:, :, :) = kpoints%cell_to_index(:, :, :)
221 ALLOCATE (index_to_cell_3c(3, 1))
222 index_to_cell_3c(:, 1) = 0
223 ALLOCATE (cell_to_index_3c(0:0, 0:0, 0:0))
224 cell_to_index_3c(0, 0, 0) = 1
228 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
230 CALL get_sub_para_kp(fm_struct_sub_kp, para_env, dimen_ri, ikp_local, first_ikp_local)
240 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
244 CALL get_cell(cell=cell, periodic=periodic)
249 CALL compute_wkp_w(qs_env, wkp_w, wkp_v, kpoints, cell%h_inv, periodic)
256 IF (do_kpoints_cubic_rpa)
THEN
257 size_p = max(num_cells_dm/2 + 1, nkp)
258 ELSE IF (do_kpoints_from_gamma)
THEN
259 size_p = max(3**(periodic(1) + periodic(2) + periodic(3)), nkp)
265 IF (my_open_shell) nspins_p_omega = 2
267 ALLOCATE (mat_p_omega(num_integ_points, size_p, nspins_p_omega))
268 DO ispin = 1, nspins_p_omega
270 DO jquad = 1, num_integ_points
271 NULLIFY (mat_p_omega(jquad, i_kp, ispin)%matrix)
272 ALLOCATE (mat_p_omega(jquad, i_kp, ispin)%matrix)
273 CALL dbcsr_create(matrix=mat_p_omega(jquad, i_kp, ispin)%matrix, &
274 template=mat_p_global%matrix)
275 CALL dbcsr_set(mat_p_omega(jquad, i_kp, ispin)%matrix, 0.0_dp)
280 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
281 CALL alloc_mat_p_omega(mat_p_omega_kp, 2, size_p, mat_p_global%matrix)
284 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
285 CALL cp_fm_create(fm_mat_ri_global_work, fm_matrix_minv_l_kpoints(1, 1)%matrix_struct, set_zero=.true.)
288 ALLOCATE (has_mat_p_blocks(num_cells_dm/2 + 1, cut_memory, cut_memory, num_3c_repl, num_3c_repl))
289 has_mat_p_blocks = .true.
291 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
292 CALL reorder_mat_l(fm_mat_minv_l_kpoints, fm_matrix_minv_l_kpoints, fm_mat_q%matrix_struct, para_env, mat_l, &
293 mat_p_global%matrix, dimen_ri, dimen_ri_red, first_ikp_local, ikp_local, fm_struct_sub_kp, &
294 allocate_mat_l=.false.)
296 CALL reorder_mat_l(fm_mat_l_kpoints, fm_matrix_l_kpoints, fm_mat_q%matrix_struct, para_env, mat_l, &
297 mat_p_global%matrix, dimen_ri, dimen_ri_red, first_ikp_local, ikp_local, fm_struct_sub_kp)
302 CALL reorder_mat_l(fm_mat_minv_l_kpoints, fm_matrix_minv_l_kpoints, fm_mat_q%matrix_struct, para_env, mat_l, &
303 mat_p_global%matrix, dimen_ri, dimen_ri_red, first_ikp_local)
307 IF (dimen_ri == dimen_ri_red)
THEN
308 CALL cp_fm_create(fm_mat_work, fm_mat_q%matrix_struct, set_zero=.true.)
311 CALL cp_fm_create(fm_mat_work, fm_mat_q%matrix_struct, nrow=dimen_ri, ncol=dimen_ri_red, &
317 IF (.NOT. (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma))
THEN
318 CALL dbcsr_get_info(mat_l%matrix, col_blk_size=col_blk_size, row_blk_size=row_blk_size)
323 CALL dbcsr_create(mat_work, template=mat_l%matrix, row_blk_size=col_blk_size, col_blk_size=row_blk_size)
326 IF (do_ri_sigma_x .OR. do_ic_model)
THEN
328 NULLIFY (mat_minvvminv%matrix)
329 ALLOCATE (mat_minvvminv%matrix)
330 CALL dbcsr_create(mat_minvvminv%matrix, template=mat_p_global%matrix)
331 CALL dbcsr_set(mat_minvvminv%matrix, 0.0_dp)
334 IF (.NOT. do_kpoints_from_gamma)
THEN
337 CALL dbcsr_multiply(
"T",
"N", 1.0_dp, mat_l%matrix, mat_l%matrix, &
338 0.0_dp, mat_minvvminv%matrix, filter_eps=eps_filter_im_time)
344 IF (do_ri_sigma_x)
THEN
346 NULLIFY (mat_dm%matrix)
347 ALLOCATE (mat_dm%matrix)
348 CALL dbcsr_create(mat_dm%matrix, template=matrix_s(1)%matrix)
352 CALL timestop(handle)
361 SUBROUTINE create_mo_coeff(cfm_mo_coeff, mo_coeff)
366 CHARACTER(LEN=*),
PARAMETER :: routinen =
'create_mo_coeff'
370 CALL timeset(routinen, handle)
373 CALL cp_fm_to_cfm(msourcer=mo_coeff, mtarget=cfm_mo_coeff)
375 CALL timestop(handle)
377 END SUBROUTINE create_mo_coeff
386 SUBROUTINE alloc_mat_p_omega(mat_P_omega, num_integ_points, size_P, template)
387 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: mat_p_omega
388 INTEGER,
INTENT(IN) :: num_integ_points, size_p
391 CHARACTER(LEN=*),
PARAMETER :: routinen =
'alloc_mat_P_omega'
393 INTEGER :: handle, i_kp, jquad
395 CALL timeset(routinen, handle)
397 NULLIFY (mat_p_omega)
400 DO jquad = 1, num_integ_points
401 ALLOCATE (mat_p_omega(jquad, i_kp)%matrix)
402 CALL dbcsr_create(matrix=mat_p_omega(jquad, i_kp)%matrix, &
404 CALL dbcsr_set(mat_p_omega(jquad, i_kp)%matrix, 0.0_dp)
408 CALL timestop(handle)
410 END SUBROUTINE alloc_mat_p_omega
427 SUBROUTINE reorder_mat_l(fm_mat_L, fm_matrix_Minv_L_kpoints, fm_struct_template, para_env, mat_L, mat_template, &
428 dimen_RI, dimen_RI_red, first_ikp_local, ikp_local, fm_struct_sub_kp, allocate_mat_L)
429 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_mat_l, fm_matrix_minv_l_kpoints
434 INTEGER,
INTENT(IN) :: dimen_ri, dimen_ri_red, first_ikp_local
435 INTEGER,
OPTIONAL :: ikp_local
437 LOGICAL,
INTENT(IN),
OPTIONAL :: allocate_mat_l
439 CHARACTER(LEN=*),
PARAMETER :: routinen =
'reorder_mat_L'
441 INTEGER :: handle, ikp, j_size, nblk
442 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
443 LOGICAL :: do_kpoints, my_allocate_mat_l
446 TYPE(
cp_fm_type) :: fm_mat_l_transposed, fmdummy
448 CALL timeset(routinen, handle)
451 IF (
PRESENT(ikp_local) .AND.
PRESENT(fm_struct_sub_kp))
THEN
457 IF (dimen_ri == dimen_ri_red)
THEN
458 fm_struct => fm_struct_template
461 CALL cp_fm_struct_create(fm_struct, nrow_global=dimen_ri_red, ncol_global=dimen_ri, template_fmstruct=fm_struct_template)
465 ALLOCATE (fm_mat_l(
SIZE(fm_matrix_minv_l_kpoints, 1),
SIZE(fm_matrix_minv_l_kpoints, 2)))
466 DO j_size = 1,
SIZE(fm_matrix_minv_l_kpoints, 2)
467 DO ikp = 1,
SIZE(fm_matrix_minv_l_kpoints, 1)
469 IF (ikp == first_ikp_local .OR. ikp_local == -1)
THEN
470 CALL cp_fm_create(fm_mat_l(ikp, j_size), fm_struct_sub_kp)
481 IF (dimen_ri == dimen_ri_red)
THEN
482 fm_struct => fm_mat_l(first_ikp_local, 1)%matrix_struct
489 template_fmstruct=fm_mat_l(first_ikp_local, 1)%matrix_struct)
503 DO j_size = 1,
SIZE(fm_matrix_minv_l_kpoints, 2)
504 DO ikp = 1,
SIZE(fm_matrix_minv_l_kpoints, 1)
506 IF (ikp_local == ikp .OR. ikp_local == -1)
THEN
507 CALL cp_fm_copy_general(fm_matrix_minv_l_kpoints(ikp, j_size), fm_mat_l_transposed, para_env)
508 CALL cp_fm_to_fm(fm_mat_l_transposed, fm_mat_l(ikp, j_size))
513 CALL cp_fm_copy_general(fm_matrix_minv_l_kpoints(ikp, j_size), fm_mat_l_transposed, blacs_env%para_env)
524 my_allocate_mat_l = .true.
525 IF (
PRESENT(allocate_mat_l)) my_allocate_mat_l = allocate_mat_l
527 IF (my_allocate_mat_l)
THEN
529 NULLIFY (mat_l%matrix)
530 ALLOCATE (mat_l%matrix)
531 IF (dimen_ri == dimen_ri_red)
THEN
534 CALL dbcsr_get_info(mat_template, nblkrows_total=nblk, col_blk_size=col_blk_size)
536 CALL calculate_equal_blk_size(row_blk_size, dimen_ri_red, nblk)
538 CALL dbcsr_create(mat_l%matrix, template=mat_template, row_blk_size=row_blk_size, col_blk_size=col_blk_size)
540 DEALLOCATE (row_blk_size)
543 IF (.NOT. (do_kpoints))
THEN
549 CALL timestop(handle)
551 END SUBROUTINE reorder_mat_l
559 SUBROUTINE calculate_equal_blk_size(blk_size_new, dimen_RI_red, nblk)
560 INTEGER,
DIMENSION(:),
POINTER :: blk_size_new
561 INTEGER,
INTENT(IN) :: dimen_ri_red, nblk
563 INTEGER :: col_per_blk, remainder
565 NULLIFY (blk_size_new)
566 ALLOCATE (blk_size_new(nblk))
568 remainder = mod(dimen_ri_red, nblk)
569 col_per_blk = dimen_ri_red/nblk
572 IF (remainder > 0) blk_size_new(1:remainder) = col_per_blk + 1
573 blk_size_new(remainder + 1:nblk) = col_per_blk
575 END SUBROUTINE calculate_equal_blk_size
600 SUBROUTINE calc_mat_q(fm_mat_S, do_ri_sos_laplace_mp2, first_cycle, virtual, &
601 Eigenval, homo, omega, omega_old, jquad, mm_style, dimen_RI, dimen_ia, alpha, fm_mat_Q, fm_mat_Q_gemm, &
602 do_bse, fm_mat_Q_static_bse_gemm, dgemm_counter, &
603 num_integ_points, count_ev_sc_GW)
605 LOGICAL,
INTENT(IN) :: do_ri_sos_laplace_mp2, first_cycle
606 INTEGER,
INTENT(IN) :: virtual
607 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: eigenval
608 INTEGER,
INTENT(IN) :: homo
609 REAL(kind=
dp),
INTENT(IN) :: omega, omega_old
610 INTEGER,
INTENT(IN) :: jquad, mm_style, dimen_ri, dimen_ia
611 REAL(kind=
dp),
INTENT(IN) :: alpha
612 TYPE(
cp_fm_type),
INTENT(IN) :: fm_mat_q, fm_mat_q_gemm
613 LOGICAL,
INTENT(IN) :: do_bse
614 TYPE(
cp_fm_type),
INTENT(IN) :: fm_mat_q_static_bse_gemm
616 INTEGER,
INTENT(IN) :: num_integ_points, count_ev_sc_gw
618 CHARACTER(LEN=*),
PARAMETER :: routinen =
'calc_mat_Q'
622 CALL timeset(routinen, handle)
624 IF (do_ri_sos_laplace_mp2)
THEN
629 homo, omega, omega_old)
632 CALL contract_s_to_q(mm_style, dimen_ri, dimen_ia, alpha, fm_mat_s, fm_mat_q_gemm, &
633 fm_mat_q, dgemm_counter)
638 IF (do_bse .AND. jquad == num_integ_points .AND. count_ev_sc_gw == 1)
THEN
639 CALL cp_fm_to_fm(fm_mat_q_gemm, fm_mat_q_static_bse_gemm)
641 CALL timestop(handle)
655 INTEGER,
INTENT(IN) :: virtual
656 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: eigenval_last
657 INTEGER,
INTENT(IN) :: homo
658 REAL(kind=
dp),
INTENT(IN) :: omega_old
660 CHARACTER(LEN=*),
PARAMETER :: routinen =
'remove_scaling_factor_rpa'
662 INTEGER :: avirt, handle, i_global, iib, iocc, &
664 INTEGER,
DIMENSION(:),
POINTER :: col_indices
665 REAL(kind=
dp) :: eigen_diff
667 CALL timeset(routinen, handle)
671 ncol_local=ncol_local, &
672 col_indices=col_indices)
676 DO iib = 1, ncol_local
677 i_global = col_indices(iib)
679 iocc = max(1, i_global - 1)/virtual + 1
680 avirt = i_global - (iocc - 1)*virtual
681 eigen_diff = eigenval_last(avirt + homo) - eigenval_last(iocc)
683 fm_mat_s%local_data(:, iib) = fm_mat_s%local_data(:, iib)/ &
684 sqrt(eigen_diff/(eigen_diff**2 + omega_old**2))
688 CALL timestop(handle)
705 LOGICAL,
INTENT(IN) :: first_cycle
706 INTEGER,
INTENT(IN) :: virtual
707 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: eigenval
708 INTEGER,
INTENT(IN) :: homo
709 REAL(kind=
dp),
INTENT(IN) :: omega, omega_old
711 CHARACTER(LEN=*),
PARAMETER :: routinen =
'calc_fm_mat_S_rpa'
713 INTEGER :: avirt, handle, i_global, iib, iocc, &
715 INTEGER,
DIMENSION(:),
POINTER :: col_indices
716 REAL(kind=
dp) :: eigen_diff
718 CALL timeset(routinen, handle)
722 ncol_local=ncol_local, &
723 col_indices=col_indices)
726 IF (first_cycle)
THEN
731 DO iib = 1, ncol_local
732 i_global = col_indices(iib)
734 iocc = max(1, i_global - 1)/virtual + 1
735 avirt = i_global - (iocc - 1)*virtual
736 eigen_diff = eigenval(avirt + homo) - eigenval(iocc)
738 fm_mat_s%local_data(:, iib) = fm_mat_s%local_data(:, iib)* &
739 sqrt(eigen_diff/(eigen_diff**2 + omega**2))
747 DO iib = 1, ncol_local
748 i_global = col_indices(iib)
750 iocc = max(1, i_global - 1)/virtual + 1
751 avirt = i_global - (iocc - 1)*virtual
752 eigen_diff = eigenval(avirt + homo) - eigenval(iocc)
754 fm_mat_s%local_data(:, iib) = fm_mat_s%local_data(:, iib)* &
755 sqrt((eigen_diff**2 + omega_old**2)/(eigen_diff**2 + omega**2))
760 CALL timestop(handle)
775 SUBROUTINE contract_s_to_q(mm_style, dimen_RI, dimen_ia, alpha, fm_mat_S, fm_mat_Q_gemm, &
776 fm_mat_Q, dgemm_counter)
778 INTEGER,
INTENT(IN) :: mm_style, dimen_ri, dimen_ia
779 REAL(kind=
dp),
INTENT(IN) :: alpha
780 TYPE(
cp_fm_type),
INTENT(IN) :: fm_mat_s, fm_mat_q_gemm, fm_mat_q
783 CHARACTER(LEN=*),
PARAMETER :: routinen =
'contract_S_to_Q'
787 CALL timeset(routinen, handle)
790 SELECT CASE (mm_style)
793 CALL parallel_gemm(transa=
"N", transb=
"T", m=dimen_ri, n=dimen_ri, k=dimen_ia, alpha=alpha, &
794 matrix_a=fm_mat_s, matrix_b=fm_mat_s, beta=0.0_dp, &
795 matrix_c=fm_mat_q_gemm)
798 CALL cp_fm_syrk(uplo=
'U', trans=
'N', k=dimen_ia, alpha=alpha, matrix_a=fm_mat_s, &
799 ia=1, ja=1, beta=0.0_dp, matrix_c=fm_mat_q_gemm)
801 cpabort(
"Unknown mm_style for contract_S_to_Q")
808 fm_mat_q_gemm%matrix_struct%context)
810 CALL timestop(handle)
812 END SUBROUTINE contract_s_to_q
822 INTEGER,
INTENT(IN) :: dimen_ri
823 REAL(kind=
dp),
DIMENSION(dimen_RI),
INTENT(OUT) :: trace_qomega
826 CHARACTER(LEN=*),
PARAMETER :: routinen =
'Q_trace_and_add_unit_matrix'
828 INTEGER :: handle, i_global, iib, j_global, jjb, &
829 ncol_local, nrow_local
830 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
833 CALL timeset(routinen, handle)
836 nrow_local=nrow_local, &
837 ncol_local=ncol_local, &
838 row_indices=row_indices, &
839 col_indices=col_indices, &
843 trace_qomega = 0.0_dp
846 DO jjb = 1, ncol_local
847 j_global = col_indices(jjb)
848 DO iib = 1, nrow_local
849 i_global = row_indices(iib)
850 IF (j_global == i_global .AND. i_global <= dimen_ri)
THEN
851 trace_qomega(i_global) = fm_mat_q%local_data(iib, jjb)
852 fm_mat_q%local_data(iib, jjb) = fm_mat_q%local_data(iib, jjb) + 1.0_dp
856 CALL para_env%sum(trace_qomega)
858 CALL timestop(handle)
873 INTEGER,
INTENT(IN) :: dimen_ri
874 REAL(kind=
dp),
DIMENSION(dimen_RI),
INTENT(IN) :: trace_qomega
877 REAL(kind=
dp),
INTENT(INOUT) :: erpa
878 REAL(kind=
dp),
INTENT(IN) :: wjquad
880 CHARACTER(LEN=*),
PARAMETER :: routinen =
'compute_Erpa_by_freq_int'
882 INTEGER :: handle, i_global, iib, info_chol, &
883 j_global, jjb, ncol_local, nrow_local
884 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
885 REAL(kind=
dp) :: fcomega
886 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: q_log
888 CALL timeset(routinen, handle)
891 nrow_local=nrow_local, &
892 ncol_local=ncol_local, &
893 row_indices=row_indices, &
894 col_indices=col_indices)
898 IF (info_chol /= 0)
THEN
899 CALL cp_warn(__location__, &
900 "The Cholesky decomposition before inverting the RPA matrix / dielectric "// &
901 "function failed. "// &
902 "In case of low-scaling RPA/GW, decreasing EPS_FILTER in the &LOW_SCALING "// &
904 "increase the overall accuracy making the matrix positive definite. "// &
908 cpassert(info_chol == 0)
910 ALLOCATE (q_log(dimen_ri))
914 DO jjb = 1, ncol_local
915 j_global = col_indices(jjb)
916 DO iib = 1, nrow_local
917 i_global = row_indices(iib)
918 IF (j_global == i_global .AND. i_global <= dimen_ri)
THEN
919 q_log(i_global) = 2.0_dp*log(fm_mat_q%local_data(iib, jjb))
923 CALL para_env_rpa%sum(q_log)
929 IF (
modulo(iib, para_env_rpa%num_pe) /= para_env_rpa%mepos) cycle
930 fcomega = fcomega + (q_log(iib) - trace_qomega(iib))/2.0_dp
932 erpa = erpa + fcomega*wjquad
936 CALL timestop(handle)
948 SUBROUTINE get_sub_para_kp(fm_struct_sub_kp, para_env, dimen_RI, &
949 ikp_local, first_ikp_local)
952 INTEGER,
INTENT(IN) :: dimen_ri
953 INTEGER,
INTENT(OUT) :: ikp_local, first_ikp_local
955 CHARACTER(len=*),
PARAMETER :: routinen =
'get_sub_para_kp'
957 INTEGER :: color_sub_kp, handle, num_proc_per_kp
961 CALL timeset(routinen, handle)
966 num_proc_per_kp = para_env%num_pe
974 color_sub_kp = para_env%mepos/num_proc_per_kp
975 ALLOCATE (para_env_sub_kp)
976 CALL para_env_sub_kp%from_split(para_env, color_sub_kp)
981 NULLIFY (blacs_env_sub_kp)
985 NULLIFY (fm_struct_sub_kp)
986 CALL cp_fm_struct_create(fm_struct_sub_kp, context=blacs_env_sub_kp, nrow_global=dimen_ri, &
987 ncol_global=dimen_ri, para_env=para_env_sub_kp)
1008 CALL timestop(handle)
1010 END SUBROUTINE get_sub_para_kp
1043 cell_to_index_3c, do_ic_model, &
1044 do_kpoints_cubic_RPA, do_kpoints_from_Gamma, do_ri_Sigma_x, &
1046 wkp_W, cfm_mat_Q, fm_mat_Minv_L_kpoints, fm_mat_L_kpoints, &
1047 fm_matrix_Minv, fm_matrix_Minv_Vtrunc_Minv, &
1048 fm_mat_RI_global_work, fm_mat_work, mat_dm, mat_L, &
1049 mat_MinvVMinv, mat_P_omega, mat_P_omega_kp, &
1050 t_3c_M, t_3c_O, t_3c_O_compressed, t_3c_O_ind, &
1053 TYPE(
cp_cfm_type),
DIMENSION(:),
INTENT(INOUT) :: cfm_mo_coeff
1054 INTEGER,
ALLOCATABLE,
DIMENSION(:, :), &
1055 INTENT(INOUT) :: index_to_cell_3c
1056 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :), &
1057 INTENT(INOUT) :: cell_to_index_3c
1058 LOGICAL,
INTENT(IN) :: do_ic_model, do_kpoints_cubic_rpa, &
1059 do_kpoints_from_gamma, do_ri_sigma_x
1060 LOGICAL,
ALLOCATABLE,
DIMENSION(:, :, :, :, :), &
1061 INTENT(INOUT) :: has_mat_p_blocks
1062 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:), &
1063 INTENT(INOUT) :: wkp_w
1065 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:, :) :: fm_mat_minv_l_kpoints, fm_mat_l_kpoints, &
1067 fm_matrix_minv_vtrunc_minv
1068 TYPE(
cp_fm_type),
INTENT(INOUT) :: fm_mat_ri_global_work, fm_mat_work
1069 TYPE(
dbcsr_p_type),
INTENT(INOUT) :: mat_dm, mat_l, mat_minvvminv
1071 DIMENSION(:, :, :),
INTENT(INOUT) :: mat_p_omega
1072 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: mat_p_omega_kp
1073 TYPE(dbt_type) :: t_3c_m
1074 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:, :) :: t_3c_o
1076 DIMENSION(:, :, :),
INTENT(INOUT) :: t_3c_o_compressed
1078 DIMENSION(:, :, :),
INTENT(INOUT) :: t_3c_o_ind
1082 CHARACTER(LEN=*),
PARAMETER :: routinen =
'dealloc_im_time'
1084 INTEGER :: cut_memory, handle, i_kp, i_mem, i_size, &
1085 ispin, j_size, jquad, nspins, unused
1086 LOGICAL :: my_open_shell
1088 CALL timeset(routinen, handle)
1090 nspins =
SIZE(cfm_mo_coeff)
1091 my_open_shell = (nspins == 2)
1093 DO ispin = 1,
SIZE(cfm_mo_coeff)
1099 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
1106 IF (.NOT. (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma))
THEN
1108 DEALLOCATE (mat_work)
1112 DEALLOCATE (mat_l%matrix)
1114 IF (do_ri_sigma_x .OR. do_ic_model)
THEN
1116 DEALLOCATE (mat_minvvminv%matrix)
1118 IF (do_ri_sigma_x)
THEN
1120 DEALLOCATE (mat_dm%matrix)
1123 DEALLOCATE (index_to_cell_3c, cell_to_index_3c)
1125 IF (
ALLOCATED(mat_p_omega))
THEN
1126 DO ispin = 1,
SIZE(mat_p_omega, 3)
1127 DO i_kp = 1,
SIZE(mat_p_omega, 2)
1128 DO jquad = 1,
SIZE(mat_p_omega, 1)
1133 DEALLOCATE (mat_p_omega)
1136 DO j_size = 1,
SIZE(t_3c_o, 2)
1137 DO i_size = 1,
SIZE(t_3c_o, 1)
1138 CALL dbt_destroy(t_3c_o(i_size, j_size))
1143 CALL dbt_destroy(t_3c_m)
1145 DEALLOCATE (has_mat_p_blocks)
1147 IF (do_kpoints_cubic_rpa .OR. do_kpoints_from_gamma)
THEN
1154 cut_memory =
SIZE(t_3c_o_compressed, 3)
1156 DEALLOCATE (t_3c_o_ind)
1157 DO i_mem = 1, cut_memory
1158 DO j_size = 1,
SIZE(t_3c_o_compressed, 2)
1159 DO i_size = 1,
SIZE(t_3c_o_compressed, 1)
1164 DEALLOCATE (t_3c_o_compressed)
1166 IF (do_kpoints_from_gamma)
THEN
1168 IF (qs_env%mp2_env%ri_g0w0%do_kpoints_Sigma)
THEN
1170 CALL kpoint_release(qs_env%mp2_env%ri_rpa_im_time%kpoints_Sigma_no_xc)
1174 CALL timestop(handle)
1191 dimen_RI_red, fm_mat_L, fm_mat_Q)
1193 TYPE(
dbcsr_type),
INTENT(IN) :: mat_p_omega, mat_l
1195 REAL(kind=
dp),
INTENT(IN) :: eps_filter_im_time
1196 TYPE(
cp_fm_type),
INTENT(INOUT) :: fm_mat_work
1197 INTEGER,
INTENT(IN) :: dimen_ri, dimen_ri_red
1198 TYPE(
cp_fm_type),
INTENT(IN) :: fm_mat_l, fm_mat_q
1200 CHARACTER(LEN=*),
PARAMETER :: routinen =
'contract_P_omega_with_mat_L'
1204 CALL timeset(routinen, handle)
1208 0.0_dp, mat_work, filter_eps=eps_filter_im_time)
1212 CALL parallel_gemm(
'N',
'N', dimen_ri_red, dimen_ri_red, dimen_ri, 1.0_dp, fm_mat_l, fm_mat_work, &
1219 CALL timestop(handle)
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Handles all functions related to the 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.
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.
subroutine, public cp_cfm_release(matrix)
Releases a full matrix.
subroutine, public cp_fm_to_cfm(msourcer, msourcei, mtarget)
Construct a complex full matrix by taking its real and imaginary parts from two separate real-value f...
subroutine, public cp_cfm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
Creates a new full matrix with the given structure.
subroutine, public cp_cfm_set_all(matrix, alpha, beta)
Set all elements of the full matrix to alpha. Besides, set all diagonal matrix elements to beta (if g...
subroutine, public dbcsr_deallocate_matrix(matrix)
...
subroutine, public dbcsr_multiply(transa, transb, alpha, matrix_a, matrix_b, beta, matrix_c, first_row, last_row, first_column, last_column, first_k, last_k, retain_sparsity, filter_eps, flop)
...
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
subroutine, public dbcsr_filter(matrix, eps)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_release(matrix)
...
DBCSR operations in CP2K.
subroutine, public copy_dbcsr_to_fm(matrix, fm, plan)
Copy a DBCSR matrix to a BLACS matrix.
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
Basic linear algebra operations for full matrices.
subroutine, public cp_fm_transpose(matrix, matrixt)
transposes a matrix matrixt = matrix ^ T
subroutine, public cp_fm_syrk(uplo, trans, k, alpha, matrix_a, ia, ja, beta, matrix_c)
performs a rank-k update of a symmetric matrix_c matrix_c = beta * matrix_c + alpha * matrix_a * tran...
various cholesky decomposition related routines
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,...
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_copy_general(source, destination, para_env)
General copy of a fm matrix to another fm matrix. Uses non-blocking MPI rather than ScaLAPACK.
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
subroutine, public cp_fm_to_fm_submat_general(source, destination, nrows, ncols, s_firstrow, s_firstcol, d_firstrow, d_firstcol, global_context)
General copy of a submatrix of fm matrix to a submatrix of another fm matrix. The two matrices can ha...
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, nrow, ncol, set_zero)
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.
subroutine, public dgemm_counter_start(dgemm_counter)
start timer of the counter
subroutine, public dgemm_counter_stop(dgemm_counter, size1, size2, size3)
stop timer of the counter and provide matrix sizes
Types and set/get functions for HFX.
subroutine, public dealloc_containers(data, memory_usage)
...
Defines the basic variable types.
integer, parameter, public dp
Types and basic routines needed for a kpoint calculation.
subroutine, public kpoint_release(kpoint)
Release a kpoint environment, deallocate all data.
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, inversion_symmetry_only, symmetry_backend, symmetry_reduction_method, gamma_centered, lattice_fft)
Retrieve information from a kpoint environment.
Definition of mathematical constants and functions.
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 MP2 energy with laplace approach.
subroutine, public calc_fm_mat_s_laplace(fm_mat_s, homo, virtual, eigenval, dajquad)
...
basic linear algebra operations for full matrixes
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, mimic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_vhxc, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
Routines treating GW and RPA calculations with kpoints.
subroutine, public compute_wkp_w(qs_env, wkp_w, wkp_v, kpoints, h_inv, periodic)
...
Utility functions for RPA calculations.
subroutine, public alloc_im_time(qs_env, para_env, dimen_ri, dimen_ri_red, num_integ_points, nspins, fm_mat_q, cfm_mo_coeff, 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, 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)
...
subroutine, public compute_erpa_by_freq_int(dimen_ri, trace_qomega, fm_mat_q, para_env_rpa, erpa, wjquad)
...
subroutine, public q_trace_and_add_unit_matrix(dimen_ri, trace_qomega, fm_mat_q)
...
subroutine, public calc_mat_q(fm_mat_s, do_ri_sos_laplace_mp2, first_cycle, virtual, eigenval, 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, count_ev_sc_gw)
...
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 dealloc_im_time(cfm_mo_coeff, 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, 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 remove_scaling_factor_rpa(fm_mat_s, virtual, eigenval_last, homo, omega_old)
...
subroutine, public calc_fm_mat_s_rpa(fm_mat_s, first_cycle, virtual, eigenval, homo, omega, omega_old)
...
Type defining parameters related to the simulation cell.
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