Merge branch 'develop' of https://github.com/deepmodeling/abacus-develop into develop

Author: pxlxingliang

source/src_pw/diago_cg.cpp

@@ -1,13 +1,15 @@
 #include "diago_cg.h"
-#include "global.h"
 #include "../src_parallel/parallel_reduce.h"
 #include "../module_base/timer.h"
+#include "module_base/constants.h"
+#include "module_base/blas_connector.h"
 
 int Diago_CG::moved = 0;
 
 
-Diago_CG::Diago_CG()
+Diago_CG::Diago_CG(Hamilt_PW* phamilt)
 {
+    this->hpw = phamilt;
     test_cg=0;
 }
 Diago_CG::~Diago_CG() {}
@@ -66,10 +68,10 @@ void Diago_CG::diag
         if (test_cg>2) GlobalV::ofs_running << "Diagonal Band : " << m << std::endl;
         for (int i=0; i<dim; i++) phi_m[i] = phi(m, i);
 
-        GlobalC::hm.hpw.s_1psi(dim, phi_m, sphi); // sphi = S|psi(m)>
+        this->hpw->s_1psi(dim, phi_m, sphi); // sphi = S|psi(m)>
         this->schmit_orth(dim, dmx, m, phi, sphi, phi_m);
 
-        GlobalC::hm.hpw.h_1psi(dim , phi_m, hphi, sphi);
+        this->hpw->h_1psi(dim , phi_m, hphi, sphi);
 
         e[m] = this->ddot_real(dim, phi_m, hphi );
 
@@ -210,7 +212,7 @@ void Diago_CG::orthogonal_gradient( const int &dim, const int &dmx,
     if (test_cg==1) ModuleBase::TITLE("Diago_CG","orthogonal_gradient");
     //ModuleBase::timer::tick("Diago_CG","orth_grad");
 
-    GlobalC::hm.hpw.s_1psi(dim , g, sg);
+    this->hpw->s_1psi(dim , g, sg);
     int inc=1;
     //<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
     //qianrui replace 2021-3-15
@@ -341,7 +343,7 @@ bool Diago_CG::update_psi(
 {
     if (test_cg==1) ModuleBase::TITLE("Diago_CG","update_psi");
     //ModuleBase::timer::tick("Diago_CG","update");
-    GlobalC::hm.hpw.h_1psi(dim, cg, hcg, scg);
+    this->hpw->h_1psi(dim, cg, hcg, scg);
     cg_norm = sqrt( this->ddot_real(dim, cg, scg) );
 
     if (cg_norm < 1.0e-10 ) return 1;
@@ -476,7 +478,7 @@ void Diago_CG::schmit_orth
     {
         psi_m[ig] /= psi_norm;
     }
-    GlobalC::hm.hpw.s_1psi(dim, psi_m, sphi); // sphi = S|psi(m)>
+    this->hpw->s_1psi(dim, psi_m, sphi); // sphi = S|psi(m)>
 
     delete [] lagrange ;
     //ModuleBase::timer::tick("Diago_CG","schmit_orth");

source/src_pw/diago_cg.h

@@ -5,11 +5,13 @@
 #include "../module_base/global_variable.h"
 #include "../module_base/complexmatrix.h"
 
+#include "src_pw/hamilt_pw.h"
+
 class Diago_CG
 {
 	public:
 
-    Diago_CG();
+    Diago_CG(Hamilt_PW* phamilt);
     ~Diago_CG();
 
     static int moved;
@@ -33,7 +35,7 @@ class Diago_CG
         int &notconv,
         double &avg_iter);
 
-    static void schmit_orth(
+    void schmit_orth(
         const int &dim,
         const int &dmx,
         const int &end,
@@ -44,6 +46,9 @@ class Diago_CG
 
 	private:
 
+    /// temp operator pointer
+    Hamilt_PW* hpw;
+
     int test_cg;
 
     void calculate_gradient(

source/src_pw/diago_david.cpp

@@ -1,11 +1,14 @@
 #include "diago_david.h"
 #include "diago_cg.h"
-#include "global.h"
 #include "../src_parallel/parallel_reduce.h"
 #include "../module_base/timer.h"
+#include "module_base/constants.h"
+#include "module_base/blas_connector.h"
+#include "module_base/lapack_connector.h"
 
-Diago_David::Diago_David()
+Diago_David::Diago_David(Hamilt_PW* phamilt)
 {
+    this->hpw = phamilt;
     test_david = 2;
     // 1: check which function is called and which step is executed
     // 2: check the eigenvalues of the result of each iteration
@@ -86,7 +89,7 @@ void Diago_David::diag
 
         this->SchmitOrth(npw, nband, m, basis, psi_m, spsi);
 
-        GlobalC::hm.hpw.h_1psi(npw, psi_m, hpsi, spsi);
+        this->hpw->h_1psi(npw, psi_m, hpsi, spsi);
 
         // basis(m) = psi_m, hp(m) = H |psi_m>, sp(m) = S |psi_m>
         for ( int ig = 0; ig < npw; ig++ )
@@ -259,7 +262,7 @@ void Diago_David::cal_grad
 */
         this->SchmitOrth(npw, nbase+notconv, nbase+m, basis, ppsi, spsi);
 
-        GlobalC::hm.hpw.h_1psi(npw, ppsi, hpsi, spsi);
+        this->hpw->h_1psi(npw, ppsi, hpsi, spsi);
 
         for ( int ig = 0; ig < npw; ig++ )
         {
@@ -534,7 +537,7 @@ void Diago_David::SchmitOrth
     assert(m >= 0);
     assert(m < n_band);
 
-    GlobalC::hm.hpw.s_1psi(npw, psi_m, spsi);
+    this->hpw->s_1psi(npw, psi_m, spsi);
 
     std::complex<double>* lagrange = new std::complex<double>[m+1];
     ModuleBase::GlobalFunc::ZEROS( lagrange, m+1 );
@@ -589,7 +592,7 @@ void Diago_David::SchmitOrth
         }
     }
 
-    GlobalC::hm.hpw.s_1psi(npw, psi_m, spsi);
+    this->hpw->s_1psi(npw, psi_m, spsi);
 
     delete[] lagrange;
     ModuleBase::timer::tick("Diago_David","SchmitOrth");

source/src_pw/diago_david.h

@@ -14,14 +14,16 @@
 #include "../module_base/global_variable.h"
 #include "../module_base/complexmatrix.h"
 
+#include "src_pw/hamilt_pw.h"
+
 class Diago_David
 {
 public:
 
-    Diago_David();
+    Diago_David(Hamilt_PW* phamilt);
     ~Diago_David();
 
-    static void SchmitOrth(
+    void SchmitOrth(
         const int& npw,
         const int n_band,
         const int m,
@@ -52,6 +54,8 @@ class Diago_David
 
 private:
 
+    Hamilt_PW* hpw;
+
     int test_david;
 
     void cal_grad(

source/src_pw/hamilt.cpp

@@ -82,7 +82,7 @@ void Hamilt::diagH_pw(
 
                     avg_iter += 1.0;
                 }
-                Diago_CG cg;
+                Diago_CG cg(&GlobalC::hm.hpw);
 
 				bool reorder = true;
 
@@ -106,7 +106,7 @@ void Hamilt::diagH_pw(
             }
 	   		else if(GlobalV::KS_SOLVER=="dav")
         	{
-				Diago_David david;
+				Diago_David david(&GlobalC::hm.hpw);
 				if(GlobalV::NPOL==1)
 				{
 					david.diag(GlobalC::wf.evc[ik0], GlobalC::wf.ekb[ik], GlobalC::kv.ngk[ik],