diff --git a/source/Makefile.Objects b/source/Makefile.Objects index 47044ce42d..c515dc426b 100644 --- a/source/Makefile.Objects +++ b/source/Makefile.Objects @@ -478,6 +478,7 @@ OBJS_IO_LCAO=cal_r_overlap_R.o\ read_dm.o\ unk_overlap_lcao.o\ read_wfc_nao.o\ + read_wfc_lcao.o\ write_wfc_nao.o\ write_HS_sparse.o\ single_R_io.o\ diff --git a/source/module_hamilt_lcao/hamilt_lcaodft/local_orbital_wfc.cpp b/source/module_hamilt_lcao/hamilt_lcaodft/local_orbital_wfc.cpp index 201c006e80..bb4171dfc5 100644 --- a/source/module_hamilt_lcao/hamilt_lcaodft/local_orbital_wfc.cpp +++ b/source/module_hamilt_lcao/hamilt_lcaodft/local_orbital_wfc.cpp @@ -1,33 +1,34 @@ #include "local_orbital_wfc.h" -#include "module_hamilt_pw/hamilt_pwdft/global.h" -#include "module_io/read_wfc_nao.h" + #include "module_base/memory.h" #include "module_base/timer.h" +#include "module_hamilt_pw/hamilt_pwdft/global.h" +#include "module_io/read_wfc_nao.h" Local_Orbital_wfc::Local_Orbital_wfc() { - allocate_flag = false; - wfck_flag = false; - complex_flag = false; + allocate_flag = false; + wfck_flag = false; + complex_flag = false; nks = 0; } Local_Orbital_wfc::~Local_Orbital_wfc() { - // used for k-points. - if(this->complex_flag) - { + // used for k-points. + if (this->complex_flag) + { delete[] this->wfc_k_grid2; } - if(this->wfck_flag) + if (this->wfck_flag) { - for(int i=0; iwfc_k_grid[i]; - } - delete[] this->wfc_k_grid; - } + for (int i = 0; i < nks; i++) + { + delete[] this->wfc_k_grid[i]; + } + delete[] this->wfc_k_grid; + } } void Local_Orbital_wfc::gamma_file(psi::Psi* psid, elecstate::ElecState* pelec) @@ -37,17 +38,15 @@ void Local_Orbital_wfc::gamma_file(psi::Psi* psid, elecstate::ElecState* double** ctot; - //allocate psi + // allocate psi int ncol = this->ParaV->ncol_bands; - if (GlobalV::KS_SOLVER == "genelpa" - || GlobalV::KS_SOLVER == "lapack_gvx" - || GlobalV::KS_SOLVER == "scalapack_gvx" - || GlobalV::KS_SOLVER == "cg_in_lcao" + if (GlobalV::KS_SOLVER == "genelpa" || GlobalV::KS_SOLVER == "lapack_gvx" || GlobalV::KS_SOLVER == "scalapack_gvx" + || GlobalV::KS_SOLVER == "cg_in_lcao" #ifdef __CUDA || GlobalV::KS_SOLVER == "cusolver" #endif - ) + ) { ncol = this->ParaV->ncol; } @@ -64,7 +63,7 @@ void Local_Orbital_wfc::gamma_file(psi::Psi* psid, elecstate::ElecState* for (int is = 0; is < GlobalV::NSPIN; ++is) { - this->error = ModuleIO::read_wfc_nao(ctot, is, GlobalV::GAMMA_ONLY_LOCAL, GlobalV::NB2D, GlobalV::NBANDS, + this->error = ModuleIO::read_wfc_nao(ctot, is, GlobalV::GAMMA_ONLY_LOCAL, GlobalV::NB2D, GlobalV::NBANDS, GlobalV::NLOCAL, GlobalV::global_readin_dir, this->ParaV, psid, pelec); #ifdef __MPI Parallel_Common::bcast_int(this->error); @@ -91,69 +90,73 @@ void Local_Orbital_wfc::gamma_file(psi::Psi* psid, elecstate::ElecState* std::cout << "WARNING: Failed to read in wavefunction, use default initialization instead." << std::endl; break; } - }//loop ispin + } // loop ispin } -void Local_Orbital_wfc::allocate_k(const int& lgd, - psi::Psi>* psi, - elecstate::ElecState* pelec, - const int& nks, - const int& nkstot, - const std::vector>& kvec_c, - const int& istep) +void Local_Orbital_wfc::allocate_k(const int& lgd, psi::Psi>* psi, elecstate::ElecState* pelec, + const int& nks, const int& nkstot, + const std::vector>& kvec_c, const int& istep) { this->nks = nks; ModuleBase::TITLE("Local_Orbital_wfc", "allocate_k"); - if(GlobalV::NLOCAL < GlobalV::NBANDS) - { - ModuleBase::WARNING_QUIT("Local_Orbital_wfc::allocate","NLOCALwfck_flag == false) - { - this->wfc_k_grid = new std::complex**[nks]; - for(int ik=0; ikwfc_k_grid[ik] = new std::complex*[GlobalV::NBANDS]; - } - this->wfck_flag = true; - } - - if(this->complex_flag) - { - delete[] this->wfc_k_grid2; - this->complex_flag = false; - } - // allocate the second part. - //if(lgd != 0) xiaohui modify 2015-02-04, fixed memory bug - if(lgd != 0) - { - const int page=GlobalV::NBANDS*lgd; - this->wfc_k_grid2=new std::complex [nks*page]; - ModuleBase::GlobalFunc::ZEROS(wfc_k_grid2, nks*page); - for(int ik=0; ikwfc_k_grid[ik][ib] = &wfc_k_grid2[ik*page+ib*lgd]; - } - ModuleBase::Memory::record("LOWF::wfc_k_grid", sizeof(std::complex) * GlobalV::NBANDS*GlobalV::NLOCAL); - this->complex_flag = true; - } - } + // mohan add the flag 2011-03-02 + // allocate the first part (only once!). + if (this->wfck_flag == false) + { + this->wfc_k_grid = new std::complex**[nks]; + for (int ik = 0; ik < nks; ik++) + { + this->wfc_k_grid[ik] = new std::complex*[GlobalV::NBANDS]; + } + this->wfck_flag = true; + } - if (INPUT.init_wfc == "atomic") + if (this->complex_flag) + { + delete[] this->wfc_k_grid2; + this->complex_flag = false; + } + // allocate the second part and initialize value as zero. + // if(lgd != 0) xiaohui modify 2015-02-04, fixed memory bug + if (lgd != 0) { + const int page = GlobalV::NBANDS * lgd; // lgd: local grid dimension + this->wfc_k_grid2 + = new std::complex[nks + * page]; // wfc_k_grid2 stores nks * nbands * lgd number of basis coefficients + ModuleBase::GlobalFunc::ZEROS(wfc_k_grid2, nks * page); + for (int ik = 0; ik < nks; ik++) + { + for (int ib = 0; ib < GlobalV::NBANDS; ib++) + { + this->wfc_k_grid[ik][ib] + = &wfc_k_grid2[ik * page + ib * lgd + + 0]; // then wfc_k_grid stores the starting address of each band + // but now there are less number of coefficients stored, if lgd < nbasis. + // this is because in grid intergration, for each grid point, only few grid points (near some + // neighboring atoms) would be involved in calculation, therefore only few basis coefficients are + // needed. + } + ModuleBase::Memory::record("LOWF::wfc_k_grid", + sizeof(std::complex) * GlobalV::NBANDS * GlobalV::NLOCAL); + this->complex_flag = true; + } } - else if (INPUT.init_wfc == "file") + + // read wavefunction from file, then divide, distribute and broadcast + if (INPUT.init_wfc == "file") // init_wfc can also be "atomic" but actually do nothing. { - if (istep > 0) - { - return; - } + // confusing, seems istep to be the index of scf step. if not the first scf, why call this function? + if (istep > 0) + { + return; + } std::cout << " Read in wave functions files: " << nkstot << std::endl; if (psi == nullptr) { @@ -166,8 +169,9 @@ void Local_Orbital_wfc::allocate_k(const int& lgd, for (int ik = 0; ik < nkstot; ++ik) { std::complex** ctot; - this->error = ModuleIO::read_wfc_nao_complex(ctot, ik, GlobalV::NB2D, GlobalV::NBANDS, GlobalV::NLOCAL, - GlobalV::global_readin_dir, kvec_c[ik], this->ParaV, psi, pelec); + this->error + = ModuleIO::read_wfc_nao_complex(ctot, ik, GlobalV::NB2D, GlobalV::NBANDS, GlobalV::NLOCAL, + GlobalV::global_readin_dir, kvec_c[ik], this->ParaV, psi, pelec); #ifdef __MPI Parallel_Common::bcast_int(this->error); #endif @@ -190,128 +194,113 @@ void Local_Orbital_wfc::allocate_k(const int& lgd, } if (this->error) { - std::cout << "WARNING: Failed to read in wavefunction, use default initialization instead." << std::endl; + std::cout << "WARNING: Failed to read in wavefunction, use default initialization instead." + << std::endl; break; } } } - else - { - ModuleBase::WARNING_QUIT("Local_Orbital_wfc","check the parameter: init_wfc"); - } - return; + return; } - int Local_Orbital_wfc::globalIndex(int localindex, int nblk, int nprocs, int myproc) { int iblock, gIndex; - iblock=localindex/nblk; - gIndex=(iblock*nprocs+myproc)*nblk+localindex%nblk; + iblock = localindex / nblk; + gIndex = (iblock * nprocs + myproc) * nblk + localindex % nblk; return gIndex; } - int Local_Orbital_wfc::localIndex(int globalindex, int nblk, int nprocs, int& myproc) { - myproc=int((globalindex%(nblk*nprocs))/nblk); - return int(globalindex/(nblk*nprocs))*nblk+globalindex%nblk; + myproc = int((globalindex % (nblk * nprocs)) / nblk); + return int(globalindex / (nblk * nprocs)) * nblk + globalindex % nblk; } #ifdef __MPI -void Local_Orbital_wfc::wfc_2d_to_grid(const double* wfc_2d, - double** wfc_grid, - const int ik, - const ModuleBase::matrix& ekb, - const ModuleBase::matrix& wg) +void Local_Orbital_wfc::wfc_2d_to_grid(const double* wfc_2d, double** wfc_grid, const int ik, + const ModuleBase::matrix& ekb, const ModuleBase::matrix& wg) { ModuleBase::TITLE(" Local_Orbital_wfc", "wfc_2d_to_grid"); - ModuleBase::timer::tick("Local_Orbital_wfc","wfc_2d_to_grid"); + ModuleBase::timer::tick("Local_Orbital_wfc", "wfc_2d_to_grid"); const Parallel_Orbitals* pv = this->ParaV; const int inc = 1; - int myid=0; + int myid = 0; MPI_Comm_rank(pv->comm_2D, &myid); - int info=0; - - //calculate maxnloc for bcasting 2d-wfc + int info = 0; + + // calculate maxnloc for bcasting 2d-wfc long maxnloc; // maximum number of elements in local matrix - info=MPI_Reduce(&pv->nloc_wfc, &maxnloc, 1, MPI_LONG, MPI_MAX, 0, pv->comm_2D); - info=MPI_Bcast(&maxnloc, 1, MPI_LONG, 0, pv->comm_2D); + info = MPI_Reduce(&pv->nloc_wfc, &maxnloc, 1, MPI_LONG, MPI_MAX, 0, pv->comm_2D); + info = MPI_Bcast(&maxnloc, 1, MPI_LONG, 0, pv->comm_2D); std::vector work(maxnloc); // work/buffer matrix int naroc[2]; // maximum number of row or column - for(int iprow=0; iprowdim0; ++iprow) + for (int iprow = 0; iprow < pv->dim0; ++iprow) { - for(int ipcol=0; ipcoldim1; ++ipcol) + for (int ipcol = 0; ipcol < pv->dim1; ++ipcol) { - const int coord[2]={iprow, ipcol}; + const int coord[2] = {iprow, ipcol}; int src_rank; - info=MPI_Cart_rank(pv->comm_2D, coord, &src_rank); - if(myid==src_rank) + info = MPI_Cart_rank(pv->comm_2D, coord, &src_rank); + if (myid == src_rank) { BlasConnector::copy(pv->nloc_wfc, wfc_2d, inc, work.data(), inc); - naroc[0]=pv->nrow; - naroc[1]=pv->ncol_bands; + naroc[0] = pv->nrow; + naroc[1] = pv->ncol_bands; } - info=MPI_Bcast(naroc, 2, MPI_INT, src_rank, pv->comm_2D); - info=MPI_Bcast(work.data(), maxnloc, MPI_DOUBLE, src_rank, pv->comm_2D); + info = MPI_Bcast(naroc, 2, MPI_INT, src_rank, pv->comm_2D); + info = MPI_Bcast(work.data(), maxnloc, MPI_DOUBLE, src_rank, pv->comm_2D); - info = this->set_wfc_grid(naroc, pv->nb, - pv->dim0, pv->dim1, iprow, ipcol, - work.data(), wfc_grid); + info = this->set_wfc_grid(naroc, pv->nb, pv->dim0, pv->dim1, iprow, ipcol, work.data(), wfc_grid); - }//loop ipcol - }//loop iprow - ModuleBase::timer::tick("Local_Orbital_wfc","wfc_2d_to_grid"); + } // loop ipcol + } // loop iprow + ModuleBase::timer::tick("Local_Orbital_wfc", "wfc_2d_to_grid"); } -void Local_Orbital_wfc::wfc_2d_to_grid(const std::complex* wfc_2d, - std::complex** wfc_grid, - int ik, - const ModuleBase::matrix& ekb, - const ModuleBase::matrix& wg, +void Local_Orbital_wfc::wfc_2d_to_grid(const std::complex* wfc_2d, std::complex** wfc_grid, int ik, + const ModuleBase::matrix& ekb, const ModuleBase::matrix& wg, const std::vector>& kvec_c) { ModuleBase::TITLE("Local_Orbital_wfc", "wfc_2d_to_grid"); - ModuleBase::timer::tick("Local_Orbital_wfc","wfc_2d_to_grid"); + ModuleBase::timer::tick("Local_Orbital_wfc", "wfc_2d_to_grid"); const Parallel_Orbitals* pv = this->ParaV; const int inc = 1; - int myid=0; + int myid = 0; MPI_Comm_rank(pv->comm_2D, &myid); - int info=0; - - //calculate maxnloc for bcasting 2d-wfc - long maxnloc=0; // maximum number of elements in local matrix - info=MPI_Reduce(&pv->nloc_wfc, &maxnloc, 1, MPI_LONG, MPI_MAX, 0, pv->comm_2D); - info=MPI_Bcast(&maxnloc, 1, MPI_LONG, 0, pv->comm_2D); + int info = 0; + + // calculate maxnloc for bcasting 2d-wfc + long maxnloc = 0; // maximum number of elements in local matrix + info = MPI_Reduce(&pv->nloc_wfc, &maxnloc, 1, MPI_LONG, MPI_MAX, 0, pv->comm_2D); + info = MPI_Bcast(&maxnloc, 1, MPI_LONG, 0, pv->comm_2D); std::vector> work(maxnloc); // work/buffer matrix - + int naroc[2] = {0}; // maximum number of row or column - for(int iprow=0; iprowdim0; ++iprow) + for (int iprow = 0; iprow < pv->dim0; ++iprow) { - for(int ipcol=0; ipcoldim1; ++ipcol) + for (int ipcol = 0; ipcol < pv->dim1; ++ipcol) { - const int coord[2]={iprow, ipcol}; + const int coord[2] = {iprow, ipcol}; int src_rank; - info=MPI_Cart_rank(pv->comm_2D, coord, &src_rank); - if(myid==src_rank) + info = MPI_Cart_rank(pv->comm_2D, coord, &src_rank); + if (myid == src_rank) { BlasConnector::copy(pv->nloc_wfc, wfc_2d, inc, work.data(), inc); - naroc[0]=pv->nrow; - naroc[1]=pv->ncol_bands; + naroc[0] = pv->nrow; + naroc[1] = pv->ncol_bands; } - info=MPI_Bcast(naroc, 2, MPI_INT, src_rank, pv->comm_2D); + info = MPI_Bcast(naroc, 2, MPI_INT, src_rank, pv->comm_2D); info = MPI_Bcast(work.data(), maxnloc, MPI_DOUBLE_COMPLEX, src_rank, pv->comm_2D); - // mohan update 2021-02-12, delte BFIELD option - info = this->set_wfc_grid(naroc, pv->nb, - pv->dim0, pv->dim1, iprow, ipcol, - work.data(), wfc_grid); - }//loop ipcol - }//loop iprow + // mohan update 2021-02-12, delte BFIELD option + info = this->set_wfc_grid(naroc, pv->nb, pv->dim0, pv->dim1, iprow, ipcol, work.data(), wfc_grid); + } // loop ipcol + } // loop iprow - ModuleBase::timer::tick("Local_Orbital_wfc","wfc_2d_to_grid"); + ModuleBase::timer::tick("Local_Orbital_wfc", "wfc_2d_to_grid"); } #endif diff --git a/source/module_io/CMakeLists.txt b/source/module_io/CMakeLists.txt index 0b40b464a1..0b621d4e91 100644 --- a/source/module_io/CMakeLists.txt +++ b/source/module_io/CMakeLists.txt @@ -54,6 +54,7 @@ if(ENABLE_LCAO) istate_envelope.cpp read_dm.cpp read_wfc_nao.cpp + read_wfc_lcao.cpp write_wfc_nao.cpp write_dm.cpp dos_nao.cpp diff --git a/source/module_io/read_wfc_lcao.cpp b/source/module_io/read_wfc_lcao.cpp new file mode 100644 index 0000000000..748c238fe1 --- /dev/null +++ b/source/module_io/read_wfc_lcao.cpp @@ -0,0 +1,355 @@ +#include "module_io/read_wfc_lcao.h" + +#include "module_base/formatter.h" +#include "module_base/tool_quit.h" + +#include +#include +#include + +#ifdef __MPI +#include "module_base/parallel_common.h" +#endif + +template +void ModuleIO::read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, int& nbands, + int& nbasis, std::vector>& lowf, std::vector& ekb, + std::vector& occ, + double& wk) //<[out] wavefunction coefficients +{ + // assert the T must be double or float + std::ifstream ifs(flowf.c_str()); + if (!ifs) + ModuleBase::WARNING_QUIT("read_abacus_lowf", "open file failed: " + flowf); + // will use line-by-line parse + std::string line; + bool read_kvec = false; + int iband = 0; + int ilocal = 0; + + wk = 1.0; + while (std::getline(ifs, line)) + { + // remove leading and trailing whitespaces + line = std::regex_replace(line, std::regex("^ +| +$|( ) +"), "$1"); + if (FmtCore::endswith(line, "(index of k points)")) + { + std::vector result = FmtCore::split(line); + ik = std::stoi(result[0]); + read_kvec = true; + continue; + } + if (read_kvec) + { + const std::vector result = FmtCore::split(line); + kvec_c.x = std::stod(result[0]); + kvec_c.y = std::stod(result[1]); + kvec_c.z = std::stod(result[2]); + read_kvec = false; + continue; + } + if (FmtCore::endswith(line, "(number of bands)")) + { + std::vector result = FmtCore::split(line); + nbands = std::stoi(result[0]); + ekb.resize(nbands, 0.0); // initialize ekb + occ.resize(nbands, 0.0); // initialize occ + } + else if (FmtCore::endswith(line, "(number of orbitals)")) + { + std::vector result = FmtCore::split(line); + nbasis = std::stoi(result[0]); + lowf.resize(nbands * nbasis, std::complex(0.0, 0.0)); // initialize lowf + } + else if (FmtCore::endswith(line, "(band)")) + { + std::vector result = FmtCore::split(line); +#ifdef __DEBUG + assert(ilocal == 0) || (ilocal == nlocal); +#endif + iband = std::stoi(result[0]) - 1; + ilocal = 0; // reset ilocal + } + else if (FmtCore::endswith(line, "(Ry)")) + { + std::vector result = FmtCore::split(line); + ekb[iband] = std::stod(result[0]); + } + else if (FmtCore::endswith(line, "(Occupations)")) + { + std::vector result = FmtCore::split(line); + occ[iband] = std::stod(result[0]); +#ifdef __DEBUG + assert(ilocal == 0); +#endif + } + else // read wavefunction coefficients + { + const std::vector result = FmtCore::split(line); + // for the case the complex number is written as a b + for (int i = 0; i < result.size(); i += 2) + { + lowf[iband * nbasis + ilocal] = std::complex(std::stod(result[i]), std::stod(result[i + 1])); + ilocal += 1; + } + } + } +#ifdef __DEBUG + assert(lowf.size() == nbands * nlocal); + assert(iband == nbands); + assert(ilocal == nlocal); +#endif +} +// instantiate the template function +template void ModuleIO::read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, + int& nbands, int& nbasis, std::vector>& lowf, + std::vector& ekb, std::vector& occ, double& wk); +template void ModuleIO::read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, + int& nbands, int& nbasis, std::vector>& lowf, + std::vector& ekb, std::vector& occ, double& wk); + +template +void ModuleIO::read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, int& nbands, + int& nbasis, std::vector& lowf, std::vector& ekb, std::vector& occ, + double& wk) +{ + std::ifstream ifs(flowf.c_str()); + if (!ifs) + ModuleBase::WARNING_QUIT("read_abacus_lowf", "open file failed: " + flowf); + // will use line-by-line parse + std::string line; + bool read_kvec = false; + int iband = 0; + int ilocal = 0; + + ik = 0; + kvec_c = ModuleBase::Vector3(0.0, 0.0, 0.0); + wk = 1.0; + while (std::getline(ifs, line)) + { + // remove leading and trailing whitespaces + line = std::regex_replace(line, std::regex("^ +| +$|( ) +"), "$1"); + if (read_kvec) + { + const std::vector result = FmtCore::split(line); + kvec_c.x = std::stod(result[0]); + kvec_c.y = std::stod(result[1]); + kvec_c.z = std::stod(result[2]); + read_kvec = false; + continue; + } + if (FmtCore::endswith(line, "(number of bands)")) + { + std::vector result = FmtCore::split(line); + nbands = std::stoi(result[0]); + ekb.resize(nbands, 0.0); // initialize ekb + occ.resize(nbands, 0.0); // initialize occ + } + else if (FmtCore::endswith(line, "(number of orbitals)")) + { + std::vector result = FmtCore::split(line); + nbasis = std::stoi(result[0]); + lowf.resize(nbands * nbasis, 0.0); // initialize lowf + } + else if (FmtCore::endswith(line, "(band)")) + { + std::vector result = FmtCore::split(line); +#ifdef __DEBUG + assert(ilocal == 0) || (ilocal == nlocal); +#endif + iband = std::stoi(result[0]) - 1; + ilocal = 0; // reset ilocal + } + else if (FmtCore::endswith(line, "(Ry)")) + { + std::vector result = FmtCore::split(line); + ekb[iband] = std::stod(result[0]); + } + else if (FmtCore::endswith(line, "(Occupations)")) + { + std::vector result = FmtCore::split(line); + occ[iband] = std::stod(result[0]); +#ifdef __DEBUG + assert(ilocal == 0); +#endif + } + else // read wavefunction coefficients + { + const std::vector result = FmtCore::split(line); + for (const auto& token: result) + { + lowf[iband * nbasis + ilocal] = static_cast(std::stod(token)); + ilocal += 1; + } + } + } +#ifdef __DEBUG + assert(lowf.size() == nbands * nlocal); + assert(iband == nbands); + assert(ilocal == nlocal); +#endif +} +// instantiate the template function +template void ModuleIO::read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, + int& nbands, int& nbasis, std::vector& lowf, std::vector& ekb, + std::vector& occ, double& wk); +template void ModuleIO::read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, + int& nbands, int& nbasis, std::vector& lowf, std::vector& ekb, + std::vector& occ, double& wk); + +#ifdef __MPI +template +void ModuleIO::restart_from_file(const std::string& out_dir, // hard-code the file name to be WFC_NAO_K*.txt? + const Parallel_2D& p2d, const int& nks, int& nbands, int& nbasis, + std::vector& lowf_loc, std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk) +{ + // reset vectors + lowf_loc.clear(); + ekb.clear(); + occ.clear(); + kvec_c.clear(); + wk.clear(); + // MPI-related variables init + int iproc; + MPI_Comm_rank(p2d.comm_2D, &iproc); + // then start + int nbands_ = -1, nbasis_ = -1; + for (int ik = 0; ik < nks; ik++) + { + // check existence of file + const std::string flowf = out_dir + "/WFC_NAO_K" + std::to_string(ik + 1) + ".txt"; + std::ifstream ifs(flowf); + if (!ifs) + ModuleBase::WARNING_QUIT("restart_from_file", "open file failed: " + flowf); + + std::vector lowf_glb; + std::vector lowf_loc_k; + std::vector ekb_; + std::vector occ_; + ModuleBase::Vector3 kvec; + double wk_; + if (iproc == 0) // only one rank is needed to read the global lowf, ekb, ... + { + int ik_; + read_abacus_lowf(flowf, ik_, kvec, nbands, nbasis, lowf_glb, ekb_, occ_, wk_); + assert(ik_ == ik + 1); // check the consistency of ik + assert(nbands == nbands_ || nbands_ == -1); // check the consistency of nbands + assert(nbasis == nbasis_ || nbasis_ == -1); // check the consistency of nbasis + nbands_ = (nbands_ == -1) ? nbands : nbands_; + nbasis_ = (nbasis_ == -1) ? nbasis : nbasis_; + ekb.insert(ekb.end(), ekb_.begin(), ekb_.end()); + occ.insert(occ.end(), occ_.begin(), occ_.end()); + wk.push_back(wk_); + kvec_c.push_back(kvec); + } + MPI_Barrier(p2d.comm_2D); // wait for finishing the reading task + // scatter the lowf_glb to lowf_loc + Parallel_2D p2d_glb; + Parallel_Common::bcast_int(nbands); + Parallel_Common::bcast_int(nbasis); + p2d_glb.init(nbasis, nbands, std::max(nbasis, nbands), p2d.comm_2D); // in the same comm world + lowf_loc_k.resize(p2d.nrow * p2d.ncol); + Cpxgemr2d(nbasis, nbands, lowf_glb.data(), 1, 1, const_cast(p2d_glb.desc), lowf_loc_k.data(), 1, 1, + const_cast(p2d.desc), p2d_glb.blacs_ctxt); + // append to the global lowf_loc + lowf_loc.insert(lowf_loc.end(), lowf_loc_k.begin(), lowf_loc_k.end()); + } + assert(lowf_loc.size() == nks * p2d.nrow * p2d.ncol); + // still something to broadcast: ekb, occ, wk and kvec. + if (iproc != 0) + { + ekb.resize(nks * nbands, 0.0); + occ.resize(nks * nbands, 0.0); + wk.resize(nks, 0.0); + kvec_c.resize(nks); + } + Parallel_Common::bcast_double(ekb.data(), nks * nbands); + Parallel_Common::bcast_double(occ.data(), nks * nbands); + Parallel_Common::bcast_double(wk.data(), nks); + // Vector3 is not a trivial datatype, need to be broadcasted element by element + for (int ik = 0; ik < nks; ik++) + { + Parallel_Common::bcast_double(kvec_c[ik].x); + Parallel_Common::bcast_double(kvec_c[ik].y); + Parallel_Common::bcast_double(kvec_c[ik].z); + } +} +// instantiate the template function +template void ModuleIO::restart_from_file(const std::string& out_dir, const Parallel_2D& p2d, const int& nks, + int& nbands, int& nbasis, std::vector& lowf_loc, + std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +template void ModuleIO::restart_from_file(const std::string& out_dir, const Parallel_2D& p2d, const int& nks, + int& nbands, int& nbasis, std::vector& lowf_loc, + std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +template void ModuleIO::restart_from_file(const std::string& out_dir, const Parallel_2D& p2d, const int& nks, + int& nbands, int& nbasis, std::vector>& lowf_loc, + std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +template void ModuleIO::restart_from_file(const std::string& out_dir, const Parallel_2D& p2d, const int& nks, + int& nbands, int& nbasis, std::vector>& lowf_loc, + std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +#endif + +template +void ModuleIO::restart_from_file(const std::string& out_dir, // hard-code the file name to be WFC_NAO_K*.txt? + const int& nks, int& nbands, int& nbasis, std::vector& lowf, + std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk) +{ + // reset vectors + lowf.clear(); + ekb.clear(); + occ.clear(); + kvec_c.clear(); + wk.clear(); + int nbands_ = -1, nbasis_ = -1; + for (int ik = 0; ik < nks; ik++) + { + // check existence of file + const std::string flowf = out_dir + "/WFC_NAO_K" + std::to_string(ik + 1) + ".txt"; + const std::ifstream ifs(flowf); + if (!ifs) + ModuleBase::WARNING_QUIT("restart_from_file", "open file failed: " + flowf); + + std::vector lowf_; + std::vector ekb_; + std::vector occ_; + ModuleBase::Vector3 kvec_; + double wk_; + int ik_; + read_abacus_lowf(flowf, ik_, kvec_, nbands, nbasis, lowf_, ekb_, occ_, wk_); + assert(nbands == nbands_ || nbands_ == -1); // check the consistency of nbands + assert(nbasis == nbasis_ || nbasis_ == -1); // check the consistency of nbasis + nbands_ = (nbands_ == -1) ? nbands : nbands_; + nbasis_ = (nbasis_ == -1) ? nbasis : nbasis_; + assert(ik_ == ik + 1); // check the consistency of ik + // append to the global lowf_loc + lowf.insert(lowf.end(), lowf_.begin(), lowf_.end()); + ekb.insert(ekb.end(), ekb_.begin(), ekb_.end()); + occ.insert(occ.end(), occ_.begin(), occ_.end()); + wk.push_back(wk_); + kvec_c.push_back(kvec_); + } + assert(ekb.size() == nks * nbands); + assert(occ.size() == nks * nbands); + assert(lowf.size() == nks * nbands * nbasis); +} +// instantiate the template function +template void ModuleIO::restart_from_file(const std::string& out_dir, const int& nks, int& nbands, int& nbasis, + std::vector& lowf, std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +template void ModuleIO::restart_from_file(const std::string& out_dir, const int& nks, int& nbands, int& nbasis, + std::vector& lowf, std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +template void ModuleIO::restart_from_file(const std::string& out_dir, const int& nks, int& nbands, int& nbasis, + std::vector>& lowf, std::vector& ekb, + std::vector& occ, std::vector>& kvec_c, + std::vector& wk); +template void ModuleIO::restart_from_file(const std::string& out_dir, const int& nks, int& nbands, int& nbasis, + std::vector>& lowf, std::vector& ekb, + std::vector& occ, std::vector>& kvec_c, + std::vector& wk); diff --git a/source/module_io/read_wfc_lcao.h b/source/module_io/read_wfc_lcao.h new file mode 100644 index 0000000000..0f96229524 --- /dev/null +++ b/source/module_io/read_wfc_lcao.h @@ -0,0 +1,113 @@ +#ifndef READ_WFC_LCAO_H +#define READ_WFC_LCAO_H + +// serial +#include "module_base/vector3.h" + +#include +#include +#include +#ifdef __MPI +// parallelization +#include "module_base/scalapack_connector.h" +#include "module_basis/module_ao/parallel_2d.h" +#endif + +/** + * @brief This class has two functions: restart psi from the previous calculation, and write psi to the disk. + * + */ +namespace ModuleIO +{ +// only when you know why you need the T, you can write function with template parameter, +// otherwise, you should overload the function for different types +// For example in this case, ONLY lowf support to be std::complex and std::complex, +// not ekb, occ, wk and kvec_c. +/** + * @brief Read the wavefunction coefficients from the file (for complex wavefunction coefficients) + * + * @tparam T + * @param flowf [in] file name like "LOWF_K_*.txt", dumped from ABACUS INPUT out_wfc_lcao 1 + * @param ik [out] the index of k points + * @param kvec_c [out] the k vector in Cartesian coordinates + * @param nbands [out] the number of bands + * @param nbasis [out] the number of orbitals + * @param lowf [out] wavefunction coefficients + * @param ekb [out] eigenvalues + * @param occ [out] occupations + * @param wk [out] weight of k points + */ +template +void read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, int& nbands, int& nbasis, + std::vector>& lowf, std::vector& ekb, std::vector& occ, + double& wk); +/** + * @brief Read the wavefunction coefficients from the file (for real wavefunction coefficients) + * + * @tparam T + * @param flowf [in] file name like "LOWF_K_*.txt", dumped from ABACUS INPUT out_wfc_lcao 1 + * @param ik [out] the index of k points + * @param kvec_c [out] the k vector in Cartesian coordinates + * @param nbands [out] the number of bands + * @param nbasis [out] the number of orbitals + * @param lowf [out] wavefunction coefficients + * @param ekb [out] eigenvalues + * @param occ [out] occupations + * @param wk [out] weight of k points + */ +template +void read_abacus_lowf(const std::string& flowf, int& ik, ModuleBase::Vector3& kvec_c, int& nbands, int& nbasis, + std::vector& lowf, std::vector& ekb, std::vector& occ, double& wk); +// the two functions above will return nbands, nbasis, lowf, ekb, occ and wk. +// the lowf is actually lowf_glb, which means the global matrix (ScaLAPACK convention), need to distribute +// to the local matrix (2D-block-cyclic parallel distribution) in the following function. + +// only-MPI-visible function, because the use of comm_world +#ifdef __MPI +/** + * @brief Restart the wavefunction coefficients from the file (MPI 2D-BCD version) + * + * @tparam T: datatype of the wavefunction coefficients, can be double, float, std::complex or + * std::complex + * @param out_dir [in] the directory where the wavefunction coefficients are stored + * @param p2d [in] the 2D parallel distribution + * @param nks [in] the number of k points + * @param nbands [out] the number of bands + * @param nbasis [out] the number of orbitals + * @param lowf_loc [out] the local wavefunction coefficients, can be used to construct psi, see constructor No.8 + * @param ekb [out] the eigenvalues + * @param occ [out] the occupations + * @param kvec_c [out] the k vectors in Cartesian coordinates + * @param wk [out] the weight of k points + * + * @warning Cpxgemr2d not implemented yet + */ +template +void restart_from_file(const std::string& out_dir, // hard-code the file name to be LOWF_K_*.txt? + const Parallel_2D& p2d, const int& nks, int& nbands, int& nbasis, std::vector& lowf_loc, + std::vector& ekb, std::vector& occ, + std::vector>& kvec_c, std::vector& wk); +#endif +// serial version, can always present +/** + * @brief Restart the wavefunction coefficients from the file (serial version) + * + * @tparam T: datatype of the wavefunction coefficients, can be double, float, std::complex or + * std::complex + * @param out_dir [in] the directory where the wavefunction coefficients are stored + * @param nks [in] the number of k points + * @param nbands [out] the number of bands + * @param nbasis [out] the number of orbitals + * @param lowf_loc [out] the local wavefunction coefficients, can be used to construct psi, see constructor No.8 + * @param ekb [out] the eigenvalues + * @param occ [out] the occupations + * @param kvec_c [out] the k vectors in Cartesian coordinates + * @param wk [out] the weight of k points + */ +template +void restart_from_file(const std::string& out_dir, // hard-code the file name to be LOWF_K_*.txt? + const int& nks, int& nbands, int& nbasis, std::vector& lowf, std::vector& ekb, + std::vector& occ, std::vector>& kvec_c, + std::vector& wk); +} // namespace ModuleIO +#endif \ No newline at end of file diff --git a/source/module_io/read_wfc_nao.cpp b/source/module_io/read_wfc_nao.cpp index 4cd0cdfeba..78146ec45d 100644 --- a/source/module_io/read_wfc_nao.cpp +++ b/source/module_io/read_wfc_nao.cpp @@ -20,16 +20,8 @@ * @param CTOT Global matrix from which the submatrix is extracted. * @return int Always returns 0 as a success indicator. */ -inline int CTOT2q(const int myid, - const int naroc[2], - const int nb, - const int dim0, - const int dim1, - const int iprow, - const int ipcol, - const int nbands, - double* work, - double** const CTOT) +inline int CTOT2q(const int myid, const int naroc[2], const int nb, const int dim0, const int dim1, const int iprow, + const int ipcol, const int nbands, double* work, double** const CTOT) { for (int j = 0; j < naroc[1]; ++j) { @@ -67,16 +59,8 @@ inline int CTOT2q(const int myid, * @param CTOT Global matrix from which the submatrix is extracted. * @return int Always returns 0 as a success indicator. */ -inline int CTOT2q_c(const int myid, - const int naroc[2], - const int nb, - const int dim0, - const int dim1, - const int iprow, - const int ipcol, - const int nbands, - std::complex* work, - std::complex** const CTOT) +inline int CTOT2q_c(const int myid, const int naroc[2], const int nb, const int dim0, const int dim1, const int iprow, + const int ipcol, const int nbands, std::complex* work, std::complex** const CTOT) { for (int j = 0; j < naroc[1]; ++j) { @@ -99,21 +83,15 @@ inline int CTOT2q_c(const int myid, } // be called in local_orbital_wfc::allocate_k -int ModuleIO::read_wfc_nao_complex(std::complex** ctot, - const int& ik, - const int& nb2d, - const int& nbands_g, - const int& nlocal_g, - const std::string& global_readin_dir, - const ModuleBase::Vector3 kvec_c, - const Parallel_Orbitals* ParaV, - psi::Psi>* psi, - elecstate::ElecState* pelec) +int ModuleIO::read_wfc_nao_complex(std::complex** ctot, const int& ik, const int& nb2d, const int& nbands_g, + const int& nlocal_g, const std::string& global_readin_dir, + const ModuleBase::Vector3 kvec_c, const Parallel_Orbitals* ParaV, + psi::Psi>* psi, elecstate::ElecState* pelec) { ModuleBase::TITLE("ModuleIO", "read_wfc_nao_complex"); ModuleBase::timer::tick("ModuleIO", "read_wfc_nao_complex"); - std::string ss = global_readin_dir + ModuleIO::wfc_nao_gen_fname(1,false,true,ik); + std::string ss = global_readin_dir + ModuleIO::wfc_nao_gen_fname(1, false, true, ik); std::ifstream ifs; int error = 0; @@ -206,6 +184,7 @@ int ModuleIO::read_wfc_nao_complex(std::complex** ctot, #ifdef __MPI Parallel_Common::bcast_int(error); + // then broadcast the "weigths" of present k-point, with length pelec->wg.nc Parallel_Common::bcast_double(&pelec->wg.c[ik * pelec->wg.nc], pelec->wg.nc); #endif @@ -257,16 +236,9 @@ int ModuleIO::read_wfc_nao_complex(std::complex** ctot, return 0; } -int ModuleIO::read_wfc_nao(double** ctot, - const int& is, - const bool& gamma_only_local, - const int& nb2d, - const int& nbands_g, - const int& nlocal_g, - const std::string& global_readin_dir, - const Parallel_Orbitals* ParaV, - psi::Psi* psid, - elecstate::ElecState* pelec) +int ModuleIO::read_wfc_nao(double** ctot, const int& is, const bool& gamma_only_local, const int& nb2d, + const int& nbands_g, const int& nlocal_g, const std::string& global_readin_dir, + const Parallel_Orbitals* ParaV, psi::Psi* psid, elecstate::ElecState* pelec) { ModuleBase::TITLE("ModuleIO", "read_wfc_nao"); ModuleBase::timer::tick("ModuleIO", "read_wfc_nao"); @@ -378,13 +350,8 @@ int ModuleIO::read_wfc_nao(double** ctot, return 0; } -void ModuleIO::distri_wfc_nao(double** ctot, - const int& is, - const int& nb2d, - const int& nbands_g, - const int& nlocal_g, - const Parallel_Orbitals* ParaV, - psi::Psi* psid) +void ModuleIO::distri_wfc_nao(double** ctot, const int& is, const int& nb2d, const int& nbands_g, const int& nlocal_g, + const Parallel_Orbitals* ParaV, psi::Psi* psid) { ModuleBase::TITLE("ModuleIO", "distri_wfc_nao"); #ifdef __MPI @@ -455,12 +422,8 @@ void ModuleIO::distri_wfc_nao(double** ctot, return; } -void ModuleIO::distri_wfc_nao_complex(std::complex** ctot, - const int& ik, - const int& nb2d, - const int& nbands_g, - const Parallel_Orbitals* ParaV, - psi::Psi>* psi) +void ModuleIO::distri_wfc_nao_complex(std::complex** ctot, const int& ik, const int& nb2d, const int& nbands_g, + const Parallel_Orbitals* ParaV, psi::Psi>* psi) { ModuleBase::TITLE("ModuleIO", "distri_wfc_nao_complex"); #ifdef __MPI diff --git a/source/module_io/test/CMakeLists.txt b/source/module_io/test/CMakeLists.txt index fb8594dea5..e10a2ad1bc 100644 --- a/source/module_io/test/CMakeLists.txt +++ b/source/module_io/test/CMakeLists.txt @@ -186,3 +186,17 @@ AddTest( LIBS base ${math_libs} device numerical_atomic_orbitals container orb SOURCES numerical_basis_test.cpp ../numerical_basis_jyjy.cpp ) + +if(ENABLE_LCAO) +AddTest( + TARGET read_wfc_lcao_test + LIBS base ${math_libs} device + SOURCES read_wfc_lcao_test.cpp ../read_wfc_lcao.cpp ../../module_basis/module_ao/parallel_2d.cpp +) + +add_test(NAME read_wfc_lcao_test_parallel + COMMAND mpirun -np 4 ./read_wfc_lcao_test + WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} +) + +endif() \ No newline at end of file diff --git a/source/module_io/test/read_wfc_lcao_test.cpp b/source/module_io/test/read_wfc_lcao_test.cpp new file mode 100644 index 0000000000..cd585e1d3c --- /dev/null +++ b/source/module_io/test/read_wfc_lcao_test.cpp @@ -0,0 +1,772 @@ +#include "module_io/read_wfc_lcao.h" + +#include "gtest/gtest.h" + +TEST(ReadWfcLcaoTest, ReadAbacusLowfComplex) +{ + + // this test should only be executed on rank 0 +#ifdef __MPI + MPI_Barrier(MPI_COMM_WORLD); + // printf("MPI environment detected, will use only rank 0\n"); + int iproc = 0; + int nprocs = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &iproc); + MPI_Comm_size(MPI_COMM_WORLD, &nprocs); + if (iproc != 0) + { + GTEST_SKIP(); + } +#endif + int ik = 0; + ModuleBase::Vector3 kvec_c; + int nbands = 0; + int nbasis = 0; + std::vector> lowf; + std::vector ekb; + std::vector occ; + double wk = -1.0; + + // first test + std::string flowf = "./support/WFC_NAO_K1.txt"; + ModuleIO::read_abacus_lowf(flowf, ik, kvec_c, nbands, nbasis, lowf, ekb, occ, wk); + EXPECT_EQ(1, ik); + EXPECT_EQ(3, nbands); + EXPECT_EQ(63, nbasis); + EXPECT_EQ(1.0, wk); // this is overwritten with 1.0 + // kvec_c + EXPECT_NEAR(-0.10540388, kvec_c.x, 1e-7); + EXPECT_NEAR(-0.060854959, kvec_c.y, 1e-7); + EXPECT_NEAR(-0.043030954, kvec_c.z, 1e-7); + // ekb + EXPECT_EQ(ekb.size(), nbands); + EXPECT_NEAR(-6.03571945e-01, ekb[0], 1e-7); + EXPECT_NEAR(-5.98035141e-01, ekb[1], 1e-7); + EXPECT_NEAR(-5.98035141e-01, ekb[2], 1e-7); + // occ + EXPECT_EQ(occ.size(), nbands); + EXPECT_NEAR(5.83090379e-03, occ[0], 1e-7); + EXPECT_NEAR(5.83090379e-03, occ[1], 1e-7); + EXPECT_NEAR(5.83090379e-03, occ[2], 1e-7); + // lowf + EXPECT_EQ(lowf.size(), nbands * nbasis); + EXPECT_NEAR(lowf[0].real(), -6.71651157e-03, 1e-7); + EXPECT_NEAR(lowf[0].imag(), 2.25946383e-02, 1e-7); + EXPECT_NEAR(lowf[1].real(), 1.43180123e-03, 1e-7); + EXPECT_NEAR(lowf[1].imag(), 1.46451488e-03, 1e-7); + EXPECT_NEAR(lowf[2].real(), 2.31452033e-03, 1e-7); + EXPECT_NEAR(lowf[2].imag(), -1.18949691e-03, 1e-7); + EXPECT_NEAR(lowf[62].real(), 1.82648757e-03, 1e-7); + EXPECT_NEAR(lowf[62].imag(), -2.11799886e-03, 1e-7); + // test reuse, expect to overwrite the previous values + flowf = "./support/WFC_NAO_K2.txt"; + ModuleIO::read_abacus_lowf(flowf, ik, kvec_c, nbands, nbasis, lowf, ekb, occ, wk); + EXPECT_EQ(2, ik); + EXPECT_EQ(3, nbands); + EXPECT_EQ(63, nbasis); + EXPECT_EQ(1.0, wk); // this is overwritten with 1.0 + // kvec_c + EXPECT_NEAR(-0.070269254, kvec_c.x, 1e-7); + EXPECT_NEAR(-0.081139946, kvec_c.y, 1e-7); + EXPECT_NEAR(-0.057374606, kvec_c.z, 1e-7); + // ekb + EXPECT_EQ(ekb.size(), nbands); + EXPECT_NEAR(-6.03667277e-01, ekb[0], 1e-7); + EXPECT_NEAR(-5.97868276e-01, ekb[1], 1e-7); + EXPECT_NEAR(-5.97662421e-01, ekb[2], 1e-7); + // occ + EXPECT_EQ(occ.size(), nbands); + EXPECT_NEAR(5.83090379e-03, occ[0], 1e-7); + EXPECT_NEAR(5.83090379e-03, occ[1], 1e-7); + EXPECT_NEAR(5.83090379e-03, occ[2], 1e-7); + // lowf + EXPECT_EQ(lowf.size(), nbands * nbasis); + EXPECT_NEAR(2.09933705e-02, lowf[0].real(), 1e-7); + EXPECT_NEAR(2.20619371e-03, lowf[0].imag(), 1e-7); + EXPECT_NEAR(1.52454416e-03, lowf[1].real(), 1e-7); + EXPECT_NEAR(-3.54139105e-04, lowf[1].imag(), 1e-7); + EXPECT_NEAR(1.31198443e-03, lowf[2].real(), 1e-7); + EXPECT_NEAR(-2.44872538e-03, lowf[2].imag(), 1e-7); + EXPECT_NEAR(lowf[62].real(), -1.15158489e-03, 1e-7); + EXPECT_NEAR(lowf[62].imag(), -1.79940038e-03, 1e-7); + // test reuse, the second time + flowf = "./support/WFC_NAO_K3.txt"; + ModuleIO::read_abacus_lowf(flowf, ik, kvec_c, nbands, nbasis, lowf, ekb, occ, wk); + EXPECT_EQ(3, ik); + EXPECT_EQ(3, nbands); + EXPECT_EQ(63, nbasis); + EXPECT_EQ(1.0, wk); // this is overwritten with 1.0 + // kvec_c + EXPECT_NEAR(-0.035134627, kvec_c.x, 1e-7); + EXPECT_NEAR(-0.10142493, kvec_c.y, 1e-7); + EXPECT_NEAR(-0.07171825, kvec_c.z, 1e-7); + // ekb + EXPECT_NEAR(-6.04664544e-01, ekb[0], 1e-7); + EXPECT_NEAR(-5.97025474e-01, ekb[1], 1e-7); + EXPECT_NEAR(-5.96870018e-01, ekb[2], 1e-7); + // occ + EXPECT_NEAR(5.83090379e-03, occ[0], 1e-7); + EXPECT_NEAR(5.83090379e-03, occ[1], 1e-7); + EXPECT_NEAR(5.83090379e-03, occ[2], 1e-7); + // lowf + EXPECT_NEAR(-6.44410072e-03, lowf[0].real(), 1e-7); + EXPECT_NEAR(2.86108252e-03, lowf[0].imag(), 1e-7); + EXPECT_NEAR(-5.81398415e-03, lowf[1].real(), 1e-7); + EXPECT_NEAR(4.01495705e-03, lowf[1].imag(), 1e-7); + EXPECT_NEAR(-2.32158666e-03, lowf[2].real(), 1e-7); + EXPECT_NEAR(2.62541166e-03, lowf[2].imag(), 1e-7); + EXPECT_NEAR(lowf[62].real(), 5.58964902e-04, 1e-7); + EXPECT_NEAR(lowf[62].imag(), 5.21866389e-04, 1e-7); +} + +TEST(ReadWfcLcaoTest, Cpzgemr2dUseTest) +{ +/* +(0,0) (0,1) (0,2) (0,3) (0,4) (0,5) (0,6) (0,7) (0,8) (0,9) (0,10) (0,11) +(1,0) (1,1) (1,2) (1,3) (1,4) (1,5) (1,6) (1,7) (1,8) (1,9) (1,10) (1,11) +(2,0) (2,1) (2,2) (2,3) (2,4) (2,5) (2,6) (2,7) (2,8) (2,9) (2,10) (2,11) +(3,0) (3,1) (3,2) (3,3) (3,4) (3,5) (3,6) (3,7) (3,8) (3,9) (3,10) (3,11) +... +(9,0) (9,1) (9,2) (9,3) (9,4) (9,5) (9,6) (9,7) (9,8) (9,9) (9,10) (9,11) +*/ +#ifdef __MPI + // this test should be run on all ranks + int iproc = 0; + int nprocs = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &iproc); + MPI_Comm_size(MPI_COMM_WORLD, &nprocs); + if (nprocs != 4) + { + if (iproc == 0) + { + printf("Please run this unittest with nprocs = 4.\n"); + } + GTEST_SKIP(); + } + // one get data on rank 0 + std::vector> lowf_glb; + const int nbands = 10; // nrow_glb + const int nbasis = 12; // ncol_glb + if (iproc == 0) + { + printf("Run unittest ScatterLowfTest::ScatterLowfComplex with MPI env:\n"); + printf("Total number of processes: %d\n\n", nprocs); + printf("Row-major processor grid is used.\n\n"); + printf("First test the \"column-major\" matrix, which means for columns their memory\n"); + printf("are consecutive. The matrix in form of (i, j), where i runs over [0, 9] and \n"); + printf("j runs over [0, 11]: (0,0), (1,0), (2,0), ..., (9,0), (0,1), ...\n"); + lowf_glb.resize(nbands * nbasis); + for (int i = 0; i < nbands * nbasis; i++) + { + const int j = i / nbands, k = i % nbands; + lowf_glb[i] = std::complex(k, j); + } + } + // the other ranks get data + MPI_Barrier(MPI_COMM_WORLD); + std::vector> lowf_loc; + Parallel_2D para2d_test; // an alternative to ParaV. But to use paraV, the desc would be desc_wfc instead of desc in + // Parallel_2D + // initialize a para2d, as if it is paraV. + // BE CAREFUL! One must be careful about defining the dimension properly. Now the original + // matrix is made column-memory-consecutive, thus the vectorized data must be cut by "the + // number of rows", then by "the number of columns". + // nbands is "the number of rows", nbasis is "the number of columns" + para2d_test.init(nbands, nbasis, 4, MPI_COMM_WORLD); + Parallel_2D para2d_glb; + para2d_glb.init(nbands, nbasis, std::max(nbands, nbasis), MPI_COMM_WORLD); + lowf_loc.resize(para2d_test.nrow * para2d_test.ncol); // the nrow and ncol are automatically + // calculated by Parallel_2D + // wait, what is the meaning of nrow here? The "nrow" again is just the number to cut + // the vectorized data into a matrix. The "nrow" is the number of rows in the matrix but + // remember the matrix is column-memory-consecutive. + + // the following function can do the scattering-gathering automatically. + // a double counterpart is pdgemr2d_, int counterpart is pigemr2d_ + // Those in C style are Cpzgemr2d, Cdgemr2d, Cigemr2d... + Cpxgemr2d(nbands, nbasis, lowf_glb.data(), 1, 1, const_cast(para2d_glb.desc), lowf_loc.data(), 1, 1, + const_cast(para2d_test.desc), para2d_glb.blacs_ctxt); + // what will happen if impose a row-major processor grid onto a column-major matrix? + // you can get correct results, expect each block is column-major. + // Have a look: + MPI_Barrier(MPI_COMM_WORLD); + std::vector sizes_loc = {4 * 4 * 3, 4 * 4 + 4 * 2, 4 * 4 * 2, 4 * 4}; + std::vector> reals = {{0, 1, 2, 3, 8, 9}, {0, 1, 2, 3, 8, 9}, {4, 5, 6, 7}, {4, 5, 6, 7}}; + std::vector> imags + = {{0, 1, 2, 3, 8, 9, 10, 11}, {4, 5, 6, 7}, {0, 1, 2, 3, 8, 9, 10, 11}, {4, 5, 6, 7}}; + for (int i = 0; i < nprocs; i++) + { + if (iproc == i) + { + EXPECT_EQ(lowf_loc.size(), sizes_loc[i]); + printf(">>> rank %d: \n", iproc); + printf("First print scattered matrix in the way that ELEMENTS WITH CONSECUTIVE\n"); + printf("MEMORIES ARE SHOWN (only shown) IN THE SAME LINE:\n"); + for (int j = 0; j < lowf_loc.size(); j++) + { + printf("(%2.0f,%2.0f)", lowf_loc[j].real(), lowf_loc[j].imag()); + if ((j + 1) % para2d_test.nrow == 0) + { + printf("\n"); + } + const int k = j % para2d_test.nrow; + EXPECT_NEAR(lowf_loc[j].real(), reals[i][k], 1e-7); + const int l = j / para2d_test.nrow; + EXPECT_NEAR(lowf_loc[j].imag(), imags[i][l], 1e-7); + } + printf("Or INDEX IT to show like \"row-major\":\n"); + // (i, j) -> (i', j') with i = j' and j = i' + // x = i*ncol + j, x' = i'*ncol' + j' with ncol' = nrow and nrow' = ncol + // i = x/ncol, j = x%ncol, x' = j*nrow + i = x%ncol*nrow + x/ncol + for (int j = 0; j < lowf_loc.size(); j++) + { + const int x = j % para2d_test.ncol * para2d_test.nrow + j / para2d_test.ncol; + printf("(%2.0f,%2.0f)", lowf_loc[x].real(), lowf_loc[x].imag()); + if ((j + 1) % para2d_test.ncol == 0) + { + printf("\n"); + } + } + printf("\n"); + usleep(10000); + } + MPI_Barrier(MPI_COMM_WORLD); + } + MPI_Barrier(MPI_COMM_WORLD); + + // test the other way around, the row-major matrix + if (iproc == 0) + { + printf("Now test the \"row-major\" matrix, which means for rows their memory\n"); + printf("are consecutive. The matrix in form of (i, j), where i runs over [0, 9] and \n"); + printf("j runs over [0, 11]: (0,0), (0,1), (0,2), ..., (0,11), (1,0), ...\n"); + for (int i = 0; i < nbands * nbasis; i++) + { + const int irow = i / nbasis, icol = i % nbasis; + lowf_glb[i] = std::complex(irow, icol); + } + } + // the other ranks get data + MPI_Barrier(MPI_COMM_WORLD); + // initialize a para2d, as if it is paraV. + Parallel_2D para2d_test_prime; + // note! this time the memory is first cut by "the number of columns", + // then by "the number of rows". Therefore the "nbasis" is put the first. + // This is how ScaLAPCK defines a matrix: the first number defines the leading dimension. + para2d_test_prime.init(nbasis, nbands, 4, MPI_COMM_WORLD); + Parallel_2D para2d_glb_prime; + para2d_glb_prime.init(nbasis, nbands, std::max(nbands, nbasis), MPI_COMM_WORLD); + lowf_loc.resize(para2d_test_prime.nrow * para2d_test_prime.ncol); + Cpxgemr2d(nbasis, nbands, lowf_glb.data(), 1, 1, const_cast(para2d_glb_prime.desc), lowf_loc.data(), 1, 1, + const_cast(para2d_test_prime.desc), para2d_glb_prime.blacs_ctxt); + MPI_Barrier(MPI_COMM_WORLD); + sizes_loc = {4 * 4 * 3, 4 * 4 * 2, 4 * 4 + 4 * 2, 4 * 4}; + reals = {{0, 1, 2, 3, 8, 9}, {4, 5, 6, 7}, {0, 1, 2, 3, 8, 9}, {4, 5, 6, 7}}; + imags = {{0, 1, 2, 3, 8, 9, 10, 11}, {0, 1, 2, 3, 8, 9, 10, 11}, {4, 5, 6, 7}, {4, 5, 6, 7}}; + for (int i = 0; i < nprocs; i++) + { + if (iproc == i) + { + EXPECT_EQ(lowf_loc.size(), sizes_loc[i]); + printf(">>> rank %d: \n", iproc); + for (int j = 0; j < lowf_loc.size(); j++) + { + printf("(%2.0f,%2.0f)", lowf_loc[j].real(), lowf_loc[j].imag()); + if ((j + 1) % para2d_test_prime.nrow == 0) + { + printf("\n"); + } + const int k = j / para2d_test_prime.nrow; + EXPECT_NEAR(lowf_loc[j].real(), reals[i][k], 1e-7); + const int l = j % para2d_test_prime.nrow; + EXPECT_NEAR(lowf_loc[j].imag(), imags[i][l], 1e-7); + } + printf("\n"); + usleep(10000); + } + MPI_Barrier(MPI_COMM_WORLD); + } + MPI_Barrier(MPI_COMM_WORLD); + if (iproc == 0) + { + printf("BE CAREFUL!\n"); + printf("You note that the PROCESSOR GRID seems to be transposed. It is because\n"); + printf("in C/C++ it is always assumed memory in the same row is consecutive, while\n"); + printf("in FORTRAN or \"what ScaLAPACK supposes\" it is column-memory-consecutive.\n"); + usleep(10000); + } + MPI_Barrier(MPI_COMM_WORLD); +#else + printf("Run unittest ScatterLowfTest::ScatterLowfComplex without MPI env:\n"); + printf("This test is not executed because MPI is not enabled.\n"); + GTEST_SKIP(); +#endif +} + +TEST(ReadWfcLcaoTest, ReadAbacusLowfReal) +{ + // this test should only be executed on rank 0 +#ifdef __MPI + MPI_Barrier(MPI_COMM_WORLD); + // printf("MPI environment detected, will use only rank 0\n"); + int iproc = 0; + int nprocs = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &iproc); + MPI_Comm_size(MPI_COMM_WORLD, &nprocs); + if (iproc != 0) + { + GTEST_SKIP(); + } +#endif + int ik = 1; // should be overwritten to 0 + ModuleBase::Vector3 kvec_c; + int nbands = 0; + int nbasis = 0; + std::vector lowf; + std::vector ekb; + std::vector occ; + double wk = -1.0; // should be overwritten to 1.0 + + // first test + const std::string flowf = "./support/WFC_NAO_GAMMA1.txt"; + ModuleIO::read_abacus_lowf(flowf, ik, kvec_c, nbands, nbasis, lowf, ekb, occ, wk); + EXPECT_EQ(0, ik); + EXPECT_EQ(3, nbands); + EXPECT_EQ(31, nbasis); + EXPECT_EQ(1.0, wk); + // kvec_c, gamma point, 0, 0, 0 + EXPECT_NEAR(0.0, kvec_c.x, 1e-7); + EXPECT_NEAR(0.0, kvec_c.y, 1e-7); + EXPECT_NEAR(0.0, kvec_c.z, 1e-7); + // ekb + EXPECT_EQ(ekb.size(), nbands); + EXPECT_NEAR(-1.22787155e+00, ekb[0], 1e-7); + EXPECT_NEAR(-3.10595658e-01, ekb[1], 1e-7); + EXPECT_NEAR(-3.00546690e-01, ekb[2], 1e-7); + // occ + EXPECT_EQ(occ.size(), nbands); + EXPECT_NEAR(2.00000000e+00, occ[0], 1e-7); + EXPECT_NEAR(2.00000000e+00, occ[1], 1e-7); + EXPECT_NEAR(2.00000000e+00, occ[2], 1e-7); + // lowf + EXPECT_EQ(lowf.size(), nbands * nbasis); + EXPECT_NEAR(-1.51728369e-02, lowf[0], 1e-7); + EXPECT_NEAR(-2.07808444e-03, lowf[1], 1e-7); + EXPECT_NEAR(1.21298954e-17, lowf[2], 1e-7); + EXPECT_NEAR(-5.44883791e-09, lowf[30], 1e-7); +} + +TEST(ReadWfcLcaoTest, Cpdgemr2dUseTest) +{ + // you can find more information in unittest Pzgemr2dUseTest, present test + // works identically to the previous one, but with real numbers. +#ifdef __MPI + // this test should be run on all ranks + int iproc = 0; + int nprocs = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &iproc); + MPI_Comm_size(MPI_COMM_WORLD, &nprocs); + if (nprocs != 4) + { + if (iproc == 0) + { + printf("Please run this unittest with nprocs = 4.\n"); + } + GTEST_SKIP(); + } + std::vector lowf_glb; + const int nbands = 10; + const int nbasis = 12; + // still, the expected matrix is organized as row-memory-consecutive but here + // just make the matrix column-memory-consecutive. + // x = i*ncol + j, x' = j*nrow + i + // i = x/ncol, j = x%ncol, x' = j*nrow + i = x%ncol*nrow + x/ncol + if (iproc == 0) + { + lowf_glb.resize(nbands * nbasis); + for (int i = 0; i < nbands * nbasis; i++) + { + lowf_glb[i] = i % nbasis * nbands + i / nbasis; + } + } + MPI_Barrier(MPI_COMM_WORLD); + std::vector lowf_loc; + Parallel_2D para2d_test; + para2d_test.init(nbasis, nbands, 4, MPI_COMM_WORLD); + Parallel_2D para2d_glb; + para2d_glb.init(nbasis, nbands, std::max(nbands, nbasis), MPI_COMM_WORLD); + lowf_loc.resize(para2d_test.nrow * para2d_test.ncol); + Cpxgemr2d(nbasis, nbands, lowf_glb.data(), 1, 1, const_cast(para2d_glb.desc), lowf_loc.data(), 1, 1, + const_cast(para2d_test.desc), para2d_glb.blacs_ctxt); + MPI_Barrier(MPI_COMM_WORLD); + std::vector sizes_loc = {4 * 4 * 3, 4 * 4 * 2, 4 * 4 + 4 * 2, 4 * 4}; + for (int i = 0; i < nprocs; i++) + { + if (iproc == i) + { + EXPECT_EQ(lowf_loc.size(), sizes_loc[i]); + } + } + MPI_Barrier(MPI_COMM_WORLD); + // test the other way around, the row-major matrix + if (iproc == 0) + { + for (int i = 0; i < nbands * nbasis; i++) + { + lowf_glb[i] = i; + } + } + MPI_Barrier(MPI_COMM_WORLD); + Parallel_2D para2d_test_prime; + para2d_test_prime.init(nbands, nbasis, 4, MPI_COMM_WORLD); + Parallel_2D para2d_glb_prime; + para2d_glb_prime.init(nbands, nbasis, std::max(nbands, nbasis), MPI_COMM_WORLD); + lowf_loc.resize(para2d_test_prime.nrow * para2d_test_prime.ncol); + Cpxgemr2d(nbands, nbasis, lowf_glb.data(), 1, 1, const_cast(para2d_glb_prime.desc), lowf_loc.data(), 1, 1, + const_cast(para2d_test_prime.desc), para2d_glb_prime.blacs_ctxt); + MPI_Barrier(MPI_COMM_WORLD); + sizes_loc = {4 * 4 * 3, 4 * 4 + 4 * 2, 4 * 4 * 2, 4 * 4}; + for (int i = 0; i < nprocs; i++) + { + if (iproc == i) + { + EXPECT_EQ(lowf_loc.size(), sizes_loc[i]); + } + } + MPI_Barrier(MPI_COMM_WORLD); +#endif +} + +TEST(ReadWfcLcaoTest, RestartFromFileParallel) +{ +#ifdef __MPI + MPI_Barrier(MPI_COMM_WORLD); + // printf("MPI environment detected, will use only rank 0\n"); + int iproc = 0; + int nprocs = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &iproc); + MPI_Comm_size(MPI_COMM_WORLD, &nprocs); + if (nprocs != 4) + { + if (iproc == 0) + { + printf("Please run this unittest with nprocs = 4.\n"); + } + GTEST_SKIP(); + } + Parallel_2D para2d_test; + // this time the nbands and nbasis are given as priori, from reading of INPUT and of orb files. + para2d_test.init(63, 3, 2, MPI_COMM_WORLD); + // therefore it is a nrows x ncols matrix, where nr = 3, nc = 63 + // scatter it with block size 2 x 2, on 4-processors grid. From testcase ReadWfcLcaoTest::Cpzgemr2dTest + // it is known the processor grid will be transposed to column-major, thus + // the global matrix will be split to 2 x 2 blocks, one 2 x 1 block, 1 x 2 blocks and 1 x 1 block + // number of elements: + // rank0: ceil(63/2)*2 = 64 + // rank1: ceil(63/2)*1 = 32 + // rank2: floor(63/2)*2 = 62 + // rank3: floor(63/2)*1 = 31 + // total size check + if (iproc == 0) + { + EXPECT_EQ(64, para2d_test.nrow * para2d_test.ncol); + } + if (iproc == 1) + { + EXPECT_EQ(32, para2d_test.nrow * para2d_test.ncol); + } + if (iproc == 2) + { + EXPECT_EQ(62, para2d_test.nrow * para2d_test.ncol); + } + if (iproc == 3) + { + EXPECT_EQ(31, para2d_test.nrow * para2d_test.ncol); + } + // nrows check + if (iproc == 0) + { + EXPECT_EQ(2, para2d_test.ncol); + } + if (iproc == 1) + { + EXPECT_EQ(1, para2d_test.ncol); + } + if (iproc == 2) + { + EXPECT_EQ(2, para2d_test.ncol); + } + if (iproc == 3) + { + EXPECT_EQ(1, para2d_test.ncol); + } + + const std::string out_dir = "./support"; + const int nks = 4; + int nbands = -1; + int nbasis = -1; + std::vector> lowf; + std::vector ekb; + std::vector occ; + std::vector> kvec_c; + std::vector wk; + ModuleIO::restart_from_file(out_dir, para2d_test, nks, nbands, nbasis, lowf, ekb, occ, kvec_c, wk); + EXPECT_EQ(3, nbands); + EXPECT_EQ(63, nbasis); + // ekb, occ, kvec_c, wk will have size irrelevant to scatter + EXPECT_EQ(nks * nbands, ekb.size()); + EXPECT_EQ(nks * nbands, occ.size()); + EXPECT_EQ(nks, kvec_c.size()); + EXPECT_EQ(nks, wk.size()); + // value test + const std::vector ekb_ref + = {-6.03571945e-01, -5.98035141e-01, -5.98035141e-01, -6.03667277e-01, -5.97868276e-01, -5.97662421e-01, + -6.04664544e-01, -5.97025474e-01, -5.96870018e-01, -6.05615293e-01, -5.96302906e-01, -5.96302906e-01}; + const std::vector occ_ref + = {5.83090379e-03, 5.83090379e-03, 5.83090379e-03, 5.83090379e-03, 5.83090379e-03, 5.83090379e-03, + 5.83090379e-03, 5.83090379e-03, 5.83090379e-03, 5.83090379e-03, 5.83090379e-03, 5.83090379e-03}; + const std::vector> kvec_c_ref + = {ModuleBase::Vector3(-0.10540388, -0.060854959, -0.043030954), + ModuleBase::Vector3(-0.070269254, -0.081139946, -0.057374606), + ModuleBase::Vector3(-0.035134627, -0.10142493, -0.07171825), + ModuleBase::Vector3(0.00000000, -0.12170991, -0.086061909)}; + const std::vector wk_ref = {1.0, 1.0, 1.0, 1.0}; + for (int i = 0; i < nprocs; i++) + { + if (iproc == i) + { + // ekb + for (int j = 0; j < nks * nbands; j++) + { + EXPECT_NEAR(ekb_ref[j], ekb[j], 1e-7); + } + // occ + for (int j = 0; j < nks * nbands; j++) + { + EXPECT_NEAR(occ_ref[j], occ[j], 1e-7); + } + // kvec_c + for (int j = 0; j < nks; j++) + { + EXPECT_NEAR(kvec_c_ref[j].x, kvec_c[j].x, 1e-7); + EXPECT_NEAR(kvec_c_ref[j].y, kvec_c[j].y, 1e-7); + EXPECT_NEAR(kvec_c_ref[j].z, kvec_c[j].z, 1e-7); + } + // wk + for (int j = 0; j < nks; j++) + { + EXPECT_NEAR(wk_ref[j], wk[j], 1e-7); + } + } + MPI_Barrier(MPI_COMM_WORLD); + } + // lowf will be the result of scattering + // rank0: 64, rank1: 32, rank2: 62, rank3: 31, each should multiply the number of k-points + if (iproc == 0) + { + EXPECT_EQ(nks * 64, lowf.size()); + } + if (iproc == 1) + { + EXPECT_EQ(nks * 32, lowf.size()); + } + if (iproc == 2) + { + EXPECT_EQ(nks * 62, lowf.size()); + } + if (iproc == 3) + { + EXPECT_EQ(nks * 31, lowf.size()); + } + // value test on lowf + // rank0 + if (iproc == 0) + { + EXPECT_NEAR(lowf[0].real(), -6.71651157e-03, 1e-7); + EXPECT_NEAR(lowf[0].imag(), 2.25946383e-02, 1e-7); + EXPECT_NEAR(lowf[1].real(), 1.43180123e-03, 1e-7); + EXPECT_NEAR(lowf[1].imag(), 1.46451488e-03, 1e-7); + EXPECT_NEAR(lowf[2].real(), -9.45760546e-03, 1e-7); + EXPECT_NEAR(lowf[2].imag(), 1.29911511e-02, 1e-7); + EXPECT_NEAR(lowf[3].real(), -5.46035106e-03, 1e-7); + EXPECT_NEAR(lowf[3].imag(), 7.50044462e-03, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].real(), -1.39799597e-03, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].imag(), -1.68192980e-03, 1e-7); + } + else if (iproc == 1) + { + EXPECT_NEAR(lowf[0].real(), 9.86470874e-13, 1e-7); + EXPECT_NEAR(lowf[0].imag(), -5.95387122e-12, 1e-7); + EXPECT_NEAR(lowf[1].real(), 4.82573453e-13, 1e-7); + EXPECT_NEAR(lowf[1].imag(), 5.54264959e-12, 1e-7); + EXPECT_NEAR(lowf[2].real(), 5.20920946e-04, 1e-7); + EXPECT_NEAR(lowf[2].imag(), -2.33076310e-03, 1e-7); + EXPECT_NEAR(lowf[3].real(), -8.35155442e-04, 1e-7); + EXPECT_NEAR(lowf[3].imag(), 3.75083842e-03, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].real(), -6.18118243e-05, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].imag(), -7.43658388e-05, 1e-7); + } + else if (iproc == 2) + { + EXPECT_NEAR(lowf[0].real(), 2.31452033e-03, 1e-7); + EXPECT_NEAR(lowf[0].imag(), -1.18949691e-03, 1e-7); + EXPECT_NEAR(lowf[1].real(), 3.86105126e-03, 1e-7); + EXPECT_NEAR(lowf[1].imag(), -5.30361525e-03, 1e-7); + EXPECT_NEAR(lowf[2].real(), 1.07440727e-03, 1e-7); + EXPECT_NEAR(lowf[2].imag(), -1.52230629e-03, 1e-7); + EXPECT_NEAR(lowf[3].real(), -2.63174959e-03, 1e-7); + EXPECT_NEAR(lowf[3].imag(), 3.72887365e-03, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].real(), 1.19038759e-04, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].imag(), 1.17824924e-04, 1e-7); + } + else if (iproc == 3) + { + EXPECT_NEAR(lowf[0].real(), 3.66087151e-13, 1e-7); + EXPECT_NEAR(lowf[0].imag(), 1.96386245e-13, 1e-7); + EXPECT_NEAR(lowf[1].real(), 9.49023673e-05, 1e-7); + EXPECT_NEAR(lowf[1].imag(), -4.04693771e-04, 1e-7); + EXPECT_NEAR(lowf[2].real(), 1.33229060e-04, 1e-7); + EXPECT_NEAR(lowf[2].imag(), 9.69176971e-04, 1e-7); + EXPECT_NEAR(lowf[3].real(), 8.23664081e-04, 1e-7); + EXPECT_NEAR(lowf[3].imag(), 5.56014508e-03, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].real(), -2.69229582e-03, 1e-7); + EXPECT_NEAR(lowf[lowf.size() - 1].imag(), -2.66484241e-03, 1e-7); + } +#else + printf("Run unittest ReadWfcLcaoTest::RestartFromFileParallel without MPI env:\n"); + printf("This test is not executed because MPI is not enabled.\n"); + GTEST_SKIP(); +#endif +} + +TEST(ReadWfcLcaoTest, RestartFromFileSerial) +{ +#ifdef __MPI + MPI_Barrier(MPI_COMM_WORLD); + // printf("MPI environment detected, will use only rank 0\n"); + int iproc = 0; + int nprocs = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &iproc); + MPI_Comm_size(MPI_COMM_WORLD, &nprocs); + if (iproc != 0) + { + GTEST_SKIP(); + } +#endif + const int nks = 4; + int nbands = -1; + int nbasis = -1; + std::vector> lowf; + std::vector ekb; + std::vector occ; + std::vector> kvec_c; + std::vector wk; + const std::string out_dir = "./support"; + ModuleIO::restart_from_file(out_dir, nks, nbands, nbasis, lowf, ekb, occ, kvec_c, wk); + EXPECT_EQ(3, nbands); + EXPECT_EQ(63, nbasis); + EXPECT_EQ(nks * nbands * nbasis, lowf.size()); + EXPECT_EQ(nks * nbands, ekb.size()); + EXPECT_EQ(nks * nbands, occ.size()); + EXPECT_EQ(nks, kvec_c.size()); + EXPECT_EQ(nks, wk.size()); + + // test the first k-point + int nbands_k0 = -1; + int nbasis_k0 = -1; + int ik_k0; + std::vector> lowf_k0; + std::vector ekb_k0; + std::vector occ_k0; + ModuleBase::Vector3 kvec_c_k0; + double wk_k0 = -1.0; + ModuleIO::read_abacus_lowf("./support/WFC_NAO_K1.txt", ik_k0, kvec_c_k0, nbands_k0, nbasis_k0, lowf_k0, ekb_k0, + occ_k0, wk_k0); + + EXPECT_EQ(1, ik_k0); + EXPECT_EQ(3, nbands_k0); + EXPECT_EQ(63, nbasis_k0); + EXPECT_EQ(1.0, wk_k0); // this is overwritten with 1.0 + // kvec_c + EXPECT_NEAR(kvec_c_k0.x, kvec_c[0].x, 1e-7); + EXPECT_NEAR(kvec_c_k0.y, kvec_c[0].y, 1e-7); + EXPECT_NEAR(kvec_c_k0.z, kvec_c[0].z, 1e-7); + // ekb + for (int i = 0; i < nbands_k0; i++) + { + EXPECT_NEAR(ekb_k0[i], ekb[i], 1e-7); + } + // occ + for (int i = 0; i < nbands_k0; i++) + { + EXPECT_NEAR(occ_k0[i], occ[i], 1e-7); + } + // lowf + for (int i = 0; i < nbands_k0 * nbasis_k0; i++) + { + EXPECT_NEAR(lowf_k0[i].real(), lowf[i].real(), 1e-7); + EXPECT_NEAR(lowf_k0[i].imag(), lowf[i].imag(), 1e-7); + } + + // test the second k-point + int nbands_k1 = -1; + int nbasis_k1 = -1; + int ik_k1; + std::vector> lowf_k1; + std::vector ekb_k1; + std::vector occ_k1; + ModuleBase::Vector3 kvec_c_k1; + double wk_k1 = -1.0; + ModuleIO::read_abacus_lowf("./support/WFC_NAO_K2.txt", ik_k1, kvec_c_k1, nbands_k1, nbasis_k1, lowf_k1, ekb_k1, + occ_k1, wk_k1); + + EXPECT_EQ(2, ik_k1); + EXPECT_EQ(3, nbands_k1); + EXPECT_EQ(63, nbasis_k1); + EXPECT_EQ(1.0, wk_k1); // this is overwritten with 1.0 + // kvec_c + EXPECT_NEAR(kvec_c_k1.x, kvec_c[1].x, 1e-7); + EXPECT_NEAR(kvec_c_k1.y, kvec_c[1].y, 1e-7); + EXPECT_NEAR(kvec_c_k1.z, kvec_c[1].z, 1e-7); + // ekb + for (int i = 0; i < nbands_k1; i++) + { + EXPECT_NEAR(ekb_k1[i], ekb[i + nbands_k0], 1e-7); + } + // occ + for (int i = 0; i < nbands_k1; i++) + { + EXPECT_NEAR(occ_k1[i], occ[i + nbands_k0], 1e-7); + } + // lowf + for (int i = 0; i < nbands_k1 * nbasis_k1; i++) + { + EXPECT_NEAR(lowf_k1[i].real(), lowf[i + nbands_k0 * nbasis_k0].real(), 1e-7); + EXPECT_NEAR(lowf_k1[i].imag(), lowf[i + nbands_k0 * nbasis_k0].imag(), 1e-7); + } +} + +int main(int argc, char** argv) +{ + ::testing::InitGoogleTest(&argc, argv); + // print cwd + char cwd[1024]; + if (getcwd(cwd, sizeof(cwd)) != nullptr) + { + printf("Current working dir: %s\n", cwd); + } + else + { + perror("getcwd() error"); + return 1; + } +#ifdef __MPI + MPI_Init(&argc, &argv); +#endif + const int status = RUN_ALL_TESTS(); + printf("Unittest read_wfc_lcao_test exits with status %d\n", status); +#ifdef __MPI + MPI_Finalize(); +#endif + return status; +} \ No newline at end of file diff --git a/source/module_io/test/support/WFC_NAO_GAMMA1.txt b/source/module_io/test/support/WFC_NAO_GAMMA1.txt new file mode 100644 index 0000000000..e05181effc --- /dev/null +++ b/source/module_io/test/support/WFC_NAO_GAMMA1.txt @@ -0,0 +1,32 @@ +3 (number of bands) +31 (number of orbitals) +1 (band) +-1.22787155e+00 (Ry) +2.00000000e+00 (Occupations) +-1.51728369e-02 -2.07808444e-03 1.21298954e-17 -7.13941907e-08 2.16402589e-07 +-2.33154682e-16 1.57981033e-07 1.37815571e-07 -4.73101688e-03 -9.50509066e-18 +-2.78717856e-18 3.97439579e-08 3.10940207e-08 5.39213060e-03 1.01734983e-18 +-6.58833740e-18 1.27090046e-07 -7.03036426e-08 -9.77585956e-01 -6.13880745e-03 +2.20316871e-15 -2.33416845e-07 -9.14787918e-07 1.07434159e-16 5.16552542e-08 +1.13539917e-07 7.73272236e-03 4.69976670e-17 1.03418318e-17 4.83739928e-07 +-5.44883791e-09 +2 (band) +-3.10595658e-01 (Ry) +2.00000000e+00 (Occupations) +8.86226706e-02 2.27254080e-02 -4.89898875e-16 4.28239622e-09 5.83161806e-08 +-4.60626214e-16 6.71493689e-09 -1.97997448e-08 -9.80282026e-01 6.43055481e-15 +2.29135711e-14 -5.33364904e-05 1.20925412e-06 -5.32138426e-03 1.33240077e-16 +5.71272071e-16 -5.49131603e-07 9.23311047e-09 1.24546190e-02 -2.68408770e-03 +-3.01869415e-16 -7.75919179e-07 -8.35043239e-07 -1.38124495e-16 -8.88948486e-08 +-3.56018509e-07 1.39706604e-03 1.17174237e-17 4.96617345e-17 -5.76417426e-09 +2.64285801e-08 +3 (band) +-3.00546690e-01 (Ry) +2.00000000e+00 (Occupations) +-8.24520852e-10 -3.26589260e-08 -1.69663394e-16 3.99552530e-04 1.58759611e-07 +2.58344137e-18 -8.99733979e-04 -3.58445734e-07 1.21571018e-06 2.25456884e-13 +-1.04214879e-13 4.17125230e-04 1.00074306e+00 -3.22028333e-09 4.79757646e-15 +-2.17277014e-15 4.39402756e-06 1.05359309e-02 1.35243111e-07 2.55990474e-09 +7.62457489e-16 -2.35220550e-02 -1.00062779e-05 6.59539788e-17 -2.44422254e-03 +-1.05746282e-06 -3.16636112e-08 2.78843989e-16 -1.97537671e-16 -8.29982457e-08 +-2.02338162e-04 diff --git a/source/module_io/test/support/WFC_NAO_K1.txt b/source/module_io/test/support/WFC_NAO_K1.txt new file mode 100644 index 0000000000..b98f72154a --- /dev/null +++ b/source/module_io/test/support/WFC_NAO_K1.txt @@ -0,0 +1,52 @@ +1 (index of k points) +-0.10540388 -0.06085495 -0.04303095 +3 (number of bands) +63 (number of orbitals) +1 (band) +-6.03571945e-01 (Ry) +5.83090379e-03 (Occupations) +-6.71651157e-03 2.25946383e-02 1.43180123e-03 1.46451488e-03 2.31452033e-03 -1.18949691e-03 3.86105126e-03 -5.30361525e-03 -9.45760546e-03 1.29911511e-02 +-5.46035106e-03 7.50044462e-03 1.07440727e-03 -1.52230629e-03 -2.63174959e-03 3.72887365e-03 -1.51944133e-03 2.15286620e-03 -1.22115353e-04 -3.37537766e-04 +2.99120310e-04 8.26795292e-04 1.72697191e-04 4.77350484e-04 1.62647420e-01 -2.84198177e-01 2.30018187e-01 -4.01916916e-01 1.32801062e-01 -2.32046840e-01 +-1.87809063e-01 3.28163787e-01 -3.25294839e-01 5.68396353e-01 -2.52010188e-03 4.07977110e-03 -3.56396226e-03 5.76966762e-03 -2.05765457e-03 3.33111915e-03 +2.90996300e-03 -4.71091388e-03 5.04020376e-03 -8.15954220e-03 -1.29707431e-03 2.57674619e-03 -1.83434009e-03 3.64406943e-03 -1.05905675e-03 2.10390446e-03 +1.49773241e-03 -2.97537023e-03 2.59414864e-03 -5.15349241e-03 -3.97550843e-04 1.44683940e-03 -6.06854639e-04 1.01883141e-03 -3.50367689e-04 5.88222591e-04 +-7.83445970e-04 1.31530570e-03 -1.35696822e-03 2.27817630e-03 -3.03660261e-14 6.17920665e-14 4.58600653e-04 -7.28132783e-05 -5.77273744e-04 -1.85618417e-04 +8.37425341e-04 -1.20858286e-03 4.83487746e-04 -6.97775642e-04 1.08111146e-03 -1.56027376e-03 1.87253999e-03 -2.70247344e-03 5.43592135e-14 1.09205667e-13 +-1.63985332e-03 1.45546190e-03 -6.90799158e-02 7.52891979e-02 -1.23303889e-02 1.18105312e-02 2.35076693e-03 -2.40917316e-03 1.73947004e-03 7.26573136e-04 +-4.26081388e-03 -1.77973354e-03 -2.45998205e-03 -1.02752963e-03 1.00801595e-03 6.43220132e-04 -2.46912478e-03 -1.57556109e-03 -1.42554985e-03 -9.09650620e-04 +8.31609942e-04 -7.66286321e-04 -2.03702004e-03 1.87701049e-03 -1.17607407e-03 1.08369251e-03 4.20973368e-04 3.56330622e-04 5.95346235e-04 5.03927605e-04 +3.43723309e-04 2.90942739e-04 -4.86098172e-04 -4.11455162e-04 -8.41946729e-04 -7.12661247e-04 -9.13243785e-04 1.05899942e-03 -1.29152174e-03 1.49765135e-03 +-7.45660429e-04 8.64669414e-04 1.05452309e-03 -1.22282721e-03 1.82648757e-03 -2.11799886e-03 +2 (band) +-5.98035141e-01 (Ry) +5.83090379e-03 (Occupations) +8.76573314e-12 -1.93189785e-11 -2.88265028e-12 2.17773729e-11 1.30377913e-12 1.08304110e-12 1.86874275e-03 -3.83200170e-03 4.77752128e-04 -9.74112406e-04 +4.93909685e-04 -1.02242215e-03 -8.69663301e-04 9.99617996e-03 -2.27023949e-04 2.54336238e-03 -2.21727797e-04 2.66313376e-03 -1.00166656e-03 -6.52953581e-04 +-2.54534300e-04 -1.66737767e-04 -2.67418877e-04 -1.72909618e-04 -2.09966472e-01 8.66452869e-01 -4.38318350e-02 8.33279277e-02 -2.53687609e-03 -3.77906886e-02 +-6.92320060e-02 2.56156425e-01 -9.70415824e-02 3.28828213e-01 7.80900421e-05 -7.50681575e-03 6.25495844e-04 1.58763114e-03 2.82147356e-04 1.37535219e-03 +2.15884927e-04 -1.47618579e-03 4.71882337e-04 -1.21702170e-03 -2.36044268e-03 -3.70876954e-03 7.57854259e-04 4.35580000e-03 5.47455951e-04 2.30034290e-03 +-3.87325116e-04 4.19702299e-04 -1.97217235e-04 1.92242675e-03 3.07908989e-04 -1.63705170e-03 -3.15643160e-04 6.71944470e-04 -2.15015255e-04 5.26439964e-04 +-3.59293523e-05 4.81481439e-04 1.08688035e-05 3.83755318e-04 -3.10906988e-05 4.97587622e-05 -3.44252715e-04 1.83027944e-03 1.26978712e-04 -1.21953671e-04 +-7.97470891e-04 -7.26654020e-05 -5.03990942e-04 -3.65931298e-05 1.70020992e-04 6.99411230e-05 3.03428152e-04 8.33737407e-05 -9.22641794e-05 -9.81936674e-06 +-1.41966518e-04 1.36348352e-04 -5.23158411e-11 9.36178642e-11 -7.61140957e-12 1.33921232e-11 1.60351645e-12 -2.80310950e-12 -3.68099832e-02 3.39093687e-02 +-9.38688683e-03 8.60833803e-03 -9.77002385e-03 9.06746575e-03 1.26048501e-02 -9.38568640e-03 3.21308226e-03 -2.38129953e-03 3.34775327e-03 -2.51215074e-03 +7.23632039e-03 -7.58482248e-03 1.84585348e-03 -1.92607628e-03 1.91973920e-03 -2.02721742e-03 -3.17187254e-03 2.40033423e-03 2.76651822e-03 -4.17266025e-03 +1.53599165e-03 -2.10378487e-03 5.44005293e-05 -7.00804085e-04 9.65945330e-04 -2.20460639e-03 -1.19968502e-03 2.33511242e-03 2.54847151e-04 4.35625490e-05 +2.17988065e-04 -1.84173316e-04 -2.34349428e-04 6.32902387e-04 -1.95343398e-04 7.57764774e-04 +3 (band) +-5.98035141e-01 (Ry) +5.83090379e-03 (Occupations) +9.86470874e-13 -5.95387122e-12 4.82573453e-13 5.54264959e-12 3.66087151e-13 1.96386245e-13 9.49023673e-05 -4.04693771e-04 5.20920946e-04 -2.33076310e-03 +-8.35155442e-04 3.75083842e-03 1.33229060e-04 9.69176971e-04 8.23664081e-04 5.56014508e-03 -1.33242086e-03 -8.94514217e-03 -1.09337525e-04 -4.04444133e-05 +-6.30599887e-04 -2.25793277e-04 1.01491773e-03 3.62486910e-04 -1.24770866e-03 8.69136907e-02 1.18610223e-02 -3.89305996e-01 -2.57019678e-02 6.86121492e-01 +-1.51258921e-02 -4.98987709e-01 6.00330846e-03 3.36374550e-01 -1.54719250e-04 -7.15401601e-04 2.04082973e-04 2.22148712e-03 -1.34550154e-04 -3.48501679e-03 +1.74768920e-03 5.61539133e-03 -9.97009876e-04 -3.45167204e-03 -3.04528154e-04 -3.01628967e-04 4.33629845e-04 -6.49120521e-04 -3.61575730e-04 2.02013352e-03 +3.23047521e-03 4.80451613e-03 -1.85836988e-03 -2.55898463e-03 -6.08492343e-06 -1.62294338e-04 -2.12474945e-05 7.08757263e-04 8.89565597e-07 -1.05044613e-03 +-1.20576279e-04 -1.99432825e-03 8.29691217e-05 1.21455445e-03 1.97898051e-04 -5.68319065e-04 6.80315141e-06 1.81450586e-04 9.43843306e-06 -1.43374981e-05 +-2.62559569e-04 7.70285858e-05 2.72016255e-04 -1.08142226e-04 1.92623971e-05 -1.71028330e-04 2.97332706e-05 1.01835109e-04 9.45113900e-04 -1.11433987e-04 +-1.05524891e-05 1.60298109e-05 -7.10509754e-12 2.44531935e-11 -1.03733385e-12 3.40512767e-12 2.18792113e-13 -7.24424287e-13 -2.76730418e-03 4.01903148e-03 +-1.56573360e-02 2.32565523e-02 2.51625219e-02 -3.74396459e-02 9.95687062e-04 -1.16460848e-03 5.65048270e-03 -6.75077464e-03 -9.08286605e-03 1.08691821e-02 +5.24168776e-04 -8.77432567e-04 2.95875368e-03 -5.07255076e-03 -4.75406840e-03 8.16547713e-03 -2.49721541e-04 2.96724144e-04 2.08318252e-04 9.39165110e-04 +3.92497623e-05 -2.65051245e-03 2.81714661e-03 -5.26537445e-03 -1.58801422e-03 2.77035026e-03 -6.36234192e-05 2.47924813e-04 1.27026963e-04 -9.86912246e-04 +-1.59672985e-04 1.68480281e-03 5.75804420e-04 -1.62364450e-03 -3.39617142e-04 1.05133952e-03 diff --git a/source/module_io/test/support/WFC_NAO_K2.txt b/source/module_io/test/support/WFC_NAO_K2.txt new file mode 100644 index 0000000000..c71b24b3ce --- /dev/null +++ b/source/module_io/test/support/WFC_NAO_K2.txt @@ -0,0 +1,52 @@ +2 (index of k points) +-0.07026925 -0.08113994 -0.05737460 +3 (number of bands) +63 (number of orbitals) +1 (band) +-6.03667277e-01 (Ry) +5.83090379e-03 (Occupations) +2.09933705e-02 2.20619371e-03 1.52454416e-03 -3.54139105e-04 1.31198443e-03 -2.44872538e-03 -9.20116330e-03 -5.45151213e-03 3.71290025e-03 1.45501785e-02 +1.30124099e-02 7.70960240e-03 3.82663106e-03 -3.84610111e-03 7.78517834e-03 -4.33511889e-03 -5.41167355e-03 5.43920833e-03 -9.15112399e-04 7.19827990e-05 +7.31599725e-04 1.31901719e-04 1.29416436e-03 -1.01799053e-04 -2.02504709e-01 -3.96665125e-01 -3.27601131e-01 -1.26268707e-01 -1.53433434e-01 -3.42508212e-01 +2.42254622e-01 4.44854270e-01 4.63297961e-01 1.78570918e-01 3.22837653e-03 5.12766840e-03 5.09299111e-03 3.04099283e-03 1.43817942e-03 4.66552064e-03 +-4.57476575e-03 -5.58236091e-03 -7.20257710e-03 -4.30061330e-03 2.71096277e-03 2.84147597e-03 4.81279492e-03 2.93105118e-03 -1.35607565e-04 2.77763433e-03 +-4.79141429e-03 -2.95749668e-03 -6.80631985e-03 -4.14513233e-03 1.38144609e-03 7.67113833e-04 8.08627314e-04 2.65802677e-04 4.36096803e-04 8.21540265e-04 +9.75142095e-04 1.83701988e-03 1.96060297e-03 6.69766736e-04 8.80212792e-05 4.35402244e-05 2.22205557e-04 -1.22295959e-03 -8.06493830e-04 3.86102824e-04 +-1.02999787e-03 -1.08410238e-03 -3.78284202e-04 -2.78671813e-04 -8.45869195e-04 -6.23129113e-04 -1.34763393e-03 -2.51513148e-03 5.51664715e-04 -5.25416787e-05 +1.57703040e-04 9.92915635e-04 5.05190321e-02 8.77355665e-02 7.71248769e-03 1.53180277e-02 -1.60787823e-03 -2.90751405e-03 4.74616908e-03 -6.51680651e-03 +-2.35857817e-03 1.35310810e-03 -6.71209668e-03 9.21615608e-03 2.60480706e-03 1.36700642e-03 -9.24054509e-04 2.36186797e-03 -3.68375348e-03 -1.93323900e-03 +-1.21442969e-04 5.99547489e-04 8.46027889e-04 2.72325634e-03 1.71746294e-04 -8.47888176e-04 1.43198828e-03 -9.83925767e-04 3.37441520e-04 -1.33906274e-04 +2.38988110e-03 7.11859012e-04 -7.90375335e-04 2.20756974e-03 -4.77214374e-04 1.89372078e-04 4.69988377e-04 7.54715897e-04 8.14293490e-04 1.27236821e-03 +-2.89004839e-05 6.29744824e-04 -8.34479479e-04 -8.61909448e-04 -1.15158489e-03 -1.79940038e-03 +2 (band) +-5.97868276e-01 (Ry) +5.83090379e-03 (Occupations) +1.44618640e-03 6.69577279e-04 1.89050182e-05 -1.73461277e-03 -1.04700276e-03 1.76062866e-03 1.85086491e-04 2.15202728e-03 2.80834863e-03 6.33731532e-03 +-2.61751869e-04 -3.04342630e-03 8.80620073e-04 -1.71701914e-04 6.96444765e-03 -6.74454643e-03 -1.24538462e-03 2.42823264e-04 5.06462111e-04 -2.72179689e-05 +-1.25714000e-03 8.17226157e-05 -7.16245627e-04 3.84920674e-05 3.14474976e-01 1.37786756e-01 -4.15349196e-01 1.42172891e-01 4.39591320e-01 2.05260088e-01 +-2.33848614e-01 -9.35128501e-02 5.87392468e-01 -2.01062827e-01 -2.47976289e-03 -3.86608243e-04 5.27243491e-03 -1.53432746e-03 -1.29394265e-03 -1.12999726e-03 +3.38011947e-03 -1.29403657e-04 -7.45634891e-03 2.16986668e-03 -1.08862164e-03 1.25888021e-03 5.36965277e-03 -1.30232731e-03 2.81236999e-03 6.22446254e-04 +3.87419368e-03 -1.74030849e-03 -7.59383564e-03 1.84176891e-03 -6.06063290e-04 -7.57571331e-05 1.04374285e-03 1.48830832e-04 -4.48554991e-04 -1.57060775e-04 +-1.00299945e-03 -3.51198586e-04 1.67129220e-03 -7.42852690e-04 -3.82545213e-04 -6.21025977e-04 5.00268162e-04 2.97217060e-04 2.96132294e-04 -2.42303693e-05 +-4.05391638e-04 1.77808661e-03 -2.33664140e-04 1.15135290e-04 -5.22488875e-04 2.57450274e-04 -1.16848478e-03 2.86716355e-04 -1.51266629e-04 -2.12996420e-03 +8.07855122e-04 -2.89222470e-04 2.09149008e-03 -3.72294710e-03 2.56202982e-03 1.96574791e-03 -1.44710850e-04 -1.20620388e-05 9.08302740e-03 1.33060713e-02 +7.53227638e-03 2.52461589e-02 -1.28453403e-02 -1.88176257e-02 -2.06589575e-03 -2.02602033e-03 6.70741575e-04 -1.53334762e-03 2.92161774e-03 2.86522525e-03 +-2.46868527e-03 -1.96217413e-03 9.59807079e-04 -4.45737658e-03 3.49124811e-03 2.77493313e-03 3.88693406e-04 1.48421465e-03 -4.06278910e-04 -1.83236707e-03 +-3.43680473e-03 -1.63094051e-03 -3.10342464e-03 -3.72398417e-03 5.74565117e-04 2.59135828e-03 1.09271660e-03 5.02722659e-04 6.08670622e-04 -1.52792331e-03 +1.32844181e-03 7.68909532e-04 -9.53290392e-04 -3.27039998e-04 -8.60790236e-04 2.16080987e-03 +3 (band) +-5.97662421e-01 (Ry) +5.83090379e-03 (Occupations) +4.19321636e-11 -2.80709607e-11 -3.33477115e-11 1.39914678e-11 1.32422571e-11 1.82935858e-11 6.72914392e-03 5.77236082e-04 5.84332768e-11 5.77396615e-12 +4.75822323e-03 4.08167546e-04 -7.62944534e-03 6.55351180e-03 -3.84311324e-11 -8.91586160e-11 -5.39483254e-03 4.63403264e-03 -1.16146065e-03 -1.99779399e-03 +-2.79410197e-12 1.19321317e-11 -8.21276705e-04 -1.41265366e-03 -6.53297129e-01 4.40450901e-01 -1.66380291e-01 1.20651972e-01 2.66707442e-01 -1.79813333e-01 +-3.77181276e-01 2.54294450e-01 -1.17648632e-01 8.53138186e-02 4.59501266e-03 -5.51617982e-03 6.77903797e-04 1.18165414e-03 -1.87590609e-03 2.25197100e-03 +2.65293181e-03 -3.18476795e-03 4.79350374e-04 8.35555768e-04 1.72516240e-03 -4.94313228e-03 -2.69530229e-04 4.32736357e-03 -7.04294587e-04 2.01802528e-03 +9.96022961e-04 -2.85391880e-03 -1.90586644e-04 3.05990822e-03 1.19138761e-03 -8.57963623e-04 -2.13194731e-04 3.98012968e-04 -1.03818250e-03 6.32572483e-05 +-9.19978000e-04 9.68574164e-04 9.53435833e-05 -1.77996832e-04 -1.65139940e-04 3.08299528e-04 -1.33201184e-03 9.59232495e-04 1.10635572e-04 6.23553988e-04 +5.15358333e-04 9.93316890e-04 -1.68295594e-03 -1.91171950e-03 6.24093859e-04 1.30442693e-04 -2.30475243e-04 -4.44224803e-04 3.99194817e-04 7.69419945e-04 +-1.23694333e-04 -6.97154565e-04 -1.07412252e-10 1.17874694e-11 1.65693418e-11 -1.29863251e-11 1.76852567e-12 -1.05920114e-14 -3.72472858e-02 -4.31938728e-03 +2.76422976e-10 1.96825832e-11 -2.63378088e-02 -3.05426793e-03 8.19050673e-03 8.13360389e-04 -5.55872600e-12 -8.50002533e-12 5.79156291e-03 5.75132596e-04 +6.68158686e-03 -1.05490337e-03 -3.93492393e-11 -2.37407840e-11 4.72459544e-03 -7.45929370e-04 -3.96307599e-03 -7.50818021e-04 3.98908133e-03 2.41422790e-05 +1.61791897e-03 3.06520202e-04 -2.28808305e-03 -4.33484978e-04 2.82070650e-03 1.70711749e-05 -2.22472624e-03 1.43644798e-03 -2.36774687e-04 4.74979602e-04 +9.08240698e-04 -5.86427447e-04 -1.28444630e-03 8.29333643e-04 -1.67424958e-04 3.35861321e-04 diff --git a/source/module_io/test/support/WFC_NAO_K3.txt b/source/module_io/test/support/WFC_NAO_K3.txt new file mode 100644 index 0000000000..e0d3852ac9 --- /dev/null +++ b/source/module_io/test/support/WFC_NAO_K3.txt @@ -0,0 +1,52 @@ +3 (index of k points) +-0.03513462 -0.10142493 -0.07171825 +3 (number of bands) +63 (number of orbitals) +1 (band) +-6.04664544e-01 (Ry) +5.83090379e-03 (Occupations) +-6.44410072e-03 2.86108252e-03 -5.81398415e-03 4.01495705e-03 -2.32158666e-03 2.62541166e-03 9.49269256e-03 6.65452048e-03 -5.02852726e-03 -6.72710070e-03 +-1.34246945e-02 -9.41091312e-03 -3.62201275e-04 3.84255737e-03 -6.19560155e-03 4.61326069e-03 5.12229964e-04 -5.43419674e-03 4.49469845e-04 -1.76093834e-04 +1.29087581e-03 4.47455297e-04 -6.35646353e-04 2.49034291e-04 3.47317658e-01 4.16938730e-01 1.57045858e-01 -3.77513068e-02 2.76511095e-01 3.55744139e-01 +-4.06048959e-01 -4.70609970e-01 -2.22096383e-01 5.33884100e-02 -4.07089162e-03 -4.74527573e-03 -3.87146974e-03 -4.12994707e-04 -2.13913520e-03 -4.44234580e-03 +5.53839409e-03 5.07784581e-03 5.47508501e-03 5.84062717e-04 -1.29777024e-03 -1.32253013e-03 -5.21386357e-03 -1.21016710e-03 1.26707763e-03 -1.53986249e-03 +3.14376317e-03 1.20184216e-03 7.37351657e-03 1.71143473e-03 -1.58981110e-03 -3.90987170e-04 -4.67410930e-04 2.49611504e-04 -7.05992818e-04 -8.28533010e-04 +-1.57864793e-03 -1.85265613e-03 -1.14424880e-03 3.84458305e-04 -5.72074381e-05 -1.00279960e-04 2.18608779e-04 1.57562717e-03 7.24875131e-04 -1.03557005e-03 +4.18740850e-04 1.00903940e-03 -1.55414326e-04 -2.16241633e-04 -3.47516996e-04 -4.83530997e-04 1.26652299e-03 1.29383203e-03 1.90635277e-04 -5.55669996e-04 +1.00949827e-03 -4.23741865e-04 -7.12585524e-02 -8.94935779e-02 -8.44214443e-03 -1.10407378e-02 2.24710652e-03 2.83840975e-03 -3.49397206e-03 5.14409284e-03 +-4.12778069e-03 -1.48295352e-03 4.94122267e-03 -7.27484581e-03 -4.40157673e-03 1.95700899e-03 3.04579120e-03 1.67939388e-03 6.22476951e-03 -2.76762867e-03 +-1.23118445e-03 5.67268289e-04 2.96994744e-04 3.21722084e-04 1.74115774e-03 -8.02238517e-04 -1.75403842e-03 1.64296557e-03 3.11170633e-04 -2.34997893e-04 +-4.42119459e-03 -2.19088074e-03 -8.81730111e-05 -4.39488646e-03 -4.40061729e-04 3.32337204e-04 2.98828807e-05 -3.88931555e-04 -3.95247872e-04 -3.69015263e-04 +4.07002218e-04 -8.70073639e-05 2.36035364e-04 6.12125722e-04 5.58964902e-04 5.21866389e-04 +2 (band) +-5.97025474e-01 (Ry) +5.83090379e-03 (Occupations) +5.69317183e-03 -3.71257283e-03 -2.47903685e-03 1.40811142e-03 -8.71555913e-04 1.11914525e-03 -1.44082193e-03 2.34277805e-03 5.05364032e-03 7.86613185e-03 +2.03762998e-03 -3.31318853e-03 2.24244925e-04 -1.03355655e-04 8.56159882e-03 -2.66562862e-03 -3.17130156e-04 1.46167038e-04 4.64810550e-04 -3.15981689e-04 +-1.45927477e-03 -2.05483986e-03 -6.57341419e-04 4.46865603e-04 8.37339722e-02 1.36725440e-01 -5.39366321e-01 5.71646184e-02 1.53973046e-01 2.50930312e-01 +-3.61561036e-02 -5.93808734e-02 7.62779164e-01 -8.08429774e-02 -4.68269214e-04 -1.03002855e-03 8.56052156e-03 1.14218437e-03 6.06313504e-04 -8.49988917e-04 +1.23979433e-03 1.18302875e-03 -1.21064056e-02 -1.61529263e-03 8.00832909e-05 -1.28365698e-04 1.15750503e-02 2.60192587e-03 2.98815079e-03 1.36987051e-03 +1.97423316e-03 1.19098051e-03 -1.63695930e-02 -3.67967885e-03 -1.74837948e-04 -1.93194996e-04 8.50833039e-04 2.46720545e-04 -6.67450729e-05 -1.24534531e-04 +-1.49246516e-04 -2.78467686e-04 2.27866920e-03 -4.65834100e-04 2.17169552e-04 -5.87464296e-04 -1.27827510e-06 1.16593270e-04 -1.15639780e-04 -1.36250437e-04 +-1.96189386e-03 9.38965027e-04 -1.49110958e-04 1.34145635e-04 -3.33422216e-04 2.99958731e-04 -3.03223768e-03 -1.46694894e-03 7.82130863e-04 -2.05914199e-03 +2.43071134e-04 -4.33592369e-04 1.71608557e-03 -1.44786774e-03 1.96406150e-03 3.69270159e-03 -1.08625036e-04 -9.16967732e-05 2.79808997e-03 6.60113769e-03 +-1.26577656e-02 1.69126755e-02 -3.95709690e-03 -9.33541820e-03 -1.74615421e-03 4.90985551e-04 8.19050687e-04 4.32530159e-03 2.46943495e-03 -6.94358459e-04 +-7.26786504e-04 -7.05484097e-04 5.87943435e-03 -4.88996984e-04 1.02783130e-03 9.97705112e-04 -3.56069910e-04 2.37424423e-04 9.58906987e-04 -1.81008349e-03 +-6.94938255e-04 -2.94861238e-03 1.25335599e-04 -2.49621492e-03 -1.35609925e-03 2.55984462e-03 3.78815412e-04 6.92269008e-04 2.44737426e-03 -6.92942400e-04 +2.13451107e-04 1.21545878e-03 -5.05194795e-04 -3.39585899e-04 -3.46110988e-03 9.79968541e-04 +3 (band) +-5.96870018e-01 (Ry) +5.83090379e-03 (Occupations) +3.61187971e-11 -3.71005110e-11 -2.89860645e-11 2.49236772e-11 1.50409985e-12 8.03323726e-12 3.32264869e-04 5.49443083e-03 2.96634123e-11 -1.14208816e-12 +2.34946713e-04 3.88514926e-03 -8.41271250e-03 5.91345914e-04 1.57108809e-11 -4.35316082e-11 -5.94868607e-03 4.18144692e-04 4.13568202e-04 -2.68426321e-03 +-1.25562752e-11 2.24235533e-12 2.92436882e-04 -1.89806070e-03 -8.08810252e-01 -2.79166381e-02 -1.05939414e-01 1.84003766e-02 3.30195407e-01 1.13969198e-02 +-4.66966817e-01 -1.61176784e-02 -7.49104766e-02 1.30110262e-02 1.09154414e-02 -2.86161677e-03 4.08886067e-05 1.49490442e-04 -4.45621032e-03 1.16825014e-03 +6.30203305e-03 -1.65215521e-03 2.89125706e-05 1.05705774e-04 1.36865755e-02 -5.28654025e-03 -8.96895600e-04 7.71771720e-04 -5.58752104e-03 2.15822100e-03 +7.90194806e-03 -3.05218544e-03 -6.34201018e-04 5.45725106e-04 1.42665481e-03 1.22665434e-04 -1.76230409e-04 1.97166175e-04 -1.07604605e-03 -6.14653363e-04 +-1.17640068e-03 1.32356987e-04 7.88126388e-05 -8.81754051e-05 -1.36507492e-04 1.52724263e-04 -1.59504857e-03 -1.37144122e-04 -1.55400391e-03 1.23842988e-03 +-1.18363694e-04 1.22314552e-03 -5.78094164e-04 -2.95838331e-03 2.06412095e-03 -1.15896247e-04 5.29338428e-05 -5.47007290e-04 -9.16841310e-05 9.47444435e-04 +1.73742918e-03 -1.38460670e-03 -9.48817639e-11 -3.30089486e-11 1.20651987e-11 7.29812260e-14 2.09180714e-12 7.72588637e-13 -1.70859295e-02 -1.77897014e-02 +4.10379087e-11 9.07776228e-11 -1.20815768e-02 -1.25792185e-02 2.43021327e-03 -4.14436032e-04 2.30488701e-11 1.29297651e-11 1.71842030e-03 -2.93050566e-04 +5.52117876e-03 9.33825523e-04 1.11334233e-11 -2.44447961e-11 3.90406296e-03 6.60314350e-04 -2.58931955e-03 -3.06573021e-03 7.63384051e-04 1.77859126e-03 +1.05708525e-03 1.25157910e-03 -1.49494434e-03 -1.77000018e-03 5.39794050e-04 1.25765395e-03 -3.23067805e-03 -4.32275378e-04 -4.95840759e-05 -7.44719598e-06 +1.31891879e-03 1.76475687e-04 -1.86523284e-03 -2.49574303e-04 -3.50612416e-05 -5.26593642e-06 diff --git a/source/module_io/test/support/WFC_NAO_K4.txt b/source/module_io/test/support/WFC_NAO_K4.txt new file mode 100644 index 0000000000..851512e24e --- /dev/null +++ b/source/module_io/test/support/WFC_NAO_K4.txt @@ -0,0 +1,52 @@ +4 (index of k points) +0.00000000 -0.12170991 -0.08606190 +3 (number of bands) +63 (number of orbitals) +1 (band) +-6.05615293e-01 (Ry) +5.83090379e-03 (Occupations) +2.63852585e-03 -1.25137325e-03 4.00871038e-03 -1.90121057e-03 1.94915250e-03 -9.24424306e-04 -5.35609014e-03 -1.12933384e-02 -2.85841773e-14 -9.28435922e-13 +7.57465531e-03 1.59711923e-02 -2.39088484e-03 -5.04119068e-03 -1.09595833e-13 -7.14903643e-13 3.38122176e-03 7.12932023e-03 4.02811230e-04 8.49329163e-04 +-4.08784780e-13 -2.32999263e-13 -5.69661103e-04 -1.20113282e-03 -2.42221221e-01 -5.10724456e-01 3.02860823e-11 6.92011964e-11 -1.91360244e-01 -4.20046054e-01 +2.84227334e-01 5.87583293e-01 -5.50191308e-11 -1.23608424e-10 2.69961427e-03 5.69214797e-03 -2.92361513e-13 -1.05855006e-12 1.16722140e-03 5.13943930e-03 +-3.85051890e-03 -6.22495709e-03 5.05196363e-13 1.67670739e-12 2.47574018e-05 5.22010776e-05 -2.44317009e-13 -1.26359525e-12 -2.20134516e-03 1.09624071e-03 +-1.59946715e-03 6.84744353e-04 4.01861793e-13 1.83366739e-12 1.02117550e-03 1.11860189e-03 -8.98144970e-14 -2.34461562e-13 5.37470485e-04 1.13325876e-03 +1.20182054e-03 2.53404362e-03 -2.19801713e-13 -3.82428306e-13 -1.02794979e-14 8.88345626e-14 -3.35601409e-04 -1.63294744e-03 -2.90360994e-04 9.75795130e-04 +1.21310342e-13 2.93409165e-13 2.81017023e-04 5.92525565e-04 6.28373167e-04 1.32492744e-03 -5.23358291e-14 5.24867534e-13 -2.84058091e-13 -8.05649475e-14 +-9.12731316e-04 -5.04127291e-04 2.77342842e-02 1.24931688e-01 2.31774737e-03 1.04405107e-02 -8.37419724e-04 -3.77223580e-03 1.44815776e-03 -3.21484659e-04 +2.59846558e-12 4.20304407e-12 -2.04800434e-03 4.54647940e-04 3.88337687e-03 -8.62092558e-04 -7.75367687e-14 -4.47249360e-13 -5.49192423e-03 1.21918299e-03 +1.07568245e-03 -2.38796764e-04 -2.67659823e-13 -8.90101405e-13 -1.52124472e-03 3.37709630e-04 1.20812666e-03 -2.68198794e-04 -3.35596973e-13 -4.60749731e-13 +2.26383042e-03 5.53518100e-03 -4.91766914e-04 4.37849795e-03 3.81998995e-13 6.66053184e-13 -8.47845044e-05 1.88217835e-05 -6.23295910e-14 -2.97589100e-13 +-1.64066226e-04 -4.11847625e-04 3.08387267e-05 -3.23820540e-04 1.14421867e-13 4.61245127e-13 +2 (band) +-5.96302906e-01 (Ry) +5.83090379e-03 (Occupations) +-2.78749601e-11 -2.03123039e-11 -4.15879440e-12 2.82749691e-12 -2.52364188e-12 1.85593727e-13 -6.86038342e-03 9.67263043e-05 -1.26665013e-04 -3.49279523e-04 +-4.85102364e-03 6.83958012e-05 -5.40427757e-03 -3.32378489e-03 -2.71962921e-04 -2.09943115e-04 -3.82140130e-03 -2.35027086e-03 1.82137797e-03 -8.76009467e-04 +-9.43394234e-06 1.09037970e-04 1.28790870e-03 -6.19432215e-04 -6.72897886e-01 -4.63959847e-01 2.10157059e-02 1.45373980e-02 2.74709415e-01 1.89410817e-01 +-3.88497773e-01 -2.67867341e-01 -2.97206969e-02 -2.05589858e-02 1.30407965e-02 6.21963679e-03 -3.26239096e-04 -3.12426357e-04 -5.32388293e-03 -2.53915610e-03 +7.52910741e-03 3.59090896e-03 4.61371761e-04 4.41837595e-04 2.00763204e-02 8.67557720e-03 -4.75943861e-04 -4.90940649e-04 -8.19612352e-03 -3.54178956e-03 +1.15910690e-02 5.00884682e-03 6.73086276e-04 6.94294933e-04 1.07782595e-03 8.55439316e-04 -2.24842293e-05 -2.08550707e-05 -7.90900486e-04 -9.47277722e-04 +-8.98619146e-04 -5.70295191e-04 -9.50927705e-05 -5.34068951e-05 -2.58748217e-05 -3.91070478e-06 -1.20504605e-03 -9.56410231e-04 -2.37008719e-03 -5.46387035e-04 +9.88056776e-05 6.07591597e-05 5.08469415e-04 -1.83046120e-03 2.52639803e-03 1.45345150e-03 7.59645012e-05 1.52833576e-04 -8.36994573e-05 9.79877086e-06 +2.64983803e-03 6.10879274e-04 -3.32930237e-11 -8.84534600e-11 -9.42799451e-12 -2.20257893e-12 1.24388548e-12 2.43183256e-12 -5.46265889e-03 -1.71890142e-02 +9.01662490e-04 3.75394865e-04 -3.86268324e-03 -1.21544686e-02 1.96040087e-03 4.69284855e-04 -1.61133855e-05 -1.07962277e-04 1.38621274e-03 3.31834481e-04 +3.96807153e-03 4.55075702e-03 -2.26874084e-04 -2.35432717e-04 2.80585028e-03 3.21787114e-03 -5.59009041e-04 -3.01411830e-03 -9.23647312e-05 -2.55809919e-05 +2.28214476e-04 1.23050861e-03 -3.22743991e-04 -1.74020200e-03 1.30623454e-04 3.61769873e-05 -2.42140006e-03 -2.91318788e-03 -8.41731141e-05 -8.33148030e-05 +9.88532434e-04 1.18930397e-03 -1.39799597e-03 -1.68192980e-03 1.19038759e-04 1.17824924e-04 +3 (band) +-5.96302906e-01 (Ry) +5.83090379e-03 (Occupations) +-6.10460263e-12 -3.29582970e-12 8.41831460e-13 -8.77218518e-13 -1.86944338e-13 -4.76566973e-13 -3.03329059e-04 4.27672185e-06 2.86477869e-03 7.89964354e-03 +-2.14486030e-04 3.02409799e-06 -2.38947936e-04 -1.46959786e-04 6.15097666e-03 4.74827645e-03 -1.68961700e-04 -1.03916266e-04 8.05314849e-05 -3.87324048e-05 +2.13367106e-04 -2.46610824e-03 5.69443559e-05 -2.73879427e-05 -2.97519063e-02 -2.05137948e-02 -4.75311553e-01 -3.28791852e-01 1.21461656e-02 8.37472217e-03 +-1.71772705e-02 -1.18436445e-02 6.72192042e-01 4.64981894e-01 5.76593513e-04 2.74998692e-04 7.37853930e-03 7.06613663e-03 -2.35393320e-04 -1.12267747e-04 +3.32896416e-04 1.58770566e-04 -1.04348303e-02 -9.99302624e-03 8.87666337e-04 3.83587088e-04 1.07644072e-02 1.11035885e-02 -3.62388267e-04 -1.56598773e-04 +5.12494397e-04 2.21464107e-04 -1.52231707e-02 -1.57028454e-02 4.76556364e-05 3.78228996e-05 5.08525091e-04 4.71678443e-04 -3.49693436e-05 -4.18834974e-05 +-3.97320797e-05 -2.52153681e-05 2.15071013e-03 1.20790195e-03 5.85210002e-04 8.84482681e-05 -5.32806212e-05 -4.22872874e-05 -1.04792440e-04 -2.41582772e-05 +-2.23468484e-03 -1.37418786e-03 2.24817659e-05 -8.09330964e-05 1.11703660e-04 6.42637636e-05 -1.71808673e-03 -3.45663200e-03 1.89302788e-03 -2.21618511e-04 +1.17161509e-04 2.70097748e-05 1.62379559e-11 6.72817010e-12 -2.24410740e-14 1.75649499e-12 -4.40615166e-13 -2.53963354e-13 -2.41529248e-04 -7.60005243e-04 +-2.03928714e-02 -8.49029295e-03 -1.70786990e-04 -5.37404884e-04 8.66783288e-05 2.07492375e-05 3.64435982e-04 2.44177930e-03 6.12908288e-05 1.46719237e-05 +1.75446671e-04 2.01209863e-04 5.13120392e-03 5.32477398e-03 1.24059528e-04 1.42276857e-04 -2.47163623e-05 -1.33268014e-04 2.08901012e-03 5.78564431e-04 +1.00904124e-05 5.44064339e-05 -1.42699917e-05 -7.69423198e-05 -2.95430645e-03 -8.18213658e-04 -1.07061220e-04 -1.28805412e-04 1.90374063e-03 1.88432814e-03 +4.37075593e-05 5.25845893e-05 -6.18118243e-05 -7.43658388e-05 -2.69229582e-03 -2.66484241e-03