comp-chem/cp2k
1
0
Fork 0
mirror of https://github.com/cp2k/cp2k.git synced 2026-07-29 06:35:28 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Regularization of minimax Fourier transform weights

Browse source
This commit is contained in:
Jan Wilhelm 2023-06-14 22:37:58 +02:00 • committed by GitHub
parent fa9b037f40
commit 186428f4e8
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
7 changed files with 142 additions and 14 deletions
Download patch file Download diff file Expand all files Collapse all files

9
src/input_cp2k_mp2.F
View file

@ -922,6 +922,15 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="REGULARIZATION_MINIMAX", &
description="Tikhonov regularization for computing weights of the Fourier transform "// &
"from imaginary time to imaginary frequency and vice versa. Needed for large minimax "// &
"grids with 20 or more points and a small range.", &
usage="REGULARIZATION_MINIMAX 1.0E-6", &
default_r_val=0.0_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
! here we generate a subsection for the periodic GW correction
CALL create_periodic_gw_correction_section(subsection)
CALL section_add_subsection(section, subsection)

48
src/mp2_grids.F
View file

@ -65,11 +65,12 @@ CONTAINS
!> \param weights_cos_tf_t_to_w ...
!> \param weights_cos_tf_w_to_t ...
!> \param weights_sin_tf_t_to_w ...
!> \param regularization ...
! **************************************************************************************************
SUBROUTINE 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)
weights_cos_tf_w_to_t, weights_sin_tf_t_to_w, regularization)
TYPE(mp_para_env_type), INTENT(IN) :: para_env
INTEGER, INTENT(IN) :: unit_nr
@ -89,6 +90,7 @@ CONTAINS
INTENT(OUT) :: weights_cos_tf_t_to_w, &
weights_cos_tf_w_to_t, &
weights_sin_tf_t_to_w
REAL(KIND=dp), INTENT(IN), OPTIONAL :: regularization
CHARACTER(LEN=*), PARAMETER :: routineN = 'get_minimax_grid'
INTEGER, PARAMETER :: num_points_per_magnitude = 200
@ -96,7 +98,7 @@ CONTAINS
INTEGER :: handle, ierr, ispin, jquad, nspins
LOGICAL :: my_do_kpoints, my_open_shell
REAL(KIND=dp) :: Emax, Emin, max_error_min, my_E_Range, &
scaling
my_regularization, scaling
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: x_tw
TYPE(section_vals_type), POINTER :: input
@ -112,6 +114,11 @@ CONTAINS
my_do_kpoints = do_kpoints_cubic_RPA
END IF
my_regularization = 0.0_dp
IF (PRESENT(regularization)) THEN
my_regularization = regularization
END IF
IF (my_do_kpoints) THEN
CALL gap_and_max_eig_diff_kpoints(qs_env, para_env, Emin, Emax, e_fermi)
my_E_Range = Emax/Emin
@ -247,14 +254,16 @@ CONTAINS
weights_cos_tf_t_to_w = 0.0_dp
CALL get_l_sq_wghts_cos_tf_t_to_w(num_integ_points, tau_tj, weights_cos_tf_t_to_w, tj, &
Emin, Emax, max_error_min, num_points_per_magnitude)
Emin, Emax, max_error_min, num_points_per_magnitude, &
my_regularization)
! get the weights for the cosine transform W^c(iw) -> W^c(it)
ALLOCATE (weights_cos_tf_w_to_t(num_integ_points, num_integ_points))
weights_cos_tf_w_to_t = 0.0_dp
CALL get_l_sq_wghts_cos_tf_w_to_t(num_integ_points, tau_tj, weights_cos_tf_w_to_t, tj, &
Emin, Emax, max_error_min, num_points_per_magnitude)
Emin, Emax, max_error_min, num_points_per_magnitude, &
my_regularization)
IF (do_gw_im_time) THEN
@ -263,7 +272,8 @@ CONTAINS
weights_sin_tf_t_to_w = 0.0_dp
CALL get_l_sq_wghts_sin_tf_t_to_w(num_integ_points, tau_tj, weights_sin_tf_t_to_w, tj, &
Emin, Emax, max_error_min, num_points_per_magnitude)
Emin, Emax, max_error_min, num_points_per_magnitude, &
my_regularization)
IF (unit_nr > 0) THEN
WRITE (UNIT=unit_nr, FMT="(T3,A,T66,ES15.2)") &
@ -721,9 +731,11 @@ CONTAINS
!> \param E_max ...
!> \param max_error ...
!> \param num_points_per_magnitude ...
!> \param regularization ...
! **************************************************************************************************
SUBROUTINE get_l_sq_wghts_cos_tf_t_to_w(num_integ_points, tau_tj, weights_cos_tf_t_to_w, omega_tj, &
E_min, E_max, max_error, num_points_per_magnitude)
E_min, E_max, max_error, num_points_per_magnitude, &
regularization)
INTEGER, INTENT(IN) :: num_integ_points
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), &
@ -735,6 +747,7 @@ CONTAINS
REAL(KIND=dp), INTENT(IN) :: E_min, E_max
REAL(KIND=dp), INTENT(INOUT) :: max_error
INTEGER, INTENT(IN) :: num_points_per_magnitude
REAL(KIND=dp), INTENT(IN) :: regularization
CHARACTER(LEN=*), PARAMETER :: routineN = 'get_l_sq_wghts_cos_tf_t_to_w'
@ -817,7 +830,9 @@ CONTAINS
! 1) V*Sigma
DO jjj = 1, num_integ_points
DO iii = 1, num_integ_points
mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)/sing_values(jjj)
! mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)/sing_values(jjj)
mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)*sing_values(jjj) &
/(regularization**2 + sing_values(jjj)**2)
END DO
END DO
@ -853,9 +868,10 @@ CONTAINS
!> \param E_max ...
!> \param max_error ...
!> \param num_points_per_magnitude ...
!> \param regularization ...
! **************************************************************************************************
SUBROUTINE get_l_sq_wghts_sin_tf_t_to_w(num_integ_points, tau_tj, weights_sin_tf_t_to_w, omega_tj, &
E_min, E_max, max_error, num_points_per_magnitude)
E_min, E_max, max_error, num_points_per_magnitude, regularization)
INTEGER, INTENT(IN) :: num_integ_points
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), &
@ -867,6 +883,7 @@ CONTAINS
REAL(KIND=dp), INTENT(IN) :: E_min, E_max
REAL(KIND=dp), INTENT(OUT) :: max_error
INTEGER, INTENT(IN) :: num_points_per_magnitude
REAL(KIND=dp), INTENT(IN) :: regularization
CHARACTER(LEN=*), PARAMETER :: routineN = 'get_l_sq_wghts_sin_tf_t_to_w'
@ -954,7 +971,9 @@ CONTAINS
! 1) V*Sigma
DO jjj = 1, num_integ_points
DO iii = 1, num_integ_points
mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)/sing_values(jjj)
! mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)/sing_values(jjj)
mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)*sing_values(jjj) &
/(regularization**2 + sing_values(jjj)**2)
END DO
END DO
@ -1187,7 +1206,8 @@ CONTAINS
END DO
CALL get_l_sq_wghts_cos_tf_t_to_w(num_integ_points, tau_tj, &
weights_cos_tf_t_to_w, tj, 1.0_dp, Rc, max_error, 200)
weights_cos_tf_t_to_w, tj, &
1.0_dp, Rc, max_error, 200, 0.0_dp)
IF (iw > 0) THEN
WRITE (iw, '(T2, I3, F12.1, ES12.3)') num_integ_points, Rc, max_error
@ -1211,9 +1231,10 @@ CONTAINS
!> \param E_max ...
!> \param max_error ...
!> \param num_points_per_magnitude ...
!> \param regularization ...
! **************************************************************************************************
SUBROUTINE get_l_sq_wghts_cos_tf_w_to_t(num_integ_points, tau_tj, weights_cos_tf_w_to_t, omega_tj, &
E_min, E_max, max_error, num_points_per_magnitude)
E_min, E_max, max_error, num_points_per_magnitude, regularization)
INTEGER, INTENT(IN) :: num_integ_points
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), &
@ -1225,6 +1246,7 @@ CONTAINS
REAL(KIND=dp), INTENT(IN) :: E_min, E_max
REAL(KIND=dp), INTENT(INOUT) :: max_error
INTEGER, INTENT(IN) :: num_points_per_magnitude
REAL(KIND=dp), INTENT(IN) :: regularization
CHARACTER(LEN=*), PARAMETER :: routineN = 'get_l_sq_wghts_cos_tf_w_to_t'
@ -1314,7 +1336,9 @@ CONTAINS
! 1) V*Sigma
DO jjj = 1, num_integ_points
DO iii = 1, num_integ_points
mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)/sing_values(jjj)
! mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)/sing_values(jjj)
mat_SinvVSinvSigma(iii, jjj) = mat_SinvVSinvT(jjj, iii)*sing_values(jjj) &
/(regularization**2 + sing_values(jjj)**2)
END DO
END DO

