Feature: enable init_chg=file for metagga

source/module_elecstate/module_charge/charge.h

@@ -145,7 +145,7 @@ class Charge
 
     // mohan add 2021-02-20
     int nrxx=0; // number of r vectors in this processor
-    int nxyz=0; // total nuber of r vectors
+    int nxyz = 0; // total number of r vectors
     int ngmc=0; // number of g vectors in this processor
     int nspin=0; // number of spins
     ModulePW::PW_Basis* rhopw = nullptr;// When double_grid is used, rhopw = rhodpw (dense grid)

source/module_elecstate/module_charge/charge_init.cpp

@@ -39,7 +39,7 @@ void Charge::init_rho(elecstate::efermi& eferm_iout,
     bool read_error = false;
     if (PARAM.inp.init_chg == "file" || PARAM.inp.init_chg == "auto")
     {
-        GlobalV::ofs_running << " try to read charge from file : " << std::endl;
+        GlobalV::ofs_running << " try to read charge from file" << std::endl;
 
         // try to read charge from binary file first, which is the same as QE
         // liuyu 2023-12-05
@@ -99,22 +99,43 @@ void Charge::init_rho(elecstate::efermi& eferm_iout,
                     break;
                 }
             }
+        }
+
+        if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
+        {
+            GlobalV::ofs_running << " try to read kinetic energy density from file" << std::endl;
+            // try to read charge from binary file first, which is the same as QE
+            std::vector<std::complex<double>> kin_g_space(PARAM.inp.nspin * this->ngmc, {0.0, 0.0});
+            std::vector<std::complex<double>*> kin_g;
+            for (int is = 0; is < PARAM.inp.nspin; is++)
+            {
+                kin_g.push_back(kin_g_space.data() + is * this->ngmc);
+            }
 
-            if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
+            std::stringstream binary;
+            binary << PARAM.globalv.global_readin_dir << PARAM.inp.suffix + "-TAU-DENSITY.restart";
+            if (ModuleIO::read_rhog(binary.str(), rhopw, kin_g.data()))
+            {
+                GlobalV::ofs_running << " Read in the kinetic energy density: " << binary.str() << std::endl;
+                for (int is = 0; is < PARAM.inp.nspin; ++is)
+                {
+                    rhopw->recip2real(kin_g[is], this->kin_r[is]);
+                }
+            }
+            else
             {
                 for (int is = 0; is < PARAM.inp.nspin; is++)
                 {
                     std::stringstream ssc;
                     ssc << PARAM.globalv.global_readin_dir << "SPIN" << is + 1 << "_TAU.cube";
-                    GlobalV::ofs_running << " try to read kinetic energy density from file : " << ssc.str()
-                                         << std::endl;
                     // mohan update 2012-02-10, sunliang update 2023-03-09
-                    if (ModuleIO::read_vdata_palgrid(pgrid,
-                        (PARAM.inp.esolver_type == "sdft" ? GlobalV::RANK_IN_STOGROUP : GlobalV::MY_RANK),
-                        GlobalV::ofs_running,
-                        ssc.str(),
-                        this->kin_r[is],
-                        ucell.nat))
+                    if (ModuleIO::read_vdata_palgrid(
+                            pgrid,
+                            (PARAM.inp.esolver_type == "sdft" ? GlobalV::RANK_IN_STOGROUP : GlobalV::MY_RANK),
+                            GlobalV::ofs_running,
+                            ssc.str(),
+                            this->kin_r[is],
+                            ucell.nat))
                     {
                         GlobalV::ofs_running << " Read in the kinetic energy density: " << ssc.str() << std::endl;
                     }

tests/integrate/101_PW_15_f_pseudopots/INPUT

@@ -1,7 +1,8 @@
 INPUT_PARAMETERS
 #Parameters (1.General)
 suffix			autotest
-calculation     scf
+calculation             scf
+dft_functional          scan
 
 nbands			8
 symmetry		1
@@ -23,3 +24,5 @@ smearing_sigma			0.002
 mixing_type		broyden
 mixing_beta		0.7
 
+init_chg                file
+read_file_dir           .

tests/integrate/101_PW_15_f_pseudopots/README

@@ -3,3 +3,5 @@ This test for:
 *SG15 pseudopotential containing f electron
 *mixing_type broyden
 *mixing_beta 0.7
+*dft_functional scan
+*init_chg file

tests/integrate/101_PW_15_f_pseudopots/STRU

@@ -2,7 +2,7 @@ ATOMIC_SPECIES
 Ce 140.115 58_Ce.UPF upf201
 
 LATTICE_CONSTANT
-8.92
+6
 
 LATTICE_VECTORS
 0.5 0.5 0.0

tests/integrate/101_PW_15_f_pseudopots/jd

@@ -1 +1 @@
-test SG15 pseudopotential containing f electron, symmetry=on
+test SG15 pseudopotential containing f electron with scan, symmetry=on, init_chg=file

tests/integrate/101_PW_15_f_pseudopots/result.ref

@@ -1,6 +1,6 @@
-etotref -3436.745619582930
-etotperatomref -3436.7456195829
+etotref -3644.833836110242
+etotperatomref -3644.8338361102
 pointgroupref O_h
 spacegroupref O_h
 nksibzref 3
-totaltimeref 0.16141
+totaltimeref 0.33