From fa7f6abd3591ba7eefcb51c70664902ba5b24b19 Mon Sep 17 00:00:00 2001 From: mohanchen Date: Tue, 23 Sep 2025 20:50:51 +0800 Subject: [PATCH 1/6] update out_freq_ion description --- docs/advanced/input_files/input-main.md | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index 5e58dcbced..f9f3d5e261 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -731,7 +731,7 @@ If only one value is set (such as `kspacing 0.5`), then kspacing values of a/b/c [back to top](#full-list-of-input-keywords) -## Variables related to input files +## Input files These variables are used to control parameters related to input files. @@ -1632,16 +1632,16 @@ These variables are used to control the geometry relaxation. [back to top](#full-list-of-input-keywords) -## Variables related to output information +## Output information These variables are used to control the output of properties. ### out_freq_ion - **Type**: Integer -- **Description**: After self-consistent-field calculations, control the interval of ionic movements for printing properties. These properties cover charge density, local potential, electrostatic potential, Hamiltonian matrix, overlap matrix, density matrix, Mulliken population analysis and so on. +- **Description**: Control the interval to print information every few ion steps. These properties cover charge density, local potential, electrostatic potential, Hamiltonian matrix, overlap matrix, density matrix, Mulliken population analysis and so on. - **Default**: 0 -- **Note**: If you want to use out_freq_elec, please set out_freq_ion to 1, otherwise out_freq_elec is useless +- **Note**: The integer indicates to print information every 'out_freq_ion' ion steps. ### out_freq_elec From 9c7d55380d29b88939872738710de382236fd270 Mon Sep 17 00:00:00 2001 From: mohanchen Date: Tue, 23 Sep 2025 21:08:01 +0800 Subject: [PATCH 2/6] update chg_init and pot_init files --- docs/advanced/input_files/input-main.md | 4 +- source/source_esolver/esolver_fp.cpp | 60 ++++++++++++++++++------- 2 files changed, 45 insertions(+), 19 deletions(-) diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index f9f3d5e261..b5a83d5b53 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -1655,11 +1655,11 @@ These variables are used to control the output of properties. - **Description**: The first integer controls whether to output the charge density on real space grids: - 1: Output the charge density (in Bohr^-3) on real space grids into the density files in the folder `OUT.${suffix}`. The files are named as: - - nspin = 1: `chgs1.cube`; + - nspin = 1: `chg.cube`; - nspin = 2: `chgs1.cube`, and `chgs2.cube`; - nspin = 4: `chgs1.cube`, `chgs2.cube`, `chgs3.cube`, and `chgs4.cube`; Note that by using the Meta-GGA functional, additional files containing the kinetic energy density will be output with the following names: - - nspin = 1: `taus1.cube`; + - nspin = 1: `tau.cube`; - nspin = 2: `taus1.cube`, and `taus2.cube`; - nspin = 4: `taus1.cube`, `taus2.cube`, `taus3.cube`, and `taus4.cube`; - 2: On top of 1, also output the initial charge density files with a suffix name as '_ini', such as `taus1_ini.cube`, etc. diff --git a/source/source_esolver/esolver_fp.cpp b/source/source_esolver/esolver_fp.cpp index 46f00cf284..3594b0e78f 100644 --- a/source/source_esolver/esolver_fp.cpp +++ b/source/source_esolver/esolver_fp.cpp @@ -184,20 +184,23 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso } // 4) write charge density + + const int nspin = PARAM.inp.nspin; + if (PARAM.inp.out_chg[0] > 0) { - for (int is = 0; is < PARAM.inp.nspin; ++is) + for (int is = 0; is < nspin; ++is) { this->pw_rhod->real2recip(this->chr.rho_save[is], this->chr.rhog_save[is]); std::string fn =PARAM.globalv.global_out_dir + "chg"; std::string spin_block; - if(PARAM.inp.nspin == 2 || PARAM.inp.nspin == 4) + if(nspin == 2 || nspin == 4) { spin_block= "s" + std::to_string(is + 1); } - else if(PARAM.inp.nspin == 1) + else if(nspin == 1) { // do nothing } @@ -207,7 +210,7 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso ModuleIO::write_vdata_palgrid(Pgrid, this->chr.rho_save[is], is, - PARAM.inp.nspin, + nspin, istep_in, fn, this->pelec->eferm.get_efval(is), @@ -224,7 +227,7 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso ModuleIO::write_vdata_palgrid(Pgrid, this->chr.kin_r_save[is], is, - PARAM.inp.nspin, + nspin, istep, fn, this->pelec->eferm.get_efval(is), @@ -236,16 +239,16 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso // 5) write potential if (PARAM.inp.out_pot == 1 || PARAM.inp.out_pot == 3) { - for (int is = 0; is < PARAM.inp.nspin; is++) + for (int is = 0; is < nspin; is++) { std::string fn =PARAM.globalv.global_out_dir + "pot"; std::string spin_block; - if(PARAM.inp.nspin == 2 || PARAM.inp.nspin == 4) + if(nspin == 2 || nspin == 4) { spin_block= "s" + std::to_string(is + 1); } - else if(PARAM.inp.nspin == 1) + else if(nspin == 1) { // do nothing } @@ -255,7 +258,7 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso ModuleIO::write_vdata_palgrid(Pgrid, this->pelec->pot->get_effective_v(is), is, - PARAM.inp.nspin, + nspin, istep_in, fn, 0.0, // efermi @@ -288,7 +291,7 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso { this->chr.cal_elf = true; Symmetry_rho srho; - for (int is = 0; is < PARAM.inp.nspin; is++) + for (int is = 0; is < nspin; is++) { srho.begin(is, this->chr, this->pw_rhod, ucell.symm); } @@ -301,7 +304,7 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso #endif out_dir, istep, - PARAM.inp.nspin, + nspin, this->chr.rho, this->chr.kin_r, this->pw_rhod, @@ -402,19 +405,32 @@ void ESolver_FP::before_scf(UnitCell& ucell, const int istep) //---------------------------------------------------------- elecstate::cal_ux(ucell); + + //---------------------------------------------------------- //! output the initial charge density //---------------------------------------------------------- + const int nspin = PARAM.inp.nspin; if (PARAM.inp.out_chg[0] == 2) { - for (int is = 0; is < PARAM.inp.nspin; is++) + for (int is = 0; is < nspin; is++) { std::stringstream ss; - ss << PARAM.globalv.global_out_dir << "chgs" << is + 1 << "_ini.cube"; + ss << PARAM.globalv.global_out_dir << "chg"; + + if(nspin==1) + { + ss << "ini.cube"; + } + else if(nspin==2 || nspin==4) + { + ss << "s" << is + 1 << "ini.cube"; + } + ModuleIO::write_vdata_palgrid(this->Pgrid, this->chr.rho[is], is, - PARAM.inp.nspin, + nspin, istep, ss.str(), this->pelec->eferm.ef, @@ -427,14 +443,24 @@ void ESolver_FP::before_scf(UnitCell& ucell, const int istep) //---------------------------------------------------------- if (PARAM.inp.out_pot == 3) { - for (int is = 0; is < PARAM.inp.nspin; is++) + for (int is = 0; is < nspin; is++) { std::stringstream ss; - ss << PARAM.globalv.global_out_dir << "pots" << is + 1 << "_ini.cube"; + ss << PARAM.globalv.global_out_dir << "pot"; + + if(nspin==1) + { + ss << "ini.cube"; + } + else if(nspin==2 || nspin==4) + { + ss << "s" << is + 1 << "ini.cube"; + } + ModuleIO::write_vdata_palgrid(this->Pgrid, this->pelec->pot->get_effective_v(is), is, - PARAM.inp.nspin, + nspin, istep, ss.str(), 0.0, // efermi From ebc48d95adfeea8fdee8517d3b937cb4e62aed13 Mon Sep 17 00:00:00 2001 From: mohanchen Date: Tue, 23 Sep 2025 21:49:50 +0800 Subject: [PATCH 3/6] cannot run, but we have include ctrl_output_fp files --- source/Makefile.Objects | 1 + source/source_esolver/esolver_fp.cpp | 174 +---------------- source/source_esolver/esolver_fp.h | 32 +-- source/source_estate/elecstate.h | 5 - source/source_io/CMakeLists.txt | 1 + source/source_io/ctrl_output_fp.cpp | 279 ++++++++++++++++----------- source/source_io/ctrl_output_fp.h | 9 +- 7 files changed, 184 insertions(+), 317 deletions(-) diff --git a/source/Makefile.Objects b/source/Makefile.Objects index ed356b4329..db5585dd3a 100644 --- a/source/Makefile.Objects +++ b/source/Makefile.Objects @@ -579,6 +579,7 @@ OBJS_IO=input_conv.o\ output_log.o\ output_mat_sparse.o\ ctrl_output_lcao.o\ + ctrl_output_fp.o\ para_json.o\ abacusjson.o\ general_info.o\ diff --git a/source/source_esolver/esolver_fp.cpp b/source/source_esolver/esolver_fp.cpp index 3594b0e78f..a0371e6616 100644 --- a/source/source_esolver/esolver_fp.cpp +++ b/source/source_esolver/esolver_fp.cpp @@ -14,14 +14,9 @@ #include "source_io/output_log.h" #include "source_io/print_info.h" #include "source_io/rhog_io.h" -#include "source_io/write_elecstat_pot.h" -#include "source_io/write_elf.h" #include "source_io/module_parameter/parameter.h" #include "source_cell/k_vector_utils.h" - -#ifdef USE_LIBXC -#include "source_io/write_libxc_r.h" -#endif +#include "source_io/ctrl_output_fp.h" namespace ModuleESolver { @@ -162,172 +157,9 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso // 3) update delta_rho for charge extrapolation CE.update_delta_rho(ucell, &(this->chr), &(this->sf)); + // 4) print out charge density, potential, elf, etc. + ModuleIO::ctrl_output_fp(ucell, this->pelec, istep); - // print out the 'g' index when istep_in!=-1 - int istep_in = -1; - if (PARAM.inp.out_freq_ion>0) // default value of out_freq_ion is 0 - { - if (istep % PARAM.inp.out_freq_ion == 0) - { - istep_in = istep; - } - } - - std::string geom_block; - if(istep_in==-1) - { - // do nothing - } - else if(istep_in>=0) - { - geom_block = "g" + std::to_string(istep + 1); - } - - // 4) write charge density - - const int nspin = PARAM.inp.nspin; - - if (PARAM.inp.out_chg[0] > 0) - { - for (int is = 0; is < nspin; ++is) - { - this->pw_rhod->real2recip(this->chr.rho_save[is], this->chr.rhog_save[is]); - - std::string fn =PARAM.globalv.global_out_dir + "chg"; - - std::string spin_block; - if(nspin == 2 || nspin == 4) - { - spin_block= "s" + std::to_string(is + 1); - } - else if(nspin == 1) - { - // do nothing - } - - fn += spin_block + geom_block + ".cube"; - - ModuleIO::write_vdata_palgrid(Pgrid, - this->chr.rho_save[is], - is, - nspin, - istep_in, - fn, - this->pelec->eferm.get_efval(is), - &(ucell), - PARAM.inp.out_chg[1], - 1); - - if (XC_Functional::get_ked_flag()) - { - fn = PARAM.globalv.global_out_dir + "tau"; - - fn += spin_block + geom_block + ".cube"; - - ModuleIO::write_vdata_palgrid(Pgrid, - this->chr.kin_r_save[is], - is, - nspin, - istep, - fn, - this->pelec->eferm.get_efval(is), - &(ucell)); - } - } - } - - // 5) write potential - if (PARAM.inp.out_pot == 1 || PARAM.inp.out_pot == 3) - { - for (int is = 0; is < nspin; is++) - { - std::string fn =PARAM.globalv.global_out_dir + "pot"; - - std::string spin_block; - if(nspin == 2 || nspin == 4) - { - spin_block= "s" + std::to_string(is + 1); - } - else if(nspin == 1) - { - // do nothing - } - - fn += spin_block + geom_block + ".cube"; - - ModuleIO::write_vdata_palgrid(Pgrid, - this->pelec->pot->get_effective_v(is), - is, - nspin, - istep_in, - fn, - 0.0, // efermi - &(ucell), - 3, // precision - 0); // out_fermi - } - } - else if (PARAM.inp.out_pot == 2) - { - std::string fn =PARAM.globalv.global_out_dir + "potes"; - fn += geom_block + ".cube"; - - ModuleIO::write_elecstat_pot( -#ifdef __MPI - this->pw_big->bz, - this->pw_big->nbz, -#endif - fn, - istep, - this->pw_rhod, - &this->chr, - &(ucell), - this->pelec->pot->get_fixed_v(), - this->solvent); - } - - // 6) write ELF - if (PARAM.inp.out_elf[0] > 0) - { - this->chr.cal_elf = true; - Symmetry_rho srho; - for (int is = 0; is < nspin; is++) - { - srho.begin(is, this->chr, this->pw_rhod, ucell.symm); - } - - std::string out_dir =PARAM.globalv.global_out_dir; - ModuleIO::write_elf( -#ifdef __MPI - this->pw_big->bz, - this->pw_big->nbz, -#endif - out_dir, - istep, - nspin, - this->chr.rho, - this->chr.kin_r, - this->pw_rhod, - this->Pgrid, - &(ucell), - PARAM.inp.out_elf[1]); - } - -#ifdef USE_LIBXC - // 7) write xc(r) - if(PARAM.inp.out_xc_r[0]>=0) - { - ModuleIO::write_libxc_r( - PARAM.inp.out_xc_r[0], - XC_Functional::get_func_id(), - this->pw_rhod->nrxx, // number of real-space grid - ucell.omega, // volume of cell - ucell.tpiba, - this->chr, - *this->pw_big, - *this->pw_rhod); - } -#endif } void ESolver_FP::before_scf(UnitCell& ucell, const int istep) diff --git a/source/source_esolver/esolver_fp.h b/source/source_esolver/esolver_fp.h index fda47313bc..7d4f34736e 100644 --- a/source/source_esolver/esolver_fp.h +++ b/source/source_esolver/esolver_fp.h @@ -7,26 +7,13 @@ #include #endif -//! plane wave basis -#include "source_basis/module_pw/pw_basis.h" - -//! symmetry analysis -#include "source_cell/module_symmetry/symmetry.h" - -//! electronic states -#include "source_estate/elecstate.h" - -//! charge extrapolation -#include "source_estate/module_charge/charge_extra.h" - -//! solvation model -#include "source_hamilt/module_surchem/surchem.h" - -//! local pseudopotential -#include "source_pw/module_pwdft/VL_in_pw.h" - -//! structure factor related to plane wave basis -#include "source_pw/module_pwdft/structure_factor.h" +#include "source_basis/module_pw/pw_basis.h" // plane wave basis +#include "source_cell/module_symmetry/symmetry.h" // symmetry analysis +#include "source_estate/elecstate.h" // electronic states +#include "source_estate/module_charge/charge_extra.h" // charge extrapolation +#include "source_hamilt/module_surchem/surchem.h" // solvation model +#include "source_pw/module_pwdft/VL_in_pw.h" // local pseudopotential +#include "source_pw/module_pwdft/structure_factor.h" // structure factor #include @@ -44,10 +31,8 @@ namespace ModuleESolver class ESolver_FP: public ESolver { public: - //! Constructor ESolver_FP(); - //! Deconstructor virtual ~ESolver_FP(); //! Initialize of the first-principels energy solver @@ -56,13 +41,10 @@ class ESolver_FP: public ESolver virtual void after_all_runners(UnitCell& ucell) override; protected: - //! Something to do before SCF iterations. virtual void before_scf(UnitCell& ucell, const int istep); - //! Something to do after SCF iterations when SCF is converged or comes to the max iter step. virtual void after_scf(UnitCell& ucell, const int istep, const bool conv_esolver); - //! Something to do after hamilt2rho function in each iter loop. virtual void iter_finish(UnitCell& ucell, const int istep, int& iter, bool &conv_esolver); //! ------------------------------------------------------------------------------ diff --git a/source/source_estate/elecstate.h b/source/source_estate/elecstate.h index bdf5f9cc77..7b1ff8fbcc 100644 --- a/source/source_estate/elecstate.h +++ b/source/source_estate/elecstate.h @@ -51,8 +51,6 @@ class ElecState { return; } - // virtual void updateRhoK(const psi::Psi> &psi) = 0; - // virtual void updateRhoK(const psi::Psi &psi)=0 virtual void cal_tau(const psi::Psi>& psi) { return; @@ -77,12 +75,9 @@ class ElecState return; } - - // use occupied weights from INPUT and skip calculate_weights // mohan updated on 2024-06-08 - // if nupdown is not 0(TWO_EFERMI case), // nelec_spin will be fixed and weights will be constrained void init_nelec_spin(); diff --git a/source/source_io/CMakeLists.txt b/source/source_io/CMakeLists.txt index 6c7b5bf349..582b7539d1 100644 --- a/source/source_io/CMakeLists.txt +++ b/source/source_io/CMakeLists.txt @@ -82,6 +82,7 @@ if(ENABLE_LCAO) cal_r_overlap_R.cpp output_mat_sparse.cpp ctrl_output_lcao.cpp + ctrl_output_fp.cpp ) endif() diff --git a/source/source_io/ctrl_output_fp.cpp b/source/source_io/ctrl_output_fp.cpp index 477612ed27..f2cc168b90 100644 --- a/source/source_io/ctrl_output_fp.cpp +++ b/source/source_io/ctrl_output_fp.cpp @@ -1,14 +1,18 @@ #include "source_io/ctrl_output_fp.h" // use ctrl_output_fp() +#include "source_io/write_elecstat_pot.h" +#include "source_io/write_elf.h" + +#ifdef USE_LIBXC +#include "source_io/write_libxc_r.h" +#endif namespace ModuleIO { -template void ctrl_output_fp(UnitCell& ucell, - elecstate::ElecStateLCAO* pelec, + elecstate::ElecState* pelec, const int istep) { -/* ModuleBase::TITLE("ModuleIO", "ctrl_output_fp"); ModuleBase::timer::tick("ModuleIO", "ctrl_output_fp"); @@ -17,121 +21,172 @@ void ctrl_output_fp(UnitCell& ucell, const int nspin = PARAM.inp.nspin; const std::string global_out_dir = PARAM.globalv.global_out_dir; - // 1) write charge density - if (PARAM.inp.out_chg[0] > 0) - { - for (int is = 0; is < PARAM.inp.nspin; is++) - { - this->pw_rhod->real2recip(this->chr.rho_save[is], this->chr.rhog_save[is]); - std::string fn =PARAM.globalv.global_out_dir + "/chgs" + std::to_string(is + 1) + ".cube"; - ModuleIO::write_vdata_palgrid(Pgrid, - this->chr.rho_save[is], - is, - PARAM.inp.nspin, - istep, - fn, - this->pelec->eferm.get_efval(is), - &(ucell), - PARAM.inp.out_chg[1], - 1); - - if (XC_Functional::get_ked_flag()) - { - fn =PARAM.globalv.global_out_dir + "/taus" + std::to_string(is + 1) + ".cube"; - ModuleIO::write_vdata_palgrid(Pgrid, - this->chr.kin_r_save[is], - is, - PARAM.inp.nspin, - istep, - fn, - this->pelec->eferm.get_efval(is), - &(ucell)); - } - } - } - - - // 2) write potential - if (PARAM.inp.out_pot == 1 || PARAM.inp.out_pot == 3) - { - for (int is = 0; is < PARAM.inp.nspin; is++) - { - std::string fn =PARAM.globalv.global_out_dir + "/pots" + std::to_string(is + 1) + ".cube"; - - ModuleIO::write_vdata_palgrid(Pgrid, - this->pelec->pot->get_effective_v(is), - is, - PARAM.inp.nspin, - istep, - fn, - 0.0, // efermi - &(ucell), - 3, // precision - 0); // out_fermi - } - } - else if (PARAM.inp.out_pot == 2) - { - std::string fn =PARAM.globalv.global_out_dir + "/pot_es.cube"; - ModuleIO::write_elecstat_pot( + + // print out the 'g' index when istep_in!=-1 + int istep_in = -1; + if (PARAM.inp.out_freq_ion>0) // default value of out_freq_ion is 0 + { + if (istep % PARAM.inp.out_freq_ion == 0) + { + istep_in = istep; + } + } + + std::string geom_block; + if(istep_in==-1) + { + // do nothing + } + else if(istep_in>=0) + { + geom_block = "g" + std::to_string(istep + 1); + } + + + // 4) write charge density + if (PARAM.inp.out_chg[0] > 0) + { + for (int is = 0; is < nspin; ++is) + { + this->pw_rhod->real2recip(this->chr.rho_save[is], this->chr.rhog_save[is]); + + std::string fn =PARAM.globalv.global_out_dir + "chg"; + + std::string spin_block; + if(nspin == 2 || nspin == 4) + { + spin_block= "s" + std::to_string(is + 1); + } + else if(nspin == 1) + { + // do nothing + } + + fn += spin_block + geom_block + ".cube"; + + ModuleIO::write_vdata_palgrid(Pgrid, + this->chr.rho_save[is], + is, + nspin, + istep_in, + fn, + pelec->eferm.get_efval(is), + &(ucell), + PARAM.inp.out_chg[1], + 1); + + if (XC_Functional::get_ked_flag()) + { + fn = PARAM.globalv.global_out_dir + "tau"; + + fn += spin_block + geom_block + ".cube"; + + ModuleIO::write_vdata_palgrid(Pgrid, + this->chr.kin_r_save[is], + is, + nspin, + istep, + fn, + pelec->eferm.get_efval(is), + &(ucell)); + } + } + } + + // 5) write potential + if (PARAM.inp.out_pot == 1 || PARAM.inp.out_pot == 3) + { + for (int is = 0; is < nspin; is++) + { + std::string fn =PARAM.globalv.global_out_dir + "pot"; + + std::string spin_block; + if(nspin == 2 || nspin == 4) + { + spin_block= "s" + std::to_string(is + 1); + } + else if(nspin == 1) + { + // do nothing + } + + fn += spin_block + geom_block + ".cube"; + + ModuleIO::write_vdata_palgrid(Pgrid, + pelec->pot->get_effective_v(is), + is, + nspin, + istep_in, + fn, + 0.0, // efermi + &(ucell), + 3, // precision + 0); // out_fermi + } + } + else if (PARAM.inp.out_pot == 2) + { + std::string fn =PARAM.globalv.global_out_dir + "potes"; + fn += geom_block + ".cube"; + + ModuleIO::write_elecstat_pot( #ifdef __MPI - this->pw_big->bz, - this->pw_big->nbz, + this->pw_big->bz, + this->pw_big->nbz, #endif - fn, - istep, - this->pw_rhod, - &this->chr, - &(ucell), - this->pelec->pot->get_fixed_v(), - this->solvent); - } - - - // 3) write ELF - if (PARAM.inp.out_elf[0] > 0) - { - this->chr.cal_elf = true; - Symmetry_rho srho; - for (int is = 0; is < PARAM.inp.nspin; is++) - { - srho.begin(is, this->chr, this->pw_rhod, ucell.symm); - } - - std::string out_dir =PARAM.globalv.global_out_dir; - ModuleIO::write_elf( + fn, + istep, + this->pw_rhod, + &this->chr, + &(ucell), + pelec->pot->get_fixed_v(), + this->solvent); + } + + // 6) write ELF + if (PARAM.inp.out_elf[0] > 0) + { + this->chr.cal_elf = true; + Symmetry_rho srho; + for (int is = 0; is < nspin; is++) + { + srho.begin(is, this->chr, this->pw_rhod, ucell.symm); + } + + std::string out_dir =PARAM.globalv.global_out_dir; + ModuleIO::write_elf( #ifdef __MPI - this->pw_big->bz, - this->pw_big->nbz, + this->pw_big->bz, + this->pw_big->nbz, +#endif + out_dir, + istep, + nspin, + this->chr.rho, + this->chr.kin_r, + this->pw_rhod, + this->Pgrid, + &(ucell), + PARAM.inp.out_elf[1]); + } + +#ifdef USE_LIBXC + // 7) write xc(r) + if(PARAM.inp.out_xc_r[0]>=0) + { + ModuleIO::write_libxc_r( + PARAM.inp.out_xc_r[0], + XC_Functional::get_func_id(), + this->pw_rhod->nrxx, // number of real-space grid + ucell.omega, // volume of cell + ucell.tpiba, + this->chr, + *this->pw_big, + *this->pw_rhod); + } #endif - out_dir, - istep, - PARAM.inp.nspin, - this->chr.rho, - this->chr.kin_r, - this->pw_rhod, - this->Pgrid, - &(ucell), - PARAM.inp.out_elf[1]); - } ModuleBase::timer::tick("ModuleIO", "ctrl_output_fp"); - -*/ - } } // End ModuleIO - - -// For gamma only -template void ModuleIO::ctrl_output_lcao(UnitCell& ucell, - const int istep); - -// For multiple k-points -template void ModuleIO::ctrl_output_lcao, double>(UnitCell& ucell, - const int istep); - -template void ModuleIO::ctrl_output_lcao, std::complex>(UnitCell& ucell, - const int istep); - diff --git a/source/source_io/ctrl_output_fp.h b/source/source_io/ctrl_output_fp.h index 98a4ad1f6d..d2686ebc7b 100644 --- a/source/source_io/ctrl_output_fp.h +++ b/source/source_io/ctrl_output_fp.h @@ -1,11 +1,12 @@ #ifndef CTRL_OUTPUT_FP_H #define CTRL_OUTPUT_FP_H +#include "source_estate/elecstate_lcao.h" + namespace ModuleIO { - template - void ctrl_output_fp(UnitCell& ucell, - elecstate::ElecStateLCAO* pelec, - const int istep); + void ctrl_output_fp(UnitCell& ucell, + elecstate::ElecState* pelec, + const int istep); } #endif From bfd1cfa33d95bb5e28e8556d90035d64f4123ea4 Mon Sep 17 00:00:00 2001 From: mohanchen Date: Tue, 23 Sep 2025 22:12:56 +0800 Subject: [PATCH 4/6] move some files from esolver_fp to ctrl_output_fp --- source/source_esolver/esolver_fp.cpp | 5 ++- source/source_esolver/esolver_of.cpp | 2 - source/source_io/ctrl_output_fp.cpp | 56 ++++++++++++++++------------ source/source_io/ctrl_output_fp.h | 9 ++++- 4 files changed, 43 insertions(+), 29 deletions(-) diff --git a/source/source_esolver/esolver_fp.cpp b/source/source_esolver/esolver_fp.cpp index a0371e6616..9d8326d2a3 100644 --- a/source/source_esolver/esolver_fp.cpp +++ b/source/source_esolver/esolver_fp.cpp @@ -8,7 +8,7 @@ #include "source_hamilt/module_vdw/vdw.h" #include "source_pw/module_pwdft/global.h" #include "source_io/cif_io.h" -#include "source_io/cube_io.h" +#include "source_io/cube_io.h" // use write_vdata_palgrid #include "source_io/json_output/init_info.h" #include "source_io/json_output/output_info.h" #include "source_io/output_log.h" @@ -158,7 +158,8 @@ void ESolver_FP::after_scf(UnitCell& ucell, const int istep, const bool conv_eso CE.update_delta_rho(ucell, &(this->chr), &(this->sf)); // 4) print out charge density, potential, elf, etc. - ModuleIO::ctrl_output_fp(ucell, this->pelec, istep); + ModuleIO::ctrl_output_fp(ucell, this->pelec, this->pw_big, this->pw_rhod, + this->chr, this->solvent, this->Pgrid, istep); } diff --git a/source/source_esolver/esolver_of.cpp b/source/source_esolver/esolver_of.cpp index 1ece91eb62..45f26babcb 100644 --- a/source/source_esolver/esolver_of.cpp +++ b/source/source_esolver/esolver_of.cpp @@ -11,9 +11,7 @@ #include "source_pw/module_pwdft/global.h" #include "source_io/print_info.h" #include "source_estate/cal_ux.h" -//-----force------------------- #include "source_pw/module_pwdft/forces.h" -//-----stress------------------ #include "source_pw/module_ofdft/of_stress_pw.h" namespace ModuleESolver diff --git a/source/source_io/ctrl_output_fp.cpp b/source/source_io/ctrl_output_fp.cpp index f2cc168b90..d79e4fbc4c 100644 --- a/source/source_io/ctrl_output_fp.cpp +++ b/source/source_io/ctrl_output_fp.cpp @@ -1,9 +1,12 @@ #include "source_io/ctrl_output_fp.h" // use ctrl_output_fp() -#include "source_io/write_elecstat_pot.h" +#include "source_estate/module_charge/symmetry_rho.h" // use Symmetry_rho +#include "source_io/write_elecstat_pot.h" // use write_elecstat_pot #include "source_io/write_elf.h" +#include "cube_io.h" // use write_vdata_palgrid #ifdef USE_LIBXC #include "source_io/write_libxc_r.h" +#include "source_hamilt/module_xc/xc_functional.h" // use XC_Functional::get_func_id() #endif namespace ModuleIO @@ -11,6 +14,11 @@ namespace ModuleIO void ctrl_output_fp(UnitCell& ucell, elecstate::ElecState* pelec, + ModulePW::PW_Basis_Big* pw_big, + ModulePW::PW_Basis* pw_rhod, + Charge &chr, + surchem &solvent, + Parallel_Grid ¶_grid, const int istep) { ModuleBase::TITLE("ModuleIO", "ctrl_output_fp"); @@ -48,7 +56,7 @@ void ctrl_output_fp(UnitCell& ucell, { for (int is = 0; is < nspin; ++is) { - this->pw_rhod->real2recip(this->chr.rho_save[is], this->chr.rhog_save[is]); + pw_rhod->real2recip(chr.rho_save[is], chr.rhog_save[is]); std::string fn =PARAM.globalv.global_out_dir + "chg"; @@ -64,8 +72,8 @@ void ctrl_output_fp(UnitCell& ucell, fn += spin_block + geom_block + ".cube"; - ModuleIO::write_vdata_palgrid(Pgrid, - this->chr.rho_save[is], + ModuleIO::write_vdata_palgrid(para_grid, + chr.rho_save[is], is, nspin, istep_in, @@ -81,8 +89,8 @@ void ctrl_output_fp(UnitCell& ucell, fn += spin_block + geom_block + ".cube"; - ModuleIO::write_vdata_palgrid(Pgrid, - this->chr.kin_r_save[is], + ModuleIO::write_vdata_palgrid(para_grid, + chr.kin_r_save[is], is, nspin, istep, @@ -112,7 +120,7 @@ void ctrl_output_fp(UnitCell& ucell, fn += spin_block + geom_block + ".cube"; - ModuleIO::write_vdata_palgrid(Pgrid, + ModuleIO::write_vdata_palgrid(para_grid, pelec->pot->get_effective_v(is), is, nspin, @@ -131,41 +139,41 @@ void ctrl_output_fp(UnitCell& ucell, ModuleIO::write_elecstat_pot( #ifdef __MPI - this->pw_big->bz, - this->pw_big->nbz, + pw_big->bz, + pw_big->nbz, #endif fn, istep, - this->pw_rhod, - &this->chr, + pw_rhod, + &chr, &(ucell), pelec->pot->get_fixed_v(), - this->solvent); + solvent); } // 6) write ELF if (PARAM.inp.out_elf[0] > 0) { - this->chr.cal_elf = true; + chr.cal_elf = true; Symmetry_rho srho; for (int is = 0; is < nspin; is++) { - srho.begin(is, this->chr, this->pw_rhod, ucell.symm); + srho.begin(is, chr, pw_rhod, ucell.symm); } std::string out_dir =PARAM.globalv.global_out_dir; ModuleIO::write_elf( #ifdef __MPI - this->pw_big->bz, - this->pw_big->nbz, + pw_big->bz, + pw_big->nbz, #endif out_dir, istep, nspin, - this->chr.rho, - this->chr.kin_r, - this->pw_rhod, - this->Pgrid, + chr.rho, + chr.kin_r, + pw_rhod, + para_grid, &(ucell), PARAM.inp.out_elf[1]); } @@ -177,12 +185,12 @@ void ctrl_output_fp(UnitCell& ucell, ModuleIO::write_libxc_r( PARAM.inp.out_xc_r[0], XC_Functional::get_func_id(), - this->pw_rhod->nrxx, // number of real-space grid + pw_rhod->nrxx, // number of real-space grid ucell.omega, // volume of cell ucell.tpiba, - this->chr, - *this->pw_big, - *this->pw_rhod); + chr, + *pw_big, + *pw_rhod); } #endif diff --git a/source/source_io/ctrl_output_fp.h b/source/source_io/ctrl_output_fp.h index d2686ebc7b..12ff1e12b0 100644 --- a/source/source_io/ctrl_output_fp.h +++ b/source/source_io/ctrl_output_fp.h @@ -5,8 +5,15 @@ namespace ModuleIO { + void ctrl_output_fp(UnitCell& ucell, - elecstate::ElecState* pelec, + elecstate::ElecState* pelec, + ModulePW::PW_Basis_Big* pw_big, + ModulePW::PW_Basis* pw_rhod, + Charge &chr, + surchem &solvent, + Parallel_Grid ¶_grid, const int istep); + } #endif From a9fda59b744b328e159afbed0a2e6bfe5755d09f Mon Sep 17 00:00:00 2001 From: mohanchen Date: Wed, 24 Sep 2025 15:27:10 +0800 Subject: [PATCH 5/6] fix a bug --- source/source_io/ctrl_output_fp.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/source/source_io/ctrl_output_fp.cpp b/source/source_io/ctrl_output_fp.cpp index d79e4fbc4c..0cde6e6011 100644 --- a/source/source_io/ctrl_output_fp.cpp +++ b/source/source_io/ctrl_output_fp.cpp @@ -3,10 +3,10 @@ #include "source_io/write_elecstat_pot.h" // use write_elecstat_pot #include "source_io/write_elf.h" #include "cube_io.h" // use write_vdata_palgrid +#include "source_hamilt/module_xc/xc_functional.h" // use XC_Functional #ifdef USE_LIBXC #include "source_io/write_libxc_r.h" -#include "source_hamilt/module_xc/xc_functional.h" // use XC_Functional::get_func_id() #endif namespace ModuleIO From 53fdb2db739cd665e1cd094033ca5dd07276194d Mon Sep 17 00:00:00 2001 From: mohanchen Date: Wed, 24 Sep 2025 22:25:48 +0800 Subject: [PATCH 6/6] move ctrl_output_fp to pw part --- source/source_io/CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/source/source_io/CMakeLists.txt b/source/source_io/CMakeLists.txt index 582b7539d1..ffa7598cd9 100644 --- a/source/source_io/CMakeLists.txt +++ b/source/source_io/CMakeLists.txt @@ -1,5 +1,6 @@ list(APPEND objects input_conv.cpp + ctrl_output_fp.cpp bessel_basis.cpp cal_test.cpp cal_dos.cpp @@ -82,7 +83,6 @@ if(ENABLE_LCAO) cal_r_overlap_R.cpp output_mat_sparse.cpp ctrl_output_lcao.cpp - ctrl_output_fp.cpp ) endif()