add write_libxc_r

docs/advanced/input_files/input-main.md

@@ -134,6 +134,7 @@
   - [Output Variables](#variables-related-to-output-information)
     - [out\_freq\_elec](#out_freq_elec)
     - [out\_chg](#out_chg)
+    - [out\_xc\_r](#out_xc_r)
     - [out\_pot](#out_pot)
     - [out\_dm](#out_dmk)
     - [out\_dm1](#out_dmr)
@@ -1652,6 +1653,25 @@ These variables are used to control the output of properties.
 - **Default**: 0 3
 - **Note**: In the 3.10-LTS version, the file names are SPIN1_CHG.cube and SPIN1_CHG_INI.cube, etc. 
 
+### out_xc_r
+
+- **Type**: Integer \[Integer\](optional)
+- **Description**: 
+  The first integer controls whether to output the exchange-correlation (in Bohr^-3) on real space grids using Libxc to folder `OUT.${suffix}`:
+  - 0: rho, amag, sigma, exc
+  - 1: vrho, vsigma
+  - 2: v2rho2, v2rhosigma, v2sigma2
+  - 3: v3rho3, v3rho2sigma, v3rhosigma2, v3sigma3
+  - 4: v4rho4, v4rho3sigma, v4rho2sigma2, v4rhosigma3, v4sigma4
+  The meaning of the files is presented in [Libxc](https://libxc.gitlab.io/manual/libxc-5.1.x/)
+
+  The second integer controls the precision of the charge density output, if not given, will use `3` as default.
+
+  ---
+  The circle order of the charge density on real space grids is: x is the outer loop, then y and finally z (z is moving fastest).
+
+- **Default**: -1 3
+
 ### out_pot
 
 - **Type**: Integer

source/Makefile.Objects

@@ -539,6 +539,7 @@ OBJS_IO=input_conv.o\
     write_dipole.o\
     td_current_io.o\
     write_wfc_r.o\
+    write_libxc_r.o\
     output_log.o\
     output_mat_sparse.o\
     ctrl_output_lcao.o\

source/module_io/CMakeLists.txt

@@ -31,6 +31,7 @@ list(APPEND objects
     write_mlkedf_descriptors.cpp
     td_current_io.cpp
     write_wfc_r.cpp
+    write_libxc_r.cpp
     output_log.cpp
     para_json.cpp
     parse_args.cpp

source/module_io/read_input_item_output.cpp

@@ -514,6 +514,27 @@ void ReadInput::item_output()
         add_intvec_bcast(input.out_wfc_re_im, para.input.out_wfc_re_im.size(), 0);
         this->add_item(item);
     }
+    {
+        Input_Item item("out_xc_r");
+        item.annotation = "if >=0, output the derivatives of exchange correlation in realspace, second parameter controls the precision";
+        item.read_value = [](const Input_Item& item, Parameter& para) {
+            size_t count = item.get_size();
+            std::vector<int> out_xc_r(count); // create a placeholder vector
+            std::transform(item.str_values.begin(), item.str_values.end(), out_xc_r.begin(), [](std::string s) { return std::stoi(s); });
+            // assign non-negative values to para.input.out_xc_r
+            std::copy(out_xc_r.begin(), out_xc_r.end(), para.input.out_xc_r.begin());
+        };
+        item.check_value = [](const Input_Item& item, const Parameter& para) {
+            if (para.input.out_xc_r[0] >= 0)
+            {
+#ifndef USE_LIBXC
+                ModuleBase::WARNING_QUIT("ReadInput", "INPUT out_xc_r is only aviailable with Libxc");
+#endif
+            }
+        };
+        sync_intvec(input.out_xc_r, 2, -1);
+        this->add_item(item);
+    }
     {
         Input_Item item("if_separate_k");
         item.annotation = "specify whether to write the partial charge densities for all k-points to individual files "