deepmodeling · Growl1234 · May 24, 2026 · May 26, 2026 · May 26, 2026 · May 26, 2026
diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml
@@ -7,6 +7,10 @@ concurrency:
   group: ${{ github.workflow }}-${{ github.ref }}
   cancel-in-progress: true
 
+defaults:
+  run:
+    shell: bash
+
 jobs:
   test:
     name: Test
@@ -32,12 +36,20 @@ jobs:
       - name: Install CI tools
         run: |
           sudo apt-get update
-          sudo apt-get install -y ccache ca-certificates python-is-python3 python3-pip
+          sudo apt-get install -y gfortran ccache ca-certificates python-is-python3 python3-pip
           sudo pip install clang-format clang-tidy
 
+      - name: Install dftd4 from toolchain
+        run: |
+          cd toolchain
+          ./install_abacus_toolchain_new.sh --with-dftd4=install --dry-run -j8
+          ./scripts/stage4/install_dftd4.sh
+          cd ..
+
       - name: Configure
         run: |
-          cmake -B build -DBUILD_TESTING=ON -DENABLE_MLALGO=ON -DENABLE_LIBXC=ON -DENABLE_LIBRI=ON -DENABLE_GOOGLEBENCH=ON -DENABLE_RAPIDJSON=ON  -DCMAKE_EXPORT_COMPILE_COMMANDS=1 -DENABLE_FLOAT_FFTW=ON
+          source toolchain/install/setup
+          cmake -B build -DBUILD_TESTING=ON -DENABLE_MLALGO=ON -DENABLE_LIBXC=ON -DENABLE_LIBRI=ON -DENABLE_GOOGLEBENCH=ON -DENABLE_RAPIDJSON=ON  -DCMAKE_EXPORT_COMPILE_COMMANDS=1 -DENABLE_FLOAT_FFTW=ON -DENABLE_DFTD4=ON
 
 # Temporarily removed because no one maintains this now.
 # And it will break the CI test workflow.

diff --git a/CMakeLists.txt b/CMakeLists.txt
@@ -34,6 +34,7 @@ option(USE_ELPA "Enable ELPA for LCAO" ON)
 option(ENABLE_LIBRI "Enable LibRI for hybrid functional" OFF)
 option(ENABLE_LIBCOMM "Enable LibComm" OFF)
 option(ENABLE_PEXSI "Enable PEXSI for LCAO" OFF)
+option(ENABLE_DFTD4 "Enable DFT-D4 dispersion correction" OFF)
 
 option(BUILD_TESTING "Build unittests" OFF)
 option(DEBUG_INFO "Print message to debug" OFF)
@@ -236,6 +237,18 @@ if(ENABLE_COVERAGE)
   add_coverage(${ABACUS_BIN_NAME})
 endif()
 
+if(ENABLE_DFTD4)
+  # DFTD4 requires enabling C and Fortran to work
+  enable_language(C)
+  enable_language(Fortran)
+  # Avoid custom-lapack fallback when resolving DFT-D4 dependencies
+  find_package(BLAS REQUIRED)
+  find_package(LAPACK REQUIRED)
+  find_package(dftd4 4.2.0 REQUIRED)
+  if(NOT TARGET dftd4::dftd4)
+    message(FATAL_ERROR "DFT-D4 was found, but target dftd4::dftd4 is missing.")
+  endif()
+endif()
 
 macro(set_if_higher VARIABLE VALUE)
   if(${VARIABLE} LESS ${VALUE})

diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md
@@ -3640,11 +3640,21 @@
   - d2: Grimme's D2 dispersion correction method
   - d3_0: Grimme's DFT-D3(0) dispersion correction method (zero-damping)
   - d3_bj: Grimme's DFTD3(BJ) dispersion correction method (BJ-damping)
+  - d4: Grimme's DFT-D4 dispersion correction method using the external DFT-D4 library
   - none: no vdW correction
 
   > Note: ABACUS supports automatic setting of DFT-D3 parameters for common functionals. To benefit from this feature, please specify the parameter dft_functional explicitly, otherwise the autoset procedure will crash. If not satisfied with the built-in parameters, any manual setting on vdw_s6, vdw_s8, vdw_a1 and vdw_a2 will overwrite the automatic values.
+
+  > Note: DFT-D4 support requires ABACUS to be configured with ENABLE_DFTD4=ON and a CMake-installed dftd4 library exporting dftd4-config.cmake. DFT-D4 damping parameters are loaded from the external library.
 - **Default**: none
 
+### vdw_d4_xc
+
+- **Type**: String
+- **Availability**: *vdw_method is set to d4*
+- **Description**: Functional name passed to the DFT-D4 library to load its internal damping parameters. If set to default, ABACUS infers the functional name from dft_functional or pseudopotential metadata.
+- **Default**: default
+
 ### vdw_s6
 
 - **Type**: String