2
src/mp2_setup.F
View file

@ -161,6 +161,8 @@ CONTAINS
l_val=mp2_env%ri_g0w0%do_gamma_only_sigma)
CALL section_vals_val_get(mp2_section, "RI_RPA%GW%UPDATE_XC_ENERGY", &
l_val=mp2_env%ri_g0w0%update_xc_energy)
CALL section_vals_val_get(mp2_section, "RI_RPA%GW%REGULARIZATION_MINIMAX", &
r_val=mp2_env%ri_g0w0%regularization_minimax)
CALL section_vals_val_get(mp2_section, "RI_RPA%GW%BSE", &
l_val=mp2_env%ri_g0w0%do_bse)

3
src/mp2_types.F
View file

@ -210,7 +210,8 @@ MODULE mp2_types
TYPE(dbcsr_p_type), DIMENSION(:), ALLOCATABLE :: &
matrix_sigma_x_minus_vxc, matrix_ks
REAL(KIND=dp) :: broadening_print_loc_bandgap, energy_window_print_loc_bandgap, &
ldos_thresh_print_loc_bandgap, energy_spacing_print_loc_bandgap
ldos_thresh_print_loc_bandgap, energy_spacing_print_loc_bandgap, &
regularization_minimax
INTEGER, DIMENSION(:), POINTER :: stride_loc_bandgap
END TYPE

