39 #include "./base/base_uses.f90"
45 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'bse_main'
69 fm_mat_Q_static_bse, fm_mat_Q_static_bse_gemm, &
70 Eigenval, homo, virtual, dimen_RI, dimen_RI_red, &
71 gd_array, color_sub, mp2_env, unit_nr)
73 TYPE(cp_fm_type),
INTENT(IN) :: fm_mat_s_ia_bse, fm_mat_s_ij_bse, &
75 TYPE(cp_fm_type),
INTENT(INOUT) :: fm_mat_q_static_bse, &
76 fm_mat_q_static_bse_gemm
77 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :), &
78 INTENT(IN) :: eigenval
79 INTEGER,
DIMENSION(:),
INTENT(IN) :: homo, virtual
80 INTEGER,
INTENT(IN) :: dimen_ri, dimen_ri_red
81 TYPE(group_dist_d1_type),
INTENT(IN) :: gd_array
82 INTEGER,
INTENT(IN) :: color_sub
83 TYPE(mp2_type) :: mp2_env
84 INTEGER,
INTENT(IN) :: unit_nr
86 CHARACTER(LEN=*),
PARAMETER :: routinen =
'start_bse_calculation'
88 INTEGER :: handle, homo_red, virtual_red
89 LOGICAL :: my_do_abba, my_do_fulldiag, &
90 my_do_iterat_diag, my_do_tda
91 REAL(kind=
dp) :: diag_runtime_est
92 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: eigenval_reduced
93 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: b_abq_bse_local, b_bar_iaq_bse_local, &
94 b_bar_ijq_bse_local, b_iaq_bse_local
95 TYPE(cp_fm_type) :: fm_a_bse, fm_b_bse, fm_c_bse, fm_inv_sqrt_a_minus_b, fm_mat_s_ab_trunc, &
96 fm_mat_s_bar_ia_bse, fm_mat_s_bar_ij_bse, fm_mat_s_ia_trunc, fm_mat_s_ij_trunc, &
98 TYPE(mp_para_env_type),
POINTER :: para_env
100 CALL timeset(routinen, handle)
102 para_env => fm_mat_s_ia_bse%matrix_struct%para_env
104 my_do_fulldiag = .false.
105 my_do_iterat_diag = .false.
109 SELECT CASE (mp2_env%ri_g0w0%bse_diag_method)
111 my_do_iterat_diag = .true.
113 cpabort(
"Iterative BSE not yet implemented")
115 my_do_fulldiag = .true.
118 SELECT CASE (mp2_env%ri_g0w0%bse_approx)
130 CALL fm_mat_s_ia_bse%matrix_struct%para_env%sync()
134 fm_mat_s_ia_trunc, fm_mat_s_ij_trunc, fm_mat_s_ab_trunc, &
135 eigenval(:, 1, 1), eigenval_reduced, &
136 homo(1), virtual(1), dimen_ri, unit_nr, &
137 homo_red, virtual_red, &
142 CALL mult_b_with_w(fm_mat_s_ij_trunc, fm_mat_s_ia_trunc, fm_mat_s_bar_ia_bse, &
143 fm_mat_s_bar_ij_bse, fm_mat_q_static_bse_gemm, &
144 dimen_ri_red, homo_red, virtual_red)
146 IF (my_do_iterat_diag)
THEN
148 fm_mat_s_bar_ia_bse, fm_mat_s_bar_ij_bse, &
149 b_bar_ijq_bse_local, b_abq_bse_local, b_bar_iaq_bse_local, &
150 b_iaq_bse_local, dimen_ri_red, homo_red, virtual_red, &
151 gd_array, color_sub, para_env)
154 IF (my_do_fulldiag)
THEN
157 para_env, diag_runtime_est)
164 CALL create_a(fm_mat_s_ia_trunc, fm_mat_s_bar_ij_bse, fm_mat_s_ab_trunc, &
165 fm_a_bse, eigenval_reduced, unit_nr, &
166 homo_red, virtual_red, dimen_ri, mp2_env, &
174 CALL create_b(fm_mat_s_ia_trunc, fm_mat_s_bar_ia_bse, fm_b_bse, &
175 homo_red, virtual_red, dimen_ri, unit_nr, mp2_env)
181 fm_sqrt_a_minus_b, fm_inv_sqrt_a_minus_b, &
182 homo_red, virtual_red, unit_nr, mp2_env, diag_runtime_est)
184 CALL cp_fm_release(fm_b_bse)
187 CALL diagonalize_a(fm_a_bse, homo_red, virtual_red, homo(1), &
188 unit_nr, diag_runtime_est, mp2_env)
191 CALL cp_fm_release(fm_a_bse)
196 CALL diagonalize_c(fm_c_bse, homo_red, virtual_red, homo(1), &
197 fm_sqrt_a_minus_b, fm_inv_sqrt_a_minus_b, &
198 unit_nr, diag_runtime_est, mp2_env)
201 CALL cp_fm_release(fm_c_bse)
205 fm_mat_s_ia_trunc, fm_mat_s_ij_trunc, fm_mat_s_ab_trunc, &
206 fm_mat_q_static_bse, fm_mat_q_static_bse_gemm)
207 DEALLOCATE (eigenval_reduced)
208 IF (my_do_iterat_diag)
THEN
211 b_iaq_bse_local, homo(1), virtual(1), mp2_env%ri_g0w0%bse_spin_config, unit_nr, &
212 eigenval(:, 1, 1), para_env, mp2_env)
214 DEALLOCATE (b_bar_ijq_bse_local, b_abq_bse_local, b_bar_iaq_bse_local, b_iaq_bse_local)
217 IF (unit_nr > 0)
THEN
218 WRITE (unit_nr,
'(T2,A4,T7,A53)')
'BSE|',
'The BSE was successfully calculated. Have a nice day!'
221 CALL timestop(handle)
Routines for the full diagonalization of GW + Bethe-Salpeter for computing electronic excitations.
subroutine, public create_a(fm_mat_S_ia_bse, fm_mat_S_bar_ij_bse, fm_mat_S_ab_bse, fm_A, Eigenval, unit_nr, homo, virtual, dimen_RI, mp2_env, para_env)
Matrix A constructed from GW energies and 3c-B-matrices (cf. subroutine mult_B_with_W) A_ia,...
subroutine, public diagonalize_c(fm_C, homo, virtual, homo_irred, fm_sqrt_A_minus_B, fm_inv_sqrt_A_minus_B, unit_nr, diag_est, mp2_env)
Solving eigenvalue equation C Z^n = (Ω^n)^2 Z^n . Here, the eigenvectors Z^n relate to X^n via Eq....
subroutine, public create_b(fm_mat_S_ia_bse, fm_mat_S_bar_ia_bse, fm_B, homo, virtual, dimen_RI, unit_nr, mp2_env)
Matrix B constructed from 3c-B-matrices (cf. subroutine mult_B_with_W) B_ia,jb = α * v_ia,...
subroutine, public create_hermitian_form_of_abba(fm_A, fm_B, fm_C, fm_sqrt_A_minus_B, fm_inv_sqrt_A_minus_B, homo, virtual, unit_nr, mp2_env, diag_est)
Construct Matrix C=(A-B)^0.5 (A+B) (A-B)^0.5 to solve full BSE matrix as a hermitian problem (cf....
subroutine, public diagonalize_a(fm_A, homo, virtual, homo_irred, unit_nr, diag_est, mp2_env)
Solving hermitian eigenvalue equation A X^n = Ω^n X^n.
Iterative routines for GW + Bethe-Salpeter for computing electronic excitations.
subroutine, public fill_local_3c_arrays(fm_mat_S_ab_bse, fm_mat_S, fm_mat_S_bar_ia_bse, fm_mat_S_bar_ij_bse, B_bar_ijQ_bse_local, B_abQ_bse_local, B_bar_iaQ_bse_local, B_iaQ_bse_local, dimen_RI, homo, virtual, gd_array, color_sub, para_env)
...
subroutine, public do_subspace_iterations(B_bar_ijQ_bse_local, B_abQ_bse_local, B_bar_iaQ_bse_local, B_iaQ_bse_local, homo, virtual, bse_spin_config, unit_nr, Eigenval, para_env, mp2_env)
...
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.
Auxiliary routines for GW + Bethe-Salpeter for computing electronic excitations.
subroutine, public estimate_bse_resources(homo_red, virtual_red, unit_nr, bse_abba, para_env, diag_runtime_est)
Roughly estimates the needed runtime and memory during the BSE run.
subroutine, public mult_b_with_w(fm_mat_S_ij_bse, fm_mat_S_ia_bse, fm_mat_S_bar_ia_bse, fm_mat_S_bar_ij_bse, fm_mat_Q_static_bse_gemm, dimen_RI, homo, virtual)
Multiplies B-matrix (RI-3c-Integrals) with W (screening) to obtain \bar{B}.
subroutine, public print_bse_start_flag(bse_tda, bse_abba, unit_nr)
...
subroutine, public deallocate_matrices_bse(fm_mat_S_bar_ia_bse, fm_mat_S_bar_ij_bse, fm_mat_S_trunc, fm_mat_S_ij_trunc, fm_mat_S_ab_trunc, fm_mat_Q_static_bse, fm_mat_Q_static_bse_gemm)
...
subroutine, public truncate_bse_matrices(fm_mat_S_ia_bse, fm_mat_S_ij_bse, fm_mat_S_ab_bse, fm_mat_S_trunc, fm_mat_S_ij_trunc, fm_mat_S_ab_trunc, Eigenval, Eigenval_reduced, homo, virtual, dimen_RI, unit_nr, homo_red, virt_red, mp2_env)
Determines indices within the given energy cutoffs and truncates Eigenvalues and matrices.
represent a full matrix distributed on many processors
Types to describe group distributions.
Defines the basic variable types.
integer, parameter, public dp
Interface to the message passing library MPI.
Types needed for MP2 calculations.