diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index 2278ddb87c..006853815e 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -1698,21 +1698,21 @@ These variables are used to control the output of properties. ### out_dmk -- **Type**: Boolean +- **Type**: Boolean \[Integer\](optional) - **Availability**: Numerical atomic orbital basis - **Description**: Whether to output the density matrix for each k-point into files in the folder `OUT.${suffix}`. The files are named as: - For gamma only case: - - nspin = 1: `dms1_nao.csr`; + - nspin = 1 and 4: `dm_nao.csr`; - nspin = 2: `dms1_nao.csr` and `dms2_nao.csr` for the two spin channels. - For multi-k points case: - - nspin = 1: `dms1k1_nao.csr`, `dms1k2_nao.csr`, ...; - - nspin = 2: `dms1k1_nao.csr`... and `dms2k1_nao.csr`... for the two spin channels. + - nspin = 1 and 4: `dmk1_nao.csr`, `dmk2_nao.csr`, ...; + - nspin = 2: `dmk1s1_nao.csr`... and `dmk1s2_nao.csr`... for the two spin channels. - **Default**: False - **Note**: In the 3.10-LTS version, the parameter is named `out_dm` and the file names are SPIN1_DM and SPIN2_DM, etc. ### out_dmr -- **Type**: Boolean +- **Type**: Boolean \[Integer\](optional) - **Availability**: Numerical atomic orbital basis (multi-k points) - **Description**: Whether to output the density matrix with Bravias lattice vector R index into files in the folder `OUT.${suffix}`. The files are named as `dmr{s}{spin index}{g}{geometry index}{_nao} + {".csr"}`. Here, 's' refers to spin, where s1 means spin up channel while s2 means spin down channel, and the sparse matrix format 'csr' is mentioned in [out_mat_hs2](#out_mat_hs2). Finally, if [out_app_flag](#out_app_flag) is set to false, the file name contains the optional 'g' index for each ionic step that may have different geometries, and if [out_app_flag](#out_app_flag) is set to true, the density matrix with respect to Bravias lattice vector R accumulates during ionic steps: - nspin = 1: `dmrs1_nao.csr`; diff --git a/examples/relax/lcao_output/INPUT b/examples/relax/lcao_output/INPUT index 5230bd6365..87261f91cb 100644 --- a/examples/relax/lcao_output/INPUT +++ b/examples/relax/lcao_output/INPUT @@ -5,10 +5,10 @@ suffix autotest nbands 8 calculation relax #cell-relax ecutwfc 10 -scf_nmax 2 +scf_nmax 20 basis_type lcao -relax_nmax 5 +relax_nmax 1 cal_stress 1 stress_thr 1e-6 @@ -25,16 +25,19 @@ relax_method bfgs pseudo_dir ../../../tests/PP_ORB orbital_dir ../../../tests/PP_ORB -nspin 4 +nspin 1 -out_chg 1 -out_pot 1 -out_dmk 1 -out_dmr 1 +out_mat_hs 0 // note +out_mat_tk 1 // note +out_mat_l 1 // note +out_chg 0 8 // note +out_pot 0 +out_dmk 0 // note +out_dmr 0 // note #out_wfc_lcao 1 // OK -out_dos 1 -out_band 1 -out_stru 1 +out_dos 0 +out_band 0 +out_stru 0 out_app_flag 0 out_interval 1 diff --git a/source/Makefile.Objects b/source/Makefile.Objects index 4ad1cbff0a..50865b5bdb 100644 --- a/source/Makefile.Objects +++ b/source/Makefile.Objects @@ -605,7 +605,7 @@ OBJS_IO_LCAO=cal_r_overlap_R.o\ write_eig_occ.o\ get_pchg_lcao.o\ get_wf_lcao.o\ - io_dmk.o\ + write_dmk.o\ unk_overlap_lcao.o\ read_wfc_nao.o\ read_wfc_lcao.o\ diff --git a/source/source_base/tool_quit.cpp b/source/source_base/tool_quit.cpp index 8798deb0e0..77c855c7a0 100644 --- a/source/source_base/tool_quit.cpp +++ b/source/source_base/tool_quit.cpp @@ -91,21 +91,21 @@ void WARNING_QUIT(const std::string &file,const std::string &description,int ret std::cout << " https://github.com/deepmodeling/abacus-develop/issues" << std::endl; #else - std::cout << " " << std::endl; - std::cout << " ---------------------------------------------------------" << std::endl; - std::cout << " !NOTICE! " << std::endl; - std::cout << " ---------------------------------------------------------" << std::endl; - std::cout << " " << std::endl; - std::cout << " " << description << std::endl; - std::cout << " CHECK IN FILE : " << PARAM.globalv.global_out_dir << "warning.log" << std::endl; - std::cout << " " << std::endl; - std::cout << " For detailed manual of ABACUS, please see the website" << std::endl; - std::cout << " https://abacus.deepmodeling.com" << std::endl; - std::cout << " For any questions, propose issues on the website" << std::endl; - std::cout << " https://github.com/deepmodeling/abacus-develop/issues" << std::endl; - std::cout << " ---------------------------------------------------------" << std::endl; - std::cout << " !NOTICE! " << std::endl; - std::cout << " ---------------------------------------------------------" << std::endl; + std::cout << " " << std::endl; + std::cout << " ---------------------------------------------------------" << std::endl; + std::cout << " !NOTICE! " << std::endl; + std::cout << " ---------------------------------------------------------" << std::endl; + std::cout << " " << std::endl; + std::cout << " " << description << std::endl; + std::cout << " CHECK IN FILE : " << PARAM.globalv.global_out_dir << "warning.log" << std::endl; + std::cout << " " << std::endl; + std::cout << " For detailed manual of ABACUS, please see the website" << std::endl; + std::cout << " https://abacus.deepmodeling.com" << std::endl; + std::cout << " For any questions, propose issues on the website" << std::endl; + std::cout << " https://github.com/deepmodeling/abacus-develop/issues" << std::endl; + std::cout << " ---------------------------------------------------------" << std::endl; + std::cout << " !NOTICE! " << std::endl; + std::cout << " ---------------------------------------------------------" << std::endl; GlobalV::ofs_running << " ---------------------------------------------------------" << std::endl; diff --git a/source/source_esolver/esolver_ks_lcao.cpp b/source/source_esolver/esolver_ks_lcao.cpp index 84844ab720..c37448c5e7 100644 --- a/source/source_esolver/esolver_ks_lcao.cpp +++ b/source/source_esolver/esolver_ks_lcao.cpp @@ -13,8 +13,7 @@ #include "source_io/berryphase.h" #include "source_io/cal_ldos.h" #include "source_io/cube_io.h" -#include "source_io/io_dmk.h" -#include "source_io/io_npz.h" +//#include "source_io/io_npz.h" #include "source_io/output_dmk.h" #include "source_io/output_log.h" #include "source_io/output_mat_sparse.h" diff --git a/source/source_hamilt/module_xc/xc_functional_libxc_tools.cpp b/source/source_hamilt/module_xc/xc_functional_libxc_tools.cpp index d7c8439d10..4d4af83a7e 100644 --- a/source/source_hamilt/module_xc/xc_functional_libxc_tools.cpp +++ b/source/source_hamilt/module_xc/xc_functional_libxc_tools.cpp @@ -16,8 +16,12 @@ std::vector XC_Functional_Libxc::convert_rho( #pragma omp parallel for collapse(2) schedule(static, 1024) #endif for( int is=0; isrho[is][ir] + 1.0/nspin*chr->rho_core[ir]; + } + } return rho; } @@ -64,7 +68,9 @@ XC_Functional_Libxc::cal_gdr( #pragma omp parallel for schedule(static, 1024) #endif for(std::size_t ir=0; ir XC_Functional_Libxc::convert_sigma( assert(nspin>0); const std::size_t nrxx = gdr[0].size(); for(std::size_t is=1; is sigma( nrxx * ((1==nspin)?1:3) ); if( 1==nspin ) @@ -126,6 +134,7 @@ std::vector XC_Functional_Libxc::cal_sgn( const std::vector &rho, const std::vector &sigma) { + //assert(nrxx>0); // adding this once will cause error in examples std::vector sgn(nrxx*nspin, 1.0); // in the case of GGA correlation for polarized case, // a cutoff for grho is required to ensure that libxc gives reasonable results @@ -137,9 +146,13 @@ std::vector XC_Functional_Libxc::cal_sgn( for( int ir=0; ir XC_Functional_Libxc::convert_vtxc_v( const double tpiba, const Charge* const chr) { + // assert(nrxx>0); // will cause error double vtxc = 0.0; ModuleBase::matrix v(nspin, nrxx); @@ -186,9 +204,10 @@ std::pair XC_Functional_Libxc::convert_vtxc_v( { for( std::size_t ir=0; ir> XC_Functional_Libxc::cal_dh( const double tpiba, const Charge* const chr) { + //assert(nrxx>0); // this line will cause bug std::vector>> h( nspin, std::vector>(nrxx) ); - if( 1==nspin ) + + if( nspin==1 ) { #ifdef _OPENMP #pragma omp parallel for schedule(static, 1024) #endif for( std::size_t ir=0; ir> XC_Functional_Libxc::cal_dh( // define two dimensional array dh [ nspin, nrxx ] std::vector> dh(nspin, std::vector(nrxx)); for( int is=0; is!=nspin; ++is ) + { XC_Functional::grad_dot( h[is].data(), dh[is].data(), chr->rhopw, tpiba); + } return dh; } @@ -267,11 +292,14 @@ ModuleBase::matrix XC_Functional_Libxc::convert_v_nspin4( const std::vector &amag, const ModuleBase::matrix &v) { + //assert(nrxx>0); assert(PARAM.inp.nspin==4); constexpr double vanishing_charge = 1.0e-10; ModuleBase::matrix v_nspin4(PARAM.inp.nspin, nrxx); for( int ir=0; irrho[ipol][ir] / amag[ir]; + } } } } return v_nspin4; } -#endif \ No newline at end of file +#endif diff --git a/source/source_io/CMakeLists.txt b/source/source_io/CMakeLists.txt index 04db6f599a..6c7b5bf349 100644 --- a/source/source_io/CMakeLists.txt +++ b/source/source_io/CMakeLists.txt @@ -60,7 +60,7 @@ if(ENABLE_LCAO) read_wfc_nao.cpp read_wfc_lcao.cpp write_wfc_nao.cpp - io_dmk.cpp + write_dmk.cpp write_dmr.cpp sparse_matrix.cpp file_reader.cpp diff --git a/source/source_io/ctrl_output_lcao.cpp b/source/source_io/ctrl_output_lcao.cpp index 6a334aa015..673e767198 100644 --- a/source/source_io/ctrl_output_lcao.cpp +++ b/source/source_io/ctrl_output_lcao.cpp @@ -8,7 +8,7 @@ // functions #include "source_io/write_dos_lcao.h" // use ModuleIO::write_dos_lcao() #include "source_io/write_dmr.h" // use ModuleIO::write_dmr() -#include "source_io/io_dmk.h" // use ModuleIO::write_dmk() +#include "source_io/write_dmk.h" // use ModuleIO::write_dmk() #include "source_io/write_HS.h" // use ModuleIO::write_hsk() #include "source_io/write_wfc_nao.h" // use ModuleIO::write_wfc_nao() #include "source_io/output_mat_sparse.h" // use ModuleIO::output_mat_sparse() @@ -92,26 +92,30 @@ void ctrl_output_lcao(UnitCell& ucell, //------------------------------------------------------------------ //! 1) Output density matrix DM(R) //------------------------------------------------------------------ - if(PARAM.inp.out_dmr) + if(PARAM.inp.out_dmr[0]) { const auto& dmr_vector = pelec->get_DM()->get_DMR_vector(); - ModuleIO::write_dmr(dmr_vector, pv, out_app_flag, + + const int precision = PARAM.inp.out_dmr[1]; + + ModuleIO::write_dmr(dmr_vector, precision, pv, out_app_flag, ucell.get_iat2iwt(), ucell.nat, istep); } //------------------------------------------------------------------ //! 2) Output density matrix DM(k) //------------------------------------------------------------------ - if (PARAM.inp.out_dmk) + if (PARAM.inp.out_dmk[0]) { std::vector efermis(nspin == 2 ? 2 : 1); for (int ispin = 0; ispin < efermis.size(); ispin++) { efermis[ispin] = pelec->eferm.get_efval(ispin); } - const int precision = 3; + const int precision = PARAM.inp.out_dmk[1]; + ModuleIO::write_dmk(pelec->get_DM()->get_DMK_vector(), - precision, efermis, &(ucell), pv); + precision, efermis, &(ucell), pv, istep); } //------------------------------------------------------------------ diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index 2e79e806a2..2286dc72ab 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -372,8 +372,8 @@ struct Input_para bool out_mul = false; ///< qifeng add 2019-9-10 bool out_proj_band = false; ///< projected band structure calculation jiyy add 2022-05-11 std::string out_level = "ie"; ///< control the output information. - bool out_dmk = false; ///< output density matrix DM(k) - bool out_dmr = false; ///< output density matrix DM(R) + std::vector out_dmr = {0, 8}; ///< output density matrix in real space DM(R) + std::vector out_dmk = {0, 8}; ///< output density matrix in reciprocal space DM(k) bool out_bandgap = false; ///< QO added for bandgap printing std::vector out_mat_hs = {0, 8}; ///< output H matrix and S matrix in local basis. std::vector out_mat_tk = {0, 8}; ///< output T(k) matrix in local basis. diff --git a/source/source_io/nscf_band.cpp b/source/source_io/nscf_band.cpp index b93b09937b..880908f57a 100644 --- a/source/source_io/nscf_band.cpp +++ b/source/source_io/nscf_band.cpp @@ -21,6 +21,8 @@ void ModuleIO::nscf_band( ModuleBase::TITLE("ModuleIO","nscf_band"); ModuleBase::timer::tick("ModuleIO", "nscf_band"); + assert(precision>0); + GlobalV::ofs_running << "\n"; GlobalV::ofs_running << " >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>" << std::endl; GlobalV::ofs_running << " | |" << std::endl; diff --git a/source/source_io/output.cpp b/source/source_io/output.cpp index cf6d5eb86f..7512e0bfe2 100644 --- a/source/source_io/output.cpp +++ b/source/source_io/output.cpp @@ -10,7 +10,6 @@ void output::printrm(std::ofstream &ofs,const std::string &s, const ModuleBase:: if (b1*b2 == 0) return; for (int i = 0;i < b1;i++) { -// ofs<<"\n row = "< limit) std::cout << std::setprecision(15) << std::setw(20) << m(i,j); if (std::abs(m(i,j)) > limit) std::cout << std::setprecision(6) << std::setw(12) << m(i,j); -// else std::cout< limit) ofs<< std::setw(12) << sqrt(norm); - else ofs< limit) std::cout<< std::setw(12) << sqrt(norm); - else std::cout<(p_ham, pv, PARAM.inp.nspin, kv.get_nks()); - auto out_dmk = ModuleIO::Output_DMK(dynamic_cast*>(pelec)->get_DM(), + auto out_s_k = ModuleIO::Output_Sk(p_ham, pv, PARAM.inp.nspin, kv.get_nks()); + auto out_dm_k = ModuleIO::Output_DMK(dynamic_cast*>(pelec)->get_DM(), pv, PARAM.inp.nspin, kv.get_nks()); - auto mulp = ModuleIO::Output_Mulliken(&(out_sk), &(out_dmk), pv, &cell_index, kv.isk, PARAM.inp.nspin); + + auto mulp = ModuleIO::Output_Mulliken(&(out_s_k), &(out_dm_k), pv, &cell_index, kv.isk, PARAM.inp.nspin); auto atom_chg = mulp.get_atom_chg(); /// used in updating mag info in STRU file ucell.atom_mulliken = mulp.get_atom_mulliken(atom_chg); diff --git a/source/source_io/read_input_item_output.cpp b/source/source_io/read_input_item_output.cpp index 5b3e880553..95b1849f62 100644 --- a/source/source_io/read_input_item_output.cpp +++ b/source/source_io/read_input_item_output.cpp @@ -43,15 +43,16 @@ void ReadInput::item_output() Input_Item item("out_chg"); item.annotation = "> 0 output charge density for selected electron steps" ", second parameter controls the precision, default is 3."; - item.read_value = [](const Input_Item& item, Parameter& para) { - const size_t count = item.get_size(); - if (count != 1 && count != 2) - { - ModuleBase::WARNING_QUIT("ReadInput", "out_chg should have 1 or 2 values"); - } - para.input.out_chg[0] = (item.str_values[0] == "-1") ? -1 : std::stoi(item.str_values[0]); - para.input.out_chg[1] = (count == 2) ? std::stoi(item.str_values[1]) : 3; - }; + item.read_value = [](const Input_Item& item, Parameter& para) { + const size_t count = item.get_size(); + if (count < 1) ModuleBase::WARNING_QUIT("ReadInput", "out_chg needs at least 1 value"); + para.input.out_chg[0] = std::stoi(item.str_values[0]); + para.input.out_chg[1] = 3; + if (count >= 2) try { para.input.out_chg[1] = std::stoi(item.str_values[1]); } + catch (const std::invalid_argument&) { /* do nothing */ } + catch (const std::out_of_range&) {/* do nothing */} + }; + // reset value in some special case item.reset_value = [](const Input_Item& item, Parameter& para) { para.input.out_chg[0] = (para.input.calculation == "get_wf" || para.input.calculation == "get_pchg") ? 1 @@ -196,31 +197,50 @@ void ReadInput::item_output() { Input_Item item("out_dmk"); item.annotation = ">0 output density matrix DM(k) for each k-point"; - item.reset_value = [](const Input_Item& item, Parameter& para) { - if (para.input.calculation == "get_pchg" || para.input.calculation == "get_wf") - { - para.input.out_dmk = false; - } - }; - read_sync_bool(input.out_dmk); + item.read_value = [](const Input_Item& item, Parameter& para) { + const size_t count = item.get_size(); + if (count < 1) ModuleBase::WARNING_QUIT("ReadInput", "out_dmk needs at least 1 value"); + para.input.out_dmk[0] = assume_as_boolean(item.str_values[0]); + para.input.out_dmk[1] = 8; + if (count >= 2) try { para.input.out_dmk[1] = std::stoi(item.str_values[1]); } + catch (const std::invalid_argument&) { /* do nothing */ } + catch (const std::out_of_range&) {/* do nothing */} + // some other case + if (para.input.calculation == "get_pchg" || para.input.calculation == "get_wf") + { + para.input.out_dmk[0] = 0; + } + }; + sync_intvec(input.out_dmk, 2, 0); this->add_item(item); } { Input_Item item("out_dmr"); - item.annotation = ">0 output density matrix DM(R) with respect to lattice vector R"; - item.reset_value = [](const Input_Item& item, Parameter& para) { - if (para.input.calculation == "get_pchg" || para.input.calculation == "get_wf") - { - para.input.out_dmr = false; - } - }; + item.annotation = "output density matrix DM(R) with respect to lattice vector R (with precision 8)"; + + item.read_value = [](const Input_Item& item, Parameter& para) { + const size_t count = item.get_size(); + if (count < 1) ModuleBase::WARNING_QUIT("ReadInput", "out_dmr needs at least 1 value"); + para.input.out_dmr[0] = assume_as_boolean(item.str_values[0]); + para.input.out_dmr[1] = 8; + if (count >= 2) try { para.input.out_dmr[1] = std::stoi(item.str_values[1]); } + catch (const std::invalid_argument&) { /* do nothing */ } + catch (const std::out_of_range&) {/* do nothing */} + // some special case + if (para.input.calculation == "get_pchg" || para.input.calculation == "get_wf") + { + para.input.out_dmr[0] = 0; + } + }; + item.check_value = [](const Input_Item& item, const Parameter& para) { - if (para.sys.gamma_only_local == true && para.input.out_dmr) + if (para.sys.gamma_only_local == true && para.input.out_dmr[0]) { ModuleBase::WARNING_QUIT("ReadInput", "out_dmr is only valid for multi-k calculation"); } }; - read_sync_bool(input.out_dmr); + + sync_intvec(input.out_dmr, 2, 0); this->add_item(item); } { @@ -232,15 +252,16 @@ void ReadInput::item_output() { Input_Item item("out_mat_hs"); item.annotation = "output H and S matrix (with precision 8)"; - item.read_value = [](const Input_Item& item, Parameter& para) { - const size_t count = item.get_size(); - if (count != 1 && count != 2) - { - ModuleBase::WARNING_QUIT("ReadInput", "out_mat_hs should have 1 or 2 values"); - } - para.input.out_mat_hs[0] = assume_as_boolean(item.str_values[0]); - para.input.out_mat_hs[1] = (count == 2) ? std::stoi(item.str_values[1]) : 8; - }; + item.read_value = [](const Input_Item& item, Parameter& para) { + const size_t count = item.get_size(); + if (count < 1) ModuleBase::WARNING_QUIT("ReadInput", "out_mat_hs needs at least 1 value"); + para.input.out_mat_hs[0] = assume_as_boolean(item.str_values[0]); + para.input.out_mat_hs[1] = 8; + if (count >= 2) try { para.input.out_mat_hs[1] = std::stoi(item.str_values[1]); } + catch (const std::invalid_argument&) { /* do nothing */ } + catch (const std::out_of_range&) {/* do nothing */} + }; + // reset value in some special case item.reset_value = [](const Input_Item& item, Parameter& para) { if (para.input.qo_switch) { @@ -252,15 +273,15 @@ void ReadInput::item_output() } { Input_Item item("out_mat_tk"); - item.annotation = "output T(k)"; - item.read_value = [](const Input_Item& item, Parameter& para) { - const size_t count = item.get_size(); - if (count != 1 && count != 2) - { - ModuleBase::WARNING_QUIT("ReadInput", "out_mat_tk should have 1 or 2 values"); - } - para.input.out_mat_tk[0] = assume_as_boolean(item.str_values[0]); - para.input.out_mat_tk[1] = (count == 2) ? std::stoi(item.str_values[1]) : 8; + item.annotation = "output kinetic matrix of electrons T(k)"; + item.read_value = [](const Input_Item& item, Parameter& para) { + const size_t count = item.get_size(); + if (count < 1) ModuleBase::WARNING_QUIT("ReadInput", "out_mat_tk needs at least 1 value"); + para.input.out_mat_tk[0] = assume_as_boolean(item.str_values[0]); + para.input.out_mat_tk[1] = 8; + if (count >= 2) try { para.input.out_mat_tk[1] = std::stoi(item.str_values[1]); } + catch (const std::invalid_argument&) { /* do nothing */ } + catch (const std::out_of_range&) {/* do nothing */} }; sync_intvec(input.out_mat_tk, 2, 0); this->add_item(item); @@ -280,14 +301,14 @@ void ReadInput::item_output() { Input_Item item("out_mat_l"); item.annotation = "output the expectation values of angular momentum operators"; - item.read_value = [](const Input_Item& item, Parameter& para) { - const size_t count = item.get_size(); - if (count != 1 && count != 2) - { - ModuleBase::WARNING_QUIT("ReadInput", "out_mat_l should have 1 or 2 values"); - } - para.input.out_mat_l[0] = assume_as_boolean(item.str_values[0]); - para.input.out_mat_l[1] = (count == 2) ? std::stoi(item.str_values[1]) : 8; + item.read_value = [](const Input_Item& item, Parameter& para) { + const size_t count = item.get_size(); + if (count < 1) ModuleBase::WARNING_QUIT("ReadInput", "out_mat_l needs at least 1 value"); + para.input.out_mat_l[0] = assume_as_boolean(item.str_values[0]); + para.input.out_mat_l[1] = 8; + if (count >= 2) try { para.input.out_mat_l[1] = std::stoi(item.str_values[1]); } + catch (const std::invalid_argument&) { /* do nothing */ } + catch (const std::out_of_range&) {/* do nothing */} }; sync_intvec(input.out_mat_l, 2, 0); this->add_item(item); @@ -502,10 +523,13 @@ void ReadInput::item_output() item.read_value = [](const Input_Item& item, Parameter& para) { size_t count = item.get_size(); std::vector out_xc_r(count); // create a placeholder vector - std::transform(item.str_values.begin(), item.str_values.end(), out_xc_r.begin(), [](std::string s) { return std::stoi(s); }); + std::transform(item.str_values.begin(), + item.str_values.end(), + out_xc_r.begin(), [](std::string s) { return std::stoi(s); }); // assign non-negative values to para.input.out_xc_r std::copy(out_xc_r.begin(), out_xc_r.end(), para.input.out_xc_r.begin()); }; + // check value item.check_value = [](const Input_Item& item, const Parameter& para) { if (para.input.out_xc_r[0] >= 0) { diff --git a/source/source_io/sparse_matrix.cpp b/source/source_io/sparse_matrix.cpp index 779760430a..10f1cf75db 100644 --- a/source/source_io/sparse_matrix.cpp +++ b/source/source_io/sparse_matrix.cpp @@ -113,4 +113,4 @@ T SparseMatrix::operator()(int row, int col) const template class SparseMatrix; template class SparseMatrix>; -} // namespace ModuleIO \ No newline at end of file +} // namespace ModuleIO diff --git a/source/source_io/test/CMakeLists.txt b/source/source_io/test/CMakeLists.txt index a5ec46d458..48b52a2054 100644 --- a/source/source_io/test/CMakeLists.txt +++ b/source/source_io/test/CMakeLists.txt @@ -231,14 +231,14 @@ add_test(NAME MODULE_IO_orb_io_test_parallel ) AddTest( - TARGET MODULE_IO_dmk_io + TARGET MODULE_IO_write_dmk LIBS parameter ${math_libs} base device cell_info - SOURCES io_dmk_test.cpp ../io_dmk.cpp ../output.cpp + SOURCES write_dmk_test.cpp ../write_dmk.cpp ../output.cpp ) add_test( - NAME MODULE_IO_dmk_io_parallel - COMMAND mpirun -np 2 ./MODULE_IO_dmk_io + NAME MODULE_IO_write_dmk_parallel + COMMAND mpirun -np 2 ./MODULE_IO_write_dmk WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} ) diff --git a/source/source_io/test/output_mulliken_test.cpp b/source/source_io/test/output_mulliken_test.cpp index 1abe3614d7..651b2a2b52 100644 --- a/source/source_io/test/output_mulliken_test.cpp +++ b/source/source_io/test/output_mulliken_test.cpp @@ -41,9 +41,9 @@ TYPED_TEST(OutputMullikenTest, nspin1) this->ncol = 13; this->paraV.init(this->nrow, this->ncol, 1, MPI_COMM_WORLD, 0); auto cell_index = CellIndex(this->atomLabels, this->atomCounts, this->lnchiCounts, 1); - auto out_sk = ModuleIO::Output_Sk(nullptr, &this->paraV, 1, 1); - auto out_dmk = ModuleIO::Output_DMK(nullptr, &this->paraV, 1, 1); - auto mulp = ModuleIO::Output_Mulliken(&(out_sk), &(out_dmk), &(this->paraV), &(cell_index), {0}, 1); + auto out_s_k = ModuleIO::Output_Sk(nullptr, &this->paraV, 1, 1); + auto out_dm_k = ModuleIO::Output_DMK(nullptr, &this->paraV, 1, 1); + auto mulp = ModuleIO::Output_Mulliken(&(out_s_k), &(out_dm_k), &(this->paraV), &(cell_index), {0}, 1); mulp.write(0, "./"); std::vector tot_chg = mulp.get_tot_chg(); EXPECT_NEAR(tot_chg[0], 4.0, 1e-5); @@ -60,9 +60,9 @@ TYPED_TEST(OutputMullikenTest, nspin2) this->ncol = 13; this->paraV.init(this->nrow, this->ncol, 1, MPI_COMM_WORLD, 0); auto cell_index = CellIndex(this->atomLabels, this->atomCounts, this->lnchiCounts, 2); - auto out_sk = ModuleIO::Output_Sk(nullptr, &this->paraV, 2, 1); - auto out_dmk = ModuleIO::Output_DMK(nullptr, &this->paraV, 2, 1); - auto mulp = ModuleIO::Output_Mulliken(&(out_sk), &(out_dmk), &(this->paraV), &(cell_index), {0, 1}, 2); + auto out_s_k = ModuleIO::Output_Sk(nullptr, &this->paraV, 2, 1); + auto out_dm_k = ModuleIO::Output_DMK(nullptr, &this->paraV, 2, 1); + auto mulp = ModuleIO::Output_Mulliken(&(out_s_k), &(out_dm_k), &(this->paraV), &(cell_index), {0, 1}, 2); mulp.write(0, "./"); std::vector tot_chg = mulp.get_tot_chg(); EXPECT_NEAR(tot_chg[0], 3.0, 1e-5); @@ -83,10 +83,10 @@ TYPED_TEST(OutputMullikenTest, nspin4) this->ncol = 26; this->paraV.init(this->nrow, this->ncol, 1, MPI_COMM_WORLD, 0); auto cell_index = CellIndex(this->atomLabels, this->atomCounts, this->lnchiCounts, 4); - auto out_sk = ModuleIO::Output_Sk>(nullptr, &this->paraV, 4, 1); - auto out_dmk = ModuleIO::Output_DMK>(nullptr, &this->paraV, 4, 1); + auto out_s_k = ModuleIO::Output_Sk>(nullptr, &this->paraV, 4, 1); + auto out_dm_k = ModuleIO::Output_DMK>(nullptr, &this->paraV, 4, 1); auto mulp - = ModuleIO::Output_Mulliken>(&(out_sk), &(out_dmk), &(this->paraV), &(cell_index), {0}, 4); + = ModuleIO::Output_Mulliken>(&(out_s_k), &(out_dm_k), &(this->paraV), &(cell_index), {0}, 4); mulp.write(0, "./"); std::vector tot_chg = mulp.get_tot_chg(); EXPECT_NEAR(tot_chg[0], 4.0, 1e-5); @@ -123,4 +123,4 @@ int main(int argc, char** argv) return result; } -#endif \ No newline at end of file +#endif diff --git a/source/source_io/test/read_input_ptest.cpp b/source/source_io/test/read_input_ptest.cpp index d6927e61b7..14b0ef5a93 100644 --- a/source/source_io/test/read_input_ptest.cpp +++ b/source/source_io/test/read_input_ptest.cpp @@ -186,8 +186,8 @@ TEST_F(InputParaTest, ParaRead) EXPECT_EQ(param.inp.out_chg[1], 3); EXPECT_EQ(param.inp.out_elf[0], 0); EXPECT_EQ(param.inp.out_elf[1], 3); - EXPECT_EQ(param.inp.out_dmk, 0); - EXPECT_EQ(param.inp.out_dmr, 0); + EXPECT_EQ(param.inp.out_dmk[0], 0); + EXPECT_EQ(param.inp.out_dmr[0], 0); EXPECT_EQ(param.inp.deepks_out_labels, 0); EXPECT_EQ(param.inp.deepks_scf, 0); EXPECT_EQ(param.inp.deepks_equiv, 0); diff --git a/source/source_io/test/support/dms1_nao.txt b/source/source_io/test/support/dm_nao.txt similarity index 99% rename from source/source_io/test/support/dms1_nao.txt rename to source/source_io/test/support/dm_nao.txt index 8f03975f1e..811fa19b2f 100644 --- a/source/source_io/test/support/dms1_nao.txt +++ b/source/source_io/test/support/dm_nao.txt @@ -1,17 +1,17 @@ -none + none 5.39761 -0.5 0 0.5 0 0.5 0.5 -0.5 0.5 0 Si 2 -Direct + Direct 0 0 0 0.75 0.75 0.75 - 1 + 1 (nspin) 0.570336288802337 (fermi energy) - 26 26 + 26 (number of basis) 3.904e-01 1.114e-02 5.918e-14 -1.032e-14 -5.650e-14 -2.278e-15 1.049e-14 7.875e-15 -1.022e-15 -9.243e-15 -4.239e-15 -4.177e-16 1.138e-14 3.904e-01 1.114e-02 -5.087e-14 @@ -116,4 +116,4 @@ Direct 9.615e-15 -3.453e-15 1.810e-01 5.280e-15 2.578e-15 -1.606e-02 -1.587e-15 9.162e-17 9.484e-17 2.460e-16 4.817e-16 -7.012e-17 3.445e-02 -1.025e-14 4.295e-15 -1.810e-01 -4.065e-16 9.526e-16 1.606e-02 2.413e-15 7.651e-16 -7.043e-18 2.049e-16 1.512e-15 - -2.685e-17 3.445e-02 \ No newline at end of file + -2.685e-17 3.445e-02 diff --git a/source/source_io/test/support/dms1k1_nao.txt b/source/source_io/test/support/dmk1_nao.txt similarity index 99% rename from source/source_io/test/support/dms1k1_nao.txt rename to source/source_io/test/support/dmk1_nao.txt index 9cacc456e0..ff611f7c7e 100644 --- a/source/source_io/test/support/dms1k1_nao.txt +++ b/source/source_io/test/support/dmk1_nao.txt @@ -1,17 +1,17 @@ -none + none 1 2.74242 2.74242 0 2.74242 0 2.74242 0 2.74242 2.74242 Si 2 -Direct + Direct 0 0 0 0.25 0.25 0.25 - 1 + 1 (nspin) 0.47370 (fermi energy) - 26 26 + 26 (number of basis) (1.411e-01,0.000e+00) (-4.479e-03,3.208e-04) (-1.892e-02,3.039e-02) (-1.892e-02,3.039e-02) (1.892e-02,-3.039e-02) (-2.049e-03,-5.032e-03) (-2.049e-03,-5.032e-03) (2.049e-03,5.032e-03) (-6.505e-17,-4.647e-18) (2.931e-03,6.477e-03) (-2.931e-03,-6.477e-03) (-2.100e-16,6.607e-17) (-2.931e-03,-6.477e-03) (-1.319e-02,3.501e-02) (9.561e-03,-5.528e-03) (-3.910e-02,-2.210e-02) @@ -116,4 +116,4 @@ Direct (-2.387e-03,7.959e-03) (1.504e-04,-6.531e-04) (-1.637e-02,-1.991e-03) (4.284e-03,-2.804e-03) (-4.284e-03,2.804e-03) (1.374e-03,8.137e-04) (4.245e-05,-5.221e-04) (-4.245e-05,5.221e-04) (-3.694e-04,2.397e-04) (4.101e-04,6.505e-04) (-4.101e-04,-6.505e-04) (-1.148e-17,-1.884e-18) (2.872e-03,1.079e-03) (-2.931e-03,-6.477e-03) (2.300e-04,7.056e-04) (1.438e-02,4.085e-03) (-3.513e-03,-6.260e-03) (3.513e-03,6.260e-03) (-1.691e-03,3.401e-04) (1.097e-04,-2.210e-04) (-1.097e-04,2.210e-04) (-2.151e-04,-3.843e-04) (5.354e-04,-4.093e-18) (-5.354e-04,-1.423e-17) - (8.904e-18,7.971e-19) (3.174e-03,0.000e+00) \ No newline at end of file + (8.904e-18,7.971e-19) (3.174e-03,0.000e+00) diff --git a/source/source_io/test/io_dmk_test.cpp b/source/source_io/test/write_dmk_test.cpp similarity index 80% rename from source/source_io/test/io_dmk_test.cpp rename to source/source_io/test/write_dmk_test.cpp index 673185648c..5c56a0bfad 100644 --- a/source/source_io/test/io_dmk_test.cpp +++ b/source/source_io/test/write_dmk_test.cpp @@ -1,4 +1,4 @@ -#include "source_io/io_dmk.h" +#include "source_io/write_dmk.h" #define private public #include "source_io/module_parameter/parameter.h" @@ -101,15 +101,30 @@ void gen_dmk(std::vector>& dmk, std::vector& efs, int ns TEST(DMKTest, GenFileName) { - std::string fname = ModuleIO::dmk_gen_fname(true, 0, 0); - EXPECT_EQ(fname, "dms1_nao.txt"); - fname = ModuleIO::dmk_gen_fname(true, 1, 1); - EXPECT_EQ(fname, "dms2_nao.txt"); - - fname = ModuleIO::dmk_gen_fname(false, 0, 0); - EXPECT_EQ(fname, "dms1k1_nao.txt"); - fname = ModuleIO::dmk_gen_fname(false, 1, 1); - EXPECT_EQ(fname, "dms2k2_nao.txt"); + bool gamma_only = true; + int ispin = 0; + int nspin = 2; + int ik = 0; + int istep = 0; + std::string fname = ModuleIO::dmk_gen_fname(gamma_only, ispin, nspin, ik, istep); + EXPECT_EQ(fname, "dms1g1_nao.txt"); + + ispin = 1; + + fname = ModuleIO::dmk_gen_fname(gamma_only, ispin, nspin, ik, istep); + EXPECT_EQ(fname, "dms2g1_nao.txt"); + + ispin = 0; + gamma_only = false; + + fname = ModuleIO::dmk_gen_fname(gamma_only, ispin, nspin, ik, istep); + EXPECT_EQ(fname, "dmk1s1g1_nao.txt"); + + ispin = 1; + ik = 1; + + fname = ModuleIO::dmk_gen_fname(gamma_only, ispin, nspin, ik, istep); + EXPECT_EQ(fname, "dmk2s2g1_nao.txt"); }; @@ -132,8 +147,9 @@ TEST(DMKTest,WriteDMK) { gen_dmk(dmk_multik, efs, nspin, nk_multik, nlocal, pv); PARAM.sys.global_out_dir = "./"; - ModuleIO::write_dmk(dmk, 3, efs, ucell, pv); - ModuleIO::write_dmk(dmk_multik, 3, efs, ucell, pv); + const int istep = -1; + ModuleIO::write_dmk(dmk, 3, efs, ucell, pv, istep); + ModuleIO::write_dmk(dmk_multik, 3, efs, ucell, pv, istep); std::ifstream ifs; int pass = 0; @@ -144,7 +160,7 @@ TEST(DMKTest,WriteDMK) { std::string str((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); EXPECT_THAT(str, testing::HasSubstr("0.00000 (fermi energy)")); - EXPECT_THAT(str, testing::HasSubstr("20 20")); + EXPECT_THAT(str, testing::HasSubstr("20 (number of basis)")); EXPECT_THAT( str, testing::HasSubstr("0.000e+00 1.000e-01 2.000e-01 3.000e-01 4.000e-01 " @@ -163,7 +179,7 @@ TEST(DMKTest,WriteDMK) { str = std::string((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); EXPECT_THAT(str, testing::HasSubstr("0.10000 (fermi energy)")); - EXPECT_THAT(str, testing::HasSubstr("20 20")); + EXPECT_THAT(str, testing::HasSubstr("20 (number of basis)")); EXPECT_THAT( str, testing::HasSubstr("1.000e+00 1.100e+00 1.200e+00 1.300e+00 1.400e+00 " @@ -177,84 +193,84 @@ TEST(DMKTest,WriteDMK) { testing::HasSubstr("2.600e+00 2.700e+00 2.800e+00 2.900e+00\n")); ifs.close(); - fn = "dms1k1_nao.txt"; + fn = "dmk1s1_nao.txt"; ifs.open(fn); str = std::string((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); EXPECT_THAT(str, testing::HasSubstr("0.00000 (fermi energy)")); - EXPECT_THAT(str, testing::HasSubstr("20 20")); + EXPECT_THAT(str, testing::HasSubstr("20 (number of basis)")); EXPECT_THAT( str, testing::HasSubstr("(0.000e+00,0.000e+00) (1.000e-01,1.000e+00) (2.000e-01,2.000e+00) " - "(3.000e-01,3.000e+00) (4.000e-01,4.000e+00) (5.000e-01,5.000e+00) " + "(3.000e-01,3.000e+00)\n (4.000e-01,4.000e+00) (5.000e-01,5.000e+00) " "(6.000e-01,6.000e+00) (7.000e-01,7.000e+00)\n")); EXPECT_THAT( str, testing::HasSubstr("(8.000e-01,8.000e+00) (9.000e-01,9.000e+00) (1.000e+00,1.000e+01) " - "(1.100e+00,1.100e+01) (1.200e+00,1.200e+01) (1.300e+00,1.300e+01) " + "(1.100e+00,1.100e+01)\n (1.200e+00,1.200e+01) (1.300e+00,1.300e+01) " "(1.400e+00,1.400e+01) (1.500e+00,1.500e+01)\n")); EXPECT_THAT( str, testing::HasSubstr("(1.600e+00,1.600e+01) (1.700e+00,1.700e+01) (1.800e+00,1.800e+01) (1.900e+00,1.900e+01)\n")); ifs.close(); - fn = "dms1k2_nao.txt"; + fn = "dmk2s1_nao.txt"; ifs.open(fn); str = std::string((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); EXPECT_THAT(str, testing::HasSubstr("0.00000 (fermi energy)")); - EXPECT_THAT(str, testing::HasSubstr("20 20")); + EXPECT_THAT(str, testing::HasSubstr("20 (number of basis)")); EXPECT_THAT( str, testing::HasSubstr("(1.000e+00,1.000e-01) (1.100e+00,1.100e+00) (1.200e+00,2.100e+00) " - "(1.300e+00,3.100e+00) (1.400e+00,4.100e+00) (1.500e+00,5.100e+00) " + "(1.300e+00,3.100e+00)\n (1.400e+00,4.100e+00) (1.500e+00,5.100e+00) " "(1.600e+00,6.100e+00) (1.700e+00,7.100e+00)\n")); EXPECT_THAT( str, testing::HasSubstr("(1.800e+00,8.100e+00) (1.900e+00,9.100e+00) (2.000e+00,1.010e+01) " - "(2.100e+00,1.110e+01) (2.200e+00,1.210e+01) (2.300e+00,1.310e+01) " + "(2.100e+00,1.110e+01)\n (2.200e+00,1.210e+01) (2.300e+00,1.310e+01) " "(2.400e+00,1.410e+01) (2.500e+00,1.510e+01)\n")); EXPECT_THAT( str, testing::HasSubstr("(2.600e+00,1.610e+01) (2.700e+00,1.710e+01) (2.800e+00,1.810e+01) (2.900e+00,1.910e+01)\n")); ifs.close(); - fn = "dms2k1_nao.txt"; + fn = "dmk1s2_nao.txt"; ifs.open(fn); str = std::string((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); EXPECT_THAT(str, testing::HasSubstr("0.10000 (fermi energy)")); - EXPECT_THAT(str, testing::HasSubstr("20 20")); + EXPECT_THAT(str, testing::HasSubstr("20 (number of basis)")); EXPECT_THAT( str, testing::HasSubstr("(2.000e+00,2.000e-01) (2.100e+00,1.200e+00) (2.200e+00,2.200e+00) " - "(2.300e+00,3.200e+00) (2.400e+00,4.200e+00) (2.500e+00,5.200e+00) " + "(2.300e+00,3.200e+00)\n (2.400e+00,4.200e+00) (2.500e+00,5.200e+00) " "(2.600e+00,6.200e+00) (2.700e+00,7.200e+00)\n")); EXPECT_THAT( str, testing::HasSubstr("(2.800e+00,8.200e+00) (2.900e+00,9.200e+00) (3.000e+00,1.020e+01) " - "(3.100e+00,1.120e+01) (3.200e+00,1.220e+01) (3.300e+00,1.320e+01) " + "(3.100e+00,1.120e+01)\n (3.200e+00,1.220e+01) (3.300e+00,1.320e+01) " "(3.400e+00,1.420e+01) (3.500e+00,1.520e+01)\n")); EXPECT_THAT( str, testing::HasSubstr("(3.600e+00,1.620e+01) (3.700e+00,1.720e+01) (3.800e+00,1.820e+01) (3.900e+00,1.920e+01)\n")); ifs.close(); - fn = "dms2k2_nao.txt"; + fn = "dmk2s2_nao.txt"; ifs.open(fn); str = std::string((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); EXPECT_THAT(str, testing::HasSubstr("0.10000 (fermi energy)")); - EXPECT_THAT(str, testing::HasSubstr("20 20")); + EXPECT_THAT(str, testing::HasSubstr("20 (number of basis)")); EXPECT_THAT( str, testing::HasSubstr("(3.000e+00,3.000e-01) (3.100e+00,1.300e+00) (3.200e+00,2.300e+00) " - "(3.300e+00,3.300e+00) (3.400e+00,4.300e+00) (3.500e+00,5.300e+00) " + "(3.300e+00,3.300e+00)\n (3.400e+00,4.300e+00) (3.500e+00,5.300e+00) " "(3.600e+00,6.300e+00) (3.700e+00,7.300e+00)\n")); EXPECT_THAT( str, testing::HasSubstr("(3.800e+00,8.300e+00) (3.900e+00,9.300e+00) (4.000e+00,1.030e+01) " - "(4.100e+00,1.130e+01) (4.200e+00,1.230e+01) (4.300e+00,1.330e+01) " + "(4.100e+00,1.130e+01)\n (4.200e+00,1.230e+01) (4.300e+00,1.330e+01) " "(4.400e+00,1.430e+01) (4.500e+00,1.530e+01)\n")); EXPECT_THAT( str, @@ -262,10 +278,10 @@ TEST(DMKTest,WriteDMK) { ifs.close(); remove("dms1_nao.txt"); remove("dms2_nao.txt"); - remove("dms1k1_nao.txt"); - remove("dms1k2_nao.txt"); - remove("dms2k1_nao.txt"); - remove("dms2k2_nao.txt"); + remove("dmk1s1_nao.txt"); + remove("dmk2s1_nao.txt"); + remove("dmk1s2_nao.txt"); + remove("dmk2s2_nao.txt"); } delete ucell; @@ -287,6 +303,8 @@ TEST(DMKTest, ReadDMK) { std::ofstream ofs_running("running_log.txt"); + GlobalV::ofs_warning.open("warning.log"); + EXPECT_TRUE(ModuleIO::read_dmk(1, 1, pv, "./support/", dmk, ofs_running)); ModuleIO::read_dmk(1, 1, pv, "./support/", dmk_multik, ofs_running); EXPECT_TRUE(ModuleIO::read_dmk(1, 1, pv, "./support/", dmk_multik, ofs_running)); @@ -302,7 +320,9 @@ TEST(DMKTest, ReadDMK) { } ofs_running.close(); + GlobalV::ofs_warning.close(); remove("running_log.txt"); + remove("warning.log"); } diff --git a/source/source_io/test_serial/read_input_item_test.cpp b/source/source_io/test_serial/read_input_item_test.cpp index d773eb9456..db842bd321 100644 --- a/source/source_io/test_serial/read_input_item_test.cpp +++ b/source/source_io/test_serial/read_input_item_test.cpp @@ -813,24 +813,27 @@ TEST_F(InputTest, Item_test) } { // out_dmk auto it = find_label("out_dmk", readinput.input_lists); - param.input.calculation = "get_wf"; - param.input.out_dmk = true; - it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.out_dmk, false); + it->second.str_values = {"1"}; + it->second.read_value(it->second, param); + EXPECT_EQ(param.input.out_dmk[0], 1); + EXPECT_EQ(param.input.out_dmk[1], 8); + + it->second.str_values = {"1", "2"}; + it->second.read_value(it->second, param); + EXPECT_EQ(param.input.out_dmk[0], 1); + EXPECT_EQ(param.input.out_dmk[1], 2); } { // out_dmr auto it = find_label("out_dmr", readinput.input_lists); - param.input.calculation = "get_wf"; - param.input.out_dmr = true; - it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.out_dmr, false); + it->second.str_values = {"1"}; + it->second.read_value(it->second, param); + EXPECT_EQ(param.input.out_dmr[0], 1); + EXPECT_EQ(param.input.out_dmr[1], 8); - param.sys.gamma_only_local = true; - param.input.out_dmr = true; - testing::internal::CaptureStdout(); - EXPECT_EXIT(it->second.check_value(it->second, param), ::testing::ExitedWithCode(1), ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("NOTICE")); + it->second.str_values = {"1", "2"}; + it->second.read_value(it->second, param); + EXPECT_EQ(param.input.out_dmr[0], 1); + EXPECT_EQ(param.input.out_dmr[1], 2); } { // method_sto auto it = find_label("method_sto", readinput.input_lists); @@ -943,12 +946,6 @@ TEST_F(InputTest, Item_test) EXPECT_EQ(param.input.out_mat_hs[0], 1); EXPECT_EQ(param.input.out_mat_hs[1], 2); - it->second.str_values = {"1", "2", "3"}; - testing::internal::CaptureStdout(); - EXPECT_EXIT(it->second.read_value(it->second, param), ::testing::ExitedWithCode(1), ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("NOTICE")); - param.input.out_mat_hs = {0}; param.input.qo_switch = true; it->second.reset_value(it->second, param); diff --git a/source/source_io/io_dmk.cpp b/source/source_io/write_dmk.cpp similarity index 76% rename from source/source_io/io_dmk.cpp rename to source/source_io/write_dmk.cpp index 01015e1e15..7a6762eeef 100644 --- a/source/source_io/io_dmk.cpp +++ b/source/source_io/write_dmk.cpp @@ -1,4 +1,4 @@ -#include "source_io/io_dmk.h" +#include "source_io/write_dmk.h" #include "source_base/parallel_common.h" #include "source_base/module_external/scalapack_connector.h" @@ -55,24 +55,39 @@ Direct ''' */ -std::string ModuleIO::dmk_gen_fname(const bool gamma_only, const int ispin, const int ik) +std::string ModuleIO::dmk_gen_fname(const bool gamma_only, const int ispin, const int nspin, const int ik, const int istep) { - if (gamma_only) - { - return std::string("dm") + "s" + std::to_string(ispin + 1) + "_nao.txt"; + // set istep = -1 if you don't want the 'g' index appears in the file name + assert(istep>=-1); + + // ik should be the correct one + + std::string fname = "dm"; + + if (!gamma_only) + { + fname += "k" + std::to_string(ik + 1); } - else + + if (nspin == 2) + { + fname += "s" + std::to_string(ispin + 1); + } + + if( istep >= 0 ) { - // mohan update 2025-05-25, the index of 'ik' should be the correct 'ik' without spin - return std::string("dm") + "s" + std::to_string(ispin + 1) - + "k" + std::to_string(ik + 1) + "_nao.txt"; + fname += "g" + std::to_string(istep + 1); } + + fname += "_nao.txt"; + + return fname; } void ModuleIO::dmk_write_ucell(std::ofstream& ofs, const UnitCell* ucell) { // write the UnitCell information - ofs << ucell->latName << std::endl; + ofs << " " << ucell->latName << std::endl; ofs << " " << ucell->lat0 * ModuleBase::BOHR_TO_A << std::endl; ofs << " " << ucell->latvec.e11 << " " << ucell->latvec.e12 << " " << ucell->latvec.e13 << std::endl; ofs << " " << ucell->latvec.e21 << " " << ucell->latvec.e22 << " " << ucell->latvec.e23 << std::endl; @@ -87,7 +102,7 @@ void ModuleIO::dmk_write_ucell(std::ofstream& ofs, const UnitCell* ucell) ofs << " " << ucell->atoms[it].na; } ofs << std::endl; - ofs << "Direct" << std::endl; + ofs << " Direct" << std::endl; for (int it = 0; it < ucell->ntype; it++) { Atom* atom = &ucell->atoms[it]; @@ -141,7 +156,7 @@ void ModuleIO::dmk_readData(std::ifstream& ifs, std::complex& data) else { ModuleBase::WARNING_QUIT("ModuleIO::dmk_readData", - "Invalid complex number format: " + complex_str); + "Invalid complex number format in dmk: " + complex_str); } } @@ -165,18 +180,6 @@ bool ModuleIO::read_dmk(const int nspin, bool gamma_only = std::is_same::value; std::vector> dmk_global(nspin * nk, std::vector(nlocal * nlocal, 0)); - // write a lambda function to check the consistency of the data - auto check_consistency - = [&](const std::string& fn, const std::string& name, const std::string& value, const int& target) { - if (std::stoi(value) != target) - { - ModuleBase::WARNING("ModuleIO::read_dmk", name + " is not consistent in file < " + fn + " >."); - std::cout << name << " = " << target << ", " << name << " in file = " << value << std::endl; - return false; - } - return true; - }; - bool read_success = true; std::string tmp; if (my_rank == 0) @@ -185,7 +188,10 @@ bool ModuleIO::read_dmk(const int nspin, { for (int ik = 0; ik < nk; ik++) { - std::string fn = dmk_dir + dmk_gen_fname(gamma_only, ispin, ik); + // to read density matrix in k space, remember to delete the step information 'g' + // set istep = -1 if you don't want the 'g' index appears in the file name + const int istep = -1; + std::string fn = dmk_dir + dmk_gen_fname(gamma_only, ispin, nspin, ik, istep); std::ifstream ifs(fn.c_str()); if (!ifs) @@ -203,37 +209,32 @@ bool ModuleIO::read_dmk(const int nspin, // read the UnitCell dmk_read_ucell(ifs); - ifs >> tmp; // nspin - if (!check_consistency(fn, "nspin", tmp, nspin)) - { - read_success = false; - ifs.close(); - break; - } - ifs >> tmp; - ifs >> tmp; - ifs >> tmp; // fermi energy - ifs >> tmp; // nlocal - if (!check_consistency(fn, "nlocal", tmp, nlocal)) - { - read_success = false; - ifs.close(); - break; - } - ifs >> tmp; // nlocal - if (!check_consistency(fn, "nlocal", tmp, nlocal)) - { - read_success = false; - ifs.close(); - break; - } + int spin_tmp = 0; + ModuleBase::GlobalFunc::READ_VALUE(ifs, spin_tmp); + + double fermi_tmp = 0.0; + ModuleBase::GlobalFunc::READ_VALUE(ifs, fermi_tmp); + + int nlocal_tmp = 0; + ModuleBase::GlobalFunc::READ_VALUE(ifs, nlocal_tmp); + + if(nlocal_tmp==nlocal) + { + ofs_running << " number of basis (nlocal) is correct: " << nlocal << std::endl; + } + else + { + ModuleBase::WARNING_QUIT("ModuleIO::read_dmk","nlocal does not match!"); + } // read the DMK data + const size_t index_k = ik + nk * ispin; for (int i = 0; i < nlocal; ++i) { + const size_t index_i = i * nlocal; for (int j = 0; j < nlocal; ++j) { - dmk_readData(ifs, dmk_global[ik + nk * ispin][i * nlocal + j]); + dmk_readData(ifs, dmk_global[index_k][index_i + j]); } } ifs.close(); @@ -283,7 +284,8 @@ void ModuleIO::write_dmk(const std::vector>& dmk, const int precision, const std::vector& efs, const UnitCell* ucell, - const Parallel_2D& pv) + const Parallel_2D& pv, + const int istep) { ModuleBase::TITLE("ModuleIO", "write_dmk"); ModuleBase::timer::tick("ModuleIO", "write_dmk"); @@ -294,13 +296,17 @@ void ModuleIO::write_dmk(const std::vector>& dmk, #endif bool gamma_only = std::is_same::value; - int nlocal = pv.get_global_row_size(); - int nspin = efs.size(); - int nk = dmk.size() / nspin; + const int nlocal = pv.get_global_row_size(); + const int nspin = efs.size(); + assert(nspin > 0); + const int nk = dmk.size() / nspin; + const double dm_thr = 1.0e-16; // mohan set 2025-09-02 + if (nk * nspin != dmk.size()) { - ModuleBase::WARNING_QUIT("write_dmk", "The size of dmk is not consistent with nspin and nk."); + ModuleBase::WARNING_QUIT("ModuleIO::write_dmk", "The size of dmk is not consistent with nspin and nk."); } + Parallel_2D pv_glb; // when nspin == 2, assume the order of K in dmk is K1_up, K2_up, ..., @@ -330,7 +336,7 @@ void ModuleIO::write_dmk(const std::vector>& dmk, if (my_rank == 0) { - std::string fn = PARAM.globalv.global_out_dir + dmk_gen_fname(gamma_only, ispin, ik); + std::string fn = PARAM.globalv.global_out_dir + dmk_gen_fname(gamma_only, ispin, nspin, ik, istep); std::ofstream ofs(fn.c_str()); if (!ofs) @@ -347,29 +353,49 @@ void ModuleIO::write_dmk(const std::vector>& dmk, dmk_write_ucell(ofs, ucell); - ofs << "\n " << nspin; // nspin + ofs << "\n " << nspin << " (nspin)"; // nspin ofs << "\n " << std::fixed << std::setprecision(5) << efs[ispin] << " (fermi energy)"; - ofs << "\n " << nlocal << " " << nlocal << std::endl; + ofs << "\n " << nlocal << " (number of basis)" << std::endl; - ofs << std::setprecision(precision); + ofs << std::fixed; ofs << std::scientific; + ofs << std::setprecision(precision); + ofs << std::right; +// ofs << std::showpos; // show positive label for (int i = 0; i < nlocal; ++i) { for (int j = 0; j < nlocal; ++j) { - if (j % 8 == 0) - { - ofs << "\n"; - } if (std::is_same::value) { + if (j % 8 == 0) + { + ofs << "\n"; + } ofs << " " << dmk_global[i * nlocal + j]; } else if (std::is_same, T>::value) { - ofs << " (" << std::real(dmk_global[i * nlocal + j]) << "," - << std::imag(dmk_global[i * nlocal + j]) << ")"; + if (j % 4 == 0) + { + ofs << "\n"; + } + + double real_v = std::real(dmk_global[i * nlocal + j]); + if(std::abs(real_v) < dm_thr) + { + real_v = 0.0; + } + double imag_v = std::imag(dmk_global[i * nlocal + j]); + if(std::abs(imag_v) < dm_thr) + { + imag_v = 0.0; + } + + ofs << " (" << real_v << "," << imag_v << ")"; + // ofs << " (" << std::real(dmk_global[i * nlocal + j]) << "," + // << std::imag(dmk_global[i * nlocal + j]) << ")"; } } } @@ -399,10 +425,12 @@ template void ModuleIO::write_dmk(const std::vector> const int precision, const std::vector& efs, const UnitCell* ucell, - const Parallel_2D& pv); + const Parallel_2D& pv, + const int istep); template void ModuleIO::write_dmk>(const std::vector>>& dmk, const int precision, const std::vector& efs, const UnitCell* ucell, - const Parallel_2D& pv); + const Parallel_2D& pv, + const int istep); diff --git a/source/source_io/io_dmk.h b/source/source_io/write_dmk.h similarity index 94% rename from source/source_io/io_dmk.h rename to source/source_io/write_dmk.h index 19a0c66323..2fe64a9554 100644 --- a/source/source_io/io_dmk.h +++ b/source/source_io/write_dmk.h @@ -1,5 +1,5 @@ -#ifndef DM_IO_H -#define DM_IO_H +#ifndef WRITE_DMK_H +#define WRITE_DMK_H #include "source_base/parallel_2d.h" #include "source_cell/unitcell.h" @@ -17,7 +17,7 @@ namespace ModuleIO { * @param ik The index of the k-point. * @return The generated filename. */ -std::string dmk_gen_fname(const bool gamma_only, const int ispin, const int ik); +std::string dmk_gen_fname(const bool gamma_only, const int ispin, const int nspin, const int ik, const int istep); /** * @brief Writes the unit cell information to a DMK file. @@ -83,8 +83,9 @@ void write_dmk(const std::vector>& dmk, const int precision, const std::vector& efs, const UnitCell* ucell, - const Parallel_2D& pv); + const Parallel_2D& pv, + const int istep); } // namespace ModuleIO -#endif // IO_DMK_H +#endif diff --git a/source/source_io/write_dmr.cpp b/source/source_io/write_dmr.cpp index fe661da5ff..ed1288cfa7 100644 --- a/source/source_io/write_dmr.cpp +++ b/source/source_io/write_dmr.cpp @@ -35,7 +35,10 @@ std::string dmr_gen_fname(const int out_type, const int ispin, const bool append return fname; } -void write_dmr_csr(std::string& fname, hamilt::HContainer* dm_serial, const int istep) +void write_dmr_csr(std::string& fname, + hamilt::HContainer* dm_serial, + const int precision, + const int istep) { // write the head: ION step number, basis number and R loop number @@ -57,22 +60,25 @@ void write_dmr_csr(std::string& fname, hamilt::HContainer* dm_serial, co // write HR_serial to ofs const double sparse_threshold = 1e-10; - const int precision = 8; - hamilt::Output_HContainer out_dmr(dm_serial, ofs, sparse_threshold, precision); - out_dmr.write(); + hamilt::Output_HContainer dmr(dm_serial, ofs, sparse_threshold, precision); + dmr.write(); ofs.close(); } void write_dmr(const std::vector*> dmr, + const int precision, const Parallel_2D& paraV, const bool append, const int* iat2iwt, const int nat, const int istep) { - for (int ispin = 0; ispin < dmr.size(); ispin++) + const int nspin = dmr.size(); + assert(nspin > 0); + for (int ispin = 0; ispin < nspin; ispin++) { const int nbasis = dmr[ispin]->get_nbasis(); + // gather the parallel matrix to serial matrix #ifdef __MPI Parallel_Orbitals serialV; @@ -86,8 +92,11 @@ void write_dmr(const std::vector*> dmr, #endif if (GlobalV::MY_RANK == 0) { - std::string fname = PARAM.globalv.global_out_dir + dmr_gen_fname(1, ispin, append, istep); - write_dmr_csr(fname, &dm_serial, istep); + // out_type = 1, csr format; + // out_type = 2, npz format (currently not support) + const int out_type = 1; + std::string fname = PARAM.globalv.global_out_dir + dmr_gen_fname(out_type, ispin, append, istep); + write_dmr_csr(fname, &dm_serial, precision, istep); } } } diff --git a/source/source_io/write_dmr.h b/source/source_io/write_dmr.h index 2619795706..4cfb4ef6af 100644 --- a/source/source_io/write_dmr.h +++ b/source/source_io/write_dmr.h @@ -24,15 +24,17 @@ std::string dmr_gen_fname(const int out_type, const int ispin, const bool append * Writes HContainer to a csr file. * * @param fname The name of the file to write the CSR representation to. + * @param precision Control the ouptut precision * @param dm_serial A pointer to the Hamiltonian container. * @param istep The current step number. */ -void write_dmr_csr(std::string& fname, hamilt::HContainer* dm_serial, const int istep); +void write_dmr_csr(std::string& fname, const int precision, hamilt::HContainer* dm_serial, const int istep); /** * Writes DMR to a file. * * @param dmr The 2D block parallel matrix representing the density matrix. The first dimension is the spin index. + * @param precision Control the output precision * @param paraV The parallel 2D object. * @param out_type The output file type. 1: csr, 2: npz. * @param sparse Whether output the sparse DM. @@ -42,6 +44,7 @@ void write_dmr_csr(std::string& fname, hamilt::HContainer* dm_serial, co * @param istep The ION step, starting from 0. */ void write_dmr(const std::vector*> dmr, + const int precision, const Parallel_2D& paraV, const bool append, const int* iat2iwt, diff --git a/source/source_lcao/module_hcontainer/output_hcontainer.cpp b/source/source_lcao/module_hcontainer/output_hcontainer.cpp index b884969815..10d21bd2ca 100644 --- a/source/source_lcao/module_hcontainer/output_hcontainer.cpp +++ b/source/source_lcao/module_hcontainer/output_hcontainer.cpp @@ -32,7 +32,9 @@ template void Output_HContainer::write(bool write_empty) { int size_for_loop_R = this->_hcontainer->size_R_loop(); - int rx, ry, rz; + int rx=0; + int ry=0; + int rz=0; int R_range[2] = {0, 0}; // find the range of R for (int iR = 0; iR < size_for_loop_R; iR++) @@ -135,4 +137,4 @@ void Output_HContainer::write_single_R(int rx, int ry, int rz) template class Output_HContainer; template class Output_HContainer>; -} // namespace hamilt \ No newline at end of file +} // namespace hamilt diff --git a/source/source_lcao/pulay_force_stress_center2.cpp b/source/source_lcao/pulay_force_stress_center2.cpp index 37f75cd5f0..9929285665 100644 --- a/source/source_lcao/pulay_force_stress_center2.cpp +++ b/source/source_lcao/pulay_force_stress_center2.cpp @@ -16,8 +16,8 @@ namespace PulayForceStress const double& factor_force, const double& factor_stress) { - ModuleBase::TITLE("Forces", "cal_pulay"); - ModuleBase::timer::tick("Forces", "cal_pulay"); + ModuleBase::TITLE("Forces", "cal_pulay_fs"); + ModuleBase::timer::tick("Forces", "cal_pulay_fs"); const int nspin = PARAM.inp.nspin; const int nlocal = PARAM.globalv.nlocal; @@ -50,9 +50,12 @@ namespace PulayForceStress } } - if (isstress) { StressTools::stress_fill(ucell.lat0, ucell.omega, s); } + if (isstress) + { + StressTools::stress_fill(ucell.lat0, ucell.omega, s); + } - ModuleBase::timer::tick("Forces", "cal_pulay"); + ModuleBase::timer::tick("Forces", "cal_pulay_fs"); } template<> //multi-k, provided xy diff --git a/source/source_main/driver_run.cpp b/source/source_main/driver_run.cpp index 8300927691..4743b95bf8 100644 --- a/source/source_main/driver_run.cpp +++ b/source/source_main/driver_run.cpp @@ -30,7 +30,6 @@ void Driver::driver_run() // this warning should not be here, mohan 2024-05-22 #ifndef __LCAO if (PARAM.inp.basis_type == "lcao_in_pw" || PARAM.inp.basis_type == "lcao") { - ModuleBase::timer::tick("Driver","driver_run"); ModuleBase::WARNING_QUIT("driver", "to use LCAO basis, compile with __LCAO"); } diff --git a/tests/02_NAO_Gamma/018_NO_GO_ODM/INPUT b/tests/02_NAO_Gamma/018_NO_GO_ODM/INPUT index df0ce8c0ca..64c4429c8e 100644 --- a/tests/02_NAO_Gamma/018_NO_GO_ODM/INPUT +++ b/tests/02_NAO_Gamma/018_NO_GO_ODM/INPUT @@ -4,7 +4,7 @@ basis_type lcao calculation scf gamma_only 1 -out_dmk 1 // file name is dms1_nao.txt +out_dmk 1 // file name is dm_nao.txt nbands 4 latname sc ecutwfc 25.0 // Rydberg diff --git a/tests/03_NAO_multik/05_NO_KP_DJ_OC2/INPUT b/tests/03_NAO_multik/05_NO_KP_DJ_OC2/INPUT index 01a4c87d0d..4aec322903 100644 --- a/tests/03_NAO_multik/05_NO_KP_DJ_OC2/INPUT +++ b/tests/03_NAO_multik/05_NO_KP_DJ_OC2/INPUT @@ -26,7 +26,7 @@ basis_type lcao ks_solver scalapack_gvx chg_extrap second-order out_dmk 0 -out_chg 2 3 +out_chg 2 3 // 2 outputs the initial charge density pw_diag_thr 0.00001 mixing_type broyden