diff --git a/source/source_base/matrix_block.h b/source/source_base/matrix_block.h new file mode 100644 index 00000000000..771aef49d97 --- /dev/null +++ b/source/source_base/matrix_block.h @@ -0,0 +1,33 @@ +#ifndef MODULE_BASE_MATRIX_BLOCK_H +#define MODULE_BASE_MATRIX_BLOCK_H + +#include + +namespace ModuleBase +{ + +/** + * @brief A non-owning description of a matrix stored in memory. + * + * It records only where the matrix lives and how it is laid out: the data + * pointer, the local number of rows and columns, and the BLACS array + * descriptor when the matrix is block-cyclically distributed. It carries no + * physical meaning, so it can be shared between the code that fills a matrix + * and the code that diagonalizes it without either side depending on the + * other. + * + * @note This is an aggregate on purpose; several call sites brace-initialize + * it as MatrixBlock{p, row, col, desc}. + */ +template struct MatrixBlock +{ + /* would change to Eigen in the future */ + T* p; + size_t row; + size_t col; + const int* desc; +}; + +} // namespace ModuleBase + +#endif // MODULE_BASE_MATRIX_BLOCK_H diff --git a/source/source_hamilt/matrixblock.h b/source/source_hamilt/matrixblock.h index 481a95b8a58..f26c585e76e 100644 --- a/source/source_hamilt/matrixblock.h +++ b/source/source_hamilt/matrixblock.h @@ -1,19 +1,19 @@ #ifndef MATRIXBLOCK_H #define MATRIXBLOCK_H -#include +#include "source_base/matrix_block.h" + namespace hamilt { -template struct MatrixBlock -{ - /* this is a simple template block of a matrix - would change to Eigen in the future */ - T* p; - size_t row; - size_t col; - const int* desc; -}; +/// MatrixBlock only describes a memory layout, so it now lives in source_base +/// and eigensolvers can use it without including the Hamiltonian interface. +/// This alias keeps the historical hamilt::MatrixBlock spelling working. +/// +/// TODO: this header is a temporary compatibility shim. Once every call site +/// spells the type as ModuleBase::MatrixBlock and includes +/// source_base/matrix_block.h directly, delete this file. +using ModuleBase::MatrixBlock; } // namespace hamilt -#endif \ No newline at end of file +#endif diff --git a/source/source_hsolver/diago_bpcg.h b/source/source_hsolver/diago_bpcg.h index 05a58b3656e..caad91e3e83 100644 --- a/source/source_hsolver/diago_bpcg.h +++ b/source/source_hsolver/diago_bpcg.h @@ -7,7 +7,6 @@ #include "source_base/module_device/memory_op.h" #include "source_base/module_device/types.h" #include "source_base/para_gemm.h" -#include "source_hamilt/hamilt.h" #include "source_hsolver/kernels/hegvd_op.h" #include "source_hsolver/para_lin_tf.h" diff --git a/source/source_hsolver/diago_cusolver.cpp b/source/source_hsolver/diago_cusolver.cpp index ed02f4cd831..b597e40ac72 100644 --- a/source/source_hsolver/diago_cusolver.cpp +++ b/source/source_hsolver/diago_cusolver.cpp @@ -32,8 +32,8 @@ DiagoCusolver::~DiagoCusolver() // Diagonalization function template void DiagoCusolver::diag( - hamilt::MatrixBlock& h_mat, - hamilt::MatrixBlock& s_mat, + ModuleBase::MatrixBlock& h_mat, + ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in) { diff --git a/source/source_hsolver/diago_cusolver.h b/source/source_hsolver/diago_cusolver.h index 50d9a393c66..d59daf87442 100644 --- a/source/source_hsolver/diago_cusolver.h +++ b/source/source_hsolver/diago_cusolver.h @@ -2,9 +2,9 @@ #define DIAGOCUSOLVER_H #include "source_base/macros.h" // GetRealType -#include "source_hamilt/hamilt.h" -#include "source_basis/module_ao/parallel_orbitals.h" +#include "source_base/matrix_block.h" #include "source_hsolver/kernels/cuda/diag_cusolver.cuh" +#include "source_psi/psi.h" namespace hsolver { @@ -26,8 +26,8 @@ class DiagoCusolver // Override the diag function for CUSOLVER diagonalization void diag( - hamilt::MatrixBlock& h_mat, - hamilt::MatrixBlock& s_mat, + ModuleBase::MatrixBlock& h_mat, + ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in); diff --git a/source/source_hsolver/diago_cusolvermp.cpp b/source/source_hsolver/diago_cusolvermp.cpp index 395add83d27..b152066c56f 100644 --- a/source/source_hsolver/diago_cusolvermp.cpp +++ b/source/source_hsolver/diago_cusolvermp.cpp @@ -2,6 +2,7 @@ #include "diago_cusolvermp.h" +#include "source_base/matrix_block.h" #include "source_base/module_external/blas_connector.h" #include "source_base/timer.h" #include "source_base/tool_title.h" @@ -14,7 +15,7 @@ template void DiagoCusolverMP::diag(hamilt::Hamilt* phm_in, psi::Psi& psi, Real* eigenvalue_in) { ModuleBase::TITLE("DiagoCusolverMP", "diag"); - hamilt::MatrixBlock h_mat, s_mat; + ModuleBase::MatrixBlock h_mat, s_mat; phm_in->matrix(h_mat, s_mat); std::vector eigen(this->nlocal, 0.0); diff --git a/source/source_hsolver/diago_elpa.cpp b/source/source_hsolver/diago_elpa.cpp index 31c2510fe2d..1673bc8091d 100644 --- a/source/source_hsolver/diago_elpa.cpp +++ b/source/source_hsolver/diago_elpa.cpp @@ -7,8 +7,8 @@ #include "source_base/tool_title.h" #include "source_base/tool_quit.h" -typedef hamilt::MatrixBlock matd; -typedef hamilt::MatrixBlock> matcd; +typedef ModuleBase::MatrixBlock matd; +typedef ModuleBase::MatrixBlock> matcd; namespace hsolver { #ifdef __MPI @@ -141,8 +141,8 @@ void DiagoElpa::diag(hamilt::Hamilt* phm_in, #ifdef __MPI template <> -void DiagoElpa>::diag_pool(hamilt::MatrixBlock>& h_mat, - hamilt::MatrixBlock>& s_mat, +void DiagoElpa>::diag_pool(ModuleBase::MatrixBlock>& h_mat, + ModuleBase::MatrixBlock>& s_mat, psi::Psi>& psi, Real* eigenvalue_in, MPI_Comm& comm) @@ -170,8 +170,8 @@ void DiagoElpa>::diag_pool(hamilt::MatrixBlock -void DiagoElpa::diag_pool(hamilt::MatrixBlock& h_mat, - hamilt::MatrixBlock& s_mat, +void DiagoElpa::diag_pool(ModuleBase::MatrixBlock& h_mat, + ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm) diff --git a/source/source_hsolver/diago_elpa.h b/source/source_hsolver/diago_elpa.h index 2bee7a14fa8..09e0dd55e4a 100644 --- a/source/source_hsolver/diago_elpa.h +++ b/source/source_hsolver/diago_elpa.h @@ -2,6 +2,7 @@ #define DIAGOELPA_H #include "source_base/macros.h" // GetRealType +#include "source_base/matrix_block.h" #include "source_hamilt/hamilt.h" #include "source_basis/module_ao/parallel_orbitals.h" @@ -22,7 +23,7 @@ class DiagoElpa void diag(hamilt::Hamilt* phm_in, psi::Psi& psi, Real* eigenvalue_in); #ifdef __MPI // diagnolization used in parallel-k case - void diag_pool(hamilt::MatrixBlock& h_mat, hamilt::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); + void diag_pool(ModuleBase::MatrixBlock& h_mat, ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); MPI_Comm setmpicomm(); // set mpi comm; static int elpa_num_thread; // need to set mpi_comm or not,-1 not,else the number of mpi needed #endif diff --git a/source/source_hsolver/diago_elpa_native.cpp b/source/source_hsolver/diago_elpa_native.cpp index 8a918377e89..a17c192137f 100644 --- a/source/source_hsolver/diago_elpa_native.cpp +++ b/source/source_hsolver/diago_elpa_native.cpp @@ -48,8 +48,8 @@ MPI_Comm DiagoElpaNative::setmpicomm() #ifdef __MPI template -void DiagoElpaNative::diag_pool(hamilt::MatrixBlock& h_mat, - hamilt::MatrixBlock& s_mat, +void DiagoElpaNative::diag_pool(ModuleBase::MatrixBlock& h_mat, + ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm) @@ -147,7 +147,7 @@ void DiagoElpaNative::diag(hamilt::Hamilt* phm_in, psi::Psi& psi, Real* { ModuleBase::TITLE("DiagoElpaNative", "diag"); #ifdef __MPI - hamilt::MatrixBlock h_mat, s_mat; + ModuleBase::MatrixBlock h_mat, s_mat; phm_in->matrix(h_mat, s_mat); MPI_Comm COMM_DIAG = setmpicomm(); // set mpi_comm needed diag_pool(h_mat, s_mat, psi, eigenvalue_in, COMM_DIAG); diff --git a/source/source_hsolver/diago_elpa_native.h b/source/source_hsolver/diago_elpa_native.h index 2cbd2c27d60..5511c8d8956 100644 --- a/source/source_hsolver/diago_elpa_native.h +++ b/source/source_hsolver/diago_elpa_native.h @@ -2,6 +2,7 @@ #define DIAGOELPANATIVE_H #include "source_base/macros.h" // GetRealType +#include "source_base/matrix_block.h" #include "source_hamilt/hamilt.h" #include "source_basis/module_ao/parallel_orbitals.h" @@ -24,7 +25,7 @@ class DiagoElpaNative void diag(hamilt::Hamilt* phm_in, psi::Psi& psi, Real* eigenvalue_in); #ifdef __MPI // diagnolization used in parallel-k case - void diag_pool(hamilt::MatrixBlock& h_mat, hamilt::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); + void diag_pool(ModuleBase::MatrixBlock& h_mat, ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); MPI_Comm setmpicomm(); // set mpi comm; static int elpa_num_thread; // need to set mpi_comm or not,-1 not,else the number of mpi needed static int lastmpinum; // last using mpi; diff --git a/source/source_hsolver/diago_lapack.cpp b/source/source_hsolver/diago_lapack.cpp index c98bac0e98e..77d42679f0a 100644 --- a/source/source_hsolver/diago_lapack.cpp +++ b/source/source_hsolver/diago_lapack.cpp @@ -8,16 +8,16 @@ #include -typedef hamilt::MatrixBlock matd; -typedef hamilt::MatrixBlock> matcd; +typedef ModuleBase::MatrixBlock matd; +typedef ModuleBase::MatrixBlock> matcd; namespace hsolver { namespace { template -void check_lapack_layout(const hamilt::MatrixBlock& h_mat, - const hamilt::MatrixBlock& s_mat, +void check_lapack_layout(const ModuleBase::MatrixBlock& h_mat, + const ModuleBase::MatrixBlock& s_mat, const std::size_t n) { if (h_mat.row != n || h_mat.col != n || s_mat.row != n || s_mat.col != n) @@ -68,8 +68,8 @@ void DiagoLapack>::diag(hamilt::Hamilt #ifdef __MPI template<> - void DiagoLapack::diag_pool(hamilt::MatrixBlock& h_mat, - hamilt::MatrixBlock& s_mat, + void DiagoLapack::diag_pool(ModuleBase::MatrixBlock& h_mat, + ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm) @@ -83,8 +83,8 @@ void DiagoLapack>::diag(hamilt::Hamilt BlasConnector::copy(this->nbands, eigen.data(), inc, eigenvalue_in, inc); } template<> - void DiagoLapack>::diag_pool(hamilt::MatrixBlock>& h_mat, - hamilt::MatrixBlock>& s_mat, + void DiagoLapack>::diag_pool(ModuleBase::MatrixBlock>& h_mat, + ModuleBase::MatrixBlock>& s_mat, psi::Psi>& psi, Real* eigenvalue_in, MPI_Comm& comm) diff --git a/source/source_hsolver/diago_lapack.h b/source/source_hsolver/diago_lapack.h index 9f674fc20fd..213d772dbcd 100644 --- a/source/source_hsolver/diago_lapack.h +++ b/source/source_hsolver/diago_lapack.h @@ -9,6 +9,7 @@ #define DIAGOLAPACK_H #include "source_base/macros.h" // GetRealType +#include "source_base/matrix_block.h" #include "source_hamilt/hamilt.h" #include "source_base/matrix.h" #include "source_basis/module_ao/parallel_orbitals.h" @@ -33,7 +34,7 @@ class DiagoLapack void diag(hamilt::Hamilt* phm_in, psi::Psi& psi, Real* eigenvalue_in); #ifdef __MPI // diagnolization used in parallel-k case - void diag_pool(hamilt::MatrixBlock& h_mat, hamilt::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); + void diag_pool(ModuleBase::MatrixBlock& h_mat, ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); #endif void dsygvx_diag(const int ncol, diff --git a/source/source_hsolver/diago_pexsi.cpp b/source/source_hsolver/diago_pexsi.cpp index 9915416a7ef..b4fac092bfc 100644 --- a/source/source_hsolver/diago_pexsi.cpp +++ b/source/source_hsolver/diago_pexsi.cpp @@ -3,13 +3,14 @@ #include #ifdef __PEXSI #include "diago_pexsi.h" +#include "source_base/matrix_block.h" #include "source_base/tool_title.h" #include "source_base/tool_quit.h" #include "source_basis/module_ao/parallel_orbitals.h" #include "module_pexsi/pexsi_solver.h" -typedef hamilt::MatrixBlock matd; -typedef hamilt::MatrixBlock> matcd; +typedef ModuleBase::MatrixBlock matd; +typedef ModuleBase::MatrixBlock> matcd; namespace hsolver { diff --git a/source/source_hsolver/diago_scalapack.cpp b/source/source_hsolver/diago_scalapack.cpp index caecd771ccd..32645855af3 100644 --- a/source/source_hsolver/diago_scalapack.cpp +++ b/source/source_hsolver/diago_scalapack.cpp @@ -8,15 +8,15 @@ #include "diago_scalapack.h" #include "source_base/global_function.h" +#include "source_base/matrix_block.h" #include "source_base/module_external/blacs_connector.h" #include "source_base/module_external/scalapack_connector.h" -#include "source_hamilt/matrixblock.h" #include #include -typedef hamilt::MatrixBlock matd; -typedef hamilt::MatrixBlock> matcd; +typedef ModuleBase::MatrixBlock matd; +typedef ModuleBase::MatrixBlock> matcd; namespace hsolver { @@ -59,8 +59,8 @@ int blacs_grid_size(const int* const desc) #ifdef __MPI template<> - void DiagoScalapack::diag_pool(hamilt::MatrixBlock& h_mat, - hamilt::MatrixBlock& s_mat, + void DiagoScalapack::diag_pool(ModuleBase::MatrixBlock& h_mat, + ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm) @@ -73,8 +73,8 @@ int blacs_grid_size(const int* const desc) BlasConnector::copy(this->nbands, eigen.data(), inc, eigenvalue_in, inc); } template<> - void DiagoScalapack>::diag_pool(hamilt::MatrixBlock>& h_mat, - hamilt::MatrixBlock>& s_mat, + void DiagoScalapack>::diag_pool(ModuleBase::MatrixBlock>& h_mat, + ModuleBase::MatrixBlock>& s_mat, psi::Psi>& psi, Real* eigenvalue_in, MPI_Comm& comm) diff --git a/source/source_hsolver/diago_scalapack.h b/source/source_hsolver/diago_scalapack.h index 7dce0306169..7cffdf8652f 100644 --- a/source/source_hsolver/diago_scalapack.h +++ b/source/source_hsolver/diago_scalapack.h @@ -13,6 +13,7 @@ #include #include "source_base/macros.h" // GetRealType +#include "source_base/matrix_block.h" #include "source_hamilt/hamilt.h" #include "source_psi/psi.h" #include "source_base/complexmatrix.h" @@ -34,7 +35,7 @@ namespace hsolver void diag(hamilt::Hamilt* phm_in, psi::Psi& psi, Real* eigenvalue_in); #ifdef __MPI // diagnolization used in parallel-k case - void diag_pool(hamilt::MatrixBlock& h_mat, hamilt::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); + void diag_pool(ModuleBase::MatrixBlock& h_mat, ModuleBase::MatrixBlock& s_mat, psi::Psi& psi, Real* eigenvalue_in, MPI_Comm& comm); #endif private: diff --git a/source/source_hsolver/hsolver_lcao.cpp b/source/source_hsolver/hsolver_lcao.cpp index 49612da3462..967ed9c153b 100644 --- a/source/source_hsolver/hsolver_lcao.cpp +++ b/source/source_hsolver/hsolver_lcao.cpp @@ -1,5 +1,7 @@ #include "hsolver_lcao.h" +#include "source_base/matrix_block.h" + #ifdef __MPI #include "diago_scalapack.h" #include "source_base/module_external/scalapack_connector.h" @@ -162,7 +164,7 @@ void HSolverLCAO::hamiltSolvePsiK(hamilt::Hamilt* hm, psi::Psi& psi, do { // Note: This branch will only be executed in the single-process case DiagoCusolver cu(this->nlocal, this->nbands); - hamilt::MatrixBlock hk, sk; + ModuleBase::MatrixBlock hk, sk; hm->matrix(hk, sk); cu.diag(hk, sk, psi, eigenvalue); } @@ -242,14 +244,14 @@ void HSolverLCAO::parakSolve(hamilt::Hamilt* pHamilt, { /// local psi in pool psi_pool.fix_k(0); - hamilt::MatrixBlock hk_pool = hamilt::MatrixBlock{k2d.hk_pool.data(), - (size_t)k2d.get_p2D_pool()->get_row_size(), - (size_t)k2d.get_p2D_pool()->get_col_size(), - k2d.get_p2D_pool()->desc}; - hamilt::MatrixBlock sk_pool = hamilt::MatrixBlock{k2d.sk_pool.data(), - (size_t)k2d.get_p2D_pool()->get_row_size(), - (size_t)k2d.get_p2D_pool()->get_col_size(), - k2d.get_p2D_pool()->desc}; + ModuleBase::MatrixBlock hk_pool = ModuleBase::MatrixBlock{k2d.hk_pool.data(), + (size_t)k2d.get_p2D_pool()->get_row_size(), + (size_t)k2d.get_p2D_pool()->get_col_size(), + k2d.get_p2D_pool()->desc}; + ModuleBase::MatrixBlock sk_pool = ModuleBase::MatrixBlock{k2d.sk_pool.data(), + (size_t)k2d.get_p2D_pool()->get_row_size(), + (size_t)k2d.get_p2D_pool()->get_col_size(), + k2d.get_p2D_pool()->desc}; /// solve eigenvector and eigenvalue for H(k) if (this->method == "scalapack_gvx") { @@ -402,7 +404,7 @@ void HSolverLCAO::parakSolve_cusolver(hamilt::Hamilt* pHamilt, sk_mat.resize(nrow * ncol); } pHamilt->updateHk(ik); - hamilt::MatrixBlock hk_2D, sk_2D; + ModuleBase::MatrixBlock hk_2D, sk_2D; pHamilt->matrix(hk_2D, sk_2D); int desc_tmp[9]; T* hk_local_ptr = hk_mat.data(); @@ -427,10 +429,10 @@ void HSolverLCAO::parakSolve_cusolver(hamilt::Hamilt* pHamilt, { psi_local.resize(1, ncol, nrow); DiagoCusolver cu(this->nlocal, this->nbands); - hamilt::MatrixBlock hk_local = hamilt::MatrixBlock{ + ModuleBase::MatrixBlock hk_local = ModuleBase::MatrixBlock{ hk_mat.data(), (size_t)nrow, (size_t)ncol, mat_para_local.desc}; - hamilt::MatrixBlock sk_local = hamilt::MatrixBlock{ + ModuleBase::MatrixBlock sk_local = ModuleBase::MatrixBlock{ sk_mat.data(), (size_t)nrow, (size_t)ncol, mat_para_local.desc}; cu.diag(hk_local, sk_local, psi_local, &(pes->ekb(kpt_assigned, 0))); diff --git a/source/source_hsolver/parallel_k2d.cpp b/source/source_hsolver/parallel_k2d.cpp index cb2a3947560..b7521fde703 100644 --- a/source/source_hsolver/parallel_k2d.cpp +++ b/source/source_hsolver/parallel_k2d.cpp @@ -1,5 +1,6 @@ #include "parallel_k2d.h" +#include "source_base/matrix_block.h" #include "source_base/parallel_global.h" #include "source_base/module_external/scalapack_connector.h" #include "source_base/timer.h" @@ -43,7 +44,7 @@ void Parallel_K2D::distribute_hsk(hamilt::Hamilt* pHamilt, for (int ipool = 0; ipool < ik_kpar.size(); ++ipool) { pHamilt->updateHk(ik_kpar[ipool]); - hamilt::MatrixBlock HK_global, SK_global; + ModuleBase::MatrixBlock HK_global, SK_global; pHamilt->matrix(HK_global, SK_global); if (this->MY_POOL == this->Pkpoints->whichpool[ik_kpar[ipool]]) { this->hk_pool.resize(this->P2D_pool->get_local_size(), 0.0); diff --git a/source/source_hsolver/parallel_k2d.h b/source/source_hsolver/parallel_k2d.h index aebd683f4bc..57e37ef278f 100644 --- a/source/source_hsolver/parallel_k2d.h +++ b/source/source_hsolver/parallel_k2d.h @@ -3,7 +3,6 @@ #include "source_base/parallel_2d.h" #include "source_cell/parallel_kpoints.h" -#include "source_hamilt/matrixblock.h" #ifdef __MPI #include "mpi.h" #endif