@@ -3737,7 +3747,7 @@
 
 - **Type**: String
 - **Availability**: *vdw_cutoff_type is set to radius*
-- **Description**: Defines the radius of the cutoff sphere when vdw_cutoff_type is set to radius. The default values depend on the chosen vdw_method.
+- **Description**: Defines the cutoff radius when vdw_cutoff_type is set to radius. The default values depend on the chosen vdw_method. For DFT-D4, this controls the two-body dispersion cutoff, while the three-body cutoff is internally limited to the DFT-D4 default value of 40 Bohr.
 - **Unit**: defined by vdw_radius_unit (default Bohr)
 
 ### vdw_radius_unit
@@ -3759,8 +3769,8 @@
 ### vdw_cn_thr
 
 - **Type**: Real
-- **Availability**: *vdw_method is set to d3_0 or d3_bj*
-- **Description**: The cutoff radius when calculating coordination numbers.
+- **Availability**: *vdw_method is set to d3_0, d3_bj, or d4*
+- **Description**: The cutoff radius when calculating coordination numbers. The default is 40 Bohr for DFT-D3 and 30 Bohr for DFT-D4.
 - **Default**: 40
 - **Unit**: defined by vdw_cn_thr_unit (default: Bohr)
 

diff --git a/docs/advanced/install.md b/docs/advanced/install.md
@@ -77,6 +77,24 @@ These two libraries are added as submodules in the [deps](https://github.com/dee
 
 If you prefer using manually downloaded libraries, provide `-DLIBRI_DIR=${path to your LibRI folder} -DLIBCOMM_DIR=${path to your LibComm folder}`.
 
+
+## Build with DFT-D4 support
+
+ABACUS can use the external [DFT-D4](https://github.com/dftd4/dftd4) library for Grimme's DFT-D4 dispersion correction. DFT-D4 support is optional and disabled by default.
+
+DFT-D4 must be built and installed with CMake, so that ABACUS can locate it through the exported `dftd4-config.cmake` package. Meson-built installations, including the Conda package of `dftd4`, are not supported unless they provide a compatible CMake package file.
+
+To build ABACUS with DFT-D4 support, pass `-DENABLE_DFTD4=ON` to CMake and provide `CMAKE_PREFIX_PATH` environment variable.
+
+In the input file, enable DFT-D4 with:
+
+```text
+vdw_method d4
+vdw_d4_xc pbe
+```
+
+If `vdw_d4_xc` is set to `default`, ABACUS will infer the functional name from `dft_functional` or pseudopotential metadata and pass it to the DFT-D4 library.
+
 ## Build Unit Tests
 
 To build tests for ABACUS, define `BUILD_TESTING` flag. You can also specify path to local installation of [Googletest](https://github.com/google/googletest) by setting `GTEST_DIR` flags. If not found in local, the configuration process will try to download it automatically.

diff --git a/docs/parameters.yaml b/docs/parameters.yaml
@@ -3738,12 +3738,23 @@ parameters:
       * d2: Grimme's D2 dispersion correction method
       * d3_0: Grimme's DFT-D3(0) dispersion correction method (zero-damping)
       * d3_bj: Grimme's DFTD3(BJ) dispersion correction method (BJ-damping)
+      * d4: Grimme's DFT-D4 dispersion correction method using the external DFT-D4 library
       * none: no vdW correction
 
       [NOTE] ABACUS supports automatic setting of DFT-D3 parameters for common functionals. To benefit from this feature, please specify the parameter dft_functional explicitly, otherwise the autoset procedure will crash. If not satisfied with the built-in parameters, any manual setting on vdw_s6, vdw_s8, vdw_a1 and vdw_a2 will overwrite the automatic values.
+
+      [NOTE] DFT-D4 support requires ABACUS to be configured with ENABLE_DFTD4=ON and a CMake-installed dftd4 library exporting `dftd4-config.cmake`. DFT-D4 damping parameters are loaded from the external library.
     default_value: none
     unit: ""
     availability: ""
+  - name: vdw_d4_xc
+    category: vdW correction
+    type: String
+    description: |
+      Functional name passed to the DFT-D4 library to load its internal damping parameters. If set to default, ABACUS infers the functional name from dft_functional or pseudopotential metadata.
+    default_value: default
+    unit: ""
+    availability: vdw_method is set to d4
   - name: vdw_s6
     category: vdW correction
     type: String
@@ -3852,7 +3863,7 @@ parameters:
     category: vdW correction
     type: String
     description: |
-      Defines the radius of the cutoff sphere when vdw_cutoff_type is set to radius. The default values depend on the chosen vdw_method.
+      Defines the cutoff radius when vdw_cutoff_type is set to radius. The default values depend on the chosen vdw_method. For DFT-D4, this controls the two-body dispersion cutoff, while the three-body cutoff is internally limited to the DFT-D4 default value of 40 Bohr.
     default_value: ""
     unit: defined by vdw_radius_unit (default Bohr)
     availability: vdw_cutoff_type is set to radius
@@ -3878,10 +3889,10 @@ parameters:
     category: vdW correction
     type: Real
     description: |
-      The cutoff radius when calculating coordination numbers.
+      The cutoff radius when calculating coordination numbers. The default is 40 Bohr for DFT-D3 and 30 Bohr for DFT-D4.
     default_value: "40"
     unit: "defined by vdw_cn_thr_unit (default: Bohr)"
-    availability: vdw_method is set to d3_0 or d3_bj
+    availability: vdw_method is set to d3_0, d3_bj, or d4
   - name: vdw_cn_thr_unit
     category: vdW correction
     type: String

diff --git a/source/source_estate/elecstate_print.cpp b/source/source_estate/elecstate_print.cpp
@@ -249,6 +249,11 @@ void print_etot(const Magnetism& magnet,
             titles.push_back("E_vdwD3");
             energies_Ry.push_back(elec.f_en.evdw);
         }
+        else if (vdw_method == "d4")
+        {
+            titles.push_back("E_vdwD4");
+            energies_Ry.push_back(elec.f_en.evdw);
+        }
 
         // mohan add 20251108
 		if (PARAM.inp.dft_plus_u)

diff --git a/source/source_hamilt/module_vdw/CMakeLists.txt b/source/source_hamilt/module_vdw/CMakeLists.txt
@@ -1,16 +1,25 @@
-add_library(
-        vdw
-        OBJECT
-        vdwd2_parameters.cpp
-        vdwd3_parameters_tab.cpp
-        vdwd3_parameters.cpp
-        vdwd2.cpp
-        vdwd3.cpp
-        vdwd3_autoset_xcname.cpp
-        vdwd3_autoset_xcparam.cpp
-        vdw.cpp
+set(vdw_sources
+    vdwd2_parameters.cpp
+    vdwd3_parameters_tab.cpp
+    vdwd3_parameters.cpp
+    vdwd2.cpp
+    vdwd3.cpp
+    vdwd3_autoset_xcname.cpp
+    vdwd3_autoset_xcparam.cpp
+    vdw.cpp
 )
 
+if(ENABLE_DFTD4)
+  list(APPEND vdw_sources vdwd4.cpp)
+endif()
+
+add_library(vdw OBJECT ${vdw_sources})
+
+if(ENABLE_DFTD4)
+  target_link_libraries(vdw PUBLIC dftd4::dftd4)
+  target_compile_definitions(vdw PUBLIC __DFTD4)
+endif()
+
 if(ENABLE_COVERAGE)
   add_coverage(vdw)
 endif()
@@ -19,4 +28,4 @@ if(BUILD_TESTING)
   if(ENABLE_MPI)
     add_subdirectory(test)
   endif()
-endif()
+endif()
diff --git a/source/source_hamilt/module_vdw/test/CMakeLists.txt b/source/source_hamilt/module_vdw/test/CMakeLists.txt
@@ -10,3 +10,8 @@ AddTest(
   LIBS parameter ${math_libs} base device vdw 
   SOURCES vdw_test.cpp 
 )
+
+if(ENABLE_DFTD4)
+  target_compile_definitions(MODULE_HAMILT_vdwTest PRIVATE __DFTD4)
+  target_link_libraries(MODULE_HAMILT_vdwTest dftd4::dftd4)
+endif()
diff --git a/source/source_hamilt/module_vdw/test/vdw_test.cpp b/source/source_hamilt/module_vdw/test/vdw_test.cpp
@@ -10,6 +10,9 @@
 #include "source_hamilt/module_vdw/vdwd3_parameters.h"
 #include "source_hamilt/module_vdw/vdwd2.h"
 #include "source_hamilt/module_vdw/vdwd3.h"
+#ifdef __DFTD4
+#include "source_hamilt/module_vdw/vdwd4.h"
+#endif
 #include "source_hamilt/module_vdw/vdw.h"
 #undef private
 
@@ -604,6 +607,73 @@ TEST_F(vdwd3abcTest, D3bjGetStress)
     EXPECT_NEAR(stress.e33, -3.4278442125590892e-05,1e-12);
 }
 
+#ifdef __DFTD4
+
+class vdwd4Test: public testing::Test
+{
+    protected:
+    UnitCell ucell;
+    Input_para input;
+
+    void SetUp(){
+        stru_ structure{std::vector<double>{0.5,0.5,0.0,0.5,0.0,0.5,0.0,0.5,0.5},
+                        std::vector<atomtype_>{atomtype_{"Si",
+                                                         std::vector<std::vector<double>>{
+                                                             {0., 0., 0.},
+                                                             {0.3, 0.25, 0.25}
+                                                         }}}};
+        construct_ucell(structure,ucell);
+
+        input.vdw_method = "d4";
+        input.vdw_d4_xc = "pbe";
+        input.vdw_cutoff_type = "radius";
+        input.vdw_radius_unit = "Bohr";
+        input.vdw_cutoff_radius = "60";
+        input.vdw_cn_thr_unit = "Bohr";
+        input.vdw_cn_thr = 30;
+    }
+
+    void TearDown(){
+        ClearUcell(ucell);
+    }
+};
+
+TEST_F(vdwd4Test, D4GetEnergy)
+{
+    auto vdw_solver = vdw::make_vdw(ucell, input);
+    double ene = vdw_solver->get_energy();
+    EXPECT_NEAR(ene, -0.04998837990336073, 1E-10);
+}
+
+TEST_F(vdwd4Test, D4GetForce)
+{
+    auto vdw_solver = vdw::make_vdw(ucell, input);
+    std::vector<ModuleBase::Vector3<double>> force = vdw_solver->get_force();
+    EXPECT_NEAR(force[0].x, -0.0023357259921368717, 1e-12);
+    EXPECT_NEAR(force[0].y, 0.0, 1e-12);
+    EXPECT_NEAR(force[0].z, 0.0, 1e-12);
+    EXPECT_NEAR(force[1].x, 0.0023357259921368730, 1e-12);
+    EXPECT_NEAR(force[1].y, 0.0, 1e-12);
+    EXPECT_NEAR(force[1].z, 0.0, 1e-12);
+}
+
+TEST_F(vdwd4Test, D4GetStress)
+{
+    auto vdw_solver = vdw::make_vdw(ucell, input);
+    ModuleBase::Matrix3 stress = vdw_solver->get_stress();
+    EXPECT_NEAR(stress.e11, 0.00015830384474877792, 1e-12);
+    EXPECT_NEAR(stress.e12, 0.0, 1e-12);
+    EXPECT_NEAR(stress.e13, 0.0, 1e-12);
+    EXPECT_NEAR(stress.e21, 0.0, 1e-12);
+    EXPECT_NEAR(stress.e22, 0.00016694998515968720, 1e-12);
+    EXPECT_NEAR(stress.e23, -1.5500973166318808e-05, 1e-12);
+    EXPECT_NEAR(stress.e31, 0.0, 1e-12);
+    EXPECT_NEAR(stress.e32, -1.5500973166318808e-05, 1e-12);
+    EXPECT_NEAR(stress.e33, 0.00016694998515968726, 1e-12);
+}
+
+#endif // __DFTD4
+
 int main(int argc, char **argv)
 {
     MPI_Init(&argc, &argv);

diff --git a/source/source_hamilt/module_vdw/vdw.cpp b/source/source_hamilt/module_vdw/vdw.cpp
@@ -6,6 +6,10 @@
 #include "vdwd3.h"
 #include "source_base/tool_quit.h"
 
+#ifdef __DFTD4
+#include "vdwd4.h"
+#endif
+
 std::string parse_xcname(const std::string &xc_input,
                          const std::vector<std::string> &xc_psp)
 {
@@ -70,6 +74,28 @@ std::unique_ptr<Vdw> make_vdw(const UnitCell &ucell,
         vdw_ptr->parameter().initial_parameters(parse_xcname(input.dft_functional, xc_psp), input, plog);
         return vdw_ptr;
     }
+    else if (input.vdw_method == "d4")
+    {
+#ifdef __DFTD4
+        std::vector<std::string> xc_psp(ucell.ntype);
+        for (int it = 0; it < ucell.ntype; it++)
+        {
+            xc_psp[it] = ucell.atoms[it].ncpp.xc_func;
+        }
+
+        std::string xc_name = input.vdw_d4_xc;
+        if (xc_name == "default")
+        {
+            xc_name = parse_xcname(input.dft_functional, xc_psp);
+        }
+
+        return vdw::make_unique<Vdwd4>(ucell, xc_name, input);
+#else
+        ModuleBase::WARNING_QUIT("ModuleHamiltGeneral::ModuleVDW::make_vdw",
+                                 "DFT-D4 support was not enabled at build time. "
+                                 "Rebuild ABACUS with -DENABLE_DFTD4=ON.");
+#endif
+    }
     else if (input.vdw_method != "none")
     {
         ModuleBase::WARNING_QUIT("ModuleHamiltGeneral::ModuleVDW::make_vdw",