Add M06-2X and PKZB Implementations (#41)

* Add EXCHCXX_ENABLE_LIBXC option to disable LibXC backend on request

* Bump to Libxc 7.0.0, update impls and scripts, add impls of M062X_{X,C}

* Added M062X as a Functional implementation

* Add PKZB implementation

* remove debug statement

* Add XC Functional impl of PKZB

Author: wavefunction91

CMakeLists.txt

@@ -57,7 +57,7 @@ endif()
 if(EXCHCXX_ENABLE_LIBXC)
 
 ## Find LibXC
-find_package( Libxc 6.2.0 CONFIG QUIET )
+find_package( Libxc 7.0.0 CONFIG QUIET )
 
 if( ${Libxc_FOUND} )
 
@@ -71,10 +71,10 @@ else()
   FetchContent_Declare(
     libxc
     GIT_REPOSITORY https://gitlab.com/libxc/libxc.git
-    GIT_TAG 6.2.0  # v6.2.0
+    GIT_TAG 7.0.0
     PATCH_COMMAND sed -i -e "s/p->info->family != XC_KINETIC/p->info->kind != XC_KINETIC/g" src/work_mgga_inc.c
   )
-  set( Libxc_VERSION 6.2.0 )
+  set( Libxc_VERSION 7.0.0 )
 
   set( OLD_BUILD_TESTING ${BUILD_TESTING} )
   set( BUILD_TESTING OFF CACHE BOOL "" FORCE )

bin/generate_primitive_kernels.py

@@ -490,7 +490,8 @@ def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff
 
 #libxc_prefix = '/Users/meji656/Projects/libxc/src/maple2c/' 
 #kernel_prefix = '/Users/meji656/Projects/ExchCXX/include/exchcxx/impl/builtin/kernels/'
-libxc_prefix = '/home/dbwy/Software/Chemistry/libxc/6.2.0/libxc-6.2.0/src/maple2c/' 
+#libxc_prefix = '/home/dbwy/Software/Chemistry/libxc/6.2.0/libxc-6.2.0/src/maple2c/' 
+libxc_prefix = '/home/davidwillia/Development/ExchCXX/libxc-7.0.0/src/maple2c/'
 kernel_prefix = 'include/exchcxx/impl/builtin/kernels/'
 gen_table = {
 
@@ -688,6 +689,71 @@ def __init__( self, local_name, libxc_file, ofname, xc_type, dens_tol, exx_coeff
     {'eta': '0.001'}
     ),
 
+  'M062X_X' : GenMetaData( 'BuiltinM062X_X',
+    libxc_prefix + 'mgga_exc/hyb_mgga_x_m05.c',
+    kernel_prefix + 'm06_2x_x.hpp',
+    'MGGA', 1e-15, 0.54,
+    {'a_0' : '0.46',
+     'a_1' : '-0.2206052',
+     'a_2' : '-9.431788e-02',  
+     'a_3' : '2.164494e+00', 
+     'a_4' : '-2.556466e+00', 
+     'a_5' : '-1.422133e+01',
+     'a_6' : '1.555044e+01',  
+     'a_7' : '3.598078e+01', 
+     'a_8' : '-2.722754e+01', 
+     'a_9' : '-3.924093e+01',  
+     'a_10': '1.522808e+01',  
+     'a_11': '1.522227e+01',
+     'csi_HF': '1.0',
+     'cx':     '0.54'}
+    ), 
+
+  'M062X_C' : GenMetaData( 'BuiltinM062X_C',
+    libxc_prefix + 'mgga_exc/mgga_c_m06l.c',
+    kernel_prefix + 'm06_2x_c.hpp',
+    'MGGA', 1e-12, 0.0,
+    {
+      'gamma_ss':  '0.06', 
+      'gamma_ab':  '0.0031', 
+      'alpha_ss':  '0.00515088', 
+      'alpha_ab':  '0.00304966',
+      'css_0':     '3.097855e-01', 
+      'css_1':    '-5.528642e+00',  
+      'css_2':     '1.347420e+01', 
+      'css_3':    '-3.213623e+01',  
+      'css_4':     '2.846742e+01',
+      'cab_0':     '8.833596e-01',  
+      'cab_1':     '3.357972e+01', 
+      'cab_2':    '-7.043548e+01',  
+      'cab_3':     '4.978271e+01', 
+      'cab_4':    '-1.852891e+01',
+      'dss_0':     '6.902145e-01',  
+      'dss_1':     '9.847204e-02',  
+      'dss_2':     '2.214797e-01', 
+      'dss_3':    '-1.968264e-03', 
+      'dss_4':    '-6.775479e-03',  
+      'dss_5':     '0.000000e+00',
+      'dab_0':     '1.166404e-01', 
+      'dab_1':    '-9.120847e-02', 
+      'dab_2':    '-6.726189e-02',  
+      'dab_3':     '6.720580e-05',  
+      'dab_4':     '8.448011e-04',  
+      'dab_5':     '0.000000e+00',
+      'Fermi_D_cnst': '1e-10'
+    }
+    ), 
+
+  'PKZB_X' : GenMetaData( 'BuiltinPKZB_X',
+    libxc_prefix + 'mgga_exc/mgga_x_pkzb.c',
+    kernel_prefix + 'pkzb_x.hpp',
+    'MGGA', 1e-15, 0., {}),
+
+  'PKZB_C' : GenMetaData( 'BuiltinPKZB_C',
+    libxc_prefix + 'mgga_exc/mgga_c_pkzb.c',
+    kernel_prefix + 'pkzb_c.hpp',
+    'MGGA', 1e-13, 0., {}),
+
   'FT98' : GenMetaData( 'BuiltinFT98_X',
     libxc_prefix + 'mgga_exc/mgga_x_ft98.c',
     kernel_prefix + 'ft98_x.hpp',

include/exchcxx/enums/functionals.hpp

@@ -60,6 +60,8 @@ enum class Functional {
   SCAN,
   R2SCAN,
   R2SCANL,
+  M062X,
+  PKZB,
   EPC17_1,
   EPC17_2,
   EPC18_1,

include/exchcxx/enums/kernels.hpp

@@ -93,6 +93,11 @@ enum class Kernel {
   R2SCANL_C,
   R2SCANL_X,
   FT98_X,
+  M062X_X,
+  M062X_C,
+  PKZB_X,
+  PKZB_C,
+  
 
   // KEDFs
   PC07_K,

include/exchcxx/impl/builtin/fwd.hpp

@@ -76,6 +76,10 @@ struct BuiltinSCAN_X;
 struct BuiltinSCAN_C;
 struct BuiltinR2SCAN_X;
 struct BuiltinR2SCAN_C;
+struct BuiltinM062X_X;
+struct BuiltinM062X_C;
+struct BuiltinPKZB_X;
+struct BuiltinPKZB_C;
 struct BuiltinFT98_X;
 
 struct BuiltinPC07_K;

include/exchcxx/impl/builtin/kernels.hpp

@@ -77,6 +77,10 @@
 #include <exchcxx/impl/builtin/kernels/scanl_x.hpp>
 #include <exchcxx/impl/builtin/kernels/r2scanl_c.hpp>
 #include <exchcxx/impl/builtin/kernels/r2scanl_x.hpp>
+#include <exchcxx/impl/builtin/kernels/m06_2x_x.hpp>
+#include <exchcxx/impl/builtin/kernels/m06_2x_c.hpp>
+#include <exchcxx/impl/builtin/kernels/pkzb_x.hpp>
+#include <exchcxx/impl/builtin/kernels/pkzb_c.hpp>
 
 #include <exchcxx/impl/builtin/kernels/pc07_k.hpp>
 #include <exchcxx/impl/builtin/kernels/pc07opt_k.hpp>

include/exchcxx/impl/builtin/kernels/ft98_x.hpp

@@ -77,7 +77,7 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t31 = safe_math::sqrt( t30 );
     const double t32 = a * t31;
     const double t33 = b1 * sigma;
-    const double t35 = t28 * t33 + 0.1e1;
+    const double t35 = t33 * t28 + 0.1e1;
     const double t36 = pow_1_4( t35 );
     const double t37 = t36 * t36;
     const double t38 = t37 * t36;
@@ -130,15 +130,15 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t109 = safe_math::sqrt( t108 );
     const double t110 = t106 + t109;
     const double t112 = piecewise_functor_5( t79, -0.4e1 * t67 + 0.4e1 * t74 + 0.2e1 * b2 - 0.1e1 / t76 / 0.2e1 + t87 / 0.8e1 - t91 / 0.16e2, t97, 0.1e1 - t68 + t75 + b2 + t85 / 0.2e1 - t89 / 0.8e1 + t100 / 0.16e2 - 0.5e1 / 0.128e3 * t102, 0.1e1 / t110 );
-    const double t114 = t112 * t60 + 0.1e1;
+    const double t114 = t60 * t112 + 0.1e1;
     const double t116 = t115 * t60;
-    const double t118 = t112 * t116 + 0.1e1;
+    const double t118 = t116 * t112 + 0.1e1;
     const double t119 = t118 * t118;
     const double t120 = t119 * t118;
     const double t121 = 0.1e1 / t120;
     const double t122 = t114 * t121;
     const double t123 = t122 * t49;
-    const double t125 = t123 * t56 + t40 * t42 + 0.1e1;
+    const double t125 = t56 * t123 + t40 * t42 + 0.1e1;
     const double t133 = b * sigma;
     const double t137 = 0.1e1 + 0.81e2 / 0.4e1 * t130 * t131 * t133 * t28;
     const double t138 = 0.1e1 / t137;
@@ -192,7 +192,7 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t31 = safe_math::sqrt( t30 );
     const double t32 = a * t31;
     const double t33 = b1 * sigma;
-    const double t35 = t28 * t33 + 0.1e1;
+    const double t35 = t33 * t28 + 0.1e1;
     const double t36 = pow_1_4( t35 );
     const double t37 = t36 * t36;
     const double t38 = t37 * t36;
@@ -245,15 +245,15 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t109 = safe_math::sqrt( t108 );
     const double t110 = t106 + t109;
     const double t112 = piecewise_functor_5( t79, -0.4e1 * t67 + 0.4e1 * t74 + 0.2e1 * b2 - 0.1e1 / t76 / 0.2e1 + t87 / 0.8e1 - t91 / 0.16e2, t97, 0.1e1 - t68 + t75 + b2 + t85 / 0.2e1 - t89 / 0.8e1 + t100 / 0.16e2 - 0.5e1 / 0.128e3 * t102, 0.1e1 / t110 );
-    const double t114 = t112 * t60 + 0.1e1;
+    const double t114 = t60 * t112 + 0.1e1;
     const double t116 = t115 * t60;
-    const double t118 = t112 * t116 + 0.1e1;
+    const double t118 = t116 * t112 + 0.1e1;
     const double t119 = t118 * t118;
     const double t120 = t119 * t118;
     const double t121 = 0.1e1 / t120;
     const double t122 = t114 * t121;
     const double t123 = t122 * t49;
-    const double t125 = t123 * t56 + t40 * t42 + 0.1e1;
+    const double t125 = t56 * t123 + t40 * t42 + 0.1e1;
     const double t133 = b * sigma;
     const double t137 = 0.1e1 + 0.81e2 / 0.4e1 * t130 * t131 * t133 * t28;
     const double t138 = 0.1e1 / t137;
@@ -296,7 +296,7 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t222 = piecewise_functor_3( t105, t201, 0.0 );
     const double t223 = 0.1e1 / t109;
     const double t224 = t223 * t106;
-    const double t226 = t222 * t224 + t222;
+    const double t226 = t224 * t222 + t222;
     const double t228 = piecewise_functor_5( t79, 0.64e2 / 0.3e1 * t192 - 0.4e2 / 0.3e1 * t196 + t198 * t201 / 0.2e1 - 0.3e1 / 0.8e1 * t204 * t201 + 0.5e1 / 0.16e2 * t207 * t201, t97, t199 - t200 + t76 * t201 - t86 * t201 / 0.2e1 + 0.3e1 / 0.8e1 * t90 * t201 - 0.5e1 / 0.16e2 * t216 * t201, -t221 * t226 );
     const double t229 = t228 * t121;
     const double t230 = t229 * t49;
@@ -315,15 +315,15 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t249 = 0.1e1 / t248;
     const double t251 = t125 * t249 * t130;
     const double t253 = t252 * t172;
-    const double t256 = t138 * t246 + 0.54e2 * t251 * t253;
+    const double t256 = t246 * t138 + 0.54e2 * t251 * t253;
     const double t261 = piecewise_functor_3( t3, 0.0, -t7 * t146 * t140 / 0.8e1 - 0.3e1 / 0.16e2 * t150 * t152 * t256 );
     const double t264 = sigma * t22;
     const double t272 = t39 * t23;
     const double t275 = t53 * t23;
     const double t276 = t275 * t27;
     const double t278 = t184 * t23;
     const double t279 = t278 * t27;
-    const double t283 = b * ( 0.2e1 * t175 * t276 - 0.2e1 * t279 * t50 );
+    const double t283 = b * ( 0.2e1 * t175 * t276 - 0.2e1 * t50 * t279 );
     const double t285 = t264 * t66;
     const double t287 = t198 * sigma;
     const double t288 = t22 * t66;
@@ -502,15 +502,15 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t153 = safe_math::cbrt( t151 );
     const double t155 = piecewise_functor_3( t152, t23, t153 * t151 );
     const double t156 = t155 * t27;
-    const double t157 = a1 * sigma_bb;
     const double t158 = rho_b * rho_b;
     const double t159 = safe_math::cbrt( rho_b );
     const double t160 = t159 * t159;
     const double t162 = 0.1e1 / t160 / t158;
-    const double t164 = t157 * t162 + 0.1e1;
+    const double t164 = a1 * sigma_bb * t162 + 0.1e1;
     const double t165 = safe_math::sqrt( t164 );
     const double t166 = a * t165;
-    const double t169 = b1 * sigma_bb * t162 + 0.1e1;
+    const double t167 = b1 * sigma_bb;
+    const double t169 = t167 * t162 + 0.1e1;
     const double t170 = pow_1_4( t169 );
     const double t171 = t170 * t170;
     const double t172 = t171 * t170;
@@ -694,15 +694,15 @@ struct kernel_traits< BuiltinFT98_X > :
     const double t153 = safe_math::cbrt( t151 );
     const double t155 = piecewise_functor_3( t152, t23, t153 * t151 );
     const double t156 = t155 * t27;
-    const double t157 = a1 * sigma_bb;
     const double t158 = rho_b * rho_b;
     const double t159 = safe_math::cbrt( rho_b );
     const double t160 = t159 * t159;
     const double t162 = 0.1e1 / t160 / t158;
-    const double t164 = t157 * t162 + 0.1e1;
+    const double t164 = a1 * sigma_bb * t162 + 0.1e1;
     const double t165 = safe_math::sqrt( t164 );
     const double t166 = a * t165;
-    const double t169 = b1 * sigma_bb * t162 + 0.1e1;
+    const double t167 = b1 * sigma_bb;
+    const double t169 = t167 * t162 + 0.1e1;
     const double t170 = pow_1_4( t169 );
     const double t171 = t170 * t170;
     const double t172 = t171 * t170;