Refactor EWALD input for TB Methods, new defaults

This commit is contained in:
Juerg Hutter 2019-03-04 16:47:55 +01:00
parent 7a86519636
commit 8d0aca42e9
53 changed files with 166 additions and 51 deletions

View file

@ -32,6 +32,7 @@ MODULE ewald_environment_types
section_vals_type,&
section_vals_val_get
USE kinds, ONLY: dp
USE mathconstants, ONLY: twopi
USE pw_poisson_types, ONLY: do_ewald_ewald,&
do_ewald_none,&
do_ewald_pme,&
@ -93,7 +94,8 @@ MODULE ewald_environment_types
ewald_env_create, &
ewald_env_retain, &
ewald_env_release, &
read_ewald_section
read_ewald_section, &
read_ewald_section_tb
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'ewald_environment_types'
INTEGER, PRIVATE, SAVE :: last_ewald_env_id_nr = 0
@ -435,6 +437,105 @@ CONTAINS
END SUBROUTINE read_ewald_section
! **************************************************************************************************
!> \brief Purpose: read the EWALD section for TB methods
!> \param ewald_env the pointer to the ewald_env
!> \param ewald_section ...
!> \param hmat ...
!> \author JGH
! **************************************************************************************************
SUBROUTINE read_ewald_section_tb(ewald_env, ewald_section, hmat)
TYPE(ewald_environment_type), POINTER :: ewald_env
TYPE(section_vals_type), POINTER :: ewald_section
REAL(KIND=dp), DIMENSION(3, 3), INTENT(IN) :: hmat
CHARACTER(len=*), PARAMETER :: routineN = 'read_ewald_section_tb', &
routineP = moduleN//':'//routineN
INTEGER :: i, iw
INTEGER, DIMENSION(:), POINTER :: gmax_read
LOGICAL :: explicit
REAL(KIND=dp) :: alat, cutoff, dummy
TYPE(cp_logger_type), POINTER :: logger
logger => cp_get_default_logger()
ewald_env%do_multipoles = .FALSE.
ewald_env%do_ipol = 0
ewald_env%eps_pol = 1.e-12_dp
ewald_env%max_multipole = 0
ewald_env%max_ipol_iter = 0
ewald_env%epsilon = 1.e-12_dp
ewald_env%ns_max = HUGE(0)
CALL section_vals_val_get(ewald_section, "EWALD_TYPE", explicit=explicit)
IF (explicit) THEN
CALL section_vals_val_get(ewald_section, "EWALD_TYPE", i_val=ewald_env%ewald_type)
IF (ewald_env%ewald_type /= do_ewald_spme) THEN
CPABORT("TB needs EWALD_TYPE SPME")
END IF
ELSE
ewald_env%ewald_type = do_ewald_spme
ENDIF
CALL section_vals_val_get(ewald_section, "ALPHA", explicit=explicit)
IF (explicit) THEN
CALL section_vals_val_get(ewald_section, "ALPHA", r_val=ewald_env%alpha)
ELSE
ewald_env%alpha = 1.0_dp
ENDIF
CALL section_vals_val_get(ewald_section, "EWALD_ACCURACY", r_val=ewald_env%precs)
CALL section_vals_val_get(ewald_section, "O_SPLINE", i_val=ewald_env%o_spline)
CALL section_vals_val_get(ewald_section, "RCUT", explicit=explicit)
IF (explicit) THEN
CALL section_vals_val_get(ewald_section, "RCUT", r_val=ewald_env%rcut)
ELSE
ewald_env%rcut = find_ewald_optimal_value(ewald_env%precs)/ewald_env%alpha
ENDIF
CALL section_vals_val_get(ewald_section, "GMAX", explicit=explicit)
IF (explicit) THEN
CALL section_vals_val_get(ewald_section, "GMAX", i_vals=gmax_read)
SELECT CASE (SIZE (gmax_read, 1))
CASE (1)
ewald_env%gmax = gmax_read(1)
CASE (3)
ewald_env%gmax = gmax_read
CASE DEFAULT
CPABORT("")
END SELECT
ELSE
! set GMAX using ECUT=alpha*45 Ry
cutoff = 45._dp*ewald_env%alpha
DO i = 1, 3
alat = SUM(hmat(:, i)**2)
CPASSERT(alat /= 0._dp)
ewald_env%gmax(i) = 2*FLOOR(SQRT(2.0_dp*cutoff*alat)/twopi)+1
ENDDO
ENDIF
iw = cp_print_key_unit_nr(logger, ewald_section, "PRINT%PROGRAM_RUN_INFO", &
extension=".log")
IF (iw > 0) THEN
WRITE (iw, '(/,T2,"EWALD| ",A,T67,A14 )') 'Summation is done by:', ADJUSTR("SPME")
dummy = cp_unit_from_cp2k(ewald_env%alpha, "angstrom^-1")
WRITE (iw, '( T2,"EWALD| ",A,A18,A,T71,F10.4 )') &
'Alpha parameter [', 'ANGSTROM^-1', ']', dummy
dummy = cp_unit_from_cp2k(ewald_env%rcut, "angstrom")
WRITE (iw, '( T2,"EWALD| ",A,A18,A,T71,F10.4 )') &
'Real Space Cutoff [', 'ANGSTROM', ']', dummy
WRITE (iw, '( T2,"EWALD| ",A,T51,3I10 )') &
'G-space max. Miller index', ewald_env%gmax
WRITE (iw, '( T2,"EWALD| ",A,T71,I10 )') &
'Spline interpolation order ', ewald_env%o_spline
END IF
CALL cp_print_key_finished_output(iw, logger, ewald_section, &
"PRINT%PROGRAM_RUN_INFO")
END SUBROUTINE read_ewald_section_tb
! **************************************************************************************************
!> \brief triggers (by bisection) the optimal value for EWALD parameter x
!> EXP(-x^2)/x^2 = EWALD_ACCURACY

