Fix: let nk_total *= nspin, in rdmft

Author: JGHan7

source/module_esolver/esolver_ks_lcao.cpp

@@ -1116,7 +1116,6 @@ void ESolver_KS_LCAO<TK, TR>::iter_finish(int& iter)
             }
 
             this->rdmft_solver.update_elec(occ_number_ks, *(this->psi));
-            std::cout << "\n******\n" << "update elec in rdmft successfully" << "\n******\n" << std::endl;
 
             //initialize the gradients of Etotal on occupation numbers and wfc, and set all elements to 0. 
             ModuleBase::matrix dE_dOccNum(this->pelec->wg.nr, this->pelec->wg.nc, true);

source/module_io/read_input_item_other.cpp

@@ -504,7 +504,7 @@ void ReadInput::item_others()
     {
         Input_Item item("rdmft_power_alpha");
         item.annotation = "the alpha parameter of power-functional, g(occ_number) = occ_number^alpha"
-                          "used in exx-type functionals such as muller and power";
+                          " used in exx-type functionals such as muller and power";
         read_sync_double(input.rdmft_power_alpha);
         item.reset_value = [](const Input_Item& item, Parameter& para) {
             if( para.input.dft_functional == "hf" || para.input.dft_functional == "pbe0" )

source/module_rdmft/rdmft.cpp

@@ -94,9 +94,10 @@ void RDMFT<TK, TR>::init(Gint_Gamma& GG_in, Gint_k& GK_in, Parallel_Orbitals& Pa
     XC_func_rdmft = XC_func_rdmft_in;
     alpha_power = alpha_power_in;
 
-    nk_total = ModuleSymmetry::Symmetry::symm_flag == -1 ? kv->get_nkstot_full(): kv->get_nks();
-    nbands_total = PARAM.inp.nbands;
     nspin = PARAM.inp.nspin;
+    nbands_total = PARAM.inp.nbands;
+    nk_total = ModuleSymmetry::Symmetry::symm_flag == -1 ? kv->get_nkstot_full(): kv->get_nks();
+    nk_total *= nspin;
     only_exx_type = ( XC_func_rdmft == "hf" || XC_func_rdmft == "muller" || XC_func_rdmft == "power" );
 
     // // create desc[] and something about MPI to Eij(nbands*nbands)

source/module_rdmft/update_state_rdmft.cpp

@@ -80,6 +80,8 @@ void RDMFT<TK, TR>::update_elec(const ModuleBase::matrix& occ_number_in, const p
     this->cal_V_hartree();
     this->cal_V_XC();
     // this->cal_Hk_Hpsi();
+
+    std::cout << "\n******\n" << "update elec in rdmft successfully" << "\n******\n" << std::endl;
 }
 
 
@@ -148,27 +150,24 @@ void RDMFT<TK, TR>::update_charge()
         charge->renormalize_rho();
     }
 
-    /*********** what's this? When we use PBE-functional, this can't be deleted *************/
+    // charge density symmetrization
     // this->pelec->calculate_weights();
     // this->pelec->calEBand();
     Symmetry_rho srho;
     for (int is = 0; is < nspin; is++)
     {
         srho.begin(is, *(this->charge), rho_basis, GlobalC::ucell.symm);
     }
-    /*********** what's this? When we use PBE-functional, this can't be deleted *************/
-
-    // // what this?
-    // if(GlobalV::VL_IN_H)
-    // {
-    //     // update Gint_K
-    //     if (!PARAM.inp.gamma_only)
-    //     {
-    //         this->UHM.GK.renew();
-    //     }
-    //     // update real space Hamiltonian
-    //     // this->p_hamilt->refresh();
-    // }
+
+    // what this? it seems that it needs to be updated at each iteration
+    if (PARAM.inp.vl_in_h)
+    {
+        // update Gint_K
+        if (!PARAM.globalv.gamma_only_local)
+        {
+            this->GK->renew();
+        }
+    }
 
 }