3
src/rpa_main.F
View file

@ -1267,7 +1267,8 @@ CONTAINS
do_ri_sos_laplace_mp2, do_print, &
tau_tj, tau_wj, qs_env, do_gw_im_time, &
do_kpoints_cubic_RPA, e_fermi(1), tj, wj, &
weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, weights_sin_tf_t_to_w)
weights_cos_tf_t_to_w, weights_cos_tf_w_to_t, weights_sin_tf_t_to_w, &
qs_env%mp2_env%ri_g0w0%regularization_minimax)
!For sos_laplace_mp2 and low-scaling RPA, potentially need to store/retrieve the initial weights
IF (qs_env%mp2_env%ri_rpa_im_time%keep_quad) THEN

90
tests/QS/regtest-gw-cubic/G0W0_H2O_PBE0_30_pts_minimax_regularization.inp Normal file
View file

@ -0,0 +1,90 @@
&GLOBAL
PROJECT G0W0_H2O_PBE0
PRINT_LEVEL MEDIUM
RUN_TYPE ENERGY
&TIMINGS
THRESHOLD 0.01
&END
&END GLOBAL
&FORCE_EVAL
METHOD Quickstep
&DFT
BASIS_SET_FILE_NAME HFX_BASIS
SORT_BASIS EXP
POTENTIAL_FILE_NAME GTH_POTENTIALS
&MGRID
CUTOFF 100
REL_CUTOFF 20
&END MGRID
&POISSON
PERIODIC NONE
POISSON_SOLVER WAVELET
&END POISSON
&QS
METHOD GPW
EPS_DEFAULT 1.0E-15
EPS_PGF_ORB 1.0E-30
&END QS
&SCF
SCF_GUESS ATOMIC
EPS_SCF 1.0E-7
MAX_SCF 100
&PRINT
&RESTART OFF
&END
&END
&END SCF
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
&WF_CORRELATION
&LOW_SCALING
EPS_FILTER 0.0
MEMORY_CUT 1
&END
&RI_RPA
&HF
FRACTION 1.0000000
&SCREENING
EPS_SCHWARZ 1.0E-8
SCREEN_ON_INITIAL_P FALSE
&END SCREENING
&END HF
! this test specifically tests a large minimax grid
RPA_NUM_QUAD_POINTS 30
&GW
CORR_MOS_OCC 10
CORR_MOS_VIRT 10
RI_SIGMA_X FALSE
! test with very large regularization parameter
REGULARIZATION_MINIMAX 0.1
&END GW
&END RI_RPA
&END
&END XC
&END DFT
&SUBSYS
&CELL
ABC [angstrom] 8.000 8.000 8.000
PERIODIC NONE
&END CELL
&KIND H
! this test specifically tests a large basis
BASIS_SET cc-TZV2P-GTH
BASIS_SET RI_AUX RI_DZVP-GTH
POTENTIAL GTH-PBE-q1
&END KIND
&KIND O
! this test specifically tests a large basis
BASIS_SET cc-TZV2P-GTH
BASIS_SET RI_AUX RI_DZVP-GTH
POTENTIAL GTH-PBE-q6
&END KIND
&TOPOLOGY
COORD_FILE_NAME H2O_gas.xyz
COORD_FILE_FORMAT xyz
&CENTER_COORDINATES
&END
&END TOPOLOGY
&END SUBSYS
&END FORCE_EVAL

1
tests/QS/regtest-gw-cubic/TEST_FILES
View file

@ -3,6 +3,7 @@ evGW_H2O_PBE_default_values.inp 78 1e-04
G0W0_H2O_PBE_GAPW.inp 78 1e-04 16.699
G0W0_H2O_PBE0.inp 78 1e-04 16.717
G0W0_H2O_PBE0_30_pts.inp 78 1e-04 15.607
G0W0_H2O_PBE0_30_pts_minimax_regularization.inp 78 1e-04 15.965
scGW0_and_evGW_H2O_PBE0_trunc_minimax.inp 11 1e-08 -17.108489005370274
scGW0_H2O_PBE0_trunc_minimax_RI_HFX.inp 11 1e-08 -13.191093644224157
G0W0_H2O_PBE_periodic.inp 78 1e-04 16.475