View file

@ -131,6 +131,9 @@ CONTAINS
ALLOCATE (pw_pools(1))
pw_pools(1)%pool => pw_big_pool
CALL pw_poisson_read_parameters(poisson_section, poisson_params)
poisson_params%ewald_type = ewald_type
poisson_params%ewald_o_spline = o_spline
poisson_params%ewald_alpha = alpha
CALL pw_poisson_set(poisson_env, cell_hmat=cell_hmat, parameters=poisson_params, &
use_level=1, pw_pools=pw_pools)
DEALLOCATE (pw_pools)

View file

@ -341,6 +341,9 @@ CONTAINS
ALLOCATE (pw_pools(1))
pw_pools(1)%pool => ewald_pw%pw_big_pool
CALL pw_poisson_read_parameters(poisson_section, poisson_params)
poisson_params%ewald_type = ewald_type
poisson_params%ewald_o_spline = o_spline
poisson_params%ewald_alpha = alpha
CALL pw_poisson_set(ewald_pw%poisson_env, cell_hmat=cell%hmat, parameters=poisson_params, &
use_level=1, pw_pools=pw_pools)
DEALLOCATE (pw_pools)

View file

@ -62,10 +62,9 @@ CONTAINS
routineP = moduleN//':'//routineN
INTEGER :: periodic
TYPE(section_vals_type), POINTER :: ewald_section, mt_section, &
wavelet_section
TYPE(section_vals_type), POINTER :: mt_section, wavelet_section
NULLIFY (ewald_section, mt_section, wavelet_section)
NULLIFY (mt_section, wavelet_section)
CALL section_vals_val_get(poisson_section, "POISSON_SOLVER", i_val=params%solver)
@ -83,14 +82,8 @@ CONTAINS
CPABORT("")
END SELECT
! parsing EWALD subsection
! Set Ewald default to NONE
params%ewald_type = do_ewald_none
ewald_section => section_vals_get_subs_vals(poisson_section, "EWALD", can_return_null=.TRUE.)
IF (ASSOCIATED(ewald_section)) THEN
CALL section_vals_val_get(ewald_section, "EWALD_TYPE", i_val=params%ewald_type)
CALL section_vals_val_get(ewald_section, "o_spline", i_val=params%ewald_o_spline)
CALL section_vals_val_get(ewald_section, "alpha", r_val=params%ewald_alpha)
ENDIF
! parsing MT subsection
mt_section => section_vals_get_subs_vals(poisson_section, "MT")

View file

