diff --git a/source/source_esolver/esolver_fp.cpp b/source/source_esolver/esolver_fp.cpp
index b0114402c6..493409716c 100644
--- a/source/source_esolver/esolver_fp.cpp
+++ b/source/source_esolver/esolver_fp.cpp
@@ -41,6 +41,21 @@ void ESolver_FP::before_all_runners(BaseCell& basecell, const Input_para& inp)
 
     this->inp_ = &inp;
 
+    SurchemParameters surchem_parameters;
+    surchem_parameters.eb_k = inp.eb_k;
+    surchem_parameters.tau = inp.tau;
+    surchem_parameters.sigma_k = inp.sigma_k;
+    surchem_parameters.nc_k = inp.nc_k;
+    this->solvent.set_parameters(surchem_parameters);
+
+    XCFunctionalParameters xc_parameters;
+    xc_parameters.xc_temperature = inp.xc_temperature;
+    xc_parameters.exx_fock_alpha = inp.exx_fock_alpha;
+    xc_parameters.exx_erfc_alpha = inp.exx_erfc_alpha;
+    xc_parameters.xc_exch_ext = inp.xc_exch_ext;
+    xc_parameters.xc_corr_ext = inp.xc_corr_ext;
+    XC_Functional::set_runtime_parameters(xc_parameters);
+
     ModuleBase::TITLE("ESolver_FP", "before_all_runners");
 
     //! 1) read pseudopotentials
@@ -211,7 +226,7 @@ void ESolver_FP::before_scf(UnitCell& ucell, const int istep)
     //! calculate ewald energy
     if (!this->inp_->test_skip_ewald)
     {
-        this->pelec->f_en.ewald_energy = H_Ewald_pw::compute_ewald(ucell, this->pw_rhod, this->sf.strucFac);
+        this->pelec->f_en.ewald_energy = H_Ewald_pw::compute_ewald(ucell, this->pw_rhod, this->sf.strucFac, this->inp_->test_energy, GlobalV::ofs_running);
     }
 
     //! set direction of magnetism, used in non-collinear case 
diff --git a/source/source_esolver/esolver_ks_lcao.cpp b/source/source_esolver/esolver_ks_lcao.cpp
index 90e1d8fa01..857f1cfd03 100644
--- a/source/source_esolver/esolver_ks_lcao.cpp
+++ b/source/source_esolver/esolver_ks_lcao.cpp
@@ -137,7 +137,12 @@ void ESolver_KS_LCAO<TK, TR>::before_scf(UnitCell& ucell, const int istep)
         this->pw_rho->nx, this->pw_rho->ny, this->pw_rho->nz,
         0, 0, this->pw_big->nbzp_start,
         this->pw_big->nbx, this->pw_big->nby, this->pw_big->nbzp,
-        orb_.Phi, ucell, this->gd));
+        orb_.Phi, ucell, this->gd,
+        this->inp_->nspin,
+        PARAM.globalv.gamma_only_local,
+        PARAM.globalv.domag,
+        this->inp_->device == "gpu",
+        this->inp_->nstream));
     ModuleGint::Gint::set_gint_info(gint_info_.get());
 
     // 7) For each atom, calculate the adjacent atoms in different cells
diff --git a/source/source_esolver/esolver_lr_lcao_bse.cpp b/source/source_esolver/esolver_lr_lcao_bse.cpp
index 409d9139fe..0f536d0ba8 100644
--- a/source/source_esolver/esolver_lr_lcao_bse.cpp
+++ b/source/source_esolver/esolver_lr_lcao_bse.cpp
@@ -123,7 +123,12 @@ void ESolver_BSE<T, TR>::before_all_runners(BaseCell& basecell, const Input_para
             this->pw_big->nbzp,
             this->orb_.Phi,
             ucell,
-            this->gd));
+            this->gd,
+            inp.nspin,
+            PARAM.globalv.gamma_only_local,
+            PARAM.globalv.domag,
+            inp.device == "gpu",
+            inp.nstream));
             ModuleGint::Gint::set_gint_info(this->gint_info_.get());
 
     this->pot.resize(this->nspin, nullptr);
diff --git a/source/source_esolver/esolver_lr_lcao_tddft.cpp b/source/source_esolver/esolver_lr_lcao_tddft.cpp
index 52bd6ad0e8..c83fa739c8 100644
--- a/source/source_esolver/esolver_lr_lcao_tddft.cpp
+++ b/source/source_esolver/esolver_lr_lcao_tddft.cpp
@@ -470,7 +470,12 @@ void ModuleESolver::ESolver_LR<T, TR>::initialize_from_unitcell_(UnitCell& ucell
         this->pw_big->nbzp,
         orb.Phi,
         ucell,
-        this->gd));
+        this->gd,
+        this->inp_->nspin,
+        PARAM.globalv.gamma_only_local,
+        PARAM.globalv.domag,
+        this->inp_->device == "gpu",
+        this->inp_->nstream));
     ModuleGint::Gint::set_gint_info(gint_info_.get());
     // if EXX from scratch, init 2-center integral and calculate Cs, Vs 
 #ifdef __EXX
diff --git a/source/source_esolver/lcao_others.cpp b/source/source_esolver/lcao_others.cpp
index 943b6632f4..4a73a20691 100644
--- a/source/source_esolver/lcao_others.cpp
+++ b/source/source_esolver/lcao_others.cpp
@@ -99,7 +99,12 @@ void ESolver_KS_LCAO<TK, TR>::others(BaseCell& basecell, const int istep)
                                               this->pw_big->nbzp,
                                               orb_.Phi,
                                               ucell,
-                                              this->gd));
+                                              this->gd,
+                                              this->inp_->nspin,
+                                              gamma_only_local,
+                                              PARAM.globalv.domag,
+                                              this->inp_->device == "gpu",
+                                              this->inp_->nstream));
     ModuleGint::Gint::set_gint_info(gint_info_.get());
 
     // (2)For each atom, calculate the adjacent atoms in different cells
diff --git a/source/source_hamilt/module_ewald/h_ewald_pw.cpp b/source/source_hamilt/module_ewald/h_ewald_pw.cpp
index 510f304b00..459b77b218 100644
--- a/source/source_hamilt/module_ewald/h_ewald_pw.cpp
+++ b/source/source_hamilt/module_ewald/h_ewald_pw.cpp
@@ -1,8 +1,7 @@
 #include "h_ewald_pw.h"
+#include "source_base/global_function.h"
 #include "source_base/parallel_comm.h"
-#include "source_io/module_parameter/parameter.h"
 #include "source_base/mymath.h" // use heapsort
-#include "source_io/module_parameter/parameter.h"
 #include "dnrm2.h"
 #include "source_base/parallel_reduce.h"
 #include "source_base/constants.h"
@@ -27,7 +26,9 @@ int H_Ewald_pw::estimate_mxr(const double &rmax, const ModuleBase::Matrix3 &bg)
 
 double H_Ewald_pw::compute_ewald(const UnitCell& cell,
                                  const ModulePW::PW_Basis* rho_basis,
-                                 const ModuleBase::ComplexMatrix& strucFac)
+                                 const ModuleBase::ComplexMatrix& strucFac,
+                                 const int test_energy,
+                                 std::ofstream& output_stream)
 {
     ModuleBase::TITLE("H_Ewald_pw","compute_ewald");
     ModuleBase::timer::start("H_Ewald_pw","compute_ewald");
@@ -73,9 +74,9 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
             charge += cell.atoms[it].na * cell.atoms[it].ncpp.zv;//mohan modify 2007-11-7
         }
     }
-    if(PARAM.inp.test_energy) 
+    if(test_energy)
     {
-        ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running,"Total ionic charge",charge);
+        ModuleBase::GlobalFunc::OUT(output_stream,"Total ionic charge",charge);
     }
 
 	// (2) calculate the converged value: alpha
@@ -94,10 +95,10 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
                      erfc(sqrt(cell.tpiba2 * rho_basis->ggecut / 4.0 / alpha));
     }
     while (upperbound > 1.0e-7);
-    if(PARAM.inp.test_energy) 
+    if(test_energy)
     {
-        ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running,"alpha",alpha);
-        ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running,"Upper bound",upperbound);
+        ModuleBase::GlobalFunc::OUT(output_stream,"alpha",alpha);
+        ModuleBase::GlobalFunc::OUT(output_stream,"Upper bound",upperbound);
     }
 
     // G-space sum here.
@@ -123,7 +124,7 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
 	// but that's not the term "gamma_only" I want to use in LCAO,  
 	fact = 1.0;
 
-    //GlobalV::ofs_running << "\n pwb.gstart = " << pwb.gstart << std::endl;
+    //output_stream << "\n pwb.gstart = " << pwb.gstart << std::endl;
     const int ig0 = rho_basis->ig_gge0;
     for (int ig = 0; ig < rho_basis->npw; ig++)
     {
@@ -165,9 +166,9 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
     rmax = 4.0 / sqrt(alpha) / cell.lat0;
     mxr = H_Ewald_pw::estimate_mxr(rmax, cell.G);
 
-    if(PARAM.inp.test_energy) 
+    if(test_energy)
     {
-        ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running,"mxr",mxr);
+        ModuleBase::GlobalFunc::OUT(output_stream,"mxr",mxr);
     }
     std::vector<ModuleBase::Vector3<double>> vec_r(mxr);
     std::vector<double> vec_r2(mxr);
@@ -177,9 +178,9 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
     double* r2 = vec_r2.data();
 
 #ifdef __MPI
-    if(PARAM.inp.test_energy) 
+    if(test_energy)
     {
-        ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running,"rmax(unit lat0)",rmax);
+        ModuleBase::GlobalFunc::OUT(output_stream,"rmax(unit lat0)",rmax);
     }
 
     int size = 0;
@@ -209,11 +210,11 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
             // calculate tau[na1]-tau[na2]
             dtau = cell.atoms[it1].tau[ia1] - cell.atoms[it2].tau[ia2];
             // generates nearest-neighbors shells
-            H_Ewald_pw::rgen(dtau, rmax, irr, cell.latvec, cell.G, r, r2, mxr, nrm);
+            H_Ewald_pw::rgen(dtau, rmax, irr, cell.latvec, cell.G, r, r2, mxr, nrm, test_energy);
             // at-->cell.latvec, bg-->G
             // and sum to the real space part
 
-            if(PARAM.inp.test_energy>1)
+            if(test_energy>1)
             {
                 ModuleBase::GlobalFunc::OUT("dtau.x",dtau.x);
                 ModuleBase::GlobalFunc::OUT("dtau.y",dtau.y);
@@ -228,7 +229,7 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
                             erfc(sqrt(alpha) * rr) / rr;
                 }
             }
-            if (PARAM.inp.test_energy>1) 
+            if (test_energy>1)
             {
                 ModuleBase::GlobalFunc::OUT("ewaldr",ewaldr);
             }
