From ef88e168249e38c526353fe7c48aab4bc697c99d Mon Sep 17 00:00:00 2001 From: Aodong Liu Date: Fri, 3 May 2024 15:42:37 -0700 Subject: [PATCH] Add Builtin LDA Electron-Proton Correlation Functionals (#38) * first commit * grabbed lda epcs from libxc * add epc functional options * enable more epc interface * added is_epc check and start generating tests * attemp to add tests for epc lda functionals * Remove unwanted changes * skip ground truth comparison for epc functionals and allow is_epc trait to automatically generate * Update generate_primitive_kernels.py * change epc checking logic in tests * check if spin is polarized before constructing epc kernel * pre libxc7 epc=false --- bin/generate_primitive_kernels.py | 58 ++++- include/exchcxx/enums/functionals.hpp | 6 +- include/exchcxx/enums/kernels.hpp | 6 + include/exchcxx/impl/builtin/fwd.hpp | 4 + include/exchcxx/impl/builtin/interface.hpp | 1 + include/exchcxx/impl/builtin/kernel_type.hpp | 4 + include/exchcxx/impl/builtin/kernels.hpp | 4 + .../exchcxx/impl/builtin/kernels/b3lyp.hpp | 1 + include/exchcxx/impl/builtin/kernels/b88.hpp | 3 +- .../impl/builtin/kernels/deorbitalized.hpp | 1 + .../exchcxx/impl/builtin/kernels/epc17_1.hpp | 211 +++++++++++++++++ .../exchcxx/impl/builtin/kernels/epc17_2.hpp | 211 +++++++++++++++++ .../exchcxx/impl/builtin/kernels/epc18_1.hpp | 216 ++++++++++++++++++ .../exchcxx/impl/builtin/kernels/epc18_2.hpp | 216 ++++++++++++++++++ .../exchcxx/impl/builtin/kernels/ft98_x.hpp | 3 +- include/exchcxx/impl/builtin/kernels/lyp.hpp | 1 + include/exchcxx/impl/builtin/kernels/pbe0.hpp | 1 + .../exchcxx/impl/builtin/kernels/pbe_c.hpp | 1 + .../exchcxx/impl/builtin/kernels/pbe_x.hpp | 3 +- .../exchcxx/impl/builtin/kernels/pc07_k.hpp | 3 +- .../impl/builtin/kernels/pc07opt_k.hpp | 3 +- .../exchcxx/impl/builtin/kernels/pw91_lda.hpp | 1 + .../impl/builtin/kernels/pw91_lda_mod.hpp | 1 + .../impl/builtin/kernels/pw91_lda_rpa.hpp | 1 + include/exchcxx/impl/builtin/kernels/pz81.hpp | 1 + .../exchcxx/impl/builtin/kernels/pz81_mod.hpp | 1 + .../exchcxx/impl/builtin/kernels/r2scan_c.hpp | 1 + .../exchcxx/impl/builtin/kernels/r2scan_x.hpp | 3 +- .../impl/builtin/kernels/rev_pbe_x.hpp | 3 +- .../exchcxx/impl/builtin/kernels/scan_c.hpp | 5 +- .../exchcxx/impl/builtin/kernels/scan_x.hpp | 3 +- .../impl/builtin/kernels/slater_exchange.hpp | 3 +- include/exchcxx/impl/builtin/kernels/vwn3.hpp | 1 + .../exchcxx/impl/builtin/kernels/vwn_rpa.hpp | 1 + include/exchcxx/impl/libxc.hpp | 1 + include/exchcxx/impl/xc_kernel.hpp | 2 + include/exchcxx/xc_functional.hpp | 8 + include/exchcxx/xc_kernel.hpp | 1 + src/builtin_interface.cxx | 22 +- src/builtin_kernel.cxx | 5 + src/cuda/builtin.cu | 6 + src/hip/builtin.hip | 5 + src/libxc.cxx | 10 + src/sycl/builtin_sycl.cxx | 5 + src/xc_functional.cxx | 22 +- src/xc_kernel.cxx | 6 +- test/reference_values.cxx | 1 - test/ut_common.hpp | 31 ++- test/xc_functional_test.cxx | 3 + test/xc_kernel_test.cxx | 89 +++++++- 50 files changed, 1161 insertions(+), 38 deletions(-) create mode 100644 include/exchcxx/impl/builtin/kernels/epc17_1.hpp create mode 100644 include/exchcxx/impl/builtin/kernels/epc17_2.hpp create mode 100644 include/exchcxx/impl/builtin/kernels/epc18_1.hpp create mode 100644 include/exchcxx/impl/builtin/kernels/epc18_2.hpp diff --git a/bin/generate_primitive_kernels.py b/bin/generate_primitive_kernels.py index 75b46a9..978a6ed 100644 --- a/bin/generate_primitive_kernels.py +++ b/bin/generate_primitive_kernels.py @@ -472,7 +472,7 @@ def get_const_lines( self, _xc_lines, xc_vars ): class GenMetaData: def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff, - params = {}, needs_laplacian = False, is_kedf = False ): + params = {}, needs_laplacian = False, is_kedf = False, is_epc = False ): self.local_name = local_name self.libxc_file = libxc_file self.ofname = ofname @@ -482,6 +482,7 @@ def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff self.params = params self.needs_laplacian = needs_laplacian self.is_kedf = is_kedf + self.is_epc = is_epc self.is_hyb = abs(float(exx_coeff)) > 1e-10 @@ -493,6 +494,46 @@ def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff kernel_prefix = 'include/exchcxx/impl/builtin/kernels/' gen_table = { + 'EPC17_1' : GenMetaData( 'BuiltinEPC17_1', + libxc_prefix + 'lda_exc/lda_c_epc17.c', + kernel_prefix + 'epc17_1.hpp', + 'LDA', 1e-24, 0., + { 'a': '2.35', + 'b': '2.40', + 'c': '3.20' }, + False, False, True + ), + + 'EPC17_2' : GenMetaData( 'BuiltinEPC17_2', + libxc_prefix + 'lda_exc/lda_c_epc17.c', + kernel_prefix + 'epc17_2.hpp', + 'LDA', 1e-24, 0., + { 'a': '2.35', + 'b': '2.40', + 'c': '6.60' }, + False, False, True + ), + + 'EPC18_1' : GenMetaData( 'BuiltinEPC18_1', + libxc_prefix + 'lda_exc/lda_c_epc18.c', + kernel_prefix + 'epc18_1.hpp', + 'LDA', 1e-24, 0., + { 'a': '1.80', + 'b': '0.10', + 'c': '0.03' }, + False, False, True + ), + + 'EPC18_2' : GenMetaData( 'BuiltinEPC18_2', + libxc_prefix + 'lda_exc/lda_c_epc18.c', + kernel_prefix + 'epc18_2.hpp', + 'LDA', 1e-24, 0., + { 'a': '3.90', + 'b': '0.50', + 'c': '0.06' }, + False, False, True + ), + 'SlaterExchange' : GenMetaData( 'BuiltinSlaterExchange', libxc_prefix + 'lda_exc/lda_x.c', kernel_prefix + 'slater_exchange.hpp', @@ -734,14 +775,15 @@ def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff static constexpr bool is_mgga = {4}; static constexpr bool needs_laplacian = {5}; static constexpr bool is_kedf = {6}; + static constexpr bool is_epc = {7}; - static constexpr double dens_tol = {7}; + static constexpr double dens_tol = {8}; static constexpr double zeta_tol = 1e-15; - static constexpr double sigma_tol = {8}; + static constexpr double sigma_tol = {9}; static constexpr double tau_tol = is_kedf ? 0.0 : 1e-20; - static constexpr bool is_hyb = {9}; - static constexpr double exx_coeff = {10}; + static constexpr bool is_hyb = {10}; + static constexpr double exx_coeff = {11}; """ lda_exc_args_unpolar = "double rho, double& eps" @@ -811,11 +853,13 @@ def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff is_mgga = xc_type == 'MGGA' needs_laplacian = (xc_type == 'MGGA') and meta_data.needs_laplacian is_kedf = meta_data.is_kedf + is_epc = meta_data.is_epc xc_struct_prefix = struct_prefix.format( meta_data.local_name, xc_type.lower(), str(is_lda).lower(), - str(is_gga).lower(), str(is_mgga).lower(), str(needs_laplacian).lower(), str(is_kedf).lower(), - meta_data.dens_tol, meta_data.dens_tol**(4.0/3.0), str(meta_data.is_hyb).lower(), meta_data.exx_coeff ) + str(is_gga).lower(), str(is_mgga).lower(), str(needs_laplacian).lower(), + str(is_kedf).lower(), str(is_epc).lower(), meta_data.dens_tol, + meta_data.dens_tol**(4.0/3.0), str(meta_data.is_hyb).lower(), meta_data.exx_coeff ) xc_param_lines = [] for pname, valstr in meta_data.params.items(): diff --git a/include/exchcxx/enums/functionals.hpp b/include/exchcxx/enums/functionals.hpp index 5da6fb5..c83928b 100644 --- a/include/exchcxx/enums/functionals.hpp +++ b/include/exchcxx/enums/functionals.hpp @@ -59,7 +59,11 @@ enum class Functional { PBE0, SCAN, R2SCAN, - R2SCANL + R2SCANL, + EPC17_1, + EPC17_2, + EPC18_1, + EPC18_2, }; extern BidirectionalMap functional_map; diff --git a/include/exchcxx/enums/kernels.hpp b/include/exchcxx/enums/kernels.hpp index 93680a1..d7eeefa 100644 --- a/include/exchcxx/enums/kernels.hpp +++ b/include/exchcxx/enums/kernels.hpp @@ -101,6 +101,12 @@ enum class Kernel { // Hybrid GGA functionals B3LYP, PBE0, + + // NEO LDA Functionals + EPC17_1, + EPC17_2, + EPC18_1, + EPC18_2, }; extern BidirectionalMap kernel_map; diff --git a/include/exchcxx/impl/builtin/fwd.hpp b/include/exchcxx/impl/builtin/fwd.hpp index 20b9a2c..fa91878 100644 --- a/include/exchcxx/impl/builtin/fwd.hpp +++ b/include/exchcxx/impl/builtin/fwd.hpp @@ -88,6 +88,10 @@ struct BuiltinSCANL_C; struct BuiltinR2SCANL_X; struct BuiltinR2SCANL_C; +struct BuiltinEPC17_1; +struct BuiltinEPC17_2; +struct BuiltinEPC18_1; +struct BuiltinEPC18_2; template struct kernel_traits; diff --git a/include/exchcxx/impl/builtin/interface.hpp b/include/exchcxx/impl/builtin/interface.hpp index cab6f5b..91cc78a 100644 --- a/include/exchcxx/impl/builtin/interface.hpp +++ b/include/exchcxx/impl/builtin/interface.hpp @@ -72,6 +72,7 @@ class BuiltinKernelInterface : public XCKernelImpl { bool is_gga_() const noexcept override; bool is_mgga_() const noexcept override; bool is_hyb_() const noexcept override; + bool is_epc_() const noexcept override; bool needs_laplacian_() const noexcept override; bool needs_tau_() const noexcept override; bool is_polarized_() const noexcept override; diff --git a/include/exchcxx/impl/builtin/kernel_type.hpp b/include/exchcxx/impl/builtin/kernel_type.hpp index 7078766..b211c71 100644 --- a/include/exchcxx/impl/builtin/kernel_type.hpp +++ b/include/exchcxx/impl/builtin/kernel_type.hpp @@ -71,6 +71,7 @@ class BuiltinKernel { virtual bool is_gga() const noexcept = 0; virtual bool is_mgga() const noexcept = 0; virtual bool is_hyb() const noexcept = 0; + virtual bool is_epc() const noexcept = 0; virtual bool needs_laplacian() const noexcept = 0; virtual bool needs_tau() const noexcept = 0; virtual double hyb_exx() const noexcept = 0; @@ -375,6 +376,7 @@ struct BuiltinKernelImpl< inline bool is_lda() const noexcept override { return traits::is_lda; } inline bool is_gga() const noexcept override { return traits::is_gga; } inline bool is_mgga() const noexcept override { return traits::is_mgga; } + inline bool is_epc() const noexcept override { return traits::is_epc; } inline double hyb_exx() const noexcept override { return traits::is_hyb ? traits::exx_coeff : 0.; @@ -525,6 +527,7 @@ struct BuiltinKernelImpl< inline bool is_lda() const noexcept override { return traits::is_lda; } inline bool is_gga() const noexcept override { return traits::is_gga; } inline bool is_mgga() const noexcept override { return traits::is_mgga; } + inline bool is_epc() const noexcept override { return traits::is_epc; } inline double hyb_exx() const noexcept override { return traits::is_hyb ? traits::exx_coeff : 0.; @@ -676,6 +679,7 @@ struct BuiltinKernelImpl< inline bool is_lda() const noexcept override { return traits::is_lda; } inline bool is_gga() const noexcept override { return traits::is_gga; } inline bool is_mgga() const noexcept override { return traits::is_mgga; } + inline bool is_epc() const noexcept override { return traits::is_epc; } inline double hyb_exx() const noexcept override { return traits::is_hyb ? traits::exx_coeff : 0.; diff --git a/include/exchcxx/impl/builtin/kernels.hpp b/include/exchcxx/impl/builtin/kernels.hpp index c1841b3..e3d5c30 100644 --- a/include/exchcxx/impl/builtin/kernels.hpp +++ b/include/exchcxx/impl/builtin/kernels.hpp @@ -81,4 +81,8 @@ #include #include +#include +#include +#include +#include diff --git a/include/exchcxx/impl/builtin/kernels/b3lyp.hpp b/include/exchcxx/impl/builtin/kernels/b3lyp.hpp index 0e19f3a..6f7b699 100644 --- a/include/exchcxx/impl/builtin/kernels/b3lyp.hpp +++ b/include/exchcxx/impl/builtin/kernels/b3lyp.hpp @@ -20,6 +20,7 @@ struct kernel_traits { static constexpr bool is_lda = false; static constexpr bool is_gga = true; static constexpr bool is_mgga = false; + static constexpr bool is_epc = false; static constexpr double exx_coeff = 0.20; static constexpr double b3lyp_ax = 0.72; // % GGA X diff --git a/include/exchcxx/impl/builtin/kernels/b88.hpp b/include/exchcxx/impl/builtin/kernels/b88.hpp index ef8bd0f..a70bfc4 100644 --- a/include/exchcxx/impl/builtin/kernels/b88.hpp +++ b/include/exchcxx/impl/builtin/kernels/b88.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinB88 > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-25; static constexpr double zeta_tol = 1e-15; @@ -391,4 +392,4 @@ struct BuiltinB88 : detail::BuiltinKernelImpl< BuiltinB88 > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/deorbitalized.hpp b/include/exchcxx/impl/builtin/kernels/deorbitalized.hpp index 7fa0991..9daab61 100644 --- a/include/exchcxx/impl/builtin/kernels/deorbitalized.hpp +++ b/include/exchcxx/impl/builtin/kernels/deorbitalized.hpp @@ -18,6 +18,7 @@ struct kernel_traits> { static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = true; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double exx_coeff = xc_traits::exx_coeff + ke_traits::exx_coeff; BUILTIN_KERNEL_EVAL_RETURN diff --git a/include/exchcxx/impl/builtin/kernels/epc17_1.hpp b/include/exchcxx/impl/builtin/kernels/epc17_1.hpp new file mode 100644 index 0000000..a33d45d --- /dev/null +++ b/include/exchcxx/impl/builtin/kernels/epc17_1.hpp @@ -0,0 +1,211 @@ +#pragma once +#include + +#include +#include +#include +#include + +#include + + + +namespace ExchCXX { + +template <> +struct kernel_traits< BuiltinEPC17_1 > : + public lda_screening_interface< BuiltinEPC17_1 > { + + static constexpr bool is_lda = true; + static constexpr bool is_gga = false; + static constexpr bool is_mgga = false; + static constexpr bool needs_laplacian = false; + static constexpr bool is_kedf = false; + static constexpr bool is_epc = true; + + static constexpr double dens_tol = 1e-24; + static constexpr double zeta_tol = 1e-15; + static constexpr double sigma_tol = 1.000000000000004e-32; + static constexpr double tau_tol = is_kedf ? 0.0 : 1e-20; + + static constexpr bool is_hyb = false; + static constexpr double exx_coeff = 0.0; + + static constexpr double a = 2.35; + static constexpr double b = 2.40; + static constexpr double c = 3.20; + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_unpolar_impl( double rho, double& eps ) { + + (void)(eps); + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t8 = square( 0.1e1 + t6 ); + const double t9 = t8 * rho; + const double t10 = rho * rho; + const double t11 = t8 * t10; + const double t12 = safe_math::sqrt( t11 ); + const double t15 = c * t8; + const double t18 = a - b * t12 / 0.2e1 + t15 * t10 / 0.4e1; + const double t19 = 0.1e1 / t18; + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t19 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_unpolar_impl( double rho, double& eps, double& vrho ) { + + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t8 = square( 0.1e1 + t6 ); + const double t9 = t8 * rho; + const double t10 = rho * rho; + const double t11 = t8 * t10; + const double t12 = safe_math::sqrt( t11 ); + const double t15 = c * t8; + const double t18 = a - b * t12 / 0.2e1 + t15 * t10 / 0.4e1; + const double t19 = 0.1e1 / t18; + const double t23 = t18 * t18; + const double t24 = 0.1e1 / t23; + const double t26 = b / t12; + const double t30 = t15 * rho / 0.2e1 - t26 * t9 / 0.2e1; + const double t35 = piecewise_functor_3( t2, 0.0, t9 * t24 * t30 / 0.4e1 - t8 * t19 / 0.4e1 ); + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t19 / 0.4e1 ); + vrho = rho * t35 + eps; + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_polar_impl( double rho_a, double rho_b, double& eps ) { + + (void)(eps); + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t23 = t4 * t4; + const double t24 = t17 * t23; + const double t25 = t24 * t22; + const double t26 = safe_math::sqrt( t25 ); + const double t29 = c * t17; + const double t30 = t23 * t22; + const double t33 = a - b * t26 / 0.2e1 + t29 * t30 / 0.4e1; + const double t34 = 0.1e1 / t33; + const double t35 = t22 * t34; + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t35 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_polar_impl( double rho_a, double rho_b, double& eps, double& vrho_a, double& vrho_b ) { + + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t23 = t4 * t4; + const double t24 = t17 * t23; + const double t25 = t24 * t22; + const double t26 = safe_math::sqrt( t25 ); + const double t29 = c * t17; + const double t30 = t23 * t22; + const double t33 = a - b * t26 / 0.2e1 + t29 * t30 / 0.4e1; + const double t34 = 0.1e1 / t33; + const double t35 = t22 * t34; + const double t38 = 0.1e1 / t23; + const double t39 = t14 * t38; + const double t41 = piecewise_functor_5( t8, 0.0, t12, 0.0, t5 - t39 ); + const double t42 = t41 * t4; + const double t44 = t17 * t22; + const double t45 = t44 * t34; + const double t46 = t19 * t38; + const double t48 = piecewise_functor_5( t12, 0.0, t8, 0.0, -t5 - t46 ); + const double t49 = t48 * t34; + const double t51 = t33 * t33; + const double t52 = 0.1e1 / t51; + const double t53 = t22 * t52; + const double t55 = b / t26; + const double t56 = t41 * t23; + const double t58 = t18 * t22; + const double t59 = 0.2e1 * t58; + const double t61 = t56 * t22 + t24 * t48 + t59; + const double t64 = c * t41; + const double t67 = t4 * t22; + const double t69 = t29 * t67 / 0.2e1; + const double t70 = t23 * t48; + const double t73 = -t55 * t61 / 0.4e1 + t64 * t30 / 0.4e1 + t69 + t29 * t70 / 0.4e1; + const double t74 = t53 * t73; + const double t78 = piecewise_functor_3( t3, 0.0, -t18 * t49 / 0.4e1 + t18 * t74 / 0.4e1 - t42 * t35 / 0.4e1 - t45 / 0.4e1 ); + const double t81 = piecewise_functor_5( t8, 0.0, t12, 0.0, -t5 - t39 ); + const double t82 = t81 * t4; + const double t85 = piecewise_functor_5( t12, 0.0, t8, 0.0, t5 - t46 ); + const double t86 = t85 * t34; + const double t88 = t81 * t23; + const double t91 = t88 * t22 + t24 * t85 + t59; + const double t94 = c * t81; + const double t97 = t23 * t85; + const double t100 = -t55 * t91 / 0.4e1 + t94 * t30 / 0.4e1 + t69 + t29 * t97 / 0.4e1; + const double t101 = t53 * t100; + const double t105 = piecewise_functor_3( t3, 0.0, t18 * t101 / 0.4e1 - t18 * t86 / 0.4e1 - t82 * t35 / 0.4e1 - t45 / 0.4e1 ); + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t35 / 0.4e1 ); + vrho_a = t4 * t78 + eps; + vrho_b = t4 * t105 + eps; + + } + + +}; + +struct BuiltinEPC17_1 : detail::BuiltinKernelImpl< BuiltinEPC17_1 > { + + BuiltinEPC17_1( Spin p ) : + detail::BuiltinKernelImpl< BuiltinEPC17_1 >(p) { } + + virtual ~BuiltinEPC17_1() = default; + +}; + + + +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/epc17_2.hpp b/include/exchcxx/impl/builtin/kernels/epc17_2.hpp new file mode 100644 index 0000000..69a8104 --- /dev/null +++ b/include/exchcxx/impl/builtin/kernels/epc17_2.hpp @@ -0,0 +1,211 @@ +#pragma once +#include + +#include +#include +#include +#include + +#include + + + +namespace ExchCXX { + +template <> +struct kernel_traits< BuiltinEPC17_2 > : + public lda_screening_interface< BuiltinEPC17_2 > { + + static constexpr bool is_lda = true; + static constexpr bool is_gga = false; + static constexpr bool is_mgga = false; + static constexpr bool needs_laplacian = false; + static constexpr bool is_kedf = false; + static constexpr bool is_epc = true; + + static constexpr double dens_tol = 1e-24; + static constexpr double zeta_tol = 1e-15; + static constexpr double sigma_tol = 1.000000000000004e-32; + static constexpr double tau_tol = is_kedf ? 0.0 : 1e-20; + + static constexpr bool is_hyb = false; + static constexpr double exx_coeff = 0.0; + + static constexpr double a = 2.35; + static constexpr double b = 2.40; + static constexpr double c = 6.60; + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_unpolar_impl( double rho, double& eps ) { + + (void)(eps); + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t8 = square( 0.1e1 + t6 ); + const double t9 = t8 * rho; + const double t10 = rho * rho; + const double t11 = t8 * t10; + const double t12 = safe_math::sqrt( t11 ); + const double t15 = c * t8; + const double t18 = a - b * t12 / 0.2e1 + t15 * t10 / 0.4e1; + const double t19 = 0.1e1 / t18; + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t19 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_unpolar_impl( double rho, double& eps, double& vrho ) { + + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t8 = square( 0.1e1 + t6 ); + const double t9 = t8 * rho; + const double t10 = rho * rho; + const double t11 = t8 * t10; + const double t12 = safe_math::sqrt( t11 ); + const double t15 = c * t8; + const double t18 = a - b * t12 / 0.2e1 + t15 * t10 / 0.4e1; + const double t19 = 0.1e1 / t18; + const double t23 = t18 * t18; + const double t24 = 0.1e1 / t23; + const double t26 = b / t12; + const double t30 = t15 * rho / 0.2e1 - t26 * t9 / 0.2e1; + const double t35 = piecewise_functor_3( t2, 0.0, t9 * t24 * t30 / 0.4e1 - t8 * t19 / 0.4e1 ); + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t19 / 0.4e1 ); + vrho = rho * t35 + eps; + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_polar_impl( double rho_a, double rho_b, double& eps ) { + + (void)(eps); + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t23 = t4 * t4; + const double t24 = t17 * t23; + const double t25 = t24 * t22; + const double t26 = safe_math::sqrt( t25 ); + const double t29 = c * t17; + const double t30 = t23 * t22; + const double t33 = a - b * t26 / 0.2e1 + t29 * t30 / 0.4e1; + const double t34 = 0.1e1 / t33; + const double t35 = t22 * t34; + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t35 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_polar_impl( double rho_a, double rho_b, double& eps, double& vrho_a, double& vrho_b ) { + + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t23 = t4 * t4; + const double t24 = t17 * t23; + const double t25 = t24 * t22; + const double t26 = safe_math::sqrt( t25 ); + const double t29 = c * t17; + const double t30 = t23 * t22; + const double t33 = a - b * t26 / 0.2e1 + t29 * t30 / 0.4e1; + const double t34 = 0.1e1 / t33; + const double t35 = t22 * t34; + const double t38 = 0.1e1 / t23; + const double t39 = t14 * t38; + const double t41 = piecewise_functor_5( t8, 0.0, t12, 0.0, t5 - t39 ); + const double t42 = t41 * t4; + const double t44 = t17 * t22; + const double t45 = t44 * t34; + const double t46 = t19 * t38; + const double t48 = piecewise_functor_5( t12, 0.0, t8, 0.0, -t5 - t46 ); + const double t49 = t48 * t34; + const double t51 = t33 * t33; + const double t52 = 0.1e1 / t51; + const double t53 = t22 * t52; + const double t55 = b / t26; + const double t56 = t41 * t23; + const double t58 = t18 * t22; + const double t59 = 0.2e1 * t58; + const double t61 = t56 * t22 + t24 * t48 + t59; + const double t64 = c * t41; + const double t67 = t4 * t22; + const double t69 = t29 * t67 / 0.2e1; + const double t70 = t23 * t48; + const double t73 = -t55 * t61 / 0.4e1 + t64 * t30 / 0.4e1 + t69 + t29 * t70 / 0.4e1; + const double t74 = t53 * t73; + const double t78 = piecewise_functor_3( t3, 0.0, -t18 * t49 / 0.4e1 + t18 * t74 / 0.4e1 - t42 * t35 / 0.4e1 - t45 / 0.4e1 ); + const double t81 = piecewise_functor_5( t8, 0.0, t12, 0.0, -t5 - t39 ); + const double t82 = t81 * t4; + const double t85 = piecewise_functor_5( t12, 0.0, t8, 0.0, t5 - t46 ); + const double t86 = t85 * t34; + const double t88 = t81 * t23; + const double t91 = t88 * t22 + t24 * t85 + t59; + const double t94 = c * t81; + const double t97 = t23 * t85; + const double t100 = -t55 * t91 / 0.4e1 + t94 * t30 / 0.4e1 + t69 + t29 * t97 / 0.4e1; + const double t101 = t53 * t100; + const double t105 = piecewise_functor_3( t3, 0.0, t18 * t101 / 0.4e1 - t18 * t86 / 0.4e1 - t82 * t35 / 0.4e1 - t45 / 0.4e1 ); + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t35 / 0.4e1 ); + vrho_a = t4 * t78 + eps; + vrho_b = t4 * t105 + eps; + + } + + +}; + +struct BuiltinEPC17_2 : detail::BuiltinKernelImpl< BuiltinEPC17_2 > { + + BuiltinEPC17_2( Spin p ) : + detail::BuiltinKernelImpl< BuiltinEPC17_2 >(p) { } + + virtual ~BuiltinEPC17_2() = default; + +}; + + + +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/epc18_1.hpp b/include/exchcxx/impl/builtin/kernels/epc18_1.hpp new file mode 100644 index 0000000..d19fe28 --- /dev/null +++ b/include/exchcxx/impl/builtin/kernels/epc18_1.hpp @@ -0,0 +1,216 @@ +#pragma once +#include + +#include +#include +#include +#include + +#include + + + +namespace ExchCXX { + +template <> +struct kernel_traits< BuiltinEPC18_1 > : + public lda_screening_interface< BuiltinEPC18_1 > { + + static constexpr bool is_lda = true; + static constexpr bool is_gga = false; + static constexpr bool is_mgga = false; + static constexpr bool needs_laplacian = false; + static constexpr bool is_kedf = false; + static constexpr bool is_epc = true; + + static constexpr double dens_tol = 1e-24; + static constexpr double zeta_tol = 1e-15; + static constexpr double sigma_tol = 1.000000000000004e-32; + static constexpr double tau_tol = is_kedf ? 0.0 : 1e-20; + + static constexpr bool is_hyb = false; + static constexpr double exx_coeff = 0.0; + + static constexpr double a = 1.80; + static constexpr double b = 0.10; + static constexpr double c = 0.03; + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_unpolar_impl( double rho, double& eps ) { + + (void)(eps); + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t7 = 0.1e1 + t6; + const double t8 = t7 * t7; + const double t9 = t8 * rho; + const double t10 = b * t7; + const double t13 = c * t8; + const double t14 = rho * rho; + const double t17 = -0.4e1 * t10 * rho + 0.16e2 * t13 * t14 + a; + const double t18 = 0.1e1 / t17; + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t18 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_unpolar_impl( double rho, double& eps, double& vrho ) { + + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t7 = 0.1e1 + t6; + const double t8 = t7 * t7; + const double t9 = t8 * rho; + const double t10 = b * t7; + const double t13 = c * t8; + const double t14 = rho * rho; + const double t17 = -0.4e1 * t10 * rho + 0.16e2 * t13 * t14 + a; + const double t18 = 0.1e1 / t17; + const double t22 = t17 * t17; + const double t23 = 0.1e1 / t22; + const double t27 = 0.32e2 * t13 * rho - 0.4e1 * t10; + const double t32 = piecewise_functor_3( t2, 0.0, t9 * t23 * t27 / 0.4e1 - t8 * t18 / 0.4e1 ); + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t18 / 0.4e1 ); + vrho = rho * t32 + eps; + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_polar_impl( double rho_a, double rho_b, double& eps ) { + + (void)(eps); + constexpr double t23 = constants::m_cbrt_2; + constexpr double t24 = t23 * t23; + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t25 = safe_math::cbrt( t18 ); + const double t27 = t22 * t4; + const double t28 = safe_math::cbrt( t27 ); + const double t31 = t24 * t25 / 0.2e1 + t24 * t28 / 0.2e1; + const double t32 = t31 * t31; + const double t33 = t32 * t31; + const double t35 = t32 * t32; + const double t38 = c * t35 * t32 - b * t33 + a; + const double t39 = 0.1e1 / t38; + const double t40 = t22 * t39; + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t40 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_polar_impl( double rho_a, double rho_b, double& eps, double& vrho_a, double& vrho_b ) { + + constexpr double t23 = constants::m_cbrt_2; + constexpr double t24 = t23 * t23; + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t25 = safe_math::cbrt( t18 ); + const double t27 = t22 * t4; + const double t28 = safe_math::cbrt( t27 ); + const double t31 = t24 * t25 / 0.2e1 + t24 * t28 / 0.2e1; + const double t32 = t31 * t31; + const double t33 = t32 * t31; + const double t35 = t32 * t32; + const double t38 = c * t35 * t32 - b * t33 + a; + const double t39 = 0.1e1 / t38; + const double t40 = t22 * t39; + const double t43 = t4 * t4; + const double t44 = 0.1e1 / t43; + const double t45 = t14 * t44; + const double t47 = piecewise_functor_5( t8, 0.0, t12, 0.0, t5 - t45 ); + const double t48 = t47 * t4; + const double t50 = t17 * t22; + const double t51 = t50 * t39; + const double t52 = t19 * t44; + const double t54 = piecewise_functor_5( t12, 0.0, t8, 0.0, -t5 - t52 ); + const double t55 = t54 * t39; + const double t57 = t38 * t38; + const double t58 = 0.1e1 / t57; + const double t59 = t22 * t58; + const double t60 = b * t32; + const double t61 = t25 * t25; + const double t63 = t24 / t61; + const double t64 = t48 + 0.1e1 + t16; + const double t66 = t28 * t28; + const double t68 = t24 / t66; + const double t70 = t54 * t4 + t21 + 0.1e1; + const double t73 = t63 * t64 / 0.6e1 + t68 * t70 / 0.6e1; + const double t77 = c * t35 * t31; + const double t80 = -0.3e1 * t60 * t73 + 0.6e1 * t77 * t73; + const double t81 = t59 * t80; + const double t85 = piecewise_functor_3( t3, 0.0, -t18 * t55 / 0.4e1 + t18 * t81 / 0.4e1 - t48 * t40 / 0.4e1 - t51 / 0.4e1 ); + const double t88 = piecewise_functor_5( t8, 0.0, t12, 0.0, -t5 - t45 ); + const double t89 = t88 * t4; + const double t92 = piecewise_functor_5( t12, 0.0, t8, 0.0, t5 - t52 ); + const double t93 = t92 * t39; + const double t95 = t89 + 0.1e1 + t16; + const double t98 = t92 * t4 + t21 + 0.1e1; + const double t101 = t63 * t95 / 0.6e1 + t68 * t98 / 0.6e1; + const double t106 = -0.3e1 * t60 * t101 + 0.6e1 * t77 * t101; + const double t107 = t59 * t106; + const double t111 = piecewise_functor_3( t3, 0.0, t18 * t107 / 0.4e1 - t18 * t93 / 0.4e1 - t89 * t40 / 0.4e1 - t51 / 0.4e1 ); + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t40 / 0.4e1 ); + vrho_a = t4 * t85 + eps; + vrho_b = t4 * t111 + eps; + + } + + +}; + +struct BuiltinEPC18_1 : detail::BuiltinKernelImpl< BuiltinEPC18_1 > { + + BuiltinEPC18_1( Spin p ) : + detail::BuiltinKernelImpl< BuiltinEPC18_1 >(p) { } + + virtual ~BuiltinEPC18_1() = default; + +}; + + + +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/epc18_2.hpp b/include/exchcxx/impl/builtin/kernels/epc18_2.hpp new file mode 100644 index 0000000..3ca66b8 --- /dev/null +++ b/include/exchcxx/impl/builtin/kernels/epc18_2.hpp @@ -0,0 +1,216 @@ +#pragma once +#include + +#include +#include +#include +#include + +#include + + + +namespace ExchCXX { + +template <> +struct kernel_traits< BuiltinEPC18_2 > : + public lda_screening_interface< BuiltinEPC18_2 > { + + static constexpr bool is_lda = true; + static constexpr bool is_gga = false; + static constexpr bool is_mgga = false; + static constexpr bool needs_laplacian = false; + static constexpr bool is_kedf = false; + static constexpr bool is_epc = true; + + static constexpr double dens_tol = 1e-24; + static constexpr double zeta_tol = 1e-15; + static constexpr double sigma_tol = 1.000000000000004e-32; + static constexpr double tau_tol = is_kedf ? 0.0 : 1e-20; + + static constexpr bool is_hyb = false; + static constexpr double exx_coeff = 0.0; + + static constexpr double a = 3.90; + static constexpr double b = 0.50; + static constexpr double c = 0.06; + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_unpolar_impl( double rho, double& eps ) { + + (void)(eps); + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t7 = 0.1e1 + t6; + const double t8 = t7 * t7; + const double t9 = t8 * rho; + const double t10 = b * t7; + const double t13 = c * t8; + const double t14 = rho * rho; + const double t17 = -0.4e1 * t10 * rho + 0.16e2 * t13 * t14 + a; + const double t18 = 0.1e1 / t17; + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t18 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_unpolar_impl( double rho, double& eps, double& vrho ) { + + + + const double t2 = rho / 0.2e1 <= dens_tol; + const double t3 = 0.1e1 <= zeta_tol; + const double t4 = zeta_tol - 0.1e1; + const double t6 = piecewise_functor_5( t3, t4, t3, -t4, 0.0 ); + const double t7 = 0.1e1 + t6; + const double t8 = t7 * t7; + const double t9 = t8 * rho; + const double t10 = b * t7; + const double t13 = c * t8; + const double t14 = rho * rho; + const double t17 = -0.4e1 * t10 * rho + 0.16e2 * t13 * t14 + a; + const double t18 = 0.1e1 / t17; + const double t22 = t17 * t17; + const double t23 = 0.1e1 / t22; + const double t27 = 0.32e2 * t13 * rho - 0.4e1 * t10; + const double t32 = piecewise_functor_3( t2, 0.0, t9 * t23 * t27 / 0.4e1 - t8 * t18 / 0.4e1 ); + + + eps = piecewise_functor_3( t2, 0.0, -t9 * t18 / 0.4e1 ); + vrho = rho * t32 + eps; + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_polar_impl( double rho_a, double rho_b, double& eps ) { + + (void)(eps); + constexpr double t23 = constants::m_cbrt_2; + constexpr double t24 = t23 * t23; + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t25 = safe_math::cbrt( t18 ); + const double t27 = t22 * t4; + const double t28 = safe_math::cbrt( t27 ); + const double t31 = t24 * t25 / 0.2e1 + t24 * t28 / 0.2e1; + const double t32 = t31 * t31; + const double t33 = t32 * t31; + const double t35 = t32 * t32; + const double t38 = c * t35 * t32 - b * t33 + a; + const double t39 = 0.1e1 / t38; + const double t40 = t22 * t39; + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t40 / 0.4e1 ); + + } + + BUILTIN_KERNEL_EVAL_RETURN + eval_exc_vxc_polar_impl( double rho_a, double rho_b, double& eps, double& vrho_a, double& vrho_b ) { + + constexpr double t23 = constants::m_cbrt_2; + constexpr double t24 = t23 * t23; + + + const double t3 = rho_a <= dens_tol && rho_b <= dens_tol; + const double t4 = rho_a + rho_b; + const double t5 = 0.1e1 / t4; + const double t8 = 0.2e1 * rho_a * t5 <= zeta_tol; + const double t9 = zeta_tol - 0.1e1; + const double t12 = 0.2e1 * rho_b * t5 <= zeta_tol; + const double t13 = -t9; + const double t14 = rho_a - rho_b; + const double t16 = piecewise_functor_5( t8, t9, t12, t13, t14 * t5 ); + const double t17 = 0.1e1 + t16; + const double t18 = t17 * t4; + const double t19 = -t14; + const double t21 = piecewise_functor_5( t12, t9, t8, t13, t19 * t5 ); + const double t22 = 0.1e1 + t21; + const double t25 = safe_math::cbrt( t18 ); + const double t27 = t22 * t4; + const double t28 = safe_math::cbrt( t27 ); + const double t31 = t24 * t25 / 0.2e1 + t24 * t28 / 0.2e1; + const double t32 = t31 * t31; + const double t33 = t32 * t31; + const double t35 = t32 * t32; + const double t38 = c * t35 * t32 - b * t33 + a; + const double t39 = 0.1e1 / t38; + const double t40 = t22 * t39; + const double t43 = t4 * t4; + const double t44 = 0.1e1 / t43; + const double t45 = t14 * t44; + const double t47 = piecewise_functor_5( t8, 0.0, t12, 0.0, t5 - t45 ); + const double t48 = t47 * t4; + const double t50 = t17 * t22; + const double t51 = t50 * t39; + const double t52 = t19 * t44; + const double t54 = piecewise_functor_5( t12, 0.0, t8, 0.0, -t5 - t52 ); + const double t55 = t54 * t39; + const double t57 = t38 * t38; + const double t58 = 0.1e1 / t57; + const double t59 = t22 * t58; + const double t60 = b * t32; + const double t61 = t25 * t25; + const double t63 = t24 / t61; + const double t64 = t48 + 0.1e1 + t16; + const double t66 = t28 * t28; + const double t68 = t24 / t66; + const double t70 = t54 * t4 + t21 + 0.1e1; + const double t73 = t63 * t64 / 0.6e1 + t68 * t70 / 0.6e1; + const double t77 = c * t35 * t31; + const double t80 = -0.3e1 * t60 * t73 + 0.6e1 * t77 * t73; + const double t81 = t59 * t80; + const double t85 = piecewise_functor_3( t3, 0.0, -t18 * t55 / 0.4e1 + t18 * t81 / 0.4e1 - t48 * t40 / 0.4e1 - t51 / 0.4e1 ); + const double t88 = piecewise_functor_5( t8, 0.0, t12, 0.0, -t5 - t45 ); + const double t89 = t88 * t4; + const double t92 = piecewise_functor_5( t12, 0.0, t8, 0.0, t5 - t52 ); + const double t93 = t92 * t39; + const double t95 = t89 + 0.1e1 + t16; + const double t98 = t92 * t4 + t21 + 0.1e1; + const double t101 = t63 * t95 / 0.6e1 + t68 * t98 / 0.6e1; + const double t106 = -0.3e1 * t60 * t101 + 0.6e1 * t77 * t101; + const double t107 = t59 * t106; + const double t111 = piecewise_functor_3( t3, 0.0, t18 * t107 / 0.4e1 - t18 * t93 / 0.4e1 - t89 * t40 / 0.4e1 - t51 / 0.4e1 ); + + + eps = piecewise_functor_3( t3, 0.0, -t18 * t40 / 0.4e1 ); + vrho_a = t4 * t85 + eps; + vrho_b = t4 * t111 + eps; + + } + + +}; + +struct BuiltinEPC18_2 : detail::BuiltinKernelImpl< BuiltinEPC18_2 > { + + BuiltinEPC18_2( Spin p ) : + detail::BuiltinKernelImpl< BuiltinEPC18_2 >(p) { } + + virtual ~BuiltinEPC18_2() = default; + +}; + + + +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/ft98_x.hpp b/include/exchcxx/impl/builtin/kernels/ft98_x.hpp index dec4949..2c49f5b 100644 --- a/include/exchcxx/impl/builtin/kernels/ft98_x.hpp +++ b/include/exchcxx/impl/builtin/kernels/ft98_x.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinFT98_X > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = true; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-12; static constexpr double zeta_tol = 1e-15; @@ -1029,4 +1030,4 @@ struct BuiltinFT98_X : detail::BuiltinKernelImpl< BuiltinFT98_X > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/lyp.hpp b/include/exchcxx/impl/builtin/kernels/lyp.hpp index b415966..ef3a028 100644 --- a/include/exchcxx/impl/builtin/kernels/lyp.hpp +++ b/include/exchcxx/impl/builtin/kernels/lyp.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinLYP > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-32; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/pbe0.hpp b/include/exchcxx/impl/builtin/kernels/pbe0.hpp index 4788d1b..80a3249 100644 --- a/include/exchcxx/impl/builtin/kernels/pbe0.hpp +++ b/include/exchcxx/impl/builtin/kernels/pbe0.hpp @@ -22,6 +22,7 @@ struct kernel_traits { static constexpr bool is_gga = true; static constexpr bool is_mgga = false; static constexpr double exx_coeff = 0.25; + static constexpr bool is_epc = false; BUILTIN_KERNEL_EVAL_RETURN diff --git a/include/exchcxx/impl/builtin/kernels/pbe_c.hpp b/include/exchcxx/impl/builtin/kernels/pbe_c.hpp index 306ed88..08049a7 100644 --- a/include/exchcxx/impl/builtin/kernels/pbe_c.hpp +++ b/include/exchcxx/impl/builtin/kernels/pbe_c.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPBE_C > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-12; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/pbe_x.hpp b/include/exchcxx/impl/builtin/kernels/pbe_x.hpp index 68794bb..a2cc279 100644 --- a/include/exchcxx/impl/builtin/kernels/pbe_x.hpp +++ b/include/exchcxx/impl/builtin/kernels/pbe_x.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPBE_X > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-32; static constexpr double zeta_tol = 1e-15; @@ -309,4 +310,4 @@ struct BuiltinPBE_X : detail::BuiltinKernelImpl< BuiltinPBE_X > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/pc07_k.hpp b/include/exchcxx/impl/builtin/kernels/pc07_k.hpp index db38b53..c991a8f 100644 --- a/include/exchcxx/impl/builtin/kernels/pc07_k.hpp +++ b/include/exchcxx/impl/builtin/kernels/pc07_k.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPC07_K > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = true; static constexpr bool is_kedf = true; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-15; static constexpr double zeta_tol = 1e-15; @@ -862,4 +863,4 @@ struct BuiltinPC07_K : detail::BuiltinKernelImpl< BuiltinPC07_K > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/pc07opt_k.hpp b/include/exchcxx/impl/builtin/kernels/pc07opt_k.hpp index 5531398..7cd3f7e 100644 --- a/include/exchcxx/impl/builtin/kernels/pc07opt_k.hpp +++ b/include/exchcxx/impl/builtin/kernels/pc07opt_k.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPC07OPT_K > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = true; static constexpr bool is_kedf = true; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-15; static constexpr double zeta_tol = 1e-15; @@ -862,4 +863,4 @@ struct BuiltinPC07OPT_K : detail::BuiltinKernelImpl< BuiltinPC07OPT_K > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/pw91_lda.hpp b/include/exchcxx/impl/builtin/kernels/pw91_lda.hpp index 5acd704..d54caa4 100644 --- a/include/exchcxx/impl/builtin/kernels/pw91_lda.hpp +++ b/include/exchcxx/impl/builtin/kernels/pw91_lda.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPW91_LDA > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/pw91_lda_mod.hpp b/include/exchcxx/impl/builtin/kernels/pw91_lda_mod.hpp index 54047fc..2c9e274 100644 --- a/include/exchcxx/impl/builtin/kernels/pw91_lda_mod.hpp +++ b/include/exchcxx/impl/builtin/kernels/pw91_lda_mod.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPW91_LDA_MOD > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/pw91_lda_rpa.hpp b/include/exchcxx/impl/builtin/kernels/pw91_lda_rpa.hpp index 3a8af70..19576d4 100644 --- a/include/exchcxx/impl/builtin/kernels/pw91_lda_rpa.hpp +++ b/include/exchcxx/impl/builtin/kernels/pw91_lda_rpa.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPW91_LDA_RPA > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/pz81.hpp b/include/exchcxx/impl/builtin/kernels/pz81.hpp index 66627ae..244098a 100644 --- a/include/exchcxx/impl/builtin/kernels/pz81.hpp +++ b/include/exchcxx/impl/builtin/kernels/pz81.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPZ81 > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/pz81_mod.hpp b/include/exchcxx/impl/builtin/kernels/pz81_mod.hpp index d55ecc8..d5a9141 100644 --- a/include/exchcxx/impl/builtin/kernels/pz81_mod.hpp +++ b/include/exchcxx/impl/builtin/kernels/pz81_mod.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinPZ81_MOD > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/r2scan_c.hpp b/include/exchcxx/impl/builtin/kernels/r2scan_c.hpp index 7a2b54b..1d14a3a 100644 --- a/include/exchcxx/impl/builtin/kernels/r2scan_c.hpp +++ b/include/exchcxx/impl/builtin/kernels/r2scan_c.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinR2SCAN_C > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-15; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/r2scan_x.hpp b/include/exchcxx/impl/builtin/kernels/r2scan_x.hpp index 4544b94..e1d0a9f 100644 --- a/include/exchcxx/impl/builtin/kernels/r2scan_x.hpp +++ b/include/exchcxx/impl/builtin/kernels/r2scan_x.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinR2SCAN_X > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-11; static constexpr double zeta_tol = 1e-15; @@ -910,4 +911,4 @@ struct BuiltinR2SCAN_X : detail::BuiltinKernelImpl< BuiltinR2SCAN_X > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/rev_pbe_x.hpp b/include/exchcxx/impl/builtin/kernels/rev_pbe_x.hpp index dcd0f8e..1a04f50 100644 --- a/include/exchcxx/impl/builtin/kernels/rev_pbe_x.hpp +++ b/include/exchcxx/impl/builtin/kernels/rev_pbe_x.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinRevPBE_X > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-32; static constexpr double zeta_tol = 1e-15; @@ -309,4 +310,4 @@ struct BuiltinRevPBE_X : detail::BuiltinKernelImpl< BuiltinRevPBE_X > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/scan_c.hpp b/include/exchcxx/impl/builtin/kernels/scan_c.hpp index 13c77db..c324e12 100644 --- a/include/exchcxx/impl/builtin/kernels/scan_c.hpp +++ b/include/exchcxx/impl/builtin/kernels/scan_c.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinSCAN_C > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-15; static constexpr double zeta_tol = 1e-15; @@ -160,7 +161,7 @@ struct kernel_traits< BuiltinSCAN_C > : constexpr double t4 = constants::m_cbrt_one_ov_pi; constexpr double t6 = constants::m_cbrt_4; constexpr double t40 = constants::m_cbrt_2; - constexpr double t61 = constants::m_pi * constants::m_pi; + constexpr double t61 = constants::m_pi_sq; constexpr double t116 = constants::m_cbrt_6; constexpr double t118 = constants::m_cbrt_pi_sq; constexpr double t5 = t2 * t4; @@ -1057,4 +1058,4 @@ struct BuiltinSCAN_C : detail::BuiltinKernelImpl< BuiltinSCAN_C > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/scan_x.hpp b/include/exchcxx/impl/builtin/kernels/scan_x.hpp index bdf0f4b..d1cb98e 100644 --- a/include/exchcxx/impl/builtin/kernels/scan_x.hpp +++ b/include/exchcxx/impl/builtin/kernels/scan_x.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinSCAN_X > : static constexpr bool is_mgga = true; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-15; static constexpr double zeta_tol = 1e-15; @@ -892,4 +893,4 @@ struct BuiltinSCAN_X : detail::BuiltinKernelImpl< BuiltinSCAN_X > { -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/slater_exchange.hpp b/include/exchcxx/impl/builtin/kernels/slater_exchange.hpp index bb354f3..fa330a3 100644 --- a/include/exchcxx/impl/builtin/kernels/slater_exchange.hpp +++ b/include/exchcxx/impl/builtin/kernels/slater_exchange.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinSlaterExchange > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; @@ -205,4 +206,4 @@ struct BuiltinSlaterExchange : detail::BuiltinKernelImpl< BuiltinSlaterExchange -} // namespace ExchCXX +} // namespace ExchCXX \ No newline at end of file diff --git a/include/exchcxx/impl/builtin/kernels/vwn3.hpp b/include/exchcxx/impl/builtin/kernels/vwn3.hpp index a14b787..f5b7bed 100644 --- a/include/exchcxx/impl/builtin/kernels/vwn3.hpp +++ b/include/exchcxx/impl/builtin/kernels/vwn3.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinVWN3 > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/builtin/kernels/vwn_rpa.hpp b/include/exchcxx/impl/builtin/kernels/vwn_rpa.hpp index cb569c1..417936c 100644 --- a/include/exchcxx/impl/builtin/kernels/vwn_rpa.hpp +++ b/include/exchcxx/impl/builtin/kernels/vwn_rpa.hpp @@ -21,6 +21,7 @@ struct kernel_traits< BuiltinVWN_RPA > : static constexpr bool is_mgga = false; static constexpr bool needs_laplacian = false; static constexpr bool is_kedf = false; + static constexpr bool is_epc = false; static constexpr double dens_tol = 1e-24; static constexpr double zeta_tol = 1e-15; diff --git a/include/exchcxx/impl/libxc.hpp b/include/exchcxx/impl/libxc.hpp index 1129ac5..bfcd17e 100644 --- a/include/exchcxx/impl/libxc.hpp +++ b/include/exchcxx/impl/libxc.hpp @@ -71,6 +71,7 @@ class LibxcKernelImpl : public XCKernelImpl { bool is_mgga_() const noexcept override; bool is_hyb_() const noexcept override; bool is_polarized_() const noexcept override; + bool is_epc_() const noexcept override; bool needs_laplacian_() const noexcept override; bool needs_tau_() const noexcept override; double hyb_exx_() const noexcept override; diff --git a/include/exchcxx/impl/xc_kernel.hpp b/include/exchcxx/impl/xc_kernel.hpp index 5b16b3c..99e283b 100644 --- a/include/exchcxx/impl/xc_kernel.hpp +++ b/include/exchcxx/impl/xc_kernel.hpp @@ -76,6 +76,7 @@ struct XCKernelImpl { bool is_mgga() const noexcept { return is_mgga_(); } bool is_hyb() const noexcept { return is_hyb_(); } bool is_polarized() const noexcept { return is_polarized_(); } + bool is_epc() const noexcept { return is_epc_(); } bool needs_laplacian() const noexcept { return needs_laplacian_(); } bool needs_tau() const noexcept { return needs_tau_(); } @@ -177,6 +178,7 @@ struct XCKernelImpl { virtual bool is_mgga_() const noexcept = 0; virtual bool is_hyb_() const noexcept = 0; virtual bool is_polarized_() const noexcept = 0; + virtual bool is_epc_() const noexcept = 0; virtual bool needs_laplacian_() const noexcept = 0; virtual bool needs_tau_() const noexcept = 0; diff --git a/include/exchcxx/xc_functional.hpp b/include/exchcxx/xc_functional.hpp index 28f3fec..470a4c3 100644 --- a/include/exchcxx/xc_functional.hpp +++ b/include/exchcxx/xc_functional.hpp @@ -155,6 +155,14 @@ class XCFunctional { ); } + inline bool is_epc() const { + throw_if_not_sane(); + return std::any_of( + kernels_.begin(), kernels_.end(), + [](const auto& x) { return x.second.is_epc(); } + ); + } + inline bool needs_laplacian() const { throw_if_not_sane(); return std::any_of( diff --git a/include/exchcxx/xc_kernel.hpp b/include/exchcxx/xc_kernel.hpp index f0c67fa..21efbd4 100644 --- a/include/exchcxx/xc_kernel.hpp +++ b/include/exchcxx/xc_kernel.hpp @@ -100,6 +100,7 @@ class XCKernel { bool is_mgga() const noexcept { return pimpl_->is_mgga(); } bool is_hyb() const noexcept { return pimpl_->is_hyb(); } bool is_polarized() const noexcept { return pimpl_->is_polarized(); } + bool is_epc() const noexcept { return pimpl_->is_epc(); } bool needs_laplacian() const noexcept { return pimpl_->needs_laplacian(); } bool needs_tau() const noexcept { return pimpl_->needs_tau(); } diff --git a/src/builtin_interface.cxx b/src/builtin_interface.cxx index 8083c86..a470cf0 100644 --- a/src/builtin_interface.cxx +++ b/src/builtin_interface.cxx @@ -46,6 +46,7 @@ #include #include #include +#include namespace ExchCXX { @@ -116,9 +117,21 @@ std::unique_ptr return std::make_unique( polar ); else if( kern == Kernel::PC07OPT_K ) return std::make_unique( polar ); - - - else + + else if( kern == Kernel::EPC17_1) { + EXCHCXX_BOOL_CHECK("EPC17_1 Needs to be Spin-Polarized!",polar==Spin::Polarized); + return std::make_unique( polar ); + } else if( kern == Kernel::EPC17_2) { + EXCHCXX_BOOL_CHECK("EPC17_2 Needs to be Spin-Polarized!",polar==Spin::Polarized); + return std::make_unique( polar ); + } else if( kern == Kernel::EPC18_1) { + EXCHCXX_BOOL_CHECK("EPC18_1 Needs to be Spin-Polarized!",polar==Spin::Polarized); + return std::make_unique( polar ); + } else if( kern == Kernel::EPC18_2) { + EXCHCXX_BOOL_CHECK("EPC18_2 Needs to be Spin-Polarized!",polar==Spin::Polarized); + return std::make_unique( polar ); + + } else throw std::runtime_error("Specified kernel does not have a builtin implementation"); @@ -168,6 +181,9 @@ bool BuiltinKernelInterface::needs_tau_() const noexcept { double BuiltinKernelInterface::hyb_exx_() const noexcept { return impl_->hyb_exx(); } +bool BuiltinKernelInterface::is_epc_() const noexcept { + return impl_->is_epc(); +} bool BuiltinKernelInterface::supports_inc_interface_() const noexcept { return true; diff --git a/src/builtin_kernel.cxx b/src/builtin_kernel.cxx index ab0e4cf..9b1c48f 100644 --- a/src/builtin_kernel.cxx +++ b/src/builtin_kernel.cxx @@ -654,6 +654,11 @@ MGGA_GENERATE_HOST_HELPERS( BuiltinSCANL_X ) MGGA_GENERATE_HOST_HELPERS( BuiltinR2SCANL_C ) MGGA_GENERATE_HOST_HELPERS( BuiltinR2SCANL_X ) +LDA_GENERATE_HOST_HELPERS( BuiltinEPC17_1 ) +LDA_GENERATE_HOST_HELPERS( BuiltinEPC17_2 ) +LDA_GENERATE_HOST_HELPERS( BuiltinEPC18_1 ) +LDA_GENERATE_HOST_HELPERS( BuiltinEPC18_2 ) + } } diff --git a/src/cuda/builtin.cu b/src/cuda/builtin.cu index b71f448..9c7d9e5 100644 --- a/src/cuda/builtin.cu +++ b/src/cuda/builtin.cu @@ -886,6 +886,12 @@ MGGA_GENERATE_DEVICE_HELPERS( BuiltinSCANL_C ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinSCANL_X ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinR2SCANL_C ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinR2SCANL_X ); + +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC17_1 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC17_2 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC18_1 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC18_2 ) + } } diff --git a/src/hip/builtin.hip b/src/hip/builtin.hip index 1eabbbe..0b71482 100644 --- a/src/hip/builtin.hip +++ b/src/hip/builtin.hip @@ -625,6 +625,11 @@ MGGA_GENERATE_DEVICE_HELPERS( BuiltinSCAN_C ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinR2SCAN_X ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinR2SCAN_C ); +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC17_1 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC17_2 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC18_1 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC18_2 ) + } } diff --git a/src/libxc.cxx b/src/libxc.cxx index 8e30825..af0c77b 100644 --- a/src/libxc.cxx +++ b/src/libxc.cxx @@ -230,6 +230,16 @@ bool LibxcKernelImpl::is_polarized_() const noexcept { return polar_ == XC_POLARIZED; } +bool LibxcKernelImpl::is_epc_() const noexcept { + int xcNumber = xc_info()->number; + return +#if XC_MAJOR_VERSION > 7 + xcNumber == XC_LDA_C_EPC17 or xcNumber == XC_LDA_C_EPC17_2 or xcNumber == XC_LDA_C_EPC18_1 or xcNumber == XC_LDA_C_EPC18_2; +#else + false; +#endif +} + double LibxcKernelImpl::hyb_exx_() const noexcept { return is_hyb_() ? xc_hyb_exx_coef(&kernel_) : 0.0; diff --git a/src/sycl/builtin_sycl.cxx b/src/sycl/builtin_sycl.cxx index 0d9f5a1..8e204b9 100644 --- a/src/sycl/builtin_sycl.cxx +++ b/src/sycl/builtin_sycl.cxx @@ -608,5 +608,10 @@ MGGA_GENERATE_DEVICE_HELPERS( BuiltinSCAN_C ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinR2SCAN_X ); MGGA_GENERATE_DEVICE_HELPERS( BuiltinR2SCAN_C ); +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC17_1 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC17_2 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC18_1 ) +LDA_GENERATE_DEVICE_HELPERS( BuiltinEPC18_2 ) + } } diff --git a/src/xc_functional.cxx b/src/xc_functional.cxx index 95bba76..4a45bfc 100644 --- a/src/xc_functional.cxx +++ b/src/xc_functional.cxx @@ -58,7 +58,11 @@ BidirectionalMap functional_map{ {"R2SCAN", Functional::R2SCAN}, {"R2SCANL", Functional::R2SCANL}, {"revPBE", Functional::revPBE}, - {"PBE0", Functional::PBE0}}}; + {"PBE0", Functional::PBE0}, + {"EPC17_1", Functional::EPC17_1}, + {"EPC17_2", Functional::EPC17_2}, + {"EPC18_1", Functional::EPC18_1}, + {"EPC18_2", Functional::EPC18_2}}}; std::ostream &operator<<(std::ostream &out, Functional functional) { out << functional_map.key(functional); @@ -121,6 +125,22 @@ std::vector< XCKernel > functional_factory( kerns = { XCKernel( backend, Kernel::PBE0, polar ) }; + else if( func == Functional::EPC17_1 ) + kerns = { + XCKernel( backend, Kernel::EPC17_1, polar ) + }; + else if( func == Functional::EPC17_2 ) + kerns = { + XCKernel( backend, Kernel::EPC17_2, polar ) + }; + else if( func == Functional::EPC18_1 ) + kerns = { + XCKernel( backend, Kernel::EPC18_1, polar ) + }; + else if( func == Functional::EPC18_2 ) + kerns = { + XCKernel( backend, Kernel::EPC18_2, polar ) + }; else { EXCHCXX_BOOL_CHECK( "Functional NYI Through Builtin Backend", false ); } diff --git a/src/xc_kernel.cxx b/src/xc_kernel.cxx index 37e530f..87bead4 100644 --- a/src/xc_kernel.cxx +++ b/src/xc_kernel.cxx @@ -75,7 +75,11 @@ BidirectionalMap kernel_map{ {"PW91_LDA", Kernel::PW91_LDA}, {"PW91_LDA_MOD", Kernel::PW91_LDA_MOD}, {"PW91_LDA_RPA", Kernel::PW91_LDA_RPA}, - {"B88", Kernel::B88}}}; + {"B88", Kernel::B88}, + {"EPC17_1", Kernel::EPC17_1}, + {"EPC17_2", Kernel::EPC17_2}, + {"EPC18_1", Kernel::EPC18_1}, + {"EPC18_2", Kernel::EPC18_2}}}; std::ostream& operator<<( std::ostream& out, Kernel kern ) { out << kernel_map.key(kern); diff --git a/test/reference_values.cxx b/test/reference_values.cxx index dbbbcbc..34c8e23 100644 --- a/test/reference_values.cxx +++ b/test/reference_values.cxx @@ -425,7 +425,6 @@ mgga_reference load_mgga_reference_values(ExchCXX::Kernel k, ExchCXX::Spin p, bo } - lda_reference gen_lda_reference_values(ExchCXX::Backend backend, ExchCXX::Kernel kern, ExchCXX::Spin polar) { diff --git a/test/ut_common.hpp b/test/ut_common.hpp index d031196..3da27e2 100644 --- a/test/ut_common.hpp +++ b/test/ut_common.hpp @@ -108,6 +108,13 @@ static std::vector mgga_kernels = { ExchCXX::Kernel::PC07OPT_K }; +static std::vector epc_lda_kernels = { + ExchCXX::Kernel::EPC17_1, + ExchCXX::Kernel::EPC17_2, + ExchCXX::Kernel::EPC18_1, + ExchCXX::Kernel::EPC18_2 +}; + static std::vector builtin_supported_kernels = { ExchCXX::Kernel::SlaterExchange, ExchCXX::Kernel::VWN3, @@ -138,7 +145,12 @@ static std::vector builtin_supported_kernels = { ExchCXX::Kernel::FT98_X, ExchCXX::Kernel::PC07_K, - ExchCXX::Kernel::PC07OPT_K + ExchCXX::Kernel::PC07OPT_K, + + ExchCXX::Kernel::EPC17_1, + ExchCXX::Kernel::EPC17_2, + ExchCXX::Kernel::EPC18_1, + ExchCXX::Kernel::EPC18_2 }; static std::vector unstable_small_kernels = { @@ -152,6 +164,11 @@ inline bool is_unstable_small(ExchCXX::Kernel kern) { kern) != unstable_small_kernels.end(); } +inline bool is_epc(ExchCXX::Kernel kern) { + return std::find(epc_lda_kernels.begin(), epc_lda_kernels.end(), + kern) != epc_lda_kernels.end(); +} + static constexpr std::array string_kernal_pairs = { std::pair("SlaterExchange", ExchCXX::Kernel::SlaterExchange), std::pair("PBE_X",ExchCXX::Kernel::PBE_X), @@ -176,7 +193,11 @@ static constexpr std::array string_kernal_pairs = { std::pair("PW91_LDA", ExchCXX::Kernel::PW91_LDA), std::pair("PW91_LDA_MOD", ExchCXX::Kernel::PW91_LDA_MOD), std::pair("PW91_LDA_RPA", ExchCXX::Kernel::PW91_LDA_RPA), - std::pair("B88", ExchCXX::Kernel::B88) + std::pair("B88", ExchCXX::Kernel::B88), + std::pair("EPC17_1", ExchCXX::Kernel::EPC17_1), + std::pair("EPC17_2", ExchCXX::Kernel::EPC17_2), + std::pair("EPC18_1", ExchCXX::Kernel::EPC18_1), + std::pair("EPC18_2", ExchCXX::Kernel::EPC18_2) }; static constexpr std::array string_functional_pairs = { @@ -187,5 +208,9 @@ static constexpr std::array string_functional_pairs = { std::pair("PBE", ExchCXX::Functional::PBE), std::pair("SCAN", ExchCXX::Functional::SCAN), std::pair("R2SCANL", ExchCXX::Functional::R2SCANL), - std::pair("PBE0", ExchCXX::Functional::PBE0) + std::pair("PBE0", ExchCXX::Functional::PBE0), + std::pair("EPC17_1", ExchCXX::Functional::EPC17_1), + std::pair("EPC17_2", ExchCXX::Functional::EPC17_2), + std::pair("EPC18_1", ExchCXX::Functional::EPC18_1), + std::pair("EPC18_2", ExchCXX::Functional::EPC18_2) }; diff --git a/test/xc_functional_test.cxx b/test/xc_functional_test.cxx index 290e2d9..b0a1288 100644 --- a/test/xc_functional_test.cxx +++ b/test/xc_functional_test.cxx @@ -191,6 +191,8 @@ void check_meta( Backend backend, Spin polar, Args&&... args ) { [](const auto& k) { return k.is_lda(); } ) and not should_be_mgga and not should_be_gga; + bool should_be_epc = std::any_of( kerns.begin(), kerns.end(), + [](const auto& k) { return k.is_epc(); } ); bool should_be_polarized = std::any_of( kerns.begin(), kerns.end(), [](const auto& k) { return k.is_polarized(); } ); @@ -204,6 +206,7 @@ void check_meta( Backend backend, Spin polar, Args&&... args ) { CHECK( func.is_lda() == should_be_lda ); CHECK( func.is_polarized() == should_be_polarized ); CHECK( func.is_hyb() == should_be_hyb ); + CHECK( func.is_epc() == should_be_epc ); CHECK( func.needs_laplacian() == should_need_lapl ); double total_exx = std::accumulate( kerns.begin(), kerns.end(), 0., diff --git a/test/xc_kernel_test.cxx b/test/xc_kernel_test.cxx index 86a4910..29d7d24 100644 --- a/test/xc_kernel_test.cxx +++ b/test/xc_kernel_test.cxx @@ -57,6 +57,7 @@ TEST_CASE( "XCKernel Metadata Validity", "[xc-kernel]" ) { auto mgga_tau_kernel_test = Kernel::SCAN_C; auto mgga_lapl_kernel_test = Kernel::R2SCANL_C; + auto epc_lda_kernel_test = Kernel::EPC17_2; Backend backend; @@ -305,6 +306,32 @@ TEST_CASE( "XCKernel Metadata Validity", "[xc-kernel]" ) { } + SECTION( "EPC LDA Polarized" ) { + + SECTION( "Builtin Backend" ) { backend = Backend::builtin; } + + XCKernel lda( backend, epc_lda_kernel_test, Spin::Polarized ); + + CHECK( lda.is_lda() ); + CHECK( lda.is_epc() ); + CHECK( lda.is_polarized() ); + CHECK( not lda.is_gga() ); + CHECK( not lda.is_mgga() ); + CHECK( not lda.is_hyb() ); + CHECK( not lda.needs_laplacian() ); + + CHECK( lda.rho_buffer_len( npts ) == 2*npts ); + CHECK( lda.sigma_buffer_len( npts ) == 0 ); + CHECK( lda.lapl_buffer_len( npts ) == 0 ); + CHECK( lda.tau_buffer_len( npts ) == 0 ); + CHECK( lda.exc_buffer_len( npts ) == npts ); + CHECK( lda.vrho_buffer_len( npts ) == 2*npts ); + CHECK( lda.vsigma_buffer_len( npts ) == 0 ); + CHECK( lda.vlapl_buffer_len( npts ) == 0 ); + CHECK( lda.vtau_buffer_len( npts ) == 0 ); + + } + } TEST_CASE( "XCKernel Metadata Correctness", "[xc-kernel]" ) { @@ -363,6 +390,22 @@ TEST_CASE( "XCKernel Metadata Correctness", "[xc-kernel]" ) { } + SECTION( "EPC LDA Kernels" ) { + + SECTION( "Builtin Backend" ) { backend = Backend::builtin; } + + for( const auto& kern : epc_lda_kernels ) { + XCKernel func( backend, kern, Spin::Polarized ); + auto exx = load_reference_exx( kern ); + + CHECK( func.is_lda() ); + CHECK( func.is_epc() ); + CHECK( exx == Approx( func.hyb_exx() ) ); + + if( std::abs(exx) > 0 ) CHECK( func.is_hyb() ); + } + + } } @@ -372,6 +415,7 @@ TEST_CASE( "XCKernel Metadata Correctness", "[xc-kernel]" ) { + void kernel_test( TestInterface interface, Backend backend, Kernel kern, Spin polar ) { @@ -728,6 +772,7 @@ void compare_libxc_builtin( TestInterface interface, EvalType evaltype, XCKernel func_libxc ( Backend::libxc, kern, polar ); XCKernel func_builtin( Backend::builtin, kern, polar ); + const int len_rho = func_libxc.rho_buffer_len( npts ); const int len_sigma = func_libxc.sigma_buffer_len( npts ); const int len_lapl = func_libxc.lapl_buffer_len( npts ); @@ -856,20 +901,25 @@ void compare_libxc_builtin( TestInterface interface, EvalType evaltype, TEST_CASE( "Builtin Corectness Test", "[xc-builtin]" ) { SECTION( "Unpolarized Regular Eval : EXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC, EvalType::Regular, kern, Spin::Unpolarized ); + } } SECTION( "Unpolarized Regular Eval : EXC + VXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC_VXC, EvalType::Regular, kern, Spin::Unpolarized ); + } } SECTION( "Unpolarized Small Eval : EXC" ) { for( auto kern : builtin_supported_kernels ) { if(is_unstable_small(kern)) continue; + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC, EvalType::Small, kern, Spin::Unpolarized ); } @@ -878,21 +928,26 @@ TEST_CASE( "Builtin Corectness Test", "[xc-builtin]" ) { SECTION( "Unpolarized Small Eval : EXC + VXC" ) { for( auto kern : builtin_supported_kernels ) { if(is_unstable_small(kern)) continue; + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC_VXC, EvalType::Small, kern, Spin::Unpolarized ); } } SECTION( "Unpolarized Zero Eval : EXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC, EvalType::Zero, kern, Spin::Unpolarized ); + } } SECTION( "Unpolarized Zero Eval : EXC + VXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC_VXC, EvalType::Zero, kern, Spin::Unpolarized ); + } } @@ -901,20 +956,25 @@ TEST_CASE( "Builtin Corectness Test", "[xc-builtin]" ) { SECTION( "Polarized Regular Eval : EXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC, EvalType::Regular, kern, Spin::Polarized ); + } } SECTION( "Polarized Regular Eval : EXC + VXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC_VXC, EvalType::Regular, kern, Spin::Polarized ); + } } SECTION( "Polarized Small Eval : EXC" ) { for( auto kern : builtin_supported_kernels ) { if(is_unstable_small(kern)) continue; + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC, EvalType::Small, kern, Spin::Polarized ); } @@ -923,21 +983,26 @@ TEST_CASE( "Builtin Corectness Test", "[xc-builtin]" ) { SECTION( "Polarized Small Eval : EXC + VXC" ) { for( auto kern : builtin_supported_kernels ) { if(is_unstable_small(kern)) continue; + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC_VXC, EvalType::Small, kern, Spin::Polarized ); } } SECTION( "Polarized Zero Eval : EXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC, EvalType::Zero, kern, Spin::Polarized ); + } } SECTION( "Polarized Zero Eval : EXC + VXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; compare_libxc_builtin( TestInterface::EXC_VXC, EvalType::Zero, kern, Spin::Polarized ); + } } } @@ -946,15 +1011,19 @@ TEST_CASE( "Builtin Corectness Test", "[xc-builtin]" ) { TEST_CASE( "Scale and Increment Interface", "[xc-inc]" ) { SECTION( "Builtin Unpolarized EXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; kernel_test( TestInterface::EXC_INC, Backend::builtin, kern, Spin::Unpolarized ); + } } SECTION( "Builtin Unpolarized EXC + VXC" ) { - for( auto kern : builtin_supported_kernels ) + for( auto kern : builtin_supported_kernels ) { + if(is_epc(kern)) continue; kernel_test( TestInterface::EXC_VXC_INC, Backend::builtin, kern, Spin::Unpolarized ); + } } SECTION( "Builtin Polarized EXC" ) {