@ -69,7 +69,8 @@ MODULE qs_environment
ewald_env_release,&
ewald_env_set,&
ewald_environment_type,&
read_ewald_section
read_ewald_section,&
read_ewald_section_tb
USE ewald_pw_methods, ONLY: ewald_pw_grid_update
USE ewald_pw_types, ONLY: ewald_pw_create,&
ewald_pw_release,&
@ -629,7 +630,7 @@ CONTAINS
ewald_section => section_vals_get_subs_vals(poisson_section, "EWALD")
print_section => section_vals_get_subs_vals(qs_env%input, "PRINT%GRID_INFORMATION")
CALL get_qs_kind_set(qs_kind_set, basis_rcut=ewald_rcut)
CALL read_ewald_section(ewald_env, ewald_section)
CALL read_ewald_section_tb(ewald_env, ewald_section, cell_ref%hmat)
CALL ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, print_section=print_section)
CALL set_qs_env(qs_env, ewald_env=ewald_env, ewald_pw=ewald_pw)
CALL ewald_env_release(ewald_env)
@ -679,7 +680,7 @@ CONTAINS
CALL ewald_env_set(ewald_env, poisson_section=poisson_section)
ewald_section => section_vals_get_subs_vals(poisson_section, "EWALD")
print_section => section_vals_get_subs_vals(qs_env%input, "PRINT%GRID_INFORMATION")
CALL read_ewald_section(ewald_env, ewald_section)
CALL read_ewald_section_tb(ewald_env, ewald_section, cell_ref%hmat)
CALL ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, print_section=print_section)
CALL set_qs_env(qs_env, ewald_env=ewald_env, ewald_pw=ewald_pw)
CALL ewald_env_release(ewald_env)

View file

@ -303,7 +303,7 @@ CONTAINS
LOGICAL, ALLOCATABLE, DIMENSION(:) :: all_present, aux_fit_present, aux_present, &
core_present, default_present, oce_present, orb_present, ppl_present, ppnl_present, &
ri_present, xb1_atom, xb2_atom
REAL(dp) :: almo_rcov, almo_rvdw, alpha, roperator, &
REAL(dp) :: almo_rcov, almo_rvdw, rcut, roperator, &
subcells
REAL(dp), ALLOCATABLE, DIMENSION(:) :: all_pot_rad, aux_fit_radius, c_radius, calpha, &
core_radius, oce_radius, orb_radius, ppl_radius, ppnl_radius, ri_radius, zeff
@ -744,8 +744,8 @@ CONTAINS
! Build the neighbor lists for the DFTB Ewald methods
IF (dft_control%qs_control%dftb_control%do_ewald) THEN
CALL get_qs_env(qs_env=qs_env, ewald_env=ewald_env)
CALL ewald_env_get(ewald_env, alpha=alpha)
c_radius = 0.5_dp*SQRT(-LOG(3.5_dp*alpha**3*1.e-12_dp))/alpha
CALL ewald_env_get(ewald_env, rcut=rcut)
c_radius = rcut
CALL pair_radius_setup(orb_present, orb_present, c_radius, c_radius, pair_radius)
CALL build_neighbor_lists(sab_tbe, particle_set, atom2d, cell, pair_radius, mic=mic, &
subcells=subcells, nlname="sab_tbe")
@ -769,11 +769,11 @@ CONTAINS
END IF
IF (xtb) THEN
! Build the neighbor lists for the DFTB Ewald methods
! Build the neighbor lists for the xTB Ewald method
IF (dft_control%qs_control%xtb_control%do_ewald) THEN
CALL get_qs_env(qs_env=qs_env, ewald_env=ewald_env)
CALL ewald_env_get(ewald_env, alpha=alpha)
c_radius = 0.5_dp*SQRT(-LOG(3.5_dp*alpha**3*1.e-12_dp))/alpha
CALL ewald_env_get(ewald_env, rcut=rcut)
c_radius = rcut
CALL pair_radius_setup(orb_present, orb_present, c_radius, c_radius, pair_radius)
CALL build_neighbor_lists(sab_tbe, particle_set, atom2d, cell, pair_radius, mic=mic, &
subcells=subcells, nlname="sab_tbe")

View file

@ -24,6 +24,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 100
&END EWALD
POISSON_SOLVER ANALYTIC

View file

