diff --git a/assets/parameters.toml b/assets/parameters.toml index eb1a340c..9d191260 100644 --- a/assets/parameters.toml +++ b/assets/parameters.toml @@ -452,4 +452,22 @@ d3.zero = {rs6=2.3933, s8=1.1129, doi="10.1039/C7CP04913G"} d3.zero = {rs6=3.3388, s8=0.0, doi="10.1039/C7CP04913G"} [parameter.pwp] -d3.zero = {rs6=2.1040, s8=0.8747, doi="10.1039/C7CP04913G"} \ No newline at end of file +d3.zero = {rs6=2.1040, s8=0.8747, doi="10.1039/C7CP04913G"} + +[parameter.pw91] +d3.bj = {a1=0.6319, s8=1.9598, a2=4.5718, doi="10.1073/pnas.1516984112"} + +[parameter.drpa75] +d3.bj = {s6=0.3754, a1=0.0, s8=0.0, a2=4.5048, doi="10.1039/C6CP00688D"} + +[parameter.scsdrpa75] +d3.bj = {s6=0.2528, a1=0.0, s8=0.0, a2=4.5050, doi="10.1021/acs.jpca.1c01295"} + +[parameter.optscsdrpa75] +d3.bj = {s6=0.2546, a1=0.0, s8=0.0, a2=4.5050, doi="10.1021/acs.jpca.1c01295"} + +[parameter.dsdpbedrpa75] +d3.bj = {s6=0.3223, a1=0.0, s8=0.0, a2=4.5050, doi="10.1021/acs.jpca.1c01295"} + +[parameter.dsdpbep86drpa75] +d3.bj = {s6=0.3012, a1=0.0, s8=0.0, a2=4.5050, doi="10.1021/acs.jpca.1c01295"} diff --git a/src/dftd3/param.f90 b/src/dftd3/param.f90 index 0e7e6425..00021773 100644 --- a/src/dftd3/param.f90 +++ b/src/dftd3/param.f90 @@ -62,8 +62,10 @@ module dftd3_param & p_m08hx_df, p_m11l_df, p_mn15l_df, p_pwp_df, p_r2scanh_df, p_r2scan0_df, & & p_r2scan50_df, p_b973c_df, p_dm21_df, p_dm21m_df, p_dm21mc_df, p_dm21mu_df, & & p_dsdpbep86_df, p_dsdpbeb95_df, p_dsdpbe_df, p_dodscan66_df, & - & p_revdsdblyp_df, p_revdsdpbep86_df, p_revdsdpbeb95_df, p_revdsdpbe_df, & - & p_revdodblyp_df, p_revdodpbep86_df, p_revdodpbeb95_df, p_revdodpbe_df + & p_revdsdblyp_df, p_revdsdpbep86_df, p_revdsdpbeb95_df, p_revdsdpbe_df, & + & p_revdodblyp_df, p_revdodpbep86_df, p_revdodpbeb95_df, p_revdodpbe_df, & + & p_pw91_df, p_drpa75_df, p_scs_drpa75_df, p_optscs_drpa75_df, & + & p_dsd_pbe_drpa75_df, p_dsd_pbep86_drpa75_df end enum contains @@ -116,10 +118,13 @@ function get_method_id(method) result(id) case("dm21m"); id = p_dm21m_df case("dm21mc"); id = p_dm21mc_df case("dm21mu"); id = p_dm21mu_df + case("drpa75"); id = p_drpa75_df case("dsdblyp"); id = p_dsdblyp_df case("dsdblypfc"); id = p_dsdblypfc_df case("dsdpbe"); id = p_dsdpbe_df + case("dsdpbedrpa75"); id = p_dsd_pbe_drpa75_df case("dsdpbep86"); id = p_dsdpbep86_df + case("dsdpbep86drpa75"); id = p_dsd_pbep86_drpa75_df case("dsdpbeb95"); id = p_dsdpbeb95_df case("dodscan66"); id = p_dodscan66_df case("hcth120"); id = p_hcth120_df @@ -164,6 +169,7 @@ function get_method_id(method) result(id) case("o3lyp"); id = p_o3lyp_df case("olyp"); id = p_olyp_df case("opbe"); id = p_opbe_df + case("optscsdrpa75"); id = p_optscs_drpa75_df case("otpss"); id = p_otpss_df case("pbe"); id = p_pbe_df case("pbe0"); id = p_pbe0_df @@ -178,6 +184,7 @@ function get_method_id(method) result(id) case("pwp", "pw91p86"); id = p_pwp_df case("pw1pw"); id = p_pw1pw_df case("pw6b95"); id = p_pw6b95_df + case("pw91"); id = p_pw91_df case("pwb6k"); id = p_pwb6k_df case("pwgga"); id = p_pwgga_df case("pwpb95"); id = p_pwpb95_df @@ -204,6 +211,7 @@ function get_method_id(method) result(id) case("rpw86pbe"); id = p_rpw86pbe_df case("rscan"); id = p_rscan_df case("scan"); id = p_scan_df + case("scsdrpa75"); id = p_scs_drpa75_df case("slaterdiracexchange"); id = p_slaterdiracexchange_df case("sogga11x"); id = p_sogga11x_df case("ssb"); id = p_ssb_df @@ -322,6 +330,8 @@ subroutine get_rational_damping(param, method, error, s9) param = d3_param(a1=0.0000_wp, s8=0.2804_wp, a2=6.5745_wp, s6=0.750_wp) case(p_pwpb95_df) param = d3_param(a1=0.0000_wp, s8=0.2904_wp, a2=7.3141_wp, s6=0.820_wp) + case(p_pw91_df) + param = d3_param(a1=0.6319_wp, s8=1.9598_wp, a2=4.5718_wp) case(p_hf_mixed_df) param = d3_param(a1=0.5607_wp, s8=3.9027_wp, a2=4.5622_wp) case(p_hf_sv_df) @@ -452,6 +462,16 @@ subroutine get_rational_damping(param, method, error, s9) param = d3_param(s6=0.4107_wp, a1=0.0_wp, s8=0.0_wp, a2=5.5_wp) case(p_revdodpbe_df) param = d3_param(s6=0.6067_wp, a1=0.0_wp, s8=0.0_wp, a2=5.5_wp) + case(p_drpa75_df) + param = d3_param(s6=0.3754_wp, a1=0.0_wp, s8=0.0_wp, a2=4.5048_wp) + case(p_scs_drpa75_df) + param = d3_param(s6=0.2528_wp, a1=0.0_wp, s8=0.0_wp, a2=4.5050_wp) + case(p_optscs_drpa75_df) + param = d3_param(s6=0.2546_wp, a1=0.0_wp, s8=0.0_wp, a2=4.5050_wp) + case(p_dsd_pbe_drpa75_df) + param = d3_param(s6=0.3223_wp, a1=0.0_wp, s8=0.0_wp, a2=4.5050_wp) + case(p_dsd_pbep86_drpa75_df) + param = d3_param(s6=0.3012_wp, a1=0.0_wp, s8=0.0_wp, a2=4.5050_wp) end select if (present(s9)) then diff --git a/test/unit/test_param.f90 b/test/unit/test_param.f90 index 4bfac922..726d3cba 100644 --- a/test/unit/test_param.f90 +++ b/test/unit/test_param.f90 @@ -81,7 +81,7 @@ subroutine test_dftd3_gen(error, mol, param, ref) call check(error, energy, ref, thr=thr) if (allocated(error)) then - print *,energy + print '(ES23.16)', energy end if end subroutine test_dftd3_gen @@ -113,7 +113,8 @@ subroutine test_d3bj_mb01(error) & "r2scanh", "r2scan0", "r2scan50", "b973c", "b3lyp3", "b3lyp5", "b3lypg", & & "dm21", "dm21m", "dm21mc", "dm21mu", "dsdpbep86", "dsdpbeb95", "dsdpbe", & & "dodscan66", "revdsdblyp", "revdsdpbep86", "revdsdpbeb95", "revdsdpbe", & - & "revdodblyp", "revdodpbep86", "revdodpbeb95", "revdodpbe"] + & "revdodblyp", "revdodpbep86", "revdodpbeb95", "revdodpbe", "pw91", "drpa75", & + & "scsdrpa75", "optscsdrpa75", "dsdpbedrpa75", "dsdpbep86drpa75"] real(wp), parameter :: ref(*) = [& &-2.9551694676908012E-2_wp,-1.6638703086788331E-2_wp,-1.6725877716130381E-2_wp, & &-3.3014429592265318E-2_wp,-2.2051435219996540E-2_wp,-3.3481565825316001E-2_wp, & @@ -152,7 +153,9 @@ subroutine test_d3bj_mb01(error) &-1.7758051042408327E-2_wp,-9.2323777698865989E-3_wp,-2.4042802546329414E-2_wp, & &-1.5411238160834008E-2_wp,-1.2978255394296128E-2_wp,-2.0231431225074755E-2_wp, & &-2.7103838130103503E-2_wp,-1.6794975103307797E-2_wp,-1.4460579192722246E-2_wp, & - &-2.1361659109385409E-2_wp] + &-2.1361659109385409E-2_wp,-9.7460845051874183E-3_wp,-2.8391262582808205E-2_wp, & + &-1.9116077511614916E-2_wp,-1.9252188823010906E-2_wp,-2.4371486479404614E-2_wp, & + &-2.2775959440262707E-02_wp] call get_structure(mol, "MB16-43", "01") do ii = 1, size(func)