Refactor: Relate to PEXSI-related changes caused by refactor #4767 (#4858)

* feat pexsi

* fix : diag not completed

* feat

* feat: pexsi hsolver

* CMake building implemented

* Works

* adapt to the new container

* Turn off USE_PEXSI

* Update LibRI to 553c91c

* modify include files

* namespace-ize

* new inputs added

* Configure Makefile Compiling, fix typos

* Fix Makefile Intel toolchains compile errors

* Fix even more PEXSI related Makefile compiling issues

* Modify inputs and update to latest version (#2)

* run INPUT.Default() in every process in InputParaTest (#3490)

Co-authored-by: kirk0830 <67682086+kirk0830@users.noreply.github.com>

* add blas support for FindLAPACK.cmake (#3497)

* more unittest of QO: towards orbital selection (#3499)

* Fix: fix bug in mulliken charge calculation (#3503)

* fix phase

* fix case test

* Refactor: namespace Conv_Coulomb_Pot_K (#3446)

* Refactor: namespace Conv_Coulomb_Pot_K

* Refactor: namespace Conv_Coulomb_Pot_K

---------

Co-authored-by: wqzhou <33364058+WHUweiqingzhou@users.noreply.github.com>

* enable the computation of all zeros in one function call (#3449)

Co-authored-by: wqzhou <33364058+WHUweiqingzhou@users.noreply.github.com>

* replace ios.eof() by ios.good() to avoid meeting badbit and failbit in reading STRU (#3506)

* Build: add ccache to accelerate the testing process (#3509)

* Build: add ccache to accelerate the testing process

* Update test.yml

* Update test.yml

* Update test.yml

* Docs: to avoid the misunderstanding in docs (#3518)

* to avoid the misunderstanding in docs

* Update docs/quick_start/hands_on.md

Co-authored-by: Chun Cai <amoycaic@gmail.com>

---------

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Docs: fix a missing depencency in conda build env (#3508)

* Feature: Add ENABLE_RAPIDJSON option to control the output of abacus.json (#3519)

Add ENABLE_RAPIDJSON option to control the output of abacus.json

* Feature: add python wrapper for math sphbes (#3475)

* recommit for review

* add python wrapper

* remove timer since performace tests add

* Feature: support segment split in kline mode in KPT file and `out_band` band output precision control, `8` as default (#3493)

* add precision control

* correct serial version of nscf_band function

* fix issue 3482

* update unit and integrated test

* update document

* correct unittest and make compatible with false and true

* fix: bug in Autotest.sh when result.ref has no totaltimeref (#3523)

* Fix : unit test of module_xc (#3524)

* Fix: omit small magnetic moments to avoid numerical instability (#3530)

* update deltalambda

* avoid numerical error in orbMulP

* add constrain on Mi

* change case reference value

* Fix: fix multiple compiler warnings (#3515)

* Fix: add noreturn attribute to warning_quit

* Add type conversion

* fix string literal

* fix small number trunctuation

* Fix system call returned value not checked

* fix missing braket

* Refactor parameter_pool.cpp and parameter_pool.h

* remove duplicated return statements

* Change WARNING_QUIT occurances in tests

* Add warning message to help debug UT

* output the default precision flag (#3496)

Co-authored-by: kirk0830 <67682086+kirk0830@users.noreply.github.com>

* Build: Improving CMake performance for finding LibXC and ELPA (#3478)

* Fix for finding LibXC and ELPA

* For compatibility to previous routines

* syntax fix for FindELPA.cmake

* Update cmake/FindELPA.cmake

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Using CMake interface as default for finding LibXC

* update docs

* fix for FindLibxc: changing imcompatible if statement

* fix for FindLibxc: changing imcompatible if statement

* fix for FindLibxc: changing imcompatible if statement

* update docs for installing pkg-config

* Update FindLibxc.cmake

* Update FindLibxc.cmake

* remove previous LibXC routine in CMakeLists.txt

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Update easy_install.md with Makefile-built LibXC supported

* Update easy_install.md to include different behavior in different version on finding ELPA

---------

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Docs: correct some docs about mp2 smearing method (#3533)

* correct some docs about mp2 smearing method

* add docs about mv method

* Feature : printing band density (#3501)

Co-authored-by: wenfei-li <liwenfei@gmail.com>
Co-authored-by: wqzhou <33364058+WHUweiqingzhou@users.noreply.github.com>

* add some docs for PR#3501 (#3537)

* Feature: enable restart charge density mixing during SCF (#3542)

* add a new parameter mixing_restart

* do not update rho if iter==mixing_restart

* do not update rho if iter==mixing_restart-1

* reset mix and rho_mdata if iter==mixing_restart

* fix SCF exit directly since drho=0 if iter=GlobalV::MIXING_RESTART

* re-set_mixing in eachiterinit for PW and LCAO

* enable SCF restarts in esolver_ks::RUN

* add some UnitTests

* add some Docs

* new inputs added

* Update input-main.md (#3551)

Solve the format problem mentioned in issue 3543

* Build: fix compatibility issue against toolchain install (#3540)

* Fix for finding LibXC and ELPA

* For compatibility to previous routines

* syntax fix for FindELPA.cmake

* Update cmake/FindELPA.cmake

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Using CMake interface as default for finding LibXC

* update docs

* fix for FindLibxc: changing imcompatible if statement

* fix for FindLibxc: changing imcompatible if statement

* fix for FindLibxc: changing imcompatible if statement

* update docs for installing pkg-config

* Update FindLibxc.cmake

* Update FindLibxc.cmake

* remove previous LibXC routine in CMakeLists.txt

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Update easy_install.md with Makefile-built LibXC supported

* Update easy_install.md to include different behavior in different version on finding ELPA

* fix compatibility issue against toolchain

* Change default ELPA install routine to old one

---------

Co-authored-by: Chun Cai <amoycaic@gmail.com>

* Test: Configure performance tests for math libraries (#3511)

* add performace test of sphbes functions.

* fix benchmark cmake errors

* add dependencies for docker

* update docs

* add performance tests for sphbes

* add google benchmark

* rewrite benchmark tests in fixtures

* disable internal testing in benchmark

* merge benchmark into integration test

---------

Co-authored-by: StarGrys <771582678@qq.com>

* Configure Makefile Compiling, fix typos

* Fix Makefile Intel toolchains compile errors

* Fix even more PEXSI related Makefile compiling issues

* Update hsolver_pw.cpp (#3556)

when use_uspp==false, overlap matrix should be E.

* Fix: cuda build target (#3276)

* Fix: cuda buid target

* Update CMakeLists.txt

---------

Co-authored-by: Denghui Lu <denghuilu@pku.edu.cn>

---------

Co-authored-by: wqzhou <33364058+WHUweiqingzhou@users.noreply.github.com>
Co-authored-by: kirk0830 <67682086+kirk0830@users.noreply.github.com>
Co-authored-by: Haozhi Han <haozhi.han@outlook.com>
Co-authored-by: Zhao Tianqi <hongriTianqi@users.noreply.github.com>
Co-authored-by: PeizeLin <78645006+PeizeLin@users.noreply.github.com>
Co-authored-by: jinzx10 <jzx016@hotmail.com>
Co-authored-by: Chun Cai <amoycaic@gmail.com>
Co-authored-by: Peng Xingliang <91927439+pxlxingliang@users.noreply.github.com>
Co-authored-by: Jie Li <76780849+jieli-matrix@users.noreply.github.com>
Co-authored-by: Wenfei Li <38569667+wenfei-li@users.noreply.github.com>
Co-authored-by: Denghui Lu <denghuilu@pku.edu.cn>
Co-authored-by: YI Zeping <18586016708@163.com>
Co-authored-by: wenfei-li <liwenfei@gmail.com>
Co-authored-by: jingan-181 <78459531+jingan-181@users.noreply.github.com>
Co-authored-by: StarGrys <771582678@qq.com>
Co-authored-by: Haozhi Han <haozhi.han@stu.pku.edu.cn>

* Revert "Modify inputs and update to latest version"

* Update FindPEXSI.cmake to fix Comments

* Fix CI errors

* Fix CI Errors and Merge with Upstream

* Resolve Pull Request Reviews

* Fix parallel communication related issue

* Fix vars in Makefile.vars, add input tests and comments for pexsi vars

* Fix nspin > 1 cases

* Improvement: take calculated mu as new initial guess, may slightly improve performance

* Fix mistakes in the last commit

* Fix: params and features
- set default pexsi_temp
- fix md in pexsi

* fix empty lines

* Fix: move params to pexsi_solver, rename USE_PEXSI to ENABLE_PEXSI

* Docs: added docs for pexsi inputs

* Fix unit test issues in input_conv

* Change default pexsi_npole from 80 to 40

* Place pexsi_EDM in DensityMatrix, set size of pexsi_dm = 1 when GlobalV::NSPIN==4, and add comments for dmToRho

* An unit test added for DiagoPexsi

* modify for changed gint interface

* correct nspin related behaviors

* add efermi passthrough

* Revert "add efermi passthrough"

This reverts commit d7b402d.

* commits to resolve conversations related to codes

* DM and EDM pointers in pexsi now handled by diagopexsi, and copying h s matrices no longer needed

* add pexsi examples

* fix pexsi unit test (original version shouldn't run)

* add building docs for pexsi

* set cxx standard to c++14, which is required in make_unique

* Fix: Fix typo related to pexsi

* update to PPEXSIDFTDriver2

* default npoints to 1, so single core pexsi will work

* Fix Compile errors

* refactor to abandon `pdiagh`

* Fix mu_buffer and nspin

---------

Co-authored-by: zhangzhihao <1900017707@pku.edu.cn>
Co-authored-by: zhangzh-pku <64026312+zhangzh-pku@users.noreply.github.com>
Co-authored-by: wqzhou <33364058+WHUweiqingzhou@users.noreply.github.com>
Co-authored-by: kirk0830 <67682086+kirk0830@users.noreply.github.com>
Co-authored-by: Haozhi Han <haozhi.han@outlook.com>
Co-authored-by: Zhao Tianqi <hongriTianqi@users.noreply.github.com>
Co-authored-by: PeizeLin <78645006+PeizeLin@users.noreply.github.com>
Co-authored-by: jinzx10 <jzx016@hotmail.com>
Co-authored-by: Chun Cai <amoycaic@gmail.com>
Co-authored-by: Peng Xingliang <91927439+pxlxingliang@users.noreply.github.com>
Co-authored-by: Jie Li <76780849+jieli-matrix@users.noreply.github.com>
Co-authored-by: Wenfei Li <38569667+wenfei-li@users.noreply.github.com>
Co-authored-by: Denghui Lu <denghuilu@pku.edu.cn>
Co-authored-by: YI Zeping <18586016708@163.com>
Co-authored-by: wenfei-li <liwenfei@gmail.com>
Co-authored-by: jingan-181 <78459531+jingan-181@users.noreply.github.com>
Co-authored-by: StarGrys <771582678@qq.com>
Co-authored-by: Haozhi Han <haozhi.han@stu.pku.edu.cn>
Co-authored-by: Mohan Chen <mohan.chen.chen.mohan@gmail.com>

Author: 18 people

source/module_hsolver/diago_pexsi.cpp

@@ -13,6 +13,9 @@ typedef hamilt::MatrixBlock<std::complex<double>> matcd;
 
 namespace hsolver
 {
+template <typename T>
+std::vector<double> DiagoPexsi<T>::mu_buffer;
+
 template <typename T>
 DiagoPexsi<T>::DiagoPexsi(const Parallel_Orbitals* ParaV_in)
 {
@@ -76,7 +79,7 @@ void DiagoPexsi<double>::diag(hamilt::Hamilt<double>* phm_in, psi::Psi<double>&
     this->totalFreeEnergy = this->ps->get_totalFreeEnergy();
     this->totalEnergyH = this->ps->get_totalEnergyH();
     this->totalEnergyS = this->ps->get_totalEnergyS();
-    this->mu_buffer[ik] = this->ps->get_mu();
+    mu_buffer[ik] = this->ps->get_mu();
 }
 
 template <>

source/module_hsolver/diago_pexsi.h

@@ -16,7 +16,7 @@ class DiagoPexsi : public DiagH<T>
 {
   private:
     using Real = typename GetTypeReal<T>::type;
-    std::vector<double> mu_buffer;
+    static std::vector<double> mu_buffer;
 
   public:
     DiagoPexsi(const Parallel_Orbitals* ParaV_in);

source/module_hsolver/hsolver_lcao.cpp

@@ -52,24 +52,23 @@ void HSolverLCAO<T, Device>::solve(hamilt::Hamilt<T>* pHamilt,
     // select the method of diagonalization
     this->method = method_in;
 
-#ifdef __PEXSI
+#ifdef __PEXSI // other purification methods should follow this routine
     if (this->method == "pexsi")
     {
-        if (this->pdiagh != nullptr)
-        {
-            if (this->pdiagh->method != this->method)
-            {
-                delete[] this->pdiagh;
-                this->pdiagh = nullptr;
-            }
-            auto tem = dynamic_cast<DiagoPexsi<T>*>(this->pdiagh);
-        }
-        if (this->pdiagh == nullptr) {
-            DiagoPexsi<T>* tem = new DiagoPexsi<T>(this->ParaV);
-            this->pdiagh = tem;
-            // this->pdiagh = dynamic_cast<DiagoPexsi<T>*>(tem);
-            this->pdiagh->method = this->method;
+        DiagoPexsi<T> pe(ParaV);
+        for (int ik = 0; ik < psi.get_nk(); ++ik) {
+            /// update H(k) for each k point
+            pHamilt->updateHk(ik);
+            psi.fix_k(ik);
+            // solve eigenvector and eigenvalue for H(k)
+            pe.diag(pHamilt, psi, nullptr);
         }
+        auto _pes = dynamic_cast<elecstate::ElecStateLCAO<T>*>(pes);
+        pes->f_en.eband = pe.totalFreeEnergy;
+        // maybe eferm could be dealt with in the future
+        _pes->dmToRho(pe.DM, pe.EDM);
+        ModuleBase::timer::tick("HSolverLCAO", "solve");
+        return;
     }
 #endif
 
@@ -100,6 +99,7 @@ void HSolverLCAO<T, Device>::solve(hamilt::Hamilt<T>* pHamilt,
             pHamilt->updateHk(ik);
 
             psi.fix_k(ik);
+            std::cout << pes->nelec_spin[0] << " " << pes->nelec_spin[1] << std::endl;
 
             // solve eigenvector and eigenvalue for H(k)
             this->hamiltSolvePsiK(pHamilt, psi, &(pes->ekb(ik, 0)));
@@ -126,22 +126,7 @@ void HSolverLCAO<T, Device>::solve(hamilt::Hamilt<T>* pHamilt,
 
     // calculate charge by psi
     // called in scf calculation
-#ifdef __PEXSI
-    if (this->method == "pexsi")
-    {
-        DiagoPexsi<T> tem = dynamic_cast<DiagoPexsi<T>*>(this->pdiagh);
-        if (tem == nullptr)
-            ModuleBase::WARNING_QUIT("HSolverLCAO", "pexsi need debug!");
-        elecstate::ElecStateLCAO<T>* _pes
-            = dynamic_cast<elecstate::ElecStateLCAO<T>*>(pes);
-        pes->f_en.eband = tem->totalFreeEnergy;
-        // maybe eferm could be dealt with in the future
-        _pes->dmToRho(tem->DM, tem->EDM);
-    } else
-#endif
-    {
-        pes->psiToRho(psi);
-    }
+    pes->psiToRho(psi);
     ModuleBase::timer::tick("HSolverLCAO", "solve");
 }
 
@@ -188,12 +173,6 @@ void HSolverLCAO<T, Device>::hamiltSolvePsiK(hamilt::Hamilt<T>* hm, psi::Psi<T>&
         la.diag(hm, psi, eigenvalue);
 #else
         ModuleBase::WARNING_QUIT("HSolverLCAO::solve", "This method of DiagH is not supported!");
-#endif
-#ifdef __PEXSI
-    else if (this->method == "pexsi")
-    {
-        this->pdiagh->diag(hm, psi, eigenvalue);
-    }
 #endif
     }
     else

source/module_hsolver/module_pexsi/simple_pexsi.cpp

@@ -110,11 +110,11 @@ int loadPEXSIOption(MPI_Comm comm,
     int_para[11] = pexsi::PEXSI_Solver::pexsi_symm;
     int_para[12] = pexsi::PEXSI_Solver::pexsi_trans;
     int_para[13] = pexsi::PEXSI_Solver::pexsi_method;
-    int_para[14] = 2;
+    int_para[14] = 1;
     int_para[15] = 0;
     int_para[16] = pexsi::PEXSI_Solver::pexsi_nproc_pole;
 
-    double_para[0] = GlobalV::NSPIN; // pexsi::PEXSI_Solver::pexsi_spin;
+    double_para[0] = 2;//GlobalV::NSPIN; // pexsi::PEXSI_Solver::pexsi_spin;
     double_para[1] = pexsi::PEXSI_Solver::pexsi_temp;
     double_para[2] = pexsi::PEXSI_Solver::pexsi_gap;
     double_para[3] = pexsi::PEXSI_Solver::pexsi_delta_e;
@@ -298,15 +298,15 @@ int simplePEXSI(MPI_Comm comm_PEXSI,
         int numTotalInertiaIter; // Number of total inertia[out]
         // LiuXh modify 2021-04-29, add DONE(ofs_running,"xx") for test
         ModuleBase::timer::tick("Diago_LCAO_Matrix", "PEXSIDFT");
-        PPEXSIDFTDriver(plan,                 // PEXSI plan[in]
-                        options,              // PEXSI Options[in]
+        PPEXSIDFTDriver2(plan,                 // PEXSI plan[in]
+                        &options,              // PEXSI Options[in]
                         numElectronExact,     // exact electron number[in]
                         &mu,                  // chemical potential[out]
                         &nelec,               // number of electrons[out]
-                        &muMinInertia,        // Lower bound for mu after the last inertia[out]
-                        &muMaxInertia,        // Upper bound for mu after the last inertia[out]
+                        // &muMinInertia,        // Lower bound for mu after the last inertia[out]
+                        // &muMaxInertia,        // Upper bound for mu after the last inertia[out]
                         &numTotalInertiaIter, // Number of total inertia[out]
-                        &numTotalPEXSIIter,   // number of total pexsi evaluation procedure[out]
+                        // &numTotalPEXSIIter,   // number of total pexsi evaluation procedure[out]
                         &info);               // 0: successful; otherwise: unsuccessful
         // LiuXh modify 2021-04-29, add DONE(ofs_running,"xx") for test
         ModuleBase::timer::tick("Diago_LCAO_Matrix", "PEXSIDFT");