diff --git a/source/source_estate/elecstate_lcao.cpp b/source/source_estate/elecstate_lcao.cpp
index 2040ee769e..0753e0338e 100644
--- a/source/source_estate/elecstate_lcao.cpp
+++ b/source/source_estate/elecstate_lcao.cpp
@@ -7,6 +7,7 @@
 #include "source_io/module_parameter/parameter.h"
 
 #include "source_lcao/module_gint/gint_interface.h"
+#include "source_lcao/rho_tau_lcao.h"
 
 #include <vector>
 
@@ -84,6 +85,15 @@ void ElecStateLCAO<std::complex<double>>::dm2rho(std::vector<std::complex<double
 }
 
 
+template <typename TK>
+void ElecStateLCAO<TK>::dmToRho(std::vector<hamilt::HContainer<double>*>& dmr,
+                                int nspin,
+                                Charge* chr,
+                                bool skip_charge)
+{
+    LCAO_domain::dm2rho(dmr, nspin, chr, skip_charge);
+}
+
 template class ElecStateLCAO<double>;               // Gamma_only case
 template class ElecStateLCAO<std::complex<double>>; // multi-k case
 
diff --git a/source/source_estate/elecstate_lcao.h b/source/source_estate/elecstate_lcao.h
index cef6653d60..fa28c5bb80 100644
--- a/source/source_estate/elecstate_lcao.h
+++ b/source/source_estate/elecstate_lcao.h
@@ -39,10 +39,23 @@ class ElecStateLCAO : public ElecState
      * @param pexsi_EDM: pointers of energy-weighed density matrix (EDMK) calculated by pexsi, needed by MD, will be
      * stored in DensityMatrix::pexsi_EDM
      */
-	void dm2rho(std::vector<TK*> pexsi_DM, 
-			std::vector<TK*> pexsi_EDM, 
+	void dm2rho(std::vector<TK*> pexsi_DM,
+			std::vector<TK*> pexsi_EDM,
 			DensityMatrix<TK, double>* dm);
 
+    /**
+     * @brief calculate electronic charge density from the density matrix (DMR)
+     *
+     * Thin wrapper over LCAO_domain::dm2rho so that HSolverLCAO delegates the
+     * charge-density calculation through the ElecState interface, mirroring the
+     * plane-wave path (ElecStatePW::psiToRho) and the pexsi branch above. This
+     * keeps the source_lcao dependency out of source_hsolver.
+     */
+    void dmToRho(std::vector<hamilt::HContainer<double>*>& dmr,
+                 int nspin,
+                 Charge* chr,
+                 bool skip_charge = false);
+
 };
 
 template <typename TK>
diff --git a/source/source_hsolver/hsolver_lcao.cpp b/source/source_hsolver/hsolver_lcao.cpp
index a210b9fe81..99a8a4d3a5 100644
--- a/source/source_hsolver/hsolver_lcao.cpp
+++ b/source/source_hsolver/hsolver_lcao.cpp
@@ -35,8 +35,6 @@
 #include "source_hsolver/parallel_k2d.h"
 #include "source_io/module_parameter/parameter.h"
 
-#include "source_lcao/rho_tau_lcao.h" // mohan add 20251024
-
 namespace hsolver
 {
 
@@ -103,7 +101,9 @@ void HSolverLCAO<TK, Device>::solve(hamilt::Hamilt<TK>* pHamilt,
         if (!skip_charge)
         {
             // compute charge density from density matrix, mohan update 20251024
-            LCAO_domain::dm2rho(dm.get_DMR_vector(), nspin, &chr);
+            // delegate to ElecStateLCAO to keep the source_lcao dependency out of
+            // source_hsolver (mirrors the pexsi branch below and the PW psiToRho path)
+            dynamic_cast<elecstate::ElecStateLCAO<TK>*>(pes)->dmToRho(dm.get_DMR_vector(), nspin, &chr);
         }
         else
         {