@ -32,6 +32,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -20,6 +20,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -20,6 +20,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -32,6 +32,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
POISSON_SOLVER ANALYTIC

View file

@ -32,6 +32,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
POISSON_SOLVER ANALYTIC

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 16
O_SPLINE 4
&END EWALD

View file

@ -33,6 +33,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 16
O_SPLINE 4
&END EWALD

View file

@ -31,6 +31,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -43,6 +43,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -41,6 +41,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -37,6 +37,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 45
O_SPLINE 8
&END EWALD

View file

@ -37,6 +37,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 45
O_SPLINE 8
&END EWALD

View file

@ -33,6 +33,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 45
O_SPLINE 8
&END EWALD

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
&END EWALD
&END POISSON

View file

@ -31,6 +31,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -33,6 +33,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -33,6 +33,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -37,6 +37,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -38,6 +38,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -38,6 +38,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -30,6 +30,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX 25
O_SPLINE 5
&END EWALD

View file

@ -46,6 +46,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -46,6 +46,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -47,6 +47,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -47,6 +47,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -48,6 +48,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -48,6 +48,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -48,6 +48,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -48,6 +48,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -48,6 +48,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -49,6 +49,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -52,6 +52,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -49,6 +49,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -49,6 +49,7 @@
&POISSON
&EWALD
EWALD_TYPE SPME
ALPHA 0.35
GMAX ${GMAXVAL} ${GMAXVAL} ${GMAXVAL}
O_SPLINE 5
&END EWALD

View file

@ -9,9 +9,9 @@ ch2o_lsd.inp 1 1.0E-12 -7.19650182
ch2o_smear.inp 1 1.0E-12 -7.19650182424155
tmol.inp 1 1.0E-12 -41.90853885813245
h2.inp 1 1.0E-12 -1.03458111733093
h2o-md.inp 1 1.0E-12 -185.15872877692459
h2o-md.inp 1 1.0E-12 -185.15881067703003
h2o_str.inp 1 1.0E-12 -5.76518308171065
h2o-atprop.inp 1 1.0E-12 -185.18224318229312
h2o-atprop.inp 1 1.0E-12 -185.18225089195914
h2o-atprop0.inp 1 1.0E-12 -187.45499307089042
si_geo.inp 1 1.0E-12 -14.55675190171168
si_kp.inp 1 1.0E-12 -14.75429924999847

View file

@ -0,0 +1 @@
##

View file

@ -25,8 +25,6 @@
&END SCF
&POISSON
&EWALD
EWALD_TYPE SPME
GMAX 45
O_SPLINE 8
&END EWALD
&END POISSON

View file

@ -18,19 +18,8 @@
&OUTER_SCF
MAX_SCF 10
&END
#&MIXING
# METHOD DIRECT_P_MIXING
# ALPHA 0.75
#&END
MAX_SCF 10
&END SCF
&POISSON
&EWALD
EWALD_TYPE SPME
GMAX 45
O_SPLINE 8
&END EWALD
&END POISSON
&END DFT
&SUBSYS
&CELL

View file

@ -4,8 +4,8 @@
# 1 compares the last total energy in the file
# for details see cp2k/tools/do_regtest
NdF3.inp 1 1.0E-12 -16.30904020314352
h2o_rtp.inp 1 1.0E-12 -5.76545575142246
h2o_emd.inp 1 1.0E-12 -5.76503346048092
h2o_rtp.inp 1 1.0E-12 -5.76545803137195
h2o_emd.inp 1 1.0E-12 -5.76503571604747
si8_wan.inp 1 1.0E-12 -14.36402454198545
si_kp.inp 1 1.0E-12 -14.75264560894033
tmol.inp 1 1.0E-12 -41.90853885813247

View file

@ -14,13 +14,6 @@
&END
MAX_SCF 20
&END SCF
&POISSON
&EWALD
EWALD_TYPE SPME
GMAX 25
O_SPLINE 5
&END EWALD
&END POISSON
&PRINT
&MULLIKEN
&END MULLIKEN

View file

@ -14,13 +14,6 @@
&END
MAX_SCF 20
&END SCF
&POISSON
&EWALD
EWALD_TYPE SPME
GMAX 25
O_SPLINE 5
&END EWALD
&END POISSON
&PRINT
&MULLIKEN
&END MULLIKEN