@@ -237,7 +238,7 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
 #else
     if (rho_basis->ig_gge0 >= 0)
     {	
-        if(PARAM.inp.test_energy) ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running,"rmax(unit lat0)",rmax);
+        if(test_energy) ModuleBase::GlobalFunc::OUT(output_stream,"rmax(unit lat0)",rmax);
         // with this choice terms up to ZiZj*erfc(4) are counted (erfc(4)=2x10^-8
         int nt1=0;
         int nt2=0;
@@ -253,11 +254,11 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
                         //calculate tau[na]-tau[nb]
                         dtau = cell.atoms[nt1].tau[na] - cell.atoms[nt2].tau[nb];
                         //generates nearest-neighbors shells
-                        H_Ewald_pw::rgen(dtau, rmax, irr, cell.latvec, cell.G, r, r2, mxr, nrm);
+                        H_Ewald_pw::rgen(dtau, rmax, irr, cell.latvec, cell.G, r, r2, mxr, nrm, test_energy);
                         // at-->cell.latvec, bg-->G
                         // and sum to the real space part
 
-                        if (PARAM.inp.test_energy>1)
+                        if (test_energy>1)
                         {
                             ModuleBase::GlobalFunc::OUT("dtau.x",dtau.x);
                             ModuleBase::GlobalFunc::OUT("dtau.y",dtau.y);
@@ -272,7 +273,7 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
                                         erfc(sqrt(alpha) * rr) / rr;
                             }
                         } // enddo
-                        if (PARAM.inp.test_energy>1) ModuleBase::GlobalFunc::OUT("ewaldr",ewaldr);
+                        if (test_energy>1) ModuleBase::GlobalFunc::OUT("ewaldr",ewaldr);
                     } // enddo
                 } // enddo
             } // nt2
@@ -285,7 +286,7 @@ double H_Ewald_pw::compute_ewald(const UnitCell& cell,
 	// mohan fix bug 2010-07-26
     Parallel_Reduce::reduce_pool(ewalds);
 
-    if (PARAM.inp.test_energy>1)
+    if (test_energy>1)
     {
         ModuleBase::GlobalFunc::OUT("ewaldg",ewaldg);
         ModuleBase::GlobalFunc::OUT("ewaldr",ewaldr);
@@ -306,7 +307,8 @@ void H_Ewald_pw::rgen(
     ModuleBase::Vector3<double> *r,
     double *r2,
     const int mxr,
-    int &nrm)
+    int &nrm,
+    const int test_energy)
 {
     //-------------------------------------------------------------------
     // generates neighbours shells (in units of alat) with length
@@ -377,7 +379,7 @@ void H_Ewald_pw::rgen(
 
     nm3 = (int)(dnrm2(3, bg1, 1) * rmax + 2);
 
-    if (PARAM.inp.test_energy>1)
+    if (test_energy>1)
     {
         ModuleBase::GlobalFunc::OUT("nm1",nm1);
         ModuleBase::GlobalFunc::OUT("nm2",nm2);
diff --git a/source/source_hamilt/module_ewald/h_ewald_pw.h b/source/source_hamilt/module_ewald/h_ewald_pw.h
index 64143c8aa8..2b942c4f5e 100644
--- a/source/source_hamilt/module_ewald/h_ewald_pw.h
+++ b/source/source_hamilt/module_ewald/h_ewald_pw.h
@@ -1,10 +1,11 @@
 #ifndef H_EWALD_PW_H
 #define H_EWALD_PW_H
 
-#include "source_base/global_function.h"
 #include "source_cell/unitcell.h"
 #include "source_basis/module_pw/pw_basis.h"
 
+#include <fstream>
+
 class H_Ewald_pw 
 {
   public:
@@ -14,7 +15,9 @@ class H_Ewald_pw
     // compute the Ewald energy
     static double compute_ewald(const UnitCell& cell,
                                 const ModulePW::PW_Basis* rho_basis,
-                                const ModuleBase::ComplexMatrix& strucFac);
+                                const ModuleBase::ComplexMatrix& strucFac,
+                                int test_energy,
+                                std::ofstream& output_stream);
 
   public:
     static int estimate_mxr(const double &rmax, const ModuleBase::Matrix3 &bg);
@@ -28,7 +31,8 @@ class H_Ewald_pw
         ModuleBase::Vector3<double> *r,
         double *r2,
       const int mxr,
-        int  &nrm
+        int  &nrm,
+        int test_energy
     );
 
 	// the coefficient of ewald method
diff --git a/source/source_hamilt/module_gint/gint.h b/source/source_hamilt/module_gint/gint.h
index 1255bae971..84ae4897c4 100644
--- a/source/source_hamilt/module_gint/gint.h
+++ b/source/source_hamilt/module_gint/gint.h
@@ -19,6 +19,11 @@ class Gint
         gint_info_ = gint_info;
     }
 
+    static const GintInfo& get_gint_info()
+    {
+        return *gint_info_;
+    }
+
     protected:
     static GintInfo* gint_info_;
 };
diff --git a/source/source_hamilt/module_gint/gint_common.cpp b/source/source_hamilt/module_gint/gint_common.cpp
index 057e92e3de..4f34657640 100644
--- a/source/source_hamilt/module_gint/gint_common.cpp
+++ b/source/source_hamilt/module_gint/gint_common.cpp
@@ -1,7 +1,6 @@
 #include "gint_common.h"
 #include "source_hamilt/module_hcontainer/hcontainer.h"
 #include "source_hamilt/module_hcontainer/hcontainer_funcs.h"
-#include "source_io/module_parameter/parameter.h"
 #include "source_base/tool_quit.h"
 #include <cassert>
 #include <type_traits>
@@ -179,7 +178,7 @@ void merge_hr_part_to_hR(const std::vector<hamilt::HContainer<double>>& hr_gint_
     std::vector<int> clx_i = {1, 0, 0, -1};
     std::vector<int> clx_j = {0, -1, 1, 0};
     for (int is = 0; is < 4; is++){
-        if(!PARAM.globalv.domag && (is==1 || is==2)) continue;
+        if(!gint_info.get_domag() && (is==1 || is==2)) continue;
         hR_tmp->set_zero();
         hamilt::HContainer<std::complex<double>>* hRGint_tmpCd = new hamilt::HContainer<std::complex<double>>(ucell_in->nat);
         hRGint_tmpCd->insert_ijrs( &(gint_info.get_ijr_info()), *(ucell_in));
@@ -311,9 +310,9 @@ void dm_2d_to_gint(
     ModuleBase::TITLE("Gint", "dm_2d_to_gint");
     ModuleBase::timer::start("Gint", "dm_2d_to_gint");
 
-    if (PARAM.inp.nspin != 4)
+    if (gint_info.get_nspin() != 4)
     {
-        // dm_gint.size() usually equals to PARAM.inp.nspin,
+        // dm_gint.size() usually equals to the configured nspin,
         // but there is exception within source_lcao/module_lr
         for (int is = 0; is < dm_gint.size(); is++)
         {
@@ -404,6 +403,7 @@ void wfc_2d_to_gint(const T* wfc_2d,
     ModuleBase::TITLE("Gint", "wfc_2d_to_gint");
     ModuleBase::timer::start("Gint", "wfc_2d_to_gint");
 
+    const int requested_nbands = nbands;
 #ifdef __MPI
     // dimension related
     nlocal = pv.desc_wfc[2];
@@ -459,7 +459,7 @@ void wfc_2d_to_gint(const T* wfc_2d,
             for (int j = 0; j < naroc[1]; ++j)
             {
                 int igcol = globalIndex(j, nb, dim1, ipcol);
-                if (igcol >= PARAM.inp.nbands)
+                if (igcol >= requested_nbands)
                 {
                     continue;
                 }
diff --git a/source/source_hamilt/module_gint/gint_info.cpp b/source/source_hamilt/module_gint/gint_info.cpp
index 572cd55287..fcc8d34185 100644
--- a/source/source_hamilt/module_gint/gint_info.cpp
+++ b/source/source_hamilt/module_gint/gint_info.cpp
@@ -1,6 +1,5 @@
 #include <cmath>
 #include <map>
-#include "source_io/module_parameter/parameter.h"
 #include "source_base/timer.h"
 #include "gint_info.h"
 #include "gint_type.h"
@@ -15,8 +14,9 @@ GintInfo::GintInfo(
     int startidx_bx, int startidx_by, int startidx_bz,
     int nbx_local, int nby_local, int nbz_local,
     const Numerical_Orbital* Phi,
-    const UnitCell& ucell, Grid_Driver& gd)
-    : ucell_(&ucell)
+    const UnitCell& ucell, Grid_Driver& gd,
+    const int nspin, const bool gamma_only, const bool domag, const bool use_gpu, const int nstream)
+    : ucell_(&ucell), nspin_(nspin), gamma_only_(gamma_only), domag_(domag), use_gpu_(use_gpu)
 {
     // initialize the unitcell information
     unitcell_info_ = std::make_shared<UnitCellInfo>(ucell_->a1 * ucell_->lat0, ucell_->a2 * ucell_->lat0, ucell_->a3 * ucell_->lat0,
@@ -46,16 +46,16 @@ GintInfo::GintInfo(
     init_atoms_(ucell_->ntype, ucell_->atoms, Phi);
 
     // initialize trace_lo_ and lgd_
-    init_trace_lo_(ucell, PARAM.inp.nspin);
+    init_trace_lo_(ucell, nspin_);
 
     // initialize the ijr_info
     // this step needs to be done after init_atoms_, because it requires the information of is_atom_on_bgrid
     init_ijr_info_(ucell, gd);
 
     #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(use_gpu_)
     {
-        streams_num_ = PARAM.inp.nstream;  // the default value of num_stream is 4
+        streams_num_ = nstream;  // the default value of num_stream is 4
         const int batch_size = nbz_local;
         init_bgrid_batches_(batch_size);
         gpu_vars_ = std::make_shared<GintGpuVars>(biggrid_info_, ucell, Phi);
@@ -151,7 +151,16 @@ void GintInfo::init_atoms_(int ntype, const Atom* atoms, const Numerical_Orbital
 
 void GintInfo::init_trace_lo_(const UnitCell& ucell, const int nspin)
 {
-    this->trace_lo_ = std::vector<int>(PARAM.globalv.nlocal, -1);
+    int nlocal = 0;
+    for (int it = 0; it < ucell.ntype; ++it)
+    {
+        nlocal += ucell.atoms[it].nw * ucell.atoms[it].na;
+    }
+    if (nspin == 4)
+    {
+        nlocal *= 2;
+    }
+    this->trace_lo_ = std::vector<int>(nlocal, -1);
     this->lgd_ = 0;
     int iat = 0;
     int iw_all = 0;
diff --git a/source/source_hamilt/module_gint/gint_info.h b/source/source_hamilt/module_gint/gint_info.h
index 4996591364..6d1736cb62 100644
--- a/source/source_hamilt/module_gint/gint_info.h
+++ b/source/source_hamilt/module_gint/gint_info.h
@@ -6,7 +6,6 @@
 #include "source_cell/unitcell.h"
 #include "source_cell/atom_spec.h"
 #include "source_hamilt/module_hcontainer/hcontainer.h"
-#include "source_io/module_parameter/parameter.h"
 #include "gint_helper.h"
 #include "big_grid.h"
 #include "gint_atom.h"
@@ -32,7 +31,8 @@ class GintInfo
         int startidx_bx, int startidx_by, int startidx_bz,
         int nbx_local, int nby_local, int nbz_local,
         const Numerical_Orbital* Phi,
-        const UnitCell& ucell, Grid_Driver& gd);
+        const UnitCell& ucell, Grid_Driver& gd,
+        int nspin, bool gamma_only, bool domag, bool use_gpu, int nstream);
 
     ~GintInfo();
 
@@ -63,6 +63,9 @@ class GintInfo
     int get_local_mgrid_num() const { return localcell_info_->get_mgrids_num(); }
     double get_mgrid_volume() const { return meshgrid_info_->get_volume(); }
     GintPrecision get_exec_precision() const { return exec_precision_; }
+    int get_nspin() const { return nspin_; }
+    bool get_domag() const { return domag_; }
+    bool use_gpu() const { return use_gpu_; }
     void set_exec_precision(const GintPrecision precision) { exec_precision_ = precision; }
 
     //=========================================
@@ -72,7 +75,7 @@ class GintInfo
     HContainer<T> get_hr(int npol = 1) const
     {
         auto hr = HContainer<T>(ucell_->nat);
-        if(PARAM.inp.gamma_only)
+        if(gamma_only_)
         {
             hr.fix_gamma();
         }
@@ -135,6 +138,10 @@ class GintInfo
     int lgd_ = 0;
 
     GintPrecision exec_precision_ = GintPrecision::fp64;
+    int nspin_ = 1;
+    bool gamma_only_ = false;
+    bool domag_ = false;
+    bool use_gpu_ = false;
 
     #ifdef __CUDA
     public:
diff --git a/source/source_hamilt/module_gint/gint_interface.cpp b/source/source_hamilt/module_gint/gint_interface.cpp
index 2ebe805f9c..bda93a17e2 100644
--- a/source/source_hamilt/module_gint/gint_interface.cpp
+++ b/source/source_hamilt/module_gint/gint_interface.cpp
@@ -1,6 +1,5 @@
 #include "gint_interface.h"
 #include "source_base/timer.h"
-#include "source_io/module_parameter/parameter.h"
 #include "gint_vl.h"
 #include "gint_vl_metagga.h"
 #include "gint_vl_nspin4.h"
@@ -31,7 +30,7 @@ void cal_gint_vl(
     HContainer<double>* hR)
 {
 #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_vl_gpu gint_vl(vr_eff, hR);
         gint_vl.cal_gint();
@@ -49,7 +48,7 @@ void cal_gint_vl(
     HContainer<std::complex<double>>* hR)
 {
     #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_vl_nspin4_gpu gint_vl_nspin4(vr_eff, hR);
         gint_vl_nspin4.cal_gint();
@@ -67,7 +66,7 @@ void cal_gint_vl_metagga(
     HContainer<double>* hR)
 {
 #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_vl_metagga_gpu gint_vl_metagga(vr_eff, vfork, hR);
         gint_vl_metagga.cal_gint();
@@ -86,7 +85,7 @@ void cal_gint_vl_metagga(
     HContainer<std::complex<double>>* hR)
 {
 #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_vl_metagga_nspin4_gpu gint_vl_metagga_nspin4(vr_eff, vofk, hR);
         gint_vl_metagga_nspin4.cal_gint();
@@ -105,7 +104,7 @@ void cal_gint_rho(
     bool is_dm_symm)
 {
     #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_rho_gpu gint_rho(dm_vec, nspin, rho, is_dm_symm);
         gint_rho.cal_gint();
@@ -135,7 +134,7 @@ void cal_gint_tau(
     double** tau)
 {
     #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_tau_gpu gint_tau(dm_vec, nspin, tau);
         gint_tau.cal_gint();
@@ -157,7 +156,7 @@ void cal_gint_fvl(
     ModuleBase::matrix* svl)
 {
 #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_fvl_gpu gint_fvl_gpu(nspin, vr_eff, dm_vec, isforce, isstress, fvl, svl);
         gint_fvl_gpu.cal_gint();
@@ -180,7 +179,7 @@ void cal_gint_fvl_meta(
     ModuleBase::matrix* svl)
 {
 #ifdef __CUDA
-    if(PARAM.inp.device == "gpu")
+    if(Gint::get_gint_info().use_gpu())
     {
         Gint_fvl_meta_gpu gint_fvl_meta(nspin, vr_eff, vofk, dm_vec, isforce, isstress, fvl, svl);
         gint_fvl_meta.cal_gint();
diff --git a/source/source_hamilt/module_hcontainer/read_hcontainer.cpp b/source/source_hamilt/module_hcontainer/read_hcontainer.cpp
index 853142467b..662d6c1ad5 100644
--- a/source/source_hamilt/module_hcontainer/read_hcontainer.cpp
+++ b/source/source_hamilt/module_hcontainer/read_hcontainer.cpp
@@ -17,8 +17,9 @@ template <typename T>
 Read_HContainer<T>::Read_HContainer(hamilt::HContainer<T>* hcontainer,
                                         const std::string& filename,
                                         const int nlocal,
-                                        const UnitCell* ucell)
-    : _hcontainer(hcontainer), _filename(filename), _nlocal(nlocal), _ucell(ucell)
+                                        const UnitCell* ucell,
+                                        const int rank)
+    : _hcontainer(hcontainer), _filename(filename), _nlocal(nlocal), _ucell(ucell), _rank(rank)
 {
 }
 
@@ -49,7 +50,7 @@ void Read_HContainer<T>::read()
     hamilt::HContainer<T> hcontainer_serial(&pv_serial);
 
 #ifdef __MPI
-    if(GlobalV::MY_RANK == 0)
+    if(this->_rank == 0)
     {
 #endif
     ModuleIO::csrFileReader<T> csr(this->_filename);
diff --git a/source/source_hamilt/module_hcontainer/read_hcontainer.h b/source/source_hamilt/module_hcontainer/read_hcontainer.h
index e8a056d7f6..f248b83b5e 100644
--- a/source/source_hamilt/module_hcontainer/read_hcontainer.h
+++ b/source/source_hamilt/module_hcontainer/read_hcontainer.h
@@ -18,7 +18,8 @@ class Read_HContainer
         hamilt::HContainer<T>* hcontainer, 
         const std::string& filename,
         const int nlocal,
-        const UnitCell* ucell
+        const UnitCell* ucell,
+        int rank
     );
     // read the matrices of all R vectors to the read stream
     void read();
@@ -40,6 +41,7 @@ class Read_HContainer
     std::string _filename;
     int _nlocal;
     const UnitCell* _ucell = nullptr;
+    int _rank = 0;
 };
 
 } // namespace hamilt
diff --git a/source/source_hamilt/module_hcontainer/test/test_hcontainer.cpp b/source/source_hamilt/module_hcontainer/test/test_hcontainer.cpp
index 9e9b5fb91e..991b074a4a 100644
--- a/source/source_hamilt/module_hcontainer/test/test_hcontainer.cpp
+++ b/source/source_hamilt/module_hcontainer/test/test_hcontainer.cpp
@@ -698,8 +698,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_hcontainer/test/test_hcontainer_complex.cpp b/source/source_hamilt/module_hcontainer/test/test_hcontainer_complex.cpp
index c5b9a906fe..c7920e202d 100644
--- a/source/source_hamilt/module_hcontainer/test/test_hcontainer_complex.cpp
+++ b/source/source_hamilt/module_hcontainer/test/test_hcontainer_complex.cpp
@@ -590,8 +590,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_hcontainer/test/test_hcontainer_time.cpp b/source/source_hamilt/module_hcontainer/test/test_hcontainer_time.cpp
index a8f85b654a..0ac84c718b 100644
--- a/source/source_hamilt/module_hcontainer/test/test_hcontainer_time.cpp
+++ b/source/source_hamilt/module_hcontainer/test/test_hcontainer_time.cpp
@@ -142,8 +142,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_surchem/cal_epsilon.cpp b/source/source_hamilt/module_surchem/cal_epsilon.cpp
index 0b7406982a..7e9b151021 100644
--- a/source/source_hamilt/module_surchem/cal_epsilon.cpp
+++ b/source/source_hamilt/module_surchem/cal_epsilon.cpp
@@ -1,6 +1,4 @@
 #include "surchem.h"
-
-#include "source_io/module_parameter/parameter.h"
 void surchem::cal_epsilon(const ModulePW::PW_Basis* rho_basis,
                           const double* PS_TOTN_real,
                           double* epsilon,
@@ -12,8 +10,9 @@ void surchem::cal_epsilon(const ModulePW::PW_Basis* rho_basis,
     double *shapefunc = new double[rho_basis->nrxx];
     for (int i = 0; i < rho_basis->nrxx; i++)
     {
-        shapefunc[i] = erfc((log(std::max(PS_TOTN_real[i], 1e-10) / PARAM.inp.nc_k)) / sqrt(2.0) / PARAM.inp.sigma_k) / 2;
-        epsilon[i] = 1 + (PARAM.inp.eb_k - 1) * shapefunc[i];
+        shapefunc[i] = erfc((log(std::max(PS_TOTN_real[i], 1e-10) / this->parameters_.nc_k))
+                            / sqrt(2.0) / this->parameters_.sigma_k) / 2;
+        epsilon[i] = 1 + (this->parameters_.eb_k - 1) * shapefunc[i];
         epsilon0[i] = 1.0;
     }
     delete[] shapefunc;
diff --git a/source/source_hamilt/module_surchem/cal_vcav.cpp b/source/source_hamilt/module_surchem/cal_vcav.cpp
index 9a56953b6a..86993b757e 100644
--- a/source/source_hamilt/module_surchem/cal_vcav.cpp
+++ b/source/source_hamilt/module_surchem/cal_vcav.cpp
@@ -1,7 +1,6 @@
 #include "source_base/timer.h"
 #include "source_base/parallel_reduce.h"
 #include "source_hamilt/module_xc/xc_functional.h"
-#include "source_io/module_parameter/parameter.h"
 #include "surchem.h"
 
 void lapl_rho(const double& tpiba2,
@@ -41,20 +40,24 @@ void lapl_rho(const double& tpiba2,
 
 // calculates first derivative of the shape function in realspace
 // exp(-(log(n/n_c))^2 /(2 sigma^2)) /(sigma * sqrt(2*pi) )/n
-void shape_gradn(const std::complex<double>* ps_totn, const ModulePW::PW_Basis* rho_basis, double* eprime)
+void shape_gradn(const std::complex<double>* ps_totn,
+                 const ModulePW::PW_Basis* rho_basis,
+                 const double nc_k,
+                 const double sigma_k,
+                 double* eprime)
 {
 
     double *ps_totn_real = new double[rho_basis->nrxx];
     ModuleBase::GlobalFunc::ZEROS(ps_totn_real, rho_basis->nrxx);
     rho_basis->recip2real(ps_totn, ps_totn_real);
 
-    double epr_c = 1.0 / sqrt(ModuleBase::TWO_PI) / PARAM.inp.sigma_k;
+    double epr_c = 1.0 / sqrt(ModuleBase::TWO_PI) / sigma_k;
     double epr_z = 0;
     double min = 1e-10;
 
     for (int ir = 0; ir < rho_basis->nrxx; ir++)
     {
-        epr_z = log(std::max(ps_totn_real[ir], min) / PARAM.inp.nc_k) / sqrt(2) / PARAM.inp.sigma_k;
+        epr_z = log(std::max(ps_totn_real[ir], min) / nc_k) / sqrt(2) / sigma_k;
         eprime[ir] = epr_c * exp(-pow(epr_z, 2)) / std::max(ps_totn_real[ir], min);
     }
 
@@ -108,7 +111,7 @@ void surchem::createcavity(const UnitCell& ucell,
     // gamma * A = exp(-(log(n/n_c))^2 /(2 sigma^2)) /(sigma * sqrt(2*pi) )
     //-------------------------------------------------------------
     double *term1 = new double[rho_basis->nrxx];
-    shape_gradn(ps_totn, rho_basis, term1);
+    shape_gradn(ps_totn, rho_basis, this->parameters_.nc_k, this->parameters_.sigma_k, term1);
 
     //-------------------------------------------------------------
     // quantum surface area, integral of (gamma*A / n) * |\nabla n|
@@ -127,7 +130,7 @@ void surchem::createcavity(const UnitCell& ucell,
     //-------------------------------------------------------------
     // cavitation energy
     //-------------------------------------------------------------
-    this->Acav = PARAM.inp.tau * qs * ucell.omega / rho_basis->nxyz;
+    this->Acav = this->parameters_.tau * qs * ucell.omega / rho_basis->nxyz;
     Parallel_Reduce::reduce_pool(this->Acav);
 
     // double Ael = cal_Acav(ucell, pwb);
@@ -153,7 +156,7 @@ void surchem::createcavity(const UnitCell& ucell,
 
     for (int ir = 0; ir < rho_basis->nrxx; ir++)
     {
-        vwork[ir] = vwork[ir] * term1[ir] * PARAM.inp.tau;
+        vwork[ir] = vwork[ir] * term1[ir] * this->parameters_.tau;
     }
 
     delete[] nablan;
diff --git a/source/source_hamilt/module_surchem/cal_vel.cpp b/source/source_hamilt/module_surchem/cal_vel.cpp
index 14b148a0f0..3736d6fc6a 100644
--- a/source/source_hamilt/module_surchem/cal_vel.cpp
+++ b/source/source_hamilt/module_surchem/cal_vel.cpp
@@ -1,18 +1,21 @@
 #include "source_base/timer.h"
 #include "source_base/parallel_reduce.h"
 #include "source_hamilt/module_xc/xc_functional.h"
-#include "source_io/module_parameter/parameter.h"
 #include "surchem.h"
 
-void shape_gradn(const double* PS_TOTN_real, const ModulePW::PW_Basis* rho_basis, double* eprime)
+void shape_gradn(const double* PS_TOTN_real,
+                 const ModulePW::PW_Basis* rho_basis,
+                 const double nc_k,
+                 const double sigma_k,
+                 double* eprime)
 {
 
-    double epr_c = 1.0 / sqrt(ModuleBase::TWO_PI) / PARAM.inp.sigma_k;
+    double epr_c = 1.0 / sqrt(ModuleBase::TWO_PI) / sigma_k;
     double epr_z = 0;
     double min = 1e-10;
     for (int ir = 0; ir < rho_basis->nrxx; ir++)
     {
-        epr_z = log(std::max(PS_TOTN_real[ir], min) / PARAM.inp.nc_k) / sqrt(2) / PARAM.inp.sigma_k;
+        epr_z = log(std::max(PS_TOTN_real[ir], min) / nc_k) / sqrt(2) / sigma_k;
         eprime[ir] = epr_c * exp(-pow(epr_z, 2)) / std::max(PS_TOTN_real[ir], min);
     }
 }
@@ -21,17 +24,20 @@ void eps_pot(const double* PS_TOTN_real,
              const double& tpiba,
              const std::complex<double>* phi,
              const ModulePW::PW_Basis* rho_basis,
+             const double nc_k,
+             const double sigma_k,
+             const double eb_k,
              double* d_eps,
              double* vwork)
 {
     double *eprime = new double[rho_basis->nrxx];
     ModuleBase::GlobalFunc::ZEROS(eprime, rho_basis->nrxx);
 
-    shape_gradn(PS_TOTN_real, rho_basis, eprime);
+    shape_gradn(PS_TOTN_real, rho_basis, nc_k, sigma_k, eprime);
 
     for (int ir = 0; ir < rho_basis->nrxx; ir++)
     {
-        eprime[ir] = eprime[ir] * (PARAM.inp.eb_k - 1);
+        eprime[ir] = eprime[ir] * (eb_k - 1);
     }
 
     ModuleBase::Vector3<double> *nabla_phi = new ModuleBase::Vector3<double>[rho_basis->nrxx];
@@ -123,7 +129,15 @@ void surchem::cal_vel(const UnitCell& cell,
     this->Ael *= cell.omega / rho_basis->nxyz;
 
     // the 2nd item of tmp_Vel
-    eps_pot(PS_TOTN_real, cell.tpiba, Sol_phi, rho_basis, epsilon, epspot);
+    eps_pot(PS_TOTN_real,
+            cell.tpiba,
+            Sol_phi,
+            rho_basis,
+            this->parameters_.nc_k,
+            this->parameters_.sigma_k,
+            this->parameters_.eb_k,
+            epsilon,
+            epspot);
 
     for (int i = 0; i < rho_basis->nrxx; i++)
     {
diff --git a/source/source_hamilt/module_surchem/sol_force.cpp b/source/source_hamilt/module_surchem/sol_force.cpp
index c998243531..48f36309ad 100644
--- a/source/source_hamilt/module_surchem/sol_force.cpp
+++ b/source/source_hamilt/module_surchem/sol_force.cpp
@@ -1,7 +1,6 @@
 #include "surchem.h"
 #include "source_base/parallel_reduce.h"
 #include "source_base/timer.h"
-#include "source_io/module_parameter/parameter.h"
 
 void surchem::force_cor_one(const UnitCell& cell,
                             const ModulePW::PW_Basis* rho_basis,
@@ -66,7 +65,10 @@ void surchem::force_cor_one(const UnitCell& cell,
 
 }
 
-void surchem::force_cor_two(const UnitCell& cell, const ModulePW::PW_Basis* rho_basis, ModuleBase::matrix& forcesol)
+void surchem::force_cor_two(const UnitCell& cell,
+                            const ModulePW::PW_Basis* rho_basis,
+                            const int nspin,
+                            ModuleBase::matrix& forcesol)
 {
    
     std::complex<double> *n_pseudo = new std::complex<double>[rho_basis->npw];
@@ -80,7 +82,7 @@ void surchem::force_cor_two(const UnitCell& cell, const ModulePW::PW_Basis* rho_
     std::complex<double> *Vel_g = new std::complex<double>[rho_basis->npw];
     ModuleBase::GlobalFunc::ZEROS(Vcav_g, rho_basis->npw);
     ModuleBase::GlobalFunc::ZEROS(Vel_g, rho_basis->npw);
-    for(int is=0; is<PARAM.inp.nspin; is++)
+    for(int is=0; is<nspin; is++)
 	{
 		for (int ir=0; ir<rho_basis->nrxx; ir++)
 		{
@@ -150,6 +152,7 @@ void surchem::force_cor_two(const UnitCell& cell, const ModulePW::PW_Basis* rho_
 void surchem::cal_force_sol(const UnitCell& cell,
                             const ModulePW::PW_Basis* rho_basis,
                             const ModuleBase::matrix& vloc,
+                            const int nspin,
                             ModuleBase::matrix& forcesol)
 {
     ModuleBase::TITLE("surchem", "cal_force_sol");
@@ -160,7 +163,7 @@ void surchem::cal_force_sol(const UnitCell& cell,
     ModuleBase::matrix force2(nat, 3);
     
     force_cor_one(cell, rho_basis, vloc, force1);
-    force_cor_two(cell, rho_basis,force2);
+    force_cor_two(cell, rho_basis, nspin, force2);
     
     int iat = 0;
     for (int it = 0;it < cell.ntype;it++)
diff --git a/source/source_hamilt/module_surchem/surchem.cpp b/source/source_hamilt/module_surchem/surchem.cpp
index 1f2a66ef86..dc175b08fb 100644
--- a/source/source_hamilt/module_surchem/surchem.cpp
+++ b/source/source_hamilt/module_surchem/surchem.cpp
@@ -13,6 +13,11 @@ surchem::surchem()
     qs = 0;
 }
 
+void surchem::set_parameters(const SurchemParameters& parameters)
+{
+    this->parameters_ = parameters;
+}
+
 void surchem::allocate(const int &nrxx, const int &nspin)
 {
     assert(nrxx >= 0);
diff --git a/source/source_hamilt/module_surchem/surchem.h b/source/source_hamilt/module_surchem/surchem.h
index ae1a7db6bb..a6a123a0cd 100644
--- a/source/source_hamilt/module_surchem/surchem.h
+++ b/source/source_hamilt/module_surchem/surchem.h
@@ -11,6 +11,14 @@
 class Parallel_Grid;
 class Structure_Factor;
 
+struct SurchemParameters
+{
+    double eb_k = 80.0;
+    double tau = 1.0798e-05;
+    double sigma_k = 0.6;
+    double nc_k = 0.00037;
+};
+
 class surchem
 {
   public:
@@ -35,6 +43,8 @@ class surchem
 
     void clear();
 
+    void set_parameters(const SurchemParameters& parameters);
+
     void cal_epsilon(const ModulePW::PW_Basis* rho_basis, const double* PS_TOTN_real, double* epsilon, double* epsilon0);
 
     void cal_pseudo(const UnitCell& cell,
@@ -120,6 +130,7 @@ class surchem
     void cal_force_sol(const UnitCell& cell,
                        const ModulePW::PW_Basis* rho_basis,
                        const ModuleBase::matrix& vloc,
+                       int nspin,
                        ModuleBase::matrix& forcesol);
 
     void force_cor_one(const UnitCell& cell,
@@ -127,13 +138,17 @@ class surchem
                        const ModuleBase::matrix& vloc,
                        ModuleBase::matrix& forcesol);
 
-    void force_cor_two(const UnitCell& cell, const ModulePW::PW_Basis* rho_basis, ModuleBase::matrix& forcesol);
+    void force_cor_two(const UnitCell& cell,
+                       const ModulePW::PW_Basis* rho_basis,
+                       int nspin,
+                       ModuleBase::matrix& forcesol);
 
     void get_totn_reci(const UnitCell& cell, const ModulePW::PW_Basis* rho_basis, std::complex<double>* totn_reci);
 
     void induced_charge(const UnitCell& cell, const ModulePW::PW_Basis* rho_basis, double* induced_rho) const;
 
   private:
+    SurchemParameters parameters_;
 };
 
 #endif
\ No newline at end of file
diff --git a/source/source_hamilt/module_surchem/test/cal_epsilon_test.cpp b/source/source_hamilt/module_surchem/test/cal_epsilon_test.cpp
index 0957dd5e4f..8c3271c4b1 100644
--- a/source/source_hamilt/module_surchem/test/cal_epsilon_test.cpp
+++ b/source/source_hamilt/module_surchem/test/cal_epsilon_test.cpp
@@ -7,7 +7,6 @@
 #include "../surchem.h"
 #include "source_base/constants.h"
 #include "source_base/global_function.h"
-#include "source_base/global_variable.h"
 #include "source_basis/module_pw/pw_basis.h"
 
 #include "gmock/gmock.h"
@@ -27,11 +26,6 @@
  *     - calculate the relative permittivity
  */
 
-namespace GlobalC
-{
-ModulePW::PW_Basis* rhopw;
-}
-
 class cal_epsilon_test : public testing::Test
 {
   protected:
@@ -62,8 +56,6 @@ TEST_F(cal_epsilon_test, cal_epsilon)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -102,6 +94,15 @@ TEST_F(cal_epsilon_test, cal_epsilon)
     EXPECT_EQ(PS_TOTN_real[0], 0.274231);
     EXPECT_EQ(epsilon[0], 1);
     EXPECT_NEAR(epsilon[12], 1.00005, doublethreshold);
+
+    SurchemParameters parameters;
+    parameters.eb_k = 40.0;
+    parameters.sigma_k = 0.8;
+    parameters.nc_k = 0.001;
+    solvent_model.set_parameters(parameters);
+    solvent_model.cal_epsilon(&pwtest, PS_TOTN_real, epsilon, epsilon0);
+    const double shape = erfc(log(PS_TOTN_real[12] / parameters.nc_k) / sqrt(2.0) / parameters.sigma_k) / 2;
+    EXPECT_NEAR(epsilon[12], 1.0 + (parameters.eb_k - 1.0) * shape, doublethreshold);
     // EXPECT_EQ(epsilon[19], 43.1009);
     // EXPECT_EQ(epsilon[26], 78.746);
     delete[] PS_TOTN_real;
@@ -113,8 +114,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_surchem/test/cal_pseudo_test.cpp b/source/source_hamilt/module_surchem/test/cal_pseudo_test.cpp
index 2bd8025f7d..e3e237547a 100644
--- a/source/source_hamilt/module_surchem/test/cal_pseudo_test.cpp
+++ b/source/source_hamilt/module_surchem/test/cal_pseudo_test.cpp
@@ -51,8 +51,6 @@ TEST_F(cal_pseudo_test, gauss_charge)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -107,8 +105,6 @@ TEST_F(cal_pseudo_test, cal_pseudo)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -151,8 +147,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_surchem/test/cal_totn_test.cpp b/source/source_hamilt/module_surchem/test/cal_totn_test.cpp
index 72e11c315c..47b5fccb4f 100644
--- a/source/source_hamilt/module_surchem/test/cal_totn_test.cpp
+++ b/source/source_hamilt/module_surchem/test/cal_totn_test.cpp
@@ -49,8 +49,6 @@ TEST_F(cal_totn_test, cal_totn)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -118,8 +116,6 @@ TEST_F(cal_totn_test, induced_charge)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -149,8 +145,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_surchem/test/cal_vcav_test.cpp b/source/source_hamilt/module_surchem/test/cal_vcav_test.cpp
index 9a0af0432d..c8c3e25c69 100644
--- a/source/source_hamilt/module_surchem/test/cal_vcav_test.cpp
+++ b/source/source_hamilt/module_surchem/test/cal_vcav_test.cpp
@@ -54,8 +54,6 @@ TEST_F(cal_vcav_test, lapl_rho)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -161,8 +159,6 @@ TEST_F(cal_vcav_test, createcavity)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -222,8 +218,6 @@ TEST_F(cal_vcav_test, cal_vcav)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -269,8 +263,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_surchem/test/cal_vel_test.cpp b/source/source_hamilt/module_surchem/test/cal_vel_test.cpp
index 9d6b120a62..eb7d3cadf5 100644
--- a/source/source_hamilt/module_surchem/test/cal_vel_test.cpp
+++ b/source/source_hamilt/module_surchem/test/cal_vel_test.cpp
@@ -11,10 +11,6 @@
 #include <fstream>
 #include <iostream>
 
-// Include parameter.h with private access for testing
-#define private public
-#include "source_io/module_parameter/parameter.h"
-#undef private
 /************************************************
  *  unit test of functions in cal_vel.cpp
  ***********************************************/
@@ -60,7 +56,7 @@ TEST_F(cal_vel_test, shape_gradn)
 
     for (int ir = 0; ir < nrxx; ir++)
     {
-        epr_z = log(std::max(PS_TOTN_real[ir], min) / PARAM.inp.nc_k) / sqrt(2) / PARAM.inp.sigma_k;
+        epr_z = log(std::max(PS_TOTN_real[ir], min) / nc_k) / sqrt(2) / sigma_k;
         eprime[ir] = epr_c * exp(-pow(epr_z, 2)) / std::max(PS_TOTN_real[ir], min);
     }
 
@@ -91,8 +87,6 @@ TEST_F(cal_vel_test, eps_pot)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -138,7 +132,7 @@ TEST_F(cal_vel_test, eps_pot)
 
     for (int ir = 0; ir < nrxx; ir++)
     {
-        eprime[ir] = eprime[ir] * (PARAM.input.eb_k - 1);
+        eprime[ir] = eprime[ir] * (80.0 - 1);
     }
 
     ModuleBase::Vector3<double>* nabla_phi = new ModuleBase::Vector3<double>[nrxx];
@@ -189,8 +183,6 @@ TEST_F(cal_vel_test, cal_vel)
 
     // init
 #ifdef __MPI
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
     MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &POOL_WORLD); // in LCAO kpar=1
 #endif
 
@@ -238,8 +230,6 @@ int main(int argc, char** argv)
 {
 #ifdef __MPI
     MPI_Init(&argc, &argv);
-    MPI_Comm_size(MPI_COMM_WORLD, &GlobalV::NPROC);
-    MPI_Comm_rank(MPI_COMM_WORLD, &GlobalV::MY_RANK);
 #endif
 
     testing::InitGoogleTest(&argc, argv);
diff --git a/source/source_hamilt/module_surchem/test/setcell.h b/source/source_hamilt/module_surchem/test/setcell.h
index f7086511cb..9136552ae5 100644
--- a/source/source_hamilt/module_surchem/test/setcell.h
+++ b/source/source_hamilt/module_surchem/test/setcell.h
@@ -10,11 +10,6 @@
 #include "source_base/parallel_grid.h"
 #include "source_pw/module_pwdft/stru_fac.h"
 
-namespace GlobalC
-{
-    ModulePW::PW_Basis* rhopw;
-}
-
 UnitCell::UnitCell(){};
 UnitCell::~UnitCell(){};
 
diff --git a/source/source_hamilt/module_vdw/test/vdw_test.cpp b/source/source_hamilt/module_vdw/test/vdw_test.cpp
index f08c92482f..23435347f6 100644
--- a/source/source_hamilt/module_vdw/test/vdw_test.cpp
+++ b/source/source_hamilt/module_vdw/test/vdw_test.cpp
@@ -203,34 +203,6 @@ TEST_F(vdwd2Test, WrongVdwType)
 }
 
 
-// mohan comment out 2025-04-05 since the original code has been removed.
-// further investigation is needed.
-/*
-TEST_F(vdwd2Test, OneAtomWarning)
-{
-    UnitCell ucell1;
-    stru_ structure1{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.}}}}};
-
-    construct_ucell(structure1,ucell1);
-
-    GlobalV::ofs_warning.open("warning.log");
-    std::ifstream ifs;
-    std::string output;
-
-    std::unique_ptr<vdw::Vdw> vdw_test = vdw::make_vdw(ucell1, input);
-
-    GlobalV::ofs_warning.close();
-	ifs.open("warning.log");
-	getline(ifs,output);
-	EXPECT_THAT(output,testing::HasSubstr("warning"));
-    EXPECT_EQ(vdw_test,nullptr);
-
-    ifs.close();
-    ClearUcell(ucell1);
-}
-*/
-
 TEST_F(vdwd2Test, D2ReadFile)
 {
     input.vdw_C6_file = "c6.txt";
diff --git a/source/source_hamilt/module_xc/libxc_abacus.h b/source/source_hamilt/module_xc/libxc_abacus.h
index 0e7b9f2df7..087a3577fc 100644
--- a/source/source_hamilt/module_xc/libxc_abacus.h
+++ b/source/source_hamilt/module_xc/libxc_abacus.h
@@ -165,7 +165,9 @@ namespace XC_Functional_Libxc
         const std::size_t nrxx,
         const Charge* const chr,
         const std::vector<double> &amag,
-        const ModuleBase::matrix &v);
+        const ModuleBase::matrix &v,
+        const bool domag,
+        const bool domag_z);
 
 
 //-------------------
diff --git a/source/source_hamilt/module_xc/libxc_pot.cpp b/source/source_hamilt/module_xc/libxc_pot.cpp
index 7e69a9a780..12351a4586 100644
--- a/source/source_hamilt/module_xc/libxc_pot.cpp
+++ b/source/source_hamilt/module_xc/libxc_pot.cpp
@@ -3,8 +3,6 @@
 #include "xc_functional.h"
 #include "libxc_abacus.h"
 #include "source_estate/module_charge/charge.h"
-#include "source_base/global_variable.h"
-#include "source_io/module_parameter/parameter.h"
 #include "source_base/parallel_reduce.h"
 #include "source_base/timer.h"
 #include "source_base/tool_title.h"
@@ -180,7 +178,7 @@ std::tuple<double,double,ModuleBase::matrix> XC_Functional_Libxc::v_xc_libxc(		/
 
     if(4==nspin_in)
     {
-        v = XC_Functional_Libxc::convert_v_nspin4(nrxx, chr, amag, v);
+        v = XC_Functional_Libxc::convert_v_nspin4(nrxx, chr, amag, v, domag, domag_z);
     }
 
     //-------------------------------------------------
diff --git a/source/source_hamilt/module_xc/libxc_setup.cpp b/source/source_hamilt/module_xc/libxc_setup.cpp
index 4ee53b9bdd..9f91c9276e 100644
--- a/source/source_hamilt/module_xc/libxc_setup.cpp
+++ b/source/source_hamilt/module_xc/libxc_setup.cpp
@@ -1,7 +1,7 @@
 #ifdef __LIBXC
 
 #include "libxc_abacus.h"
-#include "source_io/module_parameter/parameter.h"
+#include "xc_functional.h"
 #include "source_base/tool_quit.h"
 #include "source_base/formatter.h"
 
@@ -197,11 +197,11 @@ const std::vector<double> in_built_xc_func_ext_params(const int id,
     {
         // finite temperature XC functionals
         case XC_LDA_XC_KSDT:
-            return {PARAM.inp.xc_temperature * 0.5};
+            return {XC_Functional::get_runtime_parameters().xc_temperature * 0.5};
         case XC_LDA_XC_CORRKSDT:
-            return {PARAM.inp.xc_temperature * 0.5};
+            return {XC_Functional::get_runtime_parameters().xc_temperature * 0.5};
         case XC_LDA_XC_GDSMFB:
-            return {PARAM.inp.xc_temperature * 0.5};
+            return {XC_Functional::get_runtime_parameters().xc_temperature * 0.5};
 #ifdef __EXX
         // hybrid functionals
         case XC_HYB_GGA_XC_PBEH:
@@ -228,8 +228,8 @@ const std::vector<double> in_built_xc_func_ext_params(const int id,
             // This is a range-separated hybrid functional with range-separation constant  0.400,
             // and  0.0% short-range and 100.0% long-range exact exchange,
             // using the error function kernel.
-            return { std::stod(PARAM.inp.exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 1.0
-                std::stod(PARAM.inp.exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: -1.0
+            return { std::stod(XC_Functional::get_runtime_parameters().exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 1.0
+                std::stod(XC_Functional::get_runtime_parameters().exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: -1.0
                 hse_omega }; //Range separation constant: 0.4
         }
         case XC_HYB_GGA_XC_LRC_WPBE:  // Long-range corrected PBE (LRC-wPBE) by by Rohrdanz, Martins and Herbert
@@ -237,8 +237,8 @@ const std::vector<double> in_built_xc_func_ext_params(const int id,
             // This is a range-separated hybrid functional with range-separation constant  0.300,
             // and  0.0% short-range and 100.0% long-range exact exchange,
             // using the error function kernel.
-            return { std::stod(PARAM.inp.exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 1.0
-                std::stod(PARAM.inp.exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: -1.0
+            return { std::stod(XC_Functional::get_runtime_parameters().exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 1.0
+                std::stod(XC_Functional::get_runtime_parameters().exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: -1.0
                 hse_omega }; //Range separation constant: 0.3
         }
         case XC_HYB_GGA_XC_LRC_WPBEH:  // Long-range corrected short-range hybrid PBE (LRC-wPBEh) by Rohrdanz, Martins and Herbert
@@ -246,8 +246,8 @@ const std::vector<double> in_built_xc_func_ext_params(const int id,
             // This is a range-separated hybrid functional with range-separation constant  0.200,
             // and 20.0% short-range and 100.0% long-range exact exchange,
             // using the error function kernel.    
-            return { std::stod(PARAM.inp.exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 1.0
-                std::stod(PARAM.inp.exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: -0.8
+            return { std::stod(XC_Functional::get_runtime_parameters().exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 1.0
+                std::stod(XC_Functional::get_runtime_parameters().exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: -0.8
                 hse_omega }; //Range separation constant: 0.2
         }
         case XC_HYB_GGA_XC_CAM_PBEH:  // CAM hybrid screened exchange PBE version
@@ -255,8 +255,8 @@ const std::vector<double> in_built_xc_func_ext_params(const int id,
             // This is a range-separated hybrid functional with range-separation constant  0.700,
             // and 100.0% short-range and 20.0% long-range exact exchange,
             // using the error function kernel.
-            return { std::stod(PARAM.inp.exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 0.2
-                std::stod(PARAM.inp.exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: 0.8
+            return { std::stod(XC_Functional::get_runtime_parameters().exx_fock_alpha[0]),  //Fraction of Hartree-Fock exchange: 0.2
+                std::stod(XC_Functional::get_runtime_parameters().exx_erfc_alpha[0]),  //Fraction of short-range exact exchange: 0.8
                 hse_omega }; //Range separation constant: 0.7
         }
 #endif
@@ -267,13 +267,13 @@ const std::vector<double> in_built_xc_func_ext_params(const int id,
 
 const std::vector<double> external_xc_func_ext_params(const int id)
 {
-    const auto& exch_ext = PARAM.inp.xc_exch_ext;
+    const auto& exch_ext = XC_Functional::get_runtime_parameters().xc_exch_ext;
     if (!exch_ext.empty() && static_cast<int>(exch_ext.front()) == id)
     {
         return {exch_ext.begin() + 1, exch_ext.end()};
     }
 
-    const auto& corr_ext = PARAM.inp.xc_corr_ext;
+    const auto& corr_ext = XC_Functional::get_runtime_parameters().xc_corr_ext;
     if (!corr_ext.empty() && static_cast<int>(corr_ext.front()) == id)
     {
         return {corr_ext.begin() + 1, corr_ext.end()};
diff --git a/source/source_hamilt/module_xc/libxc_tools.cpp b/source/source_hamilt/module_xc/libxc_tools.cpp
index 9023f00901..5449ddbefa 100644
--- a/source/source_hamilt/module_xc/libxc_tools.cpp
+++ b/source/source_hamilt/module_xc/libxc_tools.cpp
@@ -3,7 +3,6 @@
 #include "libxc_abacus.h"
 #include "xc_functional.h"
 #include "source_estate/module_charge/charge.h"
-#include "source_io/module_parameter/parameter.h"
 
 // converting rho (abacus=>libxc)
 std::vector<double> XC_Functional_Libxc::convert_rho(
@@ -32,7 +31,7 @@ XC_Functional_Libxc::convert_rho_amag_nspin4(
 	const std::size_t nrxx,
 	const Charge* const chr)
 {
-	assert(PARAM.inp.nspin==4);
+	assert(nspin==4);
 	std::vector<double> rho(nrxx*nspin);
 	std::vector<double> amag(nrxx);
 	#ifdef _OPENMP
@@ -321,24 +320,26 @@ ModuleBase::matrix XC_Functional_Libxc::convert_v_nspin4(
 	const std::size_t nrxx,
 	const Charge* const chr,
 	const std::vector<double> &amag,
-	const ModuleBase::matrix &v)
+	const ModuleBase::matrix &v,
+    const bool domag,
+    const bool domag_z)
 {
     //assert(nrxx>0);
-	assert(PARAM.inp.nspin==4);
+	constexpr int nspin = 4;
 	constexpr double vanishing_charge = 1.0e-10;
-	ModuleBase::matrix v_nspin4(PARAM.inp.nspin, nrxx);
+	ModuleBase::matrix v_nspin4(nspin, nrxx);
 	for( int ir=0; ir<nrxx; ++ir )
 	{
 		v_nspin4(0,ir) = 0.5 * (v(0,ir)+v(1,ir));
 	}
-	if(PARAM.globalv.domag || PARAM.globalv.domag_z)
+	if(domag || domag_z)
 	{
 		for( int ir=0; ir<nrxx; ++ir )
 		{
 			if ( amag[ir] > vanishing_charge )
 			{
 				const double vs = 0.5 * (v(0,ir)-v(1,ir));
-				for(int ipol=1; ipol<PARAM.inp.nspin; ++ipol)
+				for(int ipol=1; ipol<nspin; ++ipol)
 				{
 					v_nspin4(ipol,ir) = vs * chr->rho[ipol][ir] / amag[ir];
 				}
diff --git a/source/source_hamilt/module_xc/test/CMakeLists.txt b/source/source_hamilt/module_xc/test/CMakeLists.txt
index 004901c245..51fe559c96 100644
--- a/source/source_hamilt/module_xc/test/CMakeLists.txt
+++ b/source/source_hamilt/module_xc/test/CMakeLists.txt
@@ -38,6 +38,7 @@ AddTest(
   ../xc_lda_exch.cpp ../xc_hcth.cpp
   ../../../source_base/matrix.cpp
   ../../../source_base/memory_recorder.cpp
+  ../../../source_base/global_variable.cpp
   ../../../source_base/libm/branred.cpp
   ../../../source_base/libm/sincos.cpp
   ../../../source_base/module_external/blas_connector_base.cpp ../../../source_base/module_external/blas_connector_vector.cpp ../../../source_base/module_external/blas_connector_matrix.cpp
@@ -76,6 +77,7 @@ AddTest(
     ../../../source_base/module_external/blas_connector_base.cpp ../../../source_base/module_external/blas_connector_vector.cpp ../../../source_base/module_external/blas_connector_matrix.cpp
     ../../../source_base/matrix.cpp
     ../../../source_base/memory_recorder.cpp
+  ../../../source_base/global_variable.cpp
     ../../../source_base/timer.cpp
     ../../../source_base/libm/branred.cpp
     ../../../source_base/libm/sincos.cpp
@@ -116,6 +118,7 @@ AddTest(
   ../xc_lda_exch.cpp ../xc_hcth.cpp
   ../../../source_base/matrix.cpp
   ../../../source_base/memory_recorder.cpp
+  ../../../source_base/global_variable.cpp
   ../../../source_base/libm/branred.cpp
   ../../../source_base/libm/sincos.cpp
   ../../../source_base/module_external/blas_connector_base.cpp ../../../source_base/module_external/blas_connector_vector.cpp ../../../source_base/module_external/blas_connector_matrix.cpp
diff --git a/source/source_hamilt/module_xc/test/test_xc.cpp b/source/source_hamilt/module_xc/test/test_xc.cpp
index 72b2f547bc..7085ae1851 100644
--- a/source/source_hamilt/module_xc/test/test_xc.cpp
+++ b/source/source_hamilt/module_xc/test/test_xc.cpp
@@ -1,7 +1,6 @@
 #include "gtest/gtest.h"
 #include "../xc_functional.h"
 #include "../libxc_abacus.h"
-#include "../exx_info.h"
 #include "xctest.h"
 
 /************************************************
@@ -18,19 +17,6 @@ namespace ModuleBase
     void TITLE(const std::string &class_name,const std::string &function_name,bool disable){};
 }
 
-namespace GlobalV
-{
-    std::string BASIS_TYPE = "";
-    bool CAL_STRESS = 0;
-    int CAL_FORCE = 0;
-    int NSPIN = 1;
-}
-
-namespace GlobalC
-{
-	Exx_Info exx_info;
-}
-
 class XCTest_PBE : public XCTest
 {
     protected:
diff --git a/source/source_hamilt/module_xc/test/test_xc1.cpp b/source/source_hamilt/module_xc/test/test_xc1.cpp
index b322d628c6..795868f698 100644
--- a/source/source_hamilt/module_xc/test/test_xc1.cpp
+++ b/source/source_hamilt/module_xc/test/test_xc1.cpp
@@ -1,7 +1,6 @@
 #include "gtest/gtest.h"
 #include "xctest.h"
 #include "../xc_functional.h"
-#include "../exx_info.h"
 
 /************************************************
 *  unit test of set_xc_type
@@ -18,19 +17,6 @@ namespace ModuleBase
     void TITLE(const std::string &class_name,const std::string &function_name,bool disable){};
 }
 
-namespace GlobalV
-{
-    std::string BASIS_TYPE = "";
-    bool CAL_STRESS = 0;
-    int CAL_FORCE = 0;
-    int NSPIN = 1;
-}
-
-namespace GlobalC
-{
-	Exx_Info exx_info;
-}
-
 class XCTest_HSE : public XCTest
 {
     protected:
@@ -79,6 +65,20 @@ TEST_F(XCTest_KSDT, set_xc_type)
     EXPECT_EQ(XC_Functional::get_func_type(),1);
 }
 
+TEST_F(XCTest_KSDT, runtime_parameters)
+{
+    XCFunctionalParameters parameters;
+    parameters.xc_temperature = 0.25;
+    parameters.xc_exch_ext = {101.0, 0.75};
+    parameters.xc_corr_ext = {130.0, 0.5};
+    XC_Functional::set_runtime_parameters(parameters);
+
+    const XCFunctionalParameters& stored = XC_Functional::get_runtime_parameters();
+    EXPECT_DOUBLE_EQ(stored.xc_temperature, 0.25);
+    EXPECT_EQ(stored.xc_exch_ext, parameters.xc_exch_ext);
+    EXPECT_EQ(stored.xc_corr_ext, parameters.xc_corr_ext);
+}
+
 class XCTest_KT2 : public XCTest
 {
     protected:
diff --git a/source/source_hamilt/module_xc/test/test_xc2.cpp b/source/source_hamilt/module_xc/test/test_xc2.cpp
index 5a6401e114..06c1942a4f 100644
--- a/source/source_hamilt/module_xc/test/test_xc2.cpp
+++ b/source/source_hamilt/module_xc/test/test_xc2.cpp
@@ -2,7 +2,6 @@
 #include "xctest.h"
 #include "../xc_functional.h"
 #include "../libxc_abacus.h"
-#include "../exx_info.h"
 /************************************************
 *  unit test of functionals
 ***********************************************/
@@ -17,19 +16,6 @@ namespace ModuleBase
     void TITLE(const std::string &class_name,const std::string &function_name,bool disable){};
 }
 
-namespace GlobalV
-{
-    std::string BASIS_TYPE = "";
-    bool CAL_STRESS = false;
-    int CAL_FORCE = 0;
-    int NSPIN = 2;
-}
-
-namespace GlobalC
-{
-	Exx_Info exx_info;
-}
-
 class XCTest_PBE_SPN : public XCTest
 {
     protected:
diff --git a/source/source_hamilt/module_xc/test/test_xc4.cpp b/source/source_hamilt/module_xc/test/test_xc4.cpp
index 22cb044ae2..67a4faa314 100644
--- a/source/source_hamilt/module_xc/test/test_xc4.cpp
+++ b/source/source_hamilt/module_xc/test/test_xc4.cpp
@@ -2,7 +2,6 @@
 #include "../libxc_abacus.h"
 #include "gtest/gtest.h"
 #include "xctest.h"
-#include "../exx_info.h"
 
 /************************************************
 *  unit test of functionals
@@ -18,19 +17,6 @@ namespace ModuleBase
     void TITLE(const std::string &class_name,const std::string &function_name,bool disable){};
 }
 
-namespace GlobalV
-{
-    std::string BASIS_TYPE = "";
-    bool CAL_STRESS = false;
-    int CAL_FORCE = 0;
-    int NSPIN = 1;
-}
-
-namespace GlobalC
-{
-	Exx_Info exx_info;
-}
-
 class XCTest_SCAN : public XCTest
 {
     protected:
diff --git a/source/source_hamilt/module_xc/test/test_xc6.cpp b/source/source_hamilt/module_xc/test/test_xc6.cpp
index 0e8518fadb..d4fac72332 100644
--- a/source/source_hamilt/module_xc/test/test_xc6.cpp
+++ b/source/source_hamilt/module_xc/test/test_xc6.cpp
@@ -2,7 +2,6 @@
 #include "../libxc_abacus.h"
 #include "gtest/gtest.h"
 #include "xctest.h"
-#include "../exx_info.h"
 #include <cmath>
 #include <iostream>
 #include <iomanip>
@@ -14,19 +13,6 @@ namespace ModuleBase
     void TITLE(const std::string &class_name,const std::string &function_name,bool disable){};
 }
 
-namespace GlobalV
-{
-    std::string BASIS_TYPE = "";
-    bool CAL_STRESS = false;
-    int CAL_FORCE = 0;
-    int NSPIN = 1;
-}
-
-namespace GlobalC
-{
-	Exx_Info exx_info;
-}
-
 class XCTest_SCANL_Laplacian : public XCTest
 {
     protected:
diff --git a/source/source_hamilt/module_xc/test/xc3_mock.h b/source/source_hamilt/module_xc/test/xc3_mock.h
index 83bc9ce08e..0e9e93d627 100644
--- a/source/source_hamilt/module_xc/test/xc3_mock.h
+++ b/source/source_hamilt/module_xc/test/xc3_mock.h
@@ -199,23 +199,6 @@ void TITLE(const std::string& class_name, const std::string& function_name, bool
 
 } // namespace ModuleBase
 
-namespace GlobalV
-{
-std::string BASIS_TYPE = "";
-bool CAL_STRESS = false;
-int CAL_FORCE = 0;
-int NSPIN;
-int NPOL;
-bool DOMAG;
-bool DOMAG_Z;
-std::ofstream ofs_device;
-std::ofstream ofs_running;
-} // namespace GlobalV
-
-namespace GlobalC
-{
-Exx_Info exx_info;
-}
 
 UnitCell::UnitCell() {};
 UnitCell::~UnitCell() {};
diff --git a/source/source_hamilt/module_xc/xc_functional.cpp b/source/source_hamilt/module_xc/xc_functional.cpp
index 6f290a61b4..0f2f7b2c13 100644
--- a/source/source_hamilt/module_xc/xc_functional.cpp
+++ b/source/source_hamilt/module_xc/xc_functional.cpp
@@ -1,5 +1,4 @@
 #include "xc_functional.h"
-#include "source_io/module_parameter/parameter.h"
 #include "source_base/global_function.h"
 #include "source_base/tool_title.h"
 #include "source_base/constants.h"
@@ -20,6 +19,7 @@ bool XC_Functional::need_laplacian = false;
 bool XC_Functional::use_libxc = true;
 double XC_Functional::hybrid_alpha = 0.25;
 double XC_Functional::hse_omega = 0.0;
+XCFunctionalParameters XC_Functional::runtime_parameters;
 std::map<int, double> XC_Functional::scaling_factor_xc = { {1, 1.0} }; // added by jghan, 2024-10-10
 
 void XC_Functional::set_hybrid_alpha(const double alpha_in)
@@ -32,6 +32,11 @@ void XC_Functional::set_hse_omega(const double omega_in)
     hse_omega = omega_in;
 }
 
+void XC_Functional::set_runtime_parameters(const XCFunctionalParameters& parameters)
+{
+    runtime_parameters = parameters;
+}
+
 void XC_Functional::set_xc_first_loop(const UnitCell& ucell)
 {
     ModuleBase::TITLE("XC_Functional", "set_xc_first_loop");
diff --git a/source/source_hamilt/module_xc/xc_functional.h b/source/source_hamilt/module_xc/xc_functional.h
index 66ac7adb21..01b2e4fec2 100644
--- a/source/source_hamilt/module_xc/xc_functional.h
+++ b/source/source_hamilt/module_xc/xc_functional.h
@@ -21,6 +21,16 @@
 #include <map> // added by jghan, 2024-10-10
 #include <vector>
 
+
+struct XCFunctionalParameters
+{
+    double xc_temperature = 0.0;
+    std::vector<std::string> exx_fock_alpha = {"default"};
+    std::vector<std::string> exx_erfc_alpha = {"default"};
+    std::vector<double> xc_exch_ext;
+    std::vector<double> xc_corr_ext;
+};
+
 class XC_Functional
 {
     public:
@@ -85,6 +95,13 @@ class XC_Functional
 
     static void set_hse_omega(const double omega_in);
 
+    static void set_runtime_parameters(const XCFunctionalParameters& parameters);
+
+    static const XCFunctionalParameters& get_runtime_parameters()
+    {
+        return runtime_parameters;
+    };
+
     static double get_hse_omega()
     {
         return hse_omega;
@@ -118,6 +135,8 @@ class XC_Functional
     // hse_omega for HSE functional:
     static double hse_omega;
 
+    static XCFunctionalParameters runtime_parameters;
+
     // added by jghan, 2024-07-07
     // as a scaling factor for different xc-functionals
     static std::map<int, double> scaling_factor_xc;
diff --git a/source/source_hamilt/module_xc/xc_pot.cpp b/source/source_hamilt/module_xc/xc_pot.cpp
index 1f8a3dcd32..dc6a13dabd 100644
--- a/source/source_hamilt/module_xc/xc_pot.cpp
+++ b/source/source_hamilt/module_xc/xc_pot.cpp
@@ -6,7 +6,6 @@
 
 #include "source_base/parallel_reduce.h"
 #include "source_base/timer.h"
-#include "source_io/module_parameter/parameter.h"
 #include "xc_functional.h"
 
 #ifdef __LIBXC
diff --git a/source/source_hamilt/test/rgen_test.cpp b/source/source_hamilt/test/rgen_test.cpp
index c56235fd46..d35cb9ba2b 100644
--- a/source/source_hamilt/test/rgen_test.cpp
+++ b/source/source_hamilt/test/rgen_test.cpp
@@ -48,7 +48,7 @@ TEST_F(RgenTest, ZeroRmax)
     std::vector<int> irr(mxr_test);
     int nrm = 0;
 
-    H_Ewald_pw::rgen(dtau, 0.0, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm);
+    H_Ewald_pw::rgen(dtau, 0.0, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm, 0);
 
     EXPECT_EQ(nrm, 0);
 }
@@ -64,7 +64,7 @@ TEST_F(RgenTest, SimpleCubicNearestNeighbors)
     std::vector<int> irr(mxr_test);
     int nrm = 0;
 
-    H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm);
+    H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm, 0);
 
     EXPECT_EQ(nrm, 18);
 
@@ -98,7 +98,7 @@ TEST_F(RgenTest, SimpleCubicNonZeroDtau)
     std::vector<int> irr(mxr_test);
     int nrm = 0;
 
-    H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm);
+    H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm, 0);
 
     EXPECT_EQ(nrm, 2);
     for (int i = 0; i < nrm; ++i)
@@ -130,7 +130,7 @@ TEST_F(RgenTest, LargeRmaxExceedsOriginalLimit)
     std::vector<int> irr(mxr_test);
     int nrm = 0;
 
-    H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm);
+    H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), mxr_test, nrm, 0);
 
     // Must exceed the old hard-coded limit that caused the crash
     EXPECT_GT(nrm, 200);
diff --git a/source/source_lcao/force_stress_lcao.cpp b/source/source_lcao/force_stress_lcao.cpp
index 895e03c137..165553a918 100644
--- a/source/source_lcao/force_stress_lcao.cpp
+++ b/source/source_lcao/force_stress_lcao.cpp
@@ -408,7 +408,7 @@ void Force_Stress_LCAO<T>::getForceStress(UnitCell& ucell,
     if (PARAM.inp.imp_sol && isforce)
     {
         fsol.create(nat, 3);
-        solvent.cal_force_sol(ucell, rhopw, locpp.vloc, fsol);
+        solvent.cal_force_sol(ucell, rhopw, locpp.vloc, PARAM.inp.nspin, fsol);
     }
 
     //! atomic forces from DFT+U (Quxin version)
diff --git a/source/source_lcao/lcao_set.cpp b/source/source_lcao/lcao_set.cpp
index 49b5b79699..db5c837652 100644
--- a/source/source_lcao/lcao_set.cpp
+++ b/source/source_lcao/lcao_set.cpp
@@ -128,7 +128,8 @@ void LCAO_domain::init_dm_from_file(
             dm_container,
             dmfile,
             PARAM.globalv.nlocal,
-            &ucell
+            &ucell,
+            GlobalV::MY_RANK
         );
         reader_dm.read();
     }
@@ -181,7 +182,7 @@ void LCAO_domain::init_hr_from_file(
     test_file.close();
 
     hmat->set_zero();
-    hamilt::Read_HContainer<TR> reader_hr(hmat, hrfile, PARAM.globalv.nlocal, &ucell);
+    hamilt::Read_HContainer<TR> reader_hr(hmat, hrfile, PARAM.globalv.nlocal, &ucell, GlobalV::MY_RANK);
     reader_hr.read();
     return;
 }
diff --git a/source/source_lcao/module_operator_lcao/op_exx_lcao.cpp b/source/source_lcao/module_operator_lcao/op_exx_lcao.cpp
index 5ed6d8bcf4..619e6462dd 100644
--- a/source/source_lcao/module_operator_lcao/op_exx_lcao.cpp
+++ b/source/source_lcao/module_operator_lcao/op_exx_lcao.cpp
@@ -228,7 +228,7 @@ OperatorEXX<OperatorLCAO<TK, TR>>::OperatorEXX(HS_Matrix_K<TK>* hsk_in,
                 const std::string dmfile
                     = PARAM.globalv.global_readin_dir + "/dmrs" + std::to_string(is + 1) + "_nao.csr";
                 dmR_vec[is] = new hamilt::HContainer<double>(const_cast<Parallel_Orbitals*>(pv));
-                hamilt::Read_HContainer<double> reader_dm(dmR_vec[is], dmfile, PARAM.globalv.nlocal, &ucell);
+                hamilt::Read_HContainer<double> reader_dm(dmR_vec[is], dmfile, PARAM.globalv.nlocal, &ucell, GlobalV::MY_RANK);
                 reader_dm.read();
             }
 
diff --git a/source/source_lcao/test/test_init_dm_from_file.cpp b/source/source_lcao/test/test_init_dm_from_file.cpp
index 502027b341..c0edf36773 100644
--- a/source/source_lcao/test/test_init_dm_from_file.cpp
+++ b/source/source_lcao/test/test_init_dm_from_file.cpp
@@ -149,7 +149,7 @@ TEST_F(InitDMFileTest, Nspin1_ReadSingleFile)
     ASSERT_EQ(dm->_DMR.size(), 1);
 
     hamilt::HContainer<double>* dmr0 = dm->get_DMR_vector()[0];
-    hamilt::Read_HContainer<double> reader(dmr0, "./test_dm_dir/dmrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader(dmr0, "./test_dm_dir/dmrs1_nao.csr", nlocal, &ucell, 0);
     reader.read();
 
     EXPECT_GT(dmr0->size_atom_pairs(), 0);
@@ -183,12 +183,12 @@ TEST_F(InitDMFileTest, Nspin2_ReadTwoFiles)
 
     // Read spin-up
     hamilt::HContainer<double>* dmr0 = dm->get_DMR_vector()[0];
-    hamilt::Read_HContainer<double> reader0(dmr0, "./test_dm_dir/dmrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader0(dmr0, "./test_dm_dir/dmrs1_nao.csr", nlocal, &ucell, 0);
     reader0.read();
 
     // Read spin-down
     hamilt::HContainer<double>* dmr1 = dm->get_DMR_vector()[1];
-    hamilt::Read_HContainer<double> reader1(dmr1, "./test_dm_dir/dmrs2_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader1(dmr1, "./test_dm_dir/dmrs2_nao.csr", nlocal, &ucell, 0);
     reader1.read();
 
     EXPECT_GT(dmr0->size_atom_pairs(), 0);
@@ -264,7 +264,7 @@ TEST_F(InitDMFileTest, HR_Nspin1_ReadSingleFile)
 
     // Read HR from file (same as init_hr_from_file does internally)
     hR.set_zero();
-    hamilt::Read_HContainer<double> reader(&hR, "./test_hr_dir/hrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader(&hR, "./test_hr_dir/hrs1_nao.csr", nlocal, &ucell, 0);
     reader.read();
 
     // Verify data was loaded
@@ -308,13 +308,13 @@ TEST_F(InitDMFileTest, HR_Nspin2_ReadTwoFiles)
     // Read spin-up
     auto* hR_up = create_hcontainer();
     hR_up->set_zero();
-    hamilt::Read_HContainer<double> reader_up(hR_up, "./test_hr_dir/hrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader_up(hR_up, "./test_hr_dir/hrs1_nao.csr", nlocal, &ucell, 0);
     reader_up.read();
 
     // Read spin-down
     auto* hR_down = create_hcontainer();
     hR_down->set_zero();
-    hamilt::Read_HContainer<double> reader_down(hR_down, "./test_hr_dir/hrs2_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader_down(hR_down, "./test_hr_dir/hrs2_nao.csr", nlocal, &ucell, 0);
     reader_down.read();
 
     // Verify both have data
diff --git a/source/source_lcao/test/test_output_hcontainer_consistency.cpp b/source/source_lcao/test/test_output_hcontainer_consistency.cpp
index dba4b8f0aa..f722a58069 100644
--- a/source/source_lcao/test/test_output_hcontainer_consistency.cpp
+++ b/source/source_lcao/test/test_output_hcontainer_consistency.cpp
@@ -132,7 +132,7 @@ TEST_F(OutputHContainerTest, WriteReadConsistency)
     // Read back
     auto* hc_read = create_hcontainer(0.0);
     hc_read->set_zero();
-    hamilt::Read_HContainer<double> reader(hc_read, "./test_ohc_dir/hrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader(hc_read, "./test_ohc_dir/hrs1_nao.csr", nlocal, &ucell, 0);
     reader.read();
 
     // Compare every diagonal element of every atom pair
@@ -223,7 +223,7 @@ TEST_F(OutputHContainerTest, PrecisionParameter)
 
     auto* hc_read = create_hcontainer(0.0);
     hc_read->set_zero();
-    hamilt::Read_HContainer<double> reader(hc_read, "./test_ohc_dir/hrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader(hc_read, "./test_ohc_dir/hrs1_nao.csr", nlocal, &ucell, 0);
     reader.read();
 
     int nw = ucell.atoms[0].nw;
@@ -258,13 +258,13 @@ TEST_F(OutputHContainerTest, Nspin2TwoFileConsistency)
     // Read back spin-up
     auto* hc_read_up = create_hcontainer(0.0);
     hc_read_up->set_zero();
-    hamilt::Read_HContainer<double> reader_up(hc_read_up, "./test_ohc_dir/hrs1_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader_up(hc_read_up, "./test_ohc_dir/hrs1_nao.csr", nlocal, &ucell, 0);
     reader_up.read();
 
     // Read back spin-down
     auto* hc_read_down = create_hcontainer(0.0);
     hc_read_down->set_zero();
-    hamilt::Read_HContainer<double> reader_down(hc_read_down, "./test_ohc_dir/hrs2_nao.csr", nlocal, &ucell);
+    hamilt::Read_HContainer<double> reader_down(hc_read_down, "./test_ohc_dir/hrs2_nao.csr", nlocal, &ucell, 0);
     reader_down.read();
 
     int nw = ucell.atoms[0].nw;
diff --git a/source/source_pw/module_pwdft/force_pw.cpp b/source/source_pw/module_pwdft/force_pw.cpp
index f02afc017b..787eb528bd 100644
--- a/source/source_pw/module_pwdft/force_pw.cpp
+++ b/source/source_pw/module_pwdft/force_pw.cpp
@@ -134,7 +134,7 @@ void Forces<FPTYPE, Device>::cal_force(UnitCell& ucell,
     if (PARAM.inp.imp_sol)
     {
         forcesol.create(this->nat, 3);
-        solvent.cal_force_sol(ucell, rho_basis, locpp->vloc, forcesol);
+        solvent.cal_force_sol(ucell, rho_basis, locpp->vloc, PARAM.inp.nspin, forcesol);
         if (PARAM.inp.test_force)
         {
             ModuleIO::print_force(GlobalV::ofs_running, ucell, "IMP_SOL      FORCE (Ry/Bohr)", forcesol);
@@ -686,7 +686,7 @@ void Forces<FPTYPE, Device>::cal_force_ew(const UnitCell& ucell,
                     {
                         ModuleBase::Vector3<double> d_tau
                             = ucell.atoms[T1].tau[I1] - ucell.atoms[T2].tau[I2];
-                        H_Ewald_pw::rgen(d_tau, rmax, irr.data(), ucell.latvec, ucell.G, r.data(), r2.data(), mxr, nrm);
+                        H_Ewald_pw::rgen(d_tau, rmax, irr.data(), ucell.latvec, ucell.G, r.data(), r2.data(), mxr, nrm, PARAM.inp.test_energy);
 
                         for (int n = 0; n < nrm; n++)
                         {
diff --git a/source/source_pw/module_pwdft/stress_ewa.cpp b/source/source_pw/module_pwdft/stress_ewa.cpp
index 15e966cd44..8188ccfc8d 100644
--- a/source/source_pw/module_pwdft/stress_ewa.cpp
+++ b/source/source_pw/module_pwdft/stress_ewa.cpp
@@ -144,7 +144,7 @@ void Stress_Func<FPTYPE, Device>::stress_ewa(const UnitCell& ucell,
 				//calculate tau[na]-tau[nb]
 				d_tau = ucell.atoms[it].tau[i] - ucell.atoms[jt].tau[j];
 				//generates nearest-neighbors shells 
-				H_Ewald_pw::rgen(d_tau, rmax, irr.data(), ucell.latvec, ucell.G, r.data(), r2.data(), mxr, nrm);
+				H_Ewald_pw::rgen(d_tau, rmax, irr.data(), ucell.latvec, ucell.G, r.data(), r2.data(), mxr, nrm, PARAM.inp.test_energy);
 				for(int nr=0; nr<nrm; nr++)
 				{
 					rr=sqrt(r2[nr]) * ucell.lat0;