Skip to content

Commit

Permalink
Refactor IO cube (#5362)
Browse files Browse the repository at this point in the history
* remove the redundant info in read_cube

* remove redundant code in write_cube

* make Pgrid a serial-or-parallel grid class

* Pgrid contains nxyz

* [pre-commit.ci lite] apply automatic fixes

* rename reduce

* rewrite io_cube

* rename previous read/write_cube as read/write_grid

* remove CubeInfo and xyz-zxy permutation

* minor interface changes

* revert the call in rho_sym to fix SDFT failure

* rename r/w_grid as r/w_vdata_palgrid

* [pre-commit.ci lite] apply automatic fixes

* add rank parameter to pgrid.bcast to fix SDFT

---------

Co-authored-by: pre-commit-ci-lite[bot] <117423508+pre-commit-ci-lite[bot]@users.noreply.github.com>
  • Loading branch information
maki49 and pre-commit-ci-lite[bot] authored Oct 31, 2024
1 parent 4e53045 commit 705dad7
Show file tree
Hide file tree
Showing 21 changed files with 424 additions and 898 deletions.
78 changes: 22 additions & 56 deletions source/module_elecstate/module_charge/charge_init.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -53,61 +53,37 @@ void Charge::init_rho(elecstate::efermi& eferm_iout,
{
std::stringstream ssc;
ssc << PARAM.globalv.global_readin_dir << "SPIN" << is + 1 << "_CHG.cube";
double& ef_tmp = eferm_iout.get_ef(is);
if (ModuleIO::read_cube(
#ifdef __MPI
& (GlobalC::Pgrid),
#endif
if (ModuleIO::read_vdata_palgrid(GlobalC::Pgrid,
(PARAM.inp.esolver_type == "sdft" ? GlobalV::RANK_IN_STOGROUP : GlobalV::MY_RANK),
is,
GlobalV::ofs_running,
PARAM.inp.nspin,
ssc.str(),
this->rho[is],
this->rhopw->nx,
this->rhopw->ny,
this->rhopw->nz,
ef_tmp,
& (GlobalC::ucell),
this->prenspin))
GlobalC::ucell.nat))
{
GlobalV::ofs_running << " Read in the charge density: " << ssc.str() << std::endl;
}
else if (is > 0)
else if (is > 0) // nspin=2 or 4
{
if (prenspin == 1)
if (is == 1) // failed at the second spin
{
GlobalV::ofs_running << " Didn't read in the charge density but autoset it for spin " << is + 1
<< std::endl;
for (int ir = 0; ir < this->rhopw->nrxx; ir++)
{
this->rho[is][ir] = 0.0;
}
std::cout << "Incomplete charge density file!" << std::endl;
read_error = true;
break;
}
//
else if (prenspin == 2)
{ // read up and down , then rearrange them.
if (is == 1)
{
std::cout << "Incomplete charge density file!" << std::endl;
read_error = true;
break;
}
else if (is == 2)
{
GlobalV::ofs_running << " Didn't read in the charge density but would rearrange it later. "
<< std::endl;
}
else if (is == 3)
else if (is == 2) // read 2 files when nspin=4
{
GlobalV::ofs_running << " Didn't read in the charge density but would rearrange it later. "
<< std::endl;
}
else if (is == 3) // read 2 files when nspin=4
{
GlobalV::ofs_running << " rearrange charge density " << std::endl;
for (int ir = 0; ir < this->rhopw->nrxx; ir++)
{
GlobalV::ofs_running << " rearrange charge density " << std::endl;
for (int ir = 0; ir < this->rhopw->nrxx; ir++)
{
this->rho[3][ir] = this->rho[0][ir] - this->rho[1][ir];
this->rho[0][ir] = this->rho[0][ir] + this->rho[1][ir];
this->rho[1][ir] = 0.0;
this->rho[2][ir] = 0.0;
}
this->rho[3][ir] = this->rho[0][ir] - this->rho[1][ir];
this->rho[0][ir] = this->rho[0][ir] + this->rho[1][ir];
this->rho[1][ir] = 0.0;
this->rho[2][ir] = 0.0;
}
}
}
Expand All @@ -127,22 +103,12 @@ void Charge::init_rho(elecstate::efermi& eferm_iout,
GlobalV::ofs_running << " try to read kinetic energy density from file : " << ssc.str()
<< std::endl;
// mohan update 2012-02-10, sunliang update 2023-03-09
if (ModuleIO::read_cube(
#ifdef __MPI
& (GlobalC::Pgrid),
#endif
if (ModuleIO::read_vdata_palgrid(GlobalC::Pgrid,
(PARAM.inp.esolver_type == "sdft" ? GlobalV::RANK_IN_STOGROUP : GlobalV::MY_RANK),
is,
GlobalV::ofs_running,
PARAM.inp.nspin,
ssc.str(),
this->kin_r[is],
this->rhopw->nx,
this->rhopw->ny,
this->rhopw->nz,
eferm_iout.ef,
& (GlobalC::ucell),
this->prenspin))
GlobalC::ucell.nat))
{
GlobalV::ofs_running << " Read in the kinetic energy density: " << ssc.str() << std::endl;
}
Expand Down
21 changes: 2 additions & 19 deletions source/module_elecstate/module_charge/symmetry_rho.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -97,7 +97,7 @@ void Symmetry_rho::psymm(double* rho_part,
rhotot.resize(rho_basis->nxyz);
ModuleBase::GlobalFunc::ZEROS(rhotot.data(), rho_basis->nxyz);
}
Pgrid.reduce_to_fullrho(rhotot.data(), rho_part);
Pgrid.reduce(rhotot.data(), rho_part);

// (2)
if (GlobalV::MY_RANK == 0)
Expand Down Expand Up @@ -127,24 +127,7 @@ void Symmetry_rho::psymm(double* rho_part,
}

// (3)
const int ncxy = rho_basis->nx * rho_basis->ny;
std::vector<double> zpiece(ncxy);
for (int iz = 0; iz < rho_basis->nz; iz++)
{
ModuleBase::GlobalFunc::ZEROS(zpiece.data(), ncxy);
if (GlobalV::MY_RANK == 0)
{
for (int ix = 0; ix < rho_basis->nx; ix++)
{
for (int iy = 0; iy < rho_basis->ny; iy++)
{
const int ir = ix * rho_basis->ny + iy;
zpiece[ir] = rhotot[ix * rho_basis->ny * rho_basis->nz + iy * rho_basis->nz + iz];
}
}
}
Pgrid.zpiece_to_all(zpiece.data(), iz, rho_part);
}
Pgrid.bcast(rhotot.data(), rho_part, GlobalV::MY_RANK);
#endif
return;
}
33 changes: 3 additions & 30 deletions source/module_esolver/esolver_fp.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -156,43 +156,25 @@ void ESolver_FP::after_scf(const int istep)
this->pw_rhod->real2recip(this->pelec->charge->rho_save[is], this->pelec->charge->rhog_save[is]);
}
std::string fn =PARAM.globalv.global_out_dir + "/SPIN" + std::to_string(is + 1) + "_CHG.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
data,
is,
PARAM.inp.nspin,
istep,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.get_efval(is),
&(GlobalC::ucell),
PARAM.inp.out_chg[1],
1);
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
{
fn =PARAM.globalv.global_out_dir + "/SPIN" + std::to_string(is + 1) + "_TAU.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
this->pelec->charge->kin_r_save[is],
is,
PARAM.inp.nspin,
istep,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.get_efval(is),
&(GlobalC::ucell));
}
Expand Down Expand Up @@ -224,21 +206,12 @@ void ESolver_FP::after_scf(const int istep)
{
std::string fn =PARAM.globalv.global_out_dir + "/SPIN" + std::to_string(is + 1) + "_POT.cube";

ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
this->pelec->pot->get_effective_v(is),
is,
PARAM.inp.nspin,
istep,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
0.0, // efermi
&(GlobalC::ucell),
3, // precision
Expand Down
22 changes: 2 additions & 20 deletions source/module_esolver/esolver_ks_lcao.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1029,43 +1029,25 @@ void ESolver_KS_LCAO<TK, TR>::iter_finish(int& iter)
data = this->pelec->charge->rho_save[is];
}
std::string fn = PARAM.globalv.global_out_dir + "/tmp_SPIN" + std::to_string(is + 1) + "_CHG.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
data,
is,
PARAM.inp.nspin,
0,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.get_efval(is),
&(GlobalC::ucell),
3,
1);
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
{
fn = PARAM.globalv.global_out_dir + "/tmp_SPIN" + std::to_string(is + 1) + "_TAU.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
this->pelec->charge->kin_r_save[is],
is,
PARAM.inp.nspin,
0,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.get_efval(is),
&(GlobalC::ucell));
}
Expand Down
44 changes: 4 additions & 40 deletions source/module_esolver/esolver_ks_pw.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -239,21 +239,12 @@ void ESolver_KS_PW<T, Device>::before_scf(const int istep)
{
std::stringstream ss;
ss << PARAM.globalv.global_out_dir << "SPIN" << is + 1 << "_CHG_INI.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
this->pelec->charge->rho[is],
is,
PARAM.inp.nspin,
istep,
ss.str(),
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.ef,
&(GlobalC::ucell));
}
Expand All @@ -266,21 +257,12 @@ void ESolver_KS_PW<T, Device>::before_scf(const int istep)
{
std::stringstream ss;
ss << PARAM.globalv.global_out_dir << "SPIN" << is + 1 << "_POT_INI.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
this->pelec->pot->get_effective_v(is),
is,
PARAM.inp.nspin,
istep,
ss.str(),
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
0.0, // efermi
&(GlobalC::ucell),
11, // precsion
Expand Down Expand Up @@ -480,43 +462,25 @@ void ESolver_KS_PW<T, Device>::iter_finish(int& iter)
data = this->pelec->charge->rho_save[is];
}
std::string fn = PARAM.globalv.global_out_dir + "/tmp_SPIN" + std::to_string(is + 1) + "_CHG.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
data,
is,
PARAM.inp.nspin,
0,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.get_efval(is),
&(GlobalC::ucell),
3,
1);
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
{
fn = PARAM.globalv.global_out_dir + "/tmp_SPIN" + std::to_string(is + 1) + "_TAU.cube";
ModuleIO::write_cube(
#ifdef __MPI
this->pw_big->bz,
this->pw_big->nbz,
this->pw_rhod->nplane,
this->pw_rhod->startz_current,
#endif
ModuleIO::write_vdata_palgrid(GlobalC::Pgrid,
this->pelec->charge->kin_r_save[is],
is,
PARAM.inp.nspin,
0,
fn,
this->pw_rhod->nx,
this->pw_rhod->ny,
this->pw_rhod->nz,
this->pelec->eferm.get_efval(is),
&(GlobalC::ucell));
}
Expand Down
4 changes: 2 additions & 2 deletions source/module_esolver/io_npz.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -73,7 +73,7 @@ void ESolver_KS_LCAO<TK, TR>::read_mat_npz(std::string& zipname, hamilt::HContai
const int it = GlobalC::ucell.iat2it[iat];
const int ia = GlobalC::ucell.iat2ia[iat];

//get atomic number (copied from write_cube.cpp)
//get atomic number (copied from write_vdata_palgrid.cpp)
std::string element = "";
element = GlobalC::ucell.atoms[it].label;
std::string::iterator temp = element.begin();
Expand Down Expand Up @@ -368,7 +368,7 @@ void ESolver_KS_LCAO<TK, TR>::output_mat_npz(std::string& zipname, const hamilt:
const int it = GlobalC::ucell.iat2it[iat];
const int ia = GlobalC::ucell.iat2ia[iat];

//get atomic number (copied from write_cube.cpp)
//get atomic number (copied from write_vdata_palgrid.cpp)
std::string element = "";
element = GlobalC::ucell.atoms[it].label;
std::string::iterator temp = element.begin();
Expand Down
Loading

0 comments on commit 705dad7

Please sign in to comment.