Merge pull request #1108 from deepmodeling/HSolver

Refactor: Operator for PW base, accelerated davidson diagonalization method using new Operator for a example

Author: wenfei-li

source/module_base/global_function.h

@@ -18,6 +18,8 @@
 #include "global_variable.h"
 #include "global_function-func_each_2.h"		// Peize Lin add 2016-09-07
 
+#include "blas_connector.h"
+
 namespace ModuleBase
 {
 namespace GlobalFunc
@@ -167,6 +169,23 @@ static inline void DCOPY( const T &a, T &b, const int &dim)
     for (int i=0; i<dim; ++i) b[i] = a[i];
 }
 
+template<typename T>
+inline void COPYARRAY(const T* a, T* b, const long dim);
+
+template<>
+inline void COPYARRAY(const std::complex<double>* a, std::complex<double>* b, const long dim)
+{
+    const int one = 1;
+    zcopy_(&dim, a, &one, b, &one);
+}
+
+template<>
+inline void COPYARRAY(const double* a, double* b, const long dim)
+{
+    const int one = 1;
+    dcopy_(&dim, a, &one, b, &one);
+}
+
 void BLOCK_HERE( const std::string &description );
 
 //==========================================================

source/module_hamilt/CMakeLists.txt

@@ -1,10 +1,12 @@
+add_subdirectory(ks_pw)
+
 add_library(
     hamilt
     OBJECT
     hamilt_pw.cpp
-    ekinetic_pw.cpp
-    veff_pw.cpp
-    nonlocal_pw.cpp
-    meta_pw.cpp
     hamilt_lcao.cpp
+    ks_pw/ekinetic_pw.cpp
+    ks_pw/veff_pw.cpp
+    ks_pw/nonlocal_pw.cpp
+    ks_pw/meta_pw.cpp
 )

source/module_hamilt/ekinetic_pw.cpp

@@ -33,9 +33,8 @@ class Hamilt
 
     int non_first_scf=0;
 
-  protected:
-    // array, save operations from each operators
-    std::vector<Operator*> ops;
+    // first node operator, add operations from each operators
+    Operator* ops = nullptr;
 };
 
 } // namespace hamilt

source/module_hamilt/ekinetic_pw.h

@@ -7,111 +7,98 @@
 #include "src_parallel/parallel_reduce.h"
 #include "src_pw/global.h"
 
-#include "veff_pw.h"
-#include "ekinetic_pw.h"
-#include "meta_pw.h"
-#include "nonlocal_pw.h"
+#include "ks_pw/veff_pw.h"
+#include "ks_pw/ekinetic_pw.h"
+#include "ks_pw/meta_pw.h"
+#include "ks_pw/nonlocal_pw.h"
 
 namespace hamilt
 {
 
 HamiltPW::HamiltPW()
 {
     this->classname = "HamiltPW";
-    this->ops.resize(0);
-    const int npwx = GlobalC::wf.npwx;
-    const int npol = GlobalV::NPOL;
-    const double tpiba = GlobalC::ucell.tpiba;
     const double tpiba2 = GlobalC::ucell.tpiba2;
-    const int* ngk = GlobalC::kv.ngk.data();
+    const double tpiba = GlobalC::ucell.tpiba;
     const int* isk = GlobalC::kv.isk.data();
     const double* gk2 = GlobalC::wfcpw->gk2;
 
     if (GlobalV::T_IN_H)
     {
-        Operator* ekinetic = new EkineticPW(
-            npwx, 
-            npol, 
+        Operator* ekinetic = new Ekinetic<OperatorPW>( 
             tpiba2, 
-            ngk,
-            gk2
+            gk2,
+            GlobalC::wfcpw->npwk_max
         );
-        this->ops.push_back(ekinetic);
+        if(this->ops == nullptr)
+        {
+            this->ops = ekinetic;
+        }
+        else
+        {
+            this->ops->add(ekinetic);
+        }
     }
     if (GlobalV::VL_IN_H)
     {
-        Operator* veff = new VeffPW(
-            npwx,
-            npol,
-            ngk,
+        Operator* veff = new Veff<OperatorPW>(
             isk,
             &(GlobalC::pot.vr_eff),
             GlobalC::wfcpw
         );
-        this->ops.push_back(veff);
+        if(this->ops == nullptr)
+        {
+            this->ops = veff;
+        }
+        else
+        {
+            this->ops->add(veff);
+        }
     }
     if (GlobalV::VNL_IN_H)
     {
-        Operator* nonlocal = new NonlocalPW(
-            npwx,
-            npol,
-            ngk,
+        Operator* nonlocal = new Nonlocal<OperatorPW>(
             isk,
             &GlobalC::ppcell,
             &GlobalC::ucell
         );
-        this->ops.push_back(nonlocal);
+        if(this->ops == nullptr)
+        {
+            this->ops = nonlocal;
+        }
+        else
+        {
+            this->ops->add(nonlocal);
+        }
     }
-    Operator* meta = new MetaPW(
-        npwx,
-        npol,
+    Operator* meta = new Meta<OperatorPW>(
         tpiba,
-        ngk,
         isk,
         &GlobalC::pot.vofk,
         GlobalC::wfcpw
     );
-    this->ops.push_back(meta);
+    if(this->ops == nullptr)
+    {
+        this->ops = meta;
+    }
+    else
+    {
+        this->ops->add(meta);
+    }
 }
 
 HamiltPW::~HamiltPW()
 {
     int index = 0;
-    if (GlobalV::T_IN_H)
-    {
-        delete (EkineticPW*)this->ops[index++];
-    }
-    if (GlobalV::VL_IN_H)
-    {
-        delete (VeffPW*)this->ops[index++];
-    }
-    if (GlobalV::VNL_IN_H)
-    {
-        delete (NonlocalPW*)this->ops[index++];
-    }
-    delete (MetaPW*)this->ops[index++];
+    delete this->ops;
 }
 
 void HamiltPW::updateHk(const int ik)
 {
     ModuleBase::TITLE("HamiltPW","updateHk");
 
-    for(int iter = 0; iter < this->ops.size(); ++iter)
-    {
-        this->ops[iter]->init(ik);
-    }
-
-    return;
-}
+    this->ops->init(ik);
 
-void HamiltPW::hPsi(const std::complex<double> *psi_in, std::complex<double> *hpsi, const size_t size) const
-{
-    ModuleBase::timer::tick("HamiltPW", "h_psi");
-    for(int iter = 0; iter < this->ops.size(); ++iter)
-    {
-        this->ops[iter]->act(psi_in, hpsi, size);
-    }
-    ModuleBase::timer::tick("HamiltPW", "h_psi");
     return;
 }
 
@@ -122,10 +109,7 @@ void HamiltPW::sPsi
     size_t size
 ) const
 {
-    for (size_t i=0; i<size; i++)
-    {
-        spsi[i] = psi[i];
-    }
+    ModuleBase::GlobalFunc::COPYARRAY(psi, spsi, size);
     return;
 }
 

source/module_hamilt/hamilt.h

@@ -16,7 +16,6 @@ class HamiltPW : public Hamilt
     void updateHk(const int ik) override;
 
     // core function: for solving eigenvalues of Hamiltonian with iterative method
-    virtual void hPsi(const std::complex<double> *psi_in, std::complex<double> *hpsi, const size_t size) const override;
     virtual void sPsi(const std::complex<double> *psi_in, std::complex<double> *spsi, const size_t size) const override;
 
   private:

source/module_hamilt/hamilt_pw.cpp

@@ -0,0 +1,8 @@
+add_library(
+    operator_ks_pw
+    OBJECT
+    ekinetic_pw.cpp
+    veff_pw.cpp
+    nonlocal_pw.cpp
+    meta_pw.cpp
+)