Refactor hsolver remove psi dimensions

source/module_hsolver/diago_david.h

@@ -28,6 +28,7 @@ class DiagoDavid : public DiagH<T, Device>
 
     int diag(hamilt::Hamilt<T, Device>* phm_in,   // Pointer to the Hamiltonian object for diagonalization
                       const int dim,              // Dimension of the input matrix psi to be diagonalized
+                      const int nband,            // Number of required eigenpairs
                       const int ldPsi,            // Leading dimension of the psi input
                       psi::Psi<T, Device>& psi,   // Reference to the wavefunction object for eigenvectors
                       Real* eigenvalue_in,        // Pointer to store the resulting eigenvalues
@@ -44,14 +45,10 @@ class DiagoDavid : public DiagH<T, Device>
 
     /// number of searched eigenpairs
     int n_band = 0;
-    /// dimension of the input matrix psi to be diagonalized
-    // int dim = 0;
-    // leading dimension of the matrix data
-    // const int dmx = 0;
     /// dimension of the subspace allowed in Davidson
     int david_ndim = 4;
     /// maximum dimension of the reduced basis set
-    int nbase_x = 0;
+    // int nbase_x = 0;
     /// number of unconverged eigenvalues
     int notconv = 0;
 
@@ -83,6 +80,7 @@ class DiagoDavid : public DiagH<T, Device>
     void cal_grad(hamilt::Hamilt<T, Device>* phm_in,
                   const int& dim,
                   const int& nbase,
+                  const int nbase_x,
                   const int& notconv,
                   psi::Psi<T, Device>& basis,
                   T* hphi,
@@ -93,6 +91,7 @@ class DiagoDavid : public DiagH<T, Device>
 
     void cal_elem(const int& dim,
                   int& nbase,
+                  const int nbase_x,// maximum dimension of the reduced basis set
                   const int& notconv,
                   const psi::Psi<T, Device>& basis,
                   const T* hphi,
@@ -103,6 +102,7 @@ class DiagoDavid : public DiagH<T, Device>
     void refresh(const int& dim,
                  const int& nband,
                  int& nbase,
+                 const int nbase_x,
                  const Real* eigenvalue,
                  const psi::Psi<T, Device>& psi,
                  psi::Psi<T, Device>& basis,
@@ -133,6 +133,7 @@ class DiagoDavid : public DiagH<T, Device>
 
     int diag_mock(hamilt::Hamilt<T, Device>* phm_in,
                    const int dim,
+                   const int nband,
                    const int ldPsi,
                    psi::Psi<T, Device>& psi,
                    Real* eigenvalue_in,

source/module_hsolver/hsolver_pw.cpp

@@ -549,8 +549,9 @@ void HSolverPW<T, Device>::hamiltSolvePsiK(hamilt::Hamilt<T, Device>* hm,
         const Real david_diag_thr = DiagoIterAssist<T, Device>::PW_DIAG_THR;
         const int david_maxiter = DiagoIterAssist<T, Device>::PW_DIAG_NMAX;
         // dimensions
-        const int dim = psi.get_k_first() ? psi.get_current_nbas() : psi.get_nk() * psi.get_nbasis();
-        const int ldPsi = psi.get_k_first() ? psi.get_nbasis() : psi.get_nk() * psi.get_nbasis();
+        const int dim = psi.get_current_nbas();
+        const int nband = psi.get_nbands();
+        const int ldPsi = psi.get_nbasis();
 
         DiagoDavid<T, Device> david(precondition.data(),
                                     GlobalV::PW_DIAG_NDIM,
@@ -559,6 +560,7 @@ void HSolverPW<T, Device>::hamiltSolvePsiK(hamilt::Hamilt<T, Device>* hm,
         DiagoIterAssist<T, Device>::avg_iter
             += static_cast<double>(david.diag(hm,
                                               dim,
+                                              nband,
                                               ldPsi,
                                               psi,
                                               eigenvalue,

source/module_hsolver/test/diago_david_float_test.cpp

@@ -103,10 +103,11 @@ class DiagoDavPrepare
 		start = clock();
 #endif	
 
-		const int dim = phi.get_k_first() ? phi.get_current_nbas() : phi.get_nk() * phi.get_nbasis();
-		const int ldPsi = phi.get_k_first() ? phi.get_nbasis() : phi.get_nk() * phi.get_nbasis();
+		const int dim = phi.get_current_nbas() ;
+		const int nband = phi.get_nbands();
+		const int ldPsi =phi.get_nbasis();
 		std::cout << "ldPsi = " << ldPsi << std::endl;
-		dav.diag(phm, dim, ldPsi, phi, en, eps, maxiter);
+		dav.diag(phm, dim, nband, ldPsi, phi, en, eps, maxiter);
 
 #ifdef __MPI		
 		end = MPI_Wtime();

source/module_hsolver/test/diago_david_real_test.cpp

@@ -102,9 +102,10 @@ class DiagoDavPrepare
         start = clock();
 #endif	
 
-        const int dim = phi.get_k_first() ? phi.get_current_nbas() : phi.get_nk() * phi.get_nbasis();
-        const int ldPsi = phi.get_k_first() ? phi.get_nbasis() : phi.get_nk() * phi.get_nbasis();
-        dav.diag(phm, dim, ldPsi, phi, en, eps, maxiter);
+        const int dim = phi.get_current_nbas();
+        const int nband = phi.get_nbands();
+        const int ldPsi = phi.get_nbasis();
+        dav.diag(phm, dim, nband, ldPsi, phi, en, eps, maxiter);
 
 #ifdef __MPI		
         end = MPI_Wtime();

source/module_hsolver/test/diago_david_test.cpp

@@ -105,9 +105,10 @@ class DiagoDavPrepare
 		start = clock();
 #endif	
 
-		const int dim = phi.get_k_first() ? phi.get_current_nbas() : phi.get_nk() * phi.get_nbasis();
-		const int ldPsi = phi.get_k_first() ? phi.get_nbasis() : phi.get_nk() * phi.get_nbasis();
-		dav.diag(phm, dim, ldPsi, phi, en, eps, maxiter);
+		const int dim = phi.get_current_nbas();
+		const int nband = phi.get_nbands();
+		const int ldPsi = phi.get_nbasis();
+		dav.diag(phm, dim, nband, ldPsi, phi, en, eps, maxiter);
 
 #ifdef __MPI		
 		end = MPI_Wtime();

source/module_hsolver/test/hsolver_pw_sup.h

@@ -163,6 +163,7 @@ DiagoDavid<T, Device>::~DiagoDavid() {
 template <typename T, typename Device>
 int DiagoDavid<T, Device>::diag(hamilt::Hamilt<T, Device>* phm_in,
                                 const int dim,
+                                const int nband,
                                 const int ldPsi,
                                 psi::Psi<T, Device>& psi,
                                 Real* eigenvalue_in,