Feature: Add deepks_out_base to support simultaneous output of numpy files of base and target functionals (#6483)

* add parameter deepks_out_base, when it is not "none", use a new esolver type "Esolver_DoubleXC"

* Modify Esolver_DoubleXC to calculate ebase,hbase,obase of output charge density (not charge density after charge mixing). Force to use deepks_out_freq_elec along with deepks_out_base

* some output info which helps when debuging

* output files of deepks base; add some annotation for output

* fix bug of force output related to gevdm when deepks_scf = false but deepks_out_labels = 1

* fix the bug of gamma only related to smatrix_k, add some annotations

* update annotation

* add conditions for p_hamit->refresh

* output force and stress base also in DeePKS_Labels_Elec

* revert changes to INPUT in test

* modify test for deepks_out_base

* Revert "some output info which helps when debuging" to delete output information

This reverts commit 73c3d6d.

* remove some output info for debugging in esolver_double_xc

* add some anotation for out_mat_hs used along with deepks_out_base

* add esolver_double_xc in Makefile.Objects

* Add deepks_out_base to input-main.md

Added documentation for deepks_out_base parameter.

* add ifdef for __MLALGO in esolver_double_xc, in case of not compiling with MLALGO

* Update doc for deepks_out_base

Author: xuan112358

docs/advanced/input_files/input-main.md

@@ -191,6 +191,7 @@
   - [DeePKS](#deepks)
     - [deepks\_out\_labels](#deepks_out_labels)
     - [deepks\_out\_freq\_elec](#deepks_out_freq_elec)
+    - [deepks\_out\_base](#deepks_out_base)
     - [deepks\_scf](#deepks_scf)
     - [deepks\_equiv](#deepks_equiv)
     - [deepks\_model](#deepks_model)
@@ -2173,6 +2174,13 @@ Warning: this function is not robust enough for the current version. Please try
 - **Description**: When `deepks_out_freq_elec` is greater than 0, print labels and descriptors for DeePKS in OUT.${suffix}/DeePKS_Labels_Elec per `deepks_out_freq_elec` electronic iterations, with suffix `_e*` to distinguish different steps. Often used with `deepks_out_labels` equals 1.
 - **Default**: 0
 
+### deepks_out_base
+
+- **Type**: String
+- **Availability**: Numerical atomic orbital basis and `deepks_out_freq_elec` is greater than 0
+- **Description**: Print labels and descriptors calculated by base functional ( determined by `deepks_out_base` ) and target functional ( determined by `dft_functional` ) for DeePKS in per `deepks_out_freq_elec` electronic iterations. The SCF process, labels and descriptors output of the target functional are all consistent with those when the target functional is used alone. The only additional output under this configuration is the labels of the base functional. Often used with `deepks_out_labels` equals 1.
+- **Default**: None
+
 ### deepks_scf
 
 - **Type**: Boolean

source/Makefile.Objects

@@ -275,6 +275,7 @@ OBJS_ESOLVER_LCAO=esolver_ks_lcao.o\
       esolver_gets.o\
       lcao_others.o\
       esolver_dm2rho.o\
+      esolver_double_xc.o\
 
 OBJS_GINT=gint_old.o\
       gint_gamma_env.o\

source/source_esolver/CMakeLists.txt

@@ -22,6 +22,7 @@ if(ENABLE_LCAO)
       esolver_gets.cpp
       lcao_others.cpp
       esolver_dm2rho.cpp
+      esolver_double_xc.cpp
   )
 endif()
 

source/source_esolver/esolver.cpp

@@ -6,6 +6,7 @@
 #include "source_io/module_parameter/parameter.h"
 #ifdef __LCAO
 #include "esolver_dm2rho.h"
+#include "esolver_double_xc.h"
 #include "esolver_gets.h"
 #include "esolver_ks_lcao.h"
 #include "esolver_ks_lcao_tddft.h"
@@ -197,14 +198,25 @@ ESolver* init_esolver(const Input_para& inp, UnitCell& ucell)
         }
         if (PARAM.globalv.gamma_only_local)
         {
-            return new ESolver_KS_LCAO<double, double>();
+            if (PARAM.inp.deepks_out_base != "none")
+            {
+                return new ESolver_DoubleXC<double, double>();
+            }
+            else
+            {
+                return new ESolver_KS_LCAO<double, double>();
+            }
         }
         else if (PARAM.inp.nspin < 4)
         {
             if (PARAM.inp.dm_to_rho)
             {
                 return new ESolver_DM2rho<std::complex<double>, double>();
             }
+            else if (PARAM.inp.deepks_out_base != "none")
+            {
+                return new ESolver_DoubleXC<std::complex<double>, double>();
+            }
             else
             {
                 return new ESolver_KS_LCAO<std::complex<double>, double>();
@@ -216,6 +228,10 @@ ESolver* init_esolver(const Input_para& inp, UnitCell& ucell)
             {
                 return new ESolver_DM2rho<std::complex<double>, std::complex<double>>();
             }
+            else if (PARAM.inp.deepks_out_base != "none")
+            {
+                return new ESolver_DoubleXC<std::complex<double>, std::complex<double>>();
+            }
             else
             {
                 return new ESolver_KS_LCAO<std::complex<double>, std::complex<double>>();