From 5dd75995f5915253908fc0c467ed12623063f5bf Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 06:44:51 +0000 Subject: [PATCH 01/10] merge atominput and grid together. remove delete_vector --- source/Makefile.Objects | 1 - .../module_neighbor/CMakeLists.txt | 1 - .../module_neighbor/sltk_atom_arrange.cpp | 63 +-- .../module_neighbor/sltk_atom_arrange.h | 13 - .../module_neighbor/sltk_atom_input.cpp | 326 ------------- .../module_neighbor/sltk_atom_input.h | 143 ------ .../module_cell/module_neighbor/sltk_grid.cpp | 436 ++++-------------- .../module_cell/module_neighbor/sltk_grid.h | 80 +--- .../module_neighbor/sltk_grid_driver.cpp | 11 - .../module_neighbor/sltk_grid_driver.h | 1 - .../module_neighbor/test/CMakeLists.txt | 11 +- .../test/sltk_atom_input_test.cpp | 264 ----------- .../module_neighbor/test/sltk_grid_test.cpp | 21 +- source/module_esolver/esolver_lj.cpp | 9 - .../module_deepks/test/Makefile.Objects | 1 - source/module_md/test/CMakeLists.txt | 1 - 16 files changed, 126 insertions(+), 1256 deletions(-) delete mode 100644 source/module_cell/module_neighbor/sltk_atom_input.cpp delete mode 100644 source/module_cell/module_neighbor/sltk_atom_input.h delete mode 100644 source/module_cell/module_neighbor/test/sltk_atom_input_test.cpp diff --git a/source/Makefile.Objects b/source/Makefile.Objects index 7dd7762326..a4bd83dfb4 100644 --- a/source/Makefile.Objects +++ b/source/Makefile.Objects @@ -368,7 +368,6 @@ OBJS_MD=fire.o\ OBJS_NEIGHBOR=sltk_atom.o\ sltk_atom_arrange.o\ - sltk_atom_input.o\ sltk_grid.o\ sltk_grid_driver.o\ diff --git a/source/module_cell/module_neighbor/CMakeLists.txt b/source/module_cell/module_neighbor/CMakeLists.txt index dfa2463d67..9b5fea4d14 100644 --- a/source/module_cell/module_neighbor/CMakeLists.txt +++ b/source/module_cell/module_neighbor/CMakeLists.txt @@ -3,7 +3,6 @@ add_library( OBJECT sltk_atom.cpp sltk_atom_arrange.cpp - sltk_atom_input.cpp sltk_grid.cpp sltk_grid_driver.cpp ) diff --git a/source/module_cell/module_neighbor/sltk_atom_arrange.cpp b/source/module_cell/module_neighbor/sltk_atom_arrange.cpp index d4d020c5ff..68143a5c7e 100644 --- a/source/module_cell/module_neighbor/sltk_atom_arrange.cpp +++ b/source/module_cell/module_neighbor/sltk_atom_arrange.cpp @@ -1,5 +1,4 @@ #include "sltk_atom_arrange.h" -#include "sltk_atom_input.h" #include "module_parameter/parameter.h" #include "sltk_grid.h" #include "sltk_grid_driver.h" @@ -79,52 +78,18 @@ void atom_arrange::search( { ModuleBase::TITLE("atom_arrange", "search"); ModuleBase::timer::tick("atom_arrange","search"); -/* std::cout << "pbc_flag = " << pbc_flag << std::endl; - std::cout << "search_radius_bohr = " << search_radius_bohr << std::endl; - std::cout << "test_atom_in = " << test_atom_in << std::endl; - std::cout << "test_only = " << test_only << std::endl; - */ assert( search_radius_bohr > 0.0 ); -// OUT(ofs_in,"Atom coordinates reading from",PARAM.inp.stru_file); -// OUT(ofs_in,"The coordinate type",ucell.Coordinate); -// OUT(ofs_in,"Use cartesian(unit:lat0) coordinate","TRUE"); -// if(PARAM.inp.out_level != "m") OUT(ofs_in,"searching radius is (Bohr))", search_radius_bohr); -// if(PARAM.inp.out_level != "m") OUT(ofs_in,"searching radius unit is (Bohr))",ucell.lat0); - ModuleBase::GlobalFunc::OUT(ofs_in,"searching radius is (Bohr))", search_radius_bohr); ModuleBase::GlobalFunc::OUT(ofs_in,"searching radius unit is (Bohr))",ucell.lat0); assert(ucell.nat > 0); - //============================= - // Initial Atom information - //============================= const double radius_lat0unit = search_radius_bohr / ucell.lat0; - ModuleBase::timer::tick("atom_arrange", "Atom_input"); - - Atom_input at( - ofs_in, - ucell, - ucell.nat, - ucell.ntype, - pbc_flag, - radius_lat0unit, - test_atom_in); - ModuleBase::timer::tick("atom_arrange", "Atom_input"); - - //=========================================== - // Print important information in Atom_input - //=========================================== -// at.print(std::cout); -// at.print_xyz_format("1.xyz"); - //========================================= - // Construct Grid , Cells , Adjacent atoms - //========================================= ModuleBase::timer::tick("atom_arrange", "grid_d.init"); - grid_d.init(ofs_in, ucell, at); + grid_d.init(ofs_in, ucell, radius_lat0unit, pbc_flag); ModuleBase::timer::tick("atom_arrange", "grid_d.init"); ModuleBase::timer::tick("atom_arrange", "search"); @@ -163,29 +128,3 @@ void atom_arrange::search( return; } - - -//2015-05-07 -void atom_arrange::delete_vector( - std::ofstream &ofs_in, - const bool pbc_flag, // GlobalV::SEARCH_PBC - Grid_Driver &grid_d, - const UnitCell &ucell, - const double &search_radius_bohr, - const int &test_atom_in) -{ - const double radius_lat0unit2 = search_radius_bohr / ucell.lat0; - - Atom_input at2( - ofs_in, - ucell, - ucell.nat, - ucell.ntype, - pbc_flag, - radius_lat0unit2, - test_atom_in); - - grid_d.delete_vector(at2.getGrid_layerX_minus(),at2.getGrid_layerY_minus(),at2.getGrid_layerZ_minus()); - - grid_d.delete_Cell(); -} diff --git a/source/module_cell/module_neighbor/sltk_atom_arrange.h b/source/module_cell/module_neighbor/sltk_atom_arrange.h index 14a2771487..5dd1c94f98 100644 --- a/source/module_cell/module_neighbor/sltk_atom_arrange.h +++ b/source/module_cell/module_neighbor/sltk_atom_arrange.h @@ -3,7 +3,6 @@ #include "sltk_grid.h" #include "sltk_grid_driver.h" -#include "sltk_atom_input.h" class atom_arrange @@ -29,18 +28,6 @@ class atom_arrange const double& rcutmax_Phi, const double& rcutmax_Beta, const bool gamma_only_local); - - //2015-05-07 - static void delete_vector( - std::ofstream &ofs_in, - const bool pbc_flag, - Grid_Driver &grid_d, - const UnitCell &ucell, - const double &search_radius_bohr, - const int &test_atom_in); - -private: - }; #endif diff --git a/source/module_cell/module_neighbor/sltk_atom_input.cpp b/source/module_cell/module_neighbor/sltk_atom_input.cpp deleted file mode 100644 index 2436146932..0000000000 --- a/source/module_cell/module_neighbor/sltk_atom_input.cpp +++ /dev/null @@ -1,326 +0,0 @@ -#include "sltk_atom_input.h" - -#include "module_base/memory.h" -#include "module_parameter/parameter.h" -#include "sltk_grid.h" - -//========================================================== -// define constructor and deconstructor -//========================================================== -Atom_input::Atom_input(std::ofstream& ofs_in, - const UnitCell& ucell, - const int amount, - const int ntype, - const bool boundary_in, - const double radius_in, - const int& test_atom_in) - : periodic_boundary(boundary_in), radius(radius_in), expand_flag(false), glayerX(1), glayerX_minus(0), glayerY(1), - glayerY_minus(0), glayerZ(1), glayerZ_minus(0), test_atom_input(test_atom_in) -{ - ModuleBase::TITLE("Atom_input", "Atom_input"); - - if (test_atom_input) - { - ModuleBase::GlobalFunc::OUT(ofs_in, "ntype", ntype); - ModuleBase::GlobalFunc::OUT(ofs_in, "Amount(atom number)", amount); - ModuleBase::GlobalFunc::OUT(ofs_in, "Periodic_boundary", periodic_boundary); - ModuleBase::GlobalFunc::OUT(ofs_in, "Searching radius(lat0)", radius); - } - - if (radius < 0.0) - { - ModuleBase::WARNING_QUIT("atom_arrange::init", " search radius < 0,forbidden"); - } - // random selection, in order to estimate again. - this->x_min = ucell.atoms[0].tau[0].x; - this->y_min = ucell.atoms[0].tau[0].y; - this->z_min = ucell.atoms[0].tau[0].z; - this->x_max = ucell.atoms[0].tau[0].x; - this->y_max = ucell.atoms[0].tau[0].y; - this->z_max = ucell.atoms[0].tau[0].z; - - // calculate min & max value - for (int i = 0; i < ntype; i++) - { - for (int j = 0; j < ucell.atoms[i].na; j++) - { - x_min = std::min(x_min, ucell.atoms[i].tau[j].x); - x_max = std::max(x_max, ucell.atoms[i].tau[j].x); - y_min = std::min(y_min, ucell.atoms[i].tau[j].y); - y_max = std::max(y_max, ucell.atoms[i].tau[j].y); - z_min = std::min(z_min, ucell.atoms[i].tau[j].z); - z_max = std::max(z_max, ucell.atoms[i].tau[j].z); - } - } - - if (test_atom_input) - { - ModuleBase::GlobalFunc::OUT(ofs_in, "Find the coordinate range of the input atom(unit:lat0)."); - ModuleBase::GlobalFunc::OUT(ofs_in, "min_tau", x_min, y_min, z_min); - ModuleBase::GlobalFunc::OUT(ofs_in, "max_tau", x_max, y_max, z_max); - } - - //---------------------------------------------------------- - // CALL MEMBER FUNCTION : - // NAME : Check_Expand_Condition(check if swe need to - // expand grid,and generate 6 MEMBER VARIABLE(number of - // layers for 6 dimension) - // initial value for "glayerX,Y,Z" : 1 - // (if > 2 ,expand flag = 1) - // initial value for "glayerX,Y,Z_minus" : 0 - // ( if > 1 ,expand flag = 1) - //---------------------------------------------------------- - - this->Check_Expand_Condition(ucell); - - if (test_atom_input) - { - ModuleBase::GlobalFunc::OUT(ofs_in, "glayer+", glayerX, glayerY, glayerZ); - ModuleBase::GlobalFunc::OUT(ofs_in, "glayer-", glayerX_minus, glayerY_minus, glayerZ_minus); - ModuleBase::GlobalFunc::OUT(ofs_in, "expand_flag", expand_flag); - } - - //---------------------------------------------------------- - // CALL MEMBER FUNCTION : - // NAME : calculate_cells - // Calculate how many cells we need in each direction. - //---------------------------------------------------------- - this->calculate_cells(); - if (test_atom_input) - { - ModuleBase::GlobalFunc::OUT(ofs_in, "CellDim", cell_nx, cell_ny, cell_nz); - } - return; -} - -Atom_input::~Atom_input() -{ -} - -//============================================ -// !!!! May still have bug, be very careful!! -// should use the same algorithm to generate -// dxe, dye, dze in grid_meshcell.cpp. -//============================================ -void Atom_input::Check_Expand_Condition(const UnitCell& ucell) -{ - // ModuleBase::TITLE(GlobalV::ofs_running, "Atom_input", "Check_Expand_Condition"); - - if (!periodic_boundary) - { - return; - } - - /*2016-07-19, LiuXh - // the unit of extent_1DX,Y,Z is lat0. - // means still how far can be included now. - double extent_1DX = glayerX * clength0 - dmaxX; - while (radius > extent_1DX) - { - glayerX++; - extent_1DX = glayerX * clength0 - dmaxX; - } - double extent_1DY = glayerY * clength1 - dmaxY; - while (radius > extent_1DY) - { - glayerY++; - extent_1DY = glayerY * clength1 - dmaxY; - } - double extent_1DZ = glayerZ * clength2 - dmaxZ; - while (radius > extent_1DZ) - { - glayerZ++; - extent_1DZ = glayerZ * clength2 - dmaxZ; - } - - // in case the cell is not retangle. - // mohan added 2009-10-23 - // if this is not added, it's a serious bug. - glayerX++; - glayerY++; - glayerZ++; - if(test_atom_input) - { - GlobalV::ofs_running << " Extend distance from the (maxX,maxY,maxZ) direct position in this unitcell: " << - std::endl; - } - - if(test_atom_input)OUT(GlobalV::ofs_running,"ExtentDim+",extent_1DX,extent_1DY,extent_1DZ); - - double extent_1DX_minus = glayerX_minus * clength0 + dminX; - while (radius > extent_1DX_minus) - { - glayerX_minus++; - extent_1DX_minus = glayerX_minus * clength0 + dminX; - } - double extent_1DY_minus = glayerY_minus * clength1 + dminY; - while (radius > extent_1DY_minus) - { - glayerY_minus++; - extent_1DY_minus = glayerY_minus * clength1 + dminY; - } - double extent_1DZ_minus = glayerZ_minus * clength2 + dminZ; - while (radius > extent_1DZ_minus) - { - glayerZ_minus++; - extent_1DZ_minus = glayerZ_minus * clength2 + dminZ; - } - - // in case the cell is not retangle. - // mohan added 2009-10-23 - // if this is not added, it's a serious bug. - glayerX_minus++; - glayerY_minus++; - glayerZ_minus++; - - //glayerX_minus++; - //glayerY_minus++; - //glayerZ_minus++; - 2016-07-19, LiuXh*/ - // Begin, 2016-07-19, LiuXh - double a23_1 = ucell.latvec.e22 * ucell.latvec.e33 - ucell.latvec.e23 * ucell.latvec.e32; - double a23_2 = ucell.latvec.e21 * ucell.latvec.e33 - ucell.latvec.e23 * ucell.latvec.e31; - double a23_3 = ucell.latvec.e21 * ucell.latvec.e32 - ucell.latvec.e22 * ucell.latvec.e31; - double a23_norm = sqrt(a23_1 * a23_1 + a23_2 * a23_2 + a23_3 * a23_3); - double extend_v = a23_norm * radius; - double extend_d1 = extend_v / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; - int extend_d11 = static_cast(extend_d1); - // 2016-09-05, LiuXh - if (extend_d1 - extend_d11 > 0.0) - { - extend_d11 += 1; - } - - double a31_1 = ucell.latvec.e32 * ucell.latvec.e13 - ucell.latvec.e33 * ucell.latvec.e12; - double a31_2 = ucell.latvec.e31 * ucell.latvec.e13 - ucell.latvec.e33 * ucell.latvec.e11; - double a31_3 = ucell.latvec.e31 * ucell.latvec.e12 - ucell.latvec.e32 * ucell.latvec.e11; - double a31_norm = sqrt(a31_1 * a31_1 + a31_2 * a31_2 + a31_3 * a31_3); - double extend_d2 = a31_norm * radius / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; - int extend_d22 = static_cast(extend_d2); - // 2016-09-05, LiuXh - if (extend_d2 - extend_d22 > 0.0) - { - extend_d22 += 1; - } - - double a12_1 = ucell.latvec.e12 * ucell.latvec.e23 - ucell.latvec.e13 * ucell.latvec.e22; - double a12_2 = ucell.latvec.e11 * ucell.latvec.e23 - ucell.latvec.e13 * ucell.latvec.e21; - double a12_3 = ucell.latvec.e11 * ucell.latvec.e22 - ucell.latvec.e12 * ucell.latvec.e21; - double a12_norm = sqrt(a12_1 * a12_1 + a12_2 * a12_2 + a12_3 * a12_3); - double extend_d3 = a12_norm * radius / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; - int extend_d33 = static_cast(extend_d3); - // 2016-09-05, LiuXh - if (extend_d3 - extend_d33 > 0.0) - { - extend_d33 += 1; - } - - glayerX = extend_d11 + 1; - glayerY = extend_d22 + 1; - glayerZ = extend_d33 + 1; - // Begin, 2016-09-05, LiuXh - // glayerX_minus = extend_d11 +1; - // glayerY_minus = extend_d22 +1; - // glayerZ_minus = extend_d33 +1; - glayerX_minus = extend_d11; - glayerY_minus = extend_d22; - glayerZ_minus = extend_d33; - // End, 2016-09-05, LiuXh - - if (glayerX == 1) - { - glayerX++; - } - if (glayerY == 1) - { - glayerY++; - } - if (glayerZ == 1) - { - glayerZ++; - } - if (glayerX_minus == 1) - { - glayerX_minus++; - } - if (glayerY_minus == 1) - { - glayerY_minus++; - } - if (glayerZ_minus == 1) - { - glayerZ_minus++; - } - // End, 2016-07-19, LiuXh - /* - if(test_atom_input) - { - GlobalV::ofs_running << " Extend distance from the (minX,minY,minZ) direct position in this unitcell: " << - std::endl; - } - - if(test_atom_input)OUT(GlobalV::ofs_running,"ExtentDim-",extent_1DX_minus,extent_1DY_minus,extent_1DZ_minus); - */ - //---------------------------------------------------------- - // EXPLAIN : if extent don't satisfty the searching - // requiment, we must expand one more layer - //---------------------------------------------------------- - - if (glayerX > 2 || glayerY > 2 || glayerZ > 2) - { - this->expand_flag = true; - } - else if (glayerX_minus > 1 || glayerX_minus > 1 || glayerX_minus > 1) - { - this->expand_flag = true; - } - else - { - this->expand_flag = false; - } - return; -} - -void Atom_input::calculate_cells() -{ - ModuleBase::TITLE("Atom_input", "calculate_cells"); - //---------------------------------------------------------- - // EXPLAIN : - // Expand_Case : Simple , we already know the cell numbres, - // all the trouble is only to construct adjacentset using all - // the cells. - // Not_Expand_Case : Using searching radius to construct - // the cells , trouble here,but is the convenience of searching - // time , we then only need to search 27-adjacent cell for each cell. - //---------------------------------------------------------- - if (expand_flag) - { - cell_nx = glayerX + glayerX_minus; - cell_ny = glayerY + glayerY_minus; - cell_nz = glayerZ + glayerZ_minus; - } - else - { - // maybe a bug, if we don't use direct - // coordinates, mohan note 2011-04-14 - double real_nx, real_ny, real_nz; - real_nx = (x_max - x_min) / radius; - real_ny = (y_max - y_min) / radius; - real_nz = (z_max - z_min) / radius; - cell_nx = static_cast(real_nx) + 1; - cell_ny = static_cast(real_ny) + 1; - cell_nz = static_cast(real_nz) + 1; - } - - //================ - // Wrong ! - //================ - // if(int_nx != real_nx) this->cell_nx++; - // if(int_ny != real_ny) this->cell_ny++; - // if(int_nz != real_nz) this->cell_nz++; - //======================================= - // Not need because if int_nx = real_nx, - // the position belong to the next cell - //======================================= - return; -} diff --git a/source/module_cell/module_neighbor/sltk_atom_input.h b/source/module_cell/module_neighbor/sltk_atom_input.h deleted file mode 100644 index 0b32777244..0000000000 --- a/source/module_cell/module_neighbor/sltk_atom_input.h +++ /dev/null @@ -1,143 +0,0 @@ -#ifndef ATOM_INPUT_H -#define ATOM_INPUT_H - -#include "module_cell/unitcell.h" -#include "sltk_atom.h" - -class Atom_input -{ - public: - //========================================================== - // Constructors and destructor - //========================================================== - Atom_input(std::ofstream& ofs_in, - const UnitCell& ucell, - const int amount = 0, // number of atoms - const int ntype = 0, // number of atom_types - const bool boundary = true, // 1 : periodic ocndition - const double radius_in = 0, // searching radius - const int& test_atom_in = 0 // caoyu reconst 2021-05-24 - ); - ~Atom_input(); - - public: - bool getExpandFlag() const - { - return expand_flag; - } - - int getBoundary() const - { - return periodic_boundary; - } - - double getRadius() const - { - return radius; - } - - //========================================================== - // - //========================================================== - double minX() const - { - return x_min; - } - - double minY() const - { - return y_min; - } - - double minZ() const - { - return z_min; - } - - //========================================================== - // - //========================================================== - int getCell_nX() const - { - return cell_nx; - } - - int getCell_nY() const - { - return cell_ny; - } - - int getCell_nZ() const - { - return cell_nz; - } - - //========================================================== - // - //========================================================== - int getGrid_layerX() const - { - return glayerX; - } - - int getGrid_layerX_minus() const - { - return glayerX_minus; - } - - int getGrid_layerY() const - { - return glayerY; - } - - int getGrid_layerY_minus() const - { - return glayerY_minus; - } - - int getGrid_layerZ() const - { - return glayerZ; - } - - int getGrid_layerZ_minus() const - { - return glayerZ_minus; - } - - private: - int test_atom_input; // caoyu reconst 2021-05-24 - bool periodic_boundary; - - double radius; - - double x_min; - double y_min; - double z_min; - double x_max; - double y_max; - double z_max; - //========================================================== - // MEMBRE FUNCTION : - // NAME : Check_Expand_Condition - //========================================================== - void Check_Expand_Condition(const UnitCell& ucell); - bool expand_flag; - int glayerX; - int glayerX_minus; - int glayerY; - int glayerY_minus; - int glayerZ; - int glayerZ_minus; - - //========================================================== - // MEMBRE FUNCTION : - // NAME : Expand_Grid - //========================================================== - void calculate_cells(); - int cell_nx; - int cell_ny; - int cell_nz; -}; - -#endif diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index 117a5281e1..ce7b254146 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -4,7 +4,6 @@ #include "module_base/global_variable.h" #include "module_base/memory.h" #include "module_base/timer.h" -#include "sltk_atom_input.h" //================== // Class CellSet @@ -18,70 +17,125 @@ CellSet::CellSet() Grid::Grid(const int& test_grid_in) : test_grid(test_grid_in) { - // ModuleBase::TITLE("Grid","Grid"); - //---------------------------------------------------------- - // EXPLAIN : init_cell_flag (use this flag in case memory - // leak) - //---------------------------------------------------------- - init_cell_flag = false; } Grid::~Grid() { - this->delete_Cell(); } -void Grid::init(std::ofstream& ofs_in, const UnitCell& ucell, const Atom_input& input) +void Grid::init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, bool pbc_in) { ModuleBase::TITLE("SLTK_Grid", "init"); + this->pbc = pbc_in; + this->sradius2 = radius * radius; + this->sradius = radius; + + this->x_min = ucell.atoms[0].tau[0].x; + this->y_min = ucell.atoms[0].tau[0].y; + this->z_min = ucell.atoms[0].tau[0].z; + this->x_max = ucell.atoms[0].tau[0].x; + this->y_max = ucell.atoms[0].tau[0].y; + this->z_max = ucell.atoms[0].tau[0].z; - this->setMemberVariables(ofs_in, input); - this->Build_Hash_Table(ucell, input); - this->setBoundaryAdjacent(ofs_in, input); + + // calculate min & max value + for (int i = 0; i < ucell.ntype; i++) + { + for (int j = 0; j < ucell.atoms[i].na; j++) + { + x_min = std::min(x_min, ucell.atoms[i].tau[j].x); + x_max = std::max(x_max, ucell.atoms[i].tau[j].x); + y_min = std::min(y_min, ucell.atoms[i].tau[j].y); + y_max = std::max(y_max, ucell.atoms[i].tau[j].y); + z_min = std::min(z_min, ucell.atoms[i].tau[j].z); + z_max = std::max(z_max, ucell.atoms[i].tau[j].z); + } + } + this->Check_Expand_Condition(ucell); + this->setMemberVariables(ofs_in); + this->Build_Hash_Table(ucell); + this->setBoundaryAdjacent(ofs_in); } +void Grid::Check_Expand_Condition(const UnitCell& ucell) +{ + // ModuleBase::TITLE(GlobalV::ofs_running, "Atom_input", "Check_Expand_Condition"); + + if (!pbc) + { + return; + } + double a23_1 = ucell.latvec.e22 * ucell.latvec.e33 - ucell.latvec.e23 * ucell.latvec.e32; + double a23_2 = ucell.latvec.e21 * ucell.latvec.e33 - ucell.latvec.e23 * ucell.latvec.e31; + double a23_3 = ucell.latvec.e21 * ucell.latvec.e32 - ucell.latvec.e22 * ucell.latvec.e31; + double a23_norm = sqrt(a23_1 * a23_1 + a23_2 * a23_2 + a23_3 * a23_3); + double extend_v = a23_norm * sradius; + double extend_d1 = extend_v / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; + int extend_d11 = static_cast(extend_d1); + // 2016-09-05, LiuXh + if (extend_d1 - extend_d11 > 0.0) + { + extend_d11 += 1; + } + + double a31_1 = ucell.latvec.e32 * ucell.latvec.e13 - ucell.latvec.e33 * ucell.latvec.e12; + double a31_2 = ucell.latvec.e31 * ucell.latvec.e13 - ucell.latvec.e33 * ucell.latvec.e11; + double a31_3 = ucell.latvec.e31 * ucell.latvec.e12 - ucell.latvec.e32 * ucell.latvec.e11; + double a31_norm = sqrt(a31_1 * a31_1 + a31_2 * a31_2 + a31_3 * a31_3); + double extend_d2 = a31_norm * sradius / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; + int extend_d22 = static_cast(extend_d2); + // 2016-09-05, LiuXh + if (extend_d2 - extend_d22 > 0.0) + { + extend_d22 += 1; + } + + double a12_1 = ucell.latvec.e12 * ucell.latvec.e23 - ucell.latvec.e13 * ucell.latvec.e22; + double a12_2 = ucell.latvec.e11 * ucell.latvec.e23 - ucell.latvec.e13 * ucell.latvec.e21; + double a12_3 = ucell.latvec.e11 * ucell.latvec.e22 - ucell.latvec.e12 * ucell.latvec.e21; + double a12_norm = sqrt(a12_1 * a12_1 + a12_2 * a12_2 + a12_3 * a12_3); + double extend_d3 = a12_norm * sradius / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; + int extend_d33 = static_cast(extend_d3); + // 2016-09-05, LiuXh + if (extend_d3 - extend_d33 > 0.0) + { + extend_d33 += 1; + } + glayerX = extend_d11 + 1; + glayerY = extend_d22 + 1; + glayerZ = extend_d33 + 1; + glayerX_minus = extend_d11; + glayerY_minus = extend_d22; + glayerZ_minus = extend_d33; +} //========================================================== // MEMBER FUNCTION : // NAME : setMemberVariables(read in data from Atom_input) //========================================================== -void Grid::setMemberVariables(std::ofstream& ofs_in, // output data to ofs - const Atom_input& input) +void Grid::setMemberVariables(std::ofstream& ofs_in) { ModuleBase::TITLE("SLTK_Grid", "setMemberVariables"); - this->delete_Cell(); // mohan add 2010-09-05 // AdjacentSet::call_times = 0; - this->pbc = input.getBoundary(); - this->sradius2 = input.getRadius() * input.getRadius(); - this->sradius = input.getRadius(); - this->expand_flag = input.getExpandFlag(); if (test_grid) { ModuleBase::GlobalFunc::OUT(ofs_in, "PeriodicBoundary", this->pbc); ModuleBase::GlobalFunc::OUT(ofs_in, "Radius(unit:lat0)", sradius); - ModuleBase::GlobalFunc::OUT(ofs_in, "Expand_flag", expand_flag); } //---------------------------------------------------------- // EXPLAIN : (d_minX,d_minY,d_minZ)minimal value of // x[] ,y[] , z[] //---------------------------------------------------------- - this->d_minX = input.minX(); - this->d_minY = input.minY(); - this->d_minZ = input.minZ(); - if (test_grid) - { - ModuleBase::GlobalFunc::OUT(ofs_in, "MinCoordinate", d_minX, d_minY, d_minZ); - } //---------------------------------------------------------- // set dx, dy, dz //---------------------------------------------------------- - this->cell_nx = input.getCell_nX(); - this->cell_ny = input.getCell_nY(); - this->cell_nz = input.getCell_nZ(); + this->cell_nx = glayerX + glayerX_minus; + this->cell_ny = glayerY + glayerY_minus; + this->cell_nz = glayerZ + glayerZ_minus; if (test_grid) { ModuleBase::GlobalFunc::OUT(ofs_in, "CellNumber", cell_nx, cell_ny, cell_nz); @@ -96,26 +150,18 @@ void Grid::setMemberVariables(std::ofstream& ofs_in, // output data to ofs Cell[i][j].resize(cell_nz); } } - this->init_cell_flag = true; - this->true_cell_x = input.getGrid_layerX_minus(); - this->true_cell_y = input.getGrid_layerY_minus(); - this->true_cell_z = input.getGrid_layerZ_minus(); + this->true_cell_x = glayerX_minus; + this->true_cell_y = glayerX_minus; + this->true_cell_z = glayerX_minus; } -void Grid::setBoundaryAdjacent(std::ofstream& ofs_in, const Atom_input& input) +void Grid::setBoundaryAdjacent(std::ofstream& ofs_in) { - if (expand_flag) - { - this->Construct_Adjacent_expand(true_cell_x, true_cell_y, true_cell_z); - } - else - { - this->Construct_Adjacent_begin(); - } + this->Construct_Adjacent_expand(true_cell_x, true_cell_y, true_cell_z); } -void Grid::Build_Hash_Table(const UnitCell& ucell, const Atom_input& input) +void Grid::Build_Hash_Table(const UnitCell& ucell) { ModuleBase::timer::tick("Grid", "Build_Hash_Table"); @@ -138,11 +184,11 @@ void Grid::Build_Hash_Table(const UnitCell& ucell, const Atom_input& input) ModuleBase::Vector3 vec2(ucell.latvec.e21, ucell.latvec.e22, ucell.latvec.e23); ModuleBase::Vector3 vec3(ucell.latvec.e31, ucell.latvec.e32, ucell.latvec.e33); - for (int ix = -input.getGrid_layerX_minus(); ix < input.getGrid_layerX(); ix++) + for (int ix = -glayerX_minus; ix < glayerX; ix++) { - for (int iy = -input.getGrid_layerY_minus(); iy < input.getGrid_layerY(); iy++) + for (int iy = -glayerY_minus; iy < glayerY; iy++) { - for (int iz = -input.getGrid_layerZ_minus(); iz < input.getGrid_layerZ(); iz++) + for (int iz = -glayerZ_minus; iz < glayerZ; iz++) { for (int i = 0; i < ucell.ntype; i++) { @@ -153,30 +199,9 @@ void Grid::Build_Hash_Table(const UnitCell& ucell, const Atom_input& input) double z = ucell.atoms[i].tau[j].z + vec1[2] * ix + vec2[2] * iy + vec3[2] * iz; FAtom atom(x, y, z, i, j, ix, iy, iz); int a, b, c; - if (expand_flag) - { - // EXPLAIN : In expand grid case, - // the input cell is exactly the same as input file. - a = atom.getCellX() + true_cell_x; - b = atom.getCellY() + true_cell_y; - c = atom.getCellZ() + true_cell_z; - } - else - { - //---------------------------------------------------------- - // EXPLAIN : Not expand case , the cell is 'cubic', - // the three dimension length : - // cell_x_length = |radius| - // cell_y_length = |radius| - // cell_z_length = |radius| - // - // So we don't need crystal coordinate to locate the atom. - // We use cartesian coordinate directly. - //---------------------------------------------------------- - a = static_cast(std::floor((atom.x() - this->d_minX) / this->sradius)); - b = static_cast(std::floor((atom.y() - this->d_minY) / this->sradius)); - c = static_cast(std::floor((atom.z() - this->d_minZ) / this->sradius)); - } + a = atom.getCellX() + true_cell_x; + b = atom.getCellY() + true_cell_y; + c = atom.getCellZ() + true_cell_z; this->Cell[a][b][c].atom_map[atom.getType()][atom.getNatom()] = atom; } @@ -191,40 +216,11 @@ void Grid::Construct_Adjacent_expand(const int true_i, const int true_j, const i { ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); - //----------------------------------------------------------- - // EXPLAIN : (true_i,true_j,true_k) is the cell we want - // to found AdjacentSet.And other cell save the displacement - // of center_grid in 'in_grid' - //----------------------------------------------------------- - for (int i = 0; i < this->cell_nx; i++) - { - for (int j = 0; j < this->cell_ny; j++) - { - for (int k = 0; k < this->cell_nz; k++) - { - this->Cell[i][j][k].in_grid[0] = i - true_i; - this->Cell[i][j][k].in_grid[1] = j - true_j; - this->Cell[i][j][k].in_grid[2] = k - true_k; - } - } - } - - //---------------------------------------------------------- - // EXPLAIN : Only construct AdjacentSet for 'true' cell. - //---------------------------------------------------------- for (auto& atom_vector: this->Cell[true_i][true_j][true_k].atom_map) { for (auto& fatom: atom_vector) { - if (this->pbc) - { - Construct_Adjacent_expand_periodic(true_i, true_j, true_k, fatom); - // std::cout << "fatom1 = " << fatom.getNatom() << " " << fatom.getAdjacent().size() << std::endl; - } - else - { - ModuleBase::WARNING_QUIT("Construct_Adjacent_expand", "\n Expand case, must use periodic boundary."); - } + Construct_Adjacent_expand_periodic(true_i, true_j, true_k, fatom); } } ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); @@ -254,194 +250,6 @@ void Grid::Construct_Adjacent_expand_periodic(const int true_i, const int true_j ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand_periodic"); } -void Grid::Construct_Adjacent_begin() -{ - // if (test_grid)ModuleBase::TITLE(ofs_running, "Grid", "Construct_Adjacent_begin"); - - //---------------------------------------------------------- - // EXPLAIN : Searching in all cells in this grid - //---------------------------------------------------------- - - for (int i = 0; i < this->cell_nx; i++) - { - for (int j = 0; j < this->cell_ny; j++) - { - for (int k = 0; k < this->cell_nz; k++) - { - //---------------------------------------------------------- - // EXPLAIN : Cell length == Number of atoms in this cell - //---------------------------------------------------------- - for (auto& atom_vector: this->Cell[i][j][k].atom_map) - { - for (auto& fatom2: atom_vector) - { - // pbc: periodic boundary condition - if (this->pbc) - { - Construct_Adjacent_periodic(i, j, k, fatom2); - } - else - { - Construct_Adjacent_nature(i, j, k, fatom2); - } - } - - } // ia - } // k - } // j - } // i - - return; -} - -void Grid::Construct_Adjacent_nature(const int i, const int j, const int k, FAtom& fatom1) -{ - // if(test_grid)ModuleBase::TITLE(ofs_running,"Grid","Construct_Adjacent_nature"); - for (int i2 = i - 1; i2 <= i + 1; i2++) - { - if (i2 < cell_nx && i2 >= 0) - { - for (int j2 = j - 1; j2 <= j + 1; j2++) - { - if (j2 < cell_ny && j2 >= 0) - { - for (int k2 = k - 1; k2 <= k + 1; k2++) - { - if (k2 < cell_nz && k2 >= 0) - { - for (auto& atom_vector: this->Cell[i2][j2][k2].atom_map) - { - for (auto& fatom2: atom_vector) - { - Construct_Adjacent_final(i, j, k, fatom1, i2, j2, k2, fatom2); - } // ia2 - } - } - } // k2 - } - } // j2 - } - } // 2 - - return; -} - -void Grid::Construct_Adjacent_periodic(const int i, const int j, const int k, FAtom& fatom1) -{ - // if(test_grid)ModuleBase::TITLE(ofs_running,"Grid","Construct_Adjacent_periodic"); - bool first_i = true; - - for (int i2 = i - 1; i2 <= i + 1; i2++) - { - bool first_j = true; - - for (int j2 = j - 1; j2 <= j + 1; j2++) - { - bool first_k = true; - - for (int k2 = k - 1; k2 <= k + 1; k2++) - { - int temp_i = i2; - int temp_j = j2; - int temp_k = k2; - - int g0 = 0; - int g1 = 0; - int g2 = 0; - - if (i2 < 0) - { - g0 = -1; - - if (first_i) - { - if (cell_nx >= 2) - { - i2--; - temp_i--; - } - - first_i = false; - } - - i2 += cell_nx; - } - else if (i2 >= cell_nx) - { - g0 = 1; - i2 -= cell_nx; - } - - if (j2 < 0) - { - g1 = -1; - - if (first_j) - { - if (cell_ny >= 2) - { - j2--; - temp_j--; - } - - first_j = false; - } - - j2 += cell_ny; - } - else if (j2 >= cell_ny) - { - g1 = 1; - j2 -= cell_ny; - } - - if (k2 < 0) - { - g2 = -1; - - if (first_k) - { - if (cell_nz >= 2) - { - k2--; - temp_k--; - } - - first_k = false; - } - - k2 += cell_nz; - } - else if (k2 >= cell_nz) - { - g2 = 1; - k2 -= cell_nz; - } - - Cell[i2][j2][k2].in_grid[0] = g0; - - Cell[i2][j2][k2].in_grid[1] = g1; - Cell[i2][j2][k2].in_grid[2] = g2; - - for (auto& atom_vector: this->Cell[i2][j2][k2].atom_map) - { - for (auto& fatom2: atom_vector) - { - Construct_Adjacent_final(i, j, k, fatom1, i2, j2, k2, fatom2); - } // ia2 - } - - i2 = temp_i; - - j2 = temp_j; - - k2 = temp_k; // resume i2 j2 k2 - } // k2 - } // j2 - } // i2 - - return; -} void Grid::Construct_Adjacent_final(const int i, const int j, @@ -452,51 +260,13 @@ void Grid::Construct_Adjacent_final(const int i, const int k2, FAtom& fatom2) { - //---------------------------------------------------------- - // EXPLAIN : expand_case not_expand_case - // (i,j,k,ia) only the 'true' cell only the 'true' grid - // (i2,j2,k2,ia2) all atoms in grid all atoms in 27*cell - //---------------------------------------------------------- - // (suitable for small cell periodic condition) - // Expand_Case : many 'pseudo' cells, only one true cell, - // one grid(true grid). - // Advantage : only the atoms in 'true' cell need to construct - // AdjacentSet. - // Disadvantage : must search all atoms in true grid to construct - // AdjacentSet. - // - // (suitable for large cell periodic/nature condition,here - // we discuss periodic case,once you known this case, nature - // boundary is easy to understand) - // Not_Expand_Case : 27 'pseudo' grid,only one true grid, - // many true cells. - // Advantage : (the disadvantage above is the advantage now) - // only need to search 27*cells to construct AdjacentSet - // for each cell. - // Disadvantage : (the advantave mentioned above) - // need to construct adjacent for each cell. - //---------------------------------------------------------- const double x = fatom1.x(); const double y = fatom1.y(); const double z = fatom1.z(); double x2 = fatom2.x(); double y2 = fatom2.y(); double z2 = fatom2.z(); - //---------------------------------------------------------- - // EXPLAIN : in different case 'in_grid' has different - // meaning. - //---------------------------------------------------------- - // NAME : expand_case | not_expand_case - // in_which_grid 'not available' | one of 27 adjacent grid - // in_which_cell one of all cells | 'not available' - //---------------------------------------------------------- - // The solution here is we save these datas in one structrue - // named : 'in_grid' - //---------------------------------------------------------- - //---------------------------------------------------------- - // EXPlAIN : Calculate distance between two atoms. - //---------------------------------------------------------- double delta_x = x - x2; double delta_y = y - y2; double delta_z = z - z2; @@ -508,14 +278,4 @@ void Grid::Construct_Adjacent_final(const int i, fatom1.addAdjacent(fatom2); } } -// 2015-05-07 -void Grid::delete_vector(int i, int j, int k) -{ - if (expand_flag) - { - if (this->pbc) - { - this->Cell[i][j][k].atom_map.clear(); - } - } -} + diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index 56a8b6ee3c..9f93cd40f9 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -3,7 +3,6 @@ #include "module_cell/unitcell.h" #include "sltk_atom.h" -#include "sltk_atom_input.h" #include "sltk_util.h" #include @@ -25,8 +24,6 @@ struct CellSet // Atom_input : defined elsewhere //========================================================== -class Atom_input; - //========================================================== // CLASS NAME : // Grid : @@ -41,64 +38,37 @@ class Grid Grid(const int& test_grid_in); virtual ~Grid(); - void init(std::ofstream& ofs, const UnitCell& ucell, const Atom_input& input); - - // 2015-05-07 - void delete_vector(int i, int j, int k); +void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, bool pbc_in); // Data bool pbc; // periodic boundary condition - bool expand_flag; + double sradius2; // searching radius squared double sradius; // searching radius - double d_minX; // origin of all cells - double d_minY; - double d_minZ; + int cell_nx; int cell_ny; int cell_nz; - int layer; int true_cell_x; int true_cell_y; int true_cell_z; + double x_min; + double y_min; + double z_min; + double x_max; + double y_max; + double z_max; + + void Check_Expand_Condition(const UnitCell& ucell); + int glayerX; + int glayerX_minus; + int glayerY; + int glayerY_minus; + int glayerZ; + int glayerZ_minus; std::vector>> Cell; // dx , dy ,dz is cell number in each direction,respectly. - void delete_Cell() // it will replace by container soon! - { - if (this->init_cell_flag) - { - for (int i = 0; i < this->cell_nx; i++) - { - for (int j = 0; j < this->cell_ny; j++) - { - this->Cell[i][j].clear(); - } - } - - for (int i = 0; i < this->cell_nx; i++) - { - this->Cell[i].clear(); - } - - this->Cell.clear(); - this->init_cell_flag = false; - } - } - bool init_cell_flag = false; - // LiuXh add 2019-07-15 - double getD_minX() const - { - return d_minX; - } - double getD_minY() const - { - return d_minY; - } - double getD_minZ() const - { - return d_minZ; - } int getCellX() const { @@ -127,21 +97,13 @@ class Grid private: const int test_grid; - //========================================================== - // MEMBER FUNCTIONS : - // Three Main Steps: - // NAME : setMemberVariables (read in datas from Atom_input, - // init cells.) - // NAME : setBoundaryAdjacent( Consider different situations, - // if not_expand case : nature/periodic boundary - // condition , if expand_case) - //========================================================== - void setMemberVariables(std::ofstream& ofs_in, const Atom_input& input); - void setBoundaryAdjacent(std::ofstream& ofs_in, const Atom_input& input); + void setMemberVariables(std::ofstream& ofs_in); + + void setBoundaryAdjacent(std::ofstream& ofs_in); //========================================================== - void Build_Hash_Table(const UnitCell& ucell, const Atom_input& input); + void Build_Hash_Table(const UnitCell& ucell); //========================================================== diff --git a/source/module_cell/module_neighbor/sltk_grid_driver.cpp b/source/module_cell/module_neighbor/sltk_grid_driver.cpp index 8d770ff1f9..239443f12d 100644 --- a/source/module_cell/module_neighbor/sltk_grid_driver.cpp +++ b/source/module_cell/module_neighbor/sltk_grid_driver.cpp @@ -53,18 +53,7 @@ void Grid_Driver::Find_atom( local_adjs->ntype.push_back(atom->getType()); local_adjs->natom.push_back(atom->getNatom()); local_adjs->box.push_back(ModuleBase::Vector3(atom->getCellX(), atom->getCellY(), atom->getCellZ())); - if (expand_flag) - { local_adjs->adjacent_tau.push_back(ModuleBase::Vector3(atom->x(), atom->y(), atom->z())); - } - else - { - local_adjs->adjacent_tau.push_back(Calculate_adjacent_site(atom->x(), atom->y(), atom->z(), - vec1[0], vec2[0], vec3[0], - vec1[1], vec2[1], vec3[1], - vec1[2], vec2[2], vec3[2], - atom->getCellX(), atom->getCellY(), atom->getCellZ())); - }//end if expand_flag local_adjs->adj_num++; } local_adjs->ntype.push_back(ntype); diff --git a/source/module_cell/module_neighbor/sltk_grid_driver.h b/source/module_cell/module_neighbor/sltk_grid_driver.h index fbc683cade..69e7be13fc 100644 --- a/source/module_cell/module_neighbor/sltk_grid_driver.h +++ b/source/module_cell/module_neighbor/sltk_grid_driver.h @@ -7,7 +7,6 @@ #include "module_cell/unitcell.h" #include "module_hamilt_pw/hamilt_pwdft/structure_factor.h" #include "sltk_atom.h" -#include "sltk_atom_input.h" #include "sltk_grid.h" #include diff --git a/source/module_cell/module_neighbor/test/CMakeLists.txt b/source/module_cell/module_neighbor/test/CMakeLists.txt index 72742a4ced..933c00f8bb 100644 --- a/source/module_cell/module_neighbor/test/CMakeLists.txt +++ b/source/module_cell/module_neighbor/test/CMakeLists.txt @@ -9,17 +9,10 @@ AddTest( SOURCES sltk_atom_test.cpp ../sltk_atom.cpp ) -AddTest( - TARGET cell_neighbor_sltk_atom_input - LIBS parameter ${math_libs} base device cell_info - SOURCES sltk_atom_input_test.cpp ../sltk_atom_input.cpp ../sltk_atom.cpp - ../../../module_io/output.cpp -) - AddTest( TARGET cell_neighbor_sltk_grid LIBS parameter ${math_libs} base device cell_info - SOURCES sltk_grid_test.cpp ../sltk_grid.cpp ../sltk_atom_input.cpp ../sltk_atom.cpp + SOURCES sltk_grid_test.cpp ../sltk_grid.cpp ../sltk_atom.cpp ../../../module_io/output.cpp ) @@ -27,6 +20,6 @@ AddTest( TARGET cell_neighbor_sltk_atom_arrange LIBS parameter ${math_libs} base device cell_info SOURCES sltk_atom_arrange_test.cpp ../sltk_atom_arrange.cpp ../sltk_grid_driver.cpp ../sltk_grid.cpp - ../sltk_atom_input.cpp ../sltk_atom.cpp + ../sltk_atom.cpp ../../../module_io/output.cpp ) \ No newline at end of file diff --git a/source/module_cell/module_neighbor/test/sltk_atom_input_test.cpp b/source/module_cell/module_neighbor/test/sltk_atom_input_test.cpp deleted file mode 100644 index bed6207560..0000000000 --- a/source/module_cell/module_neighbor/test/sltk_atom_input_test.cpp +++ /dev/null @@ -1,264 +0,0 @@ -#include "../sltk_atom_input.h" - -#define private public -#include "module_parameter/parameter.h" -#undef private -#include - -#include "gmock/gmock.h" -#include "gtest/gtest.h" -#include "prepare_unitcell.h" - -#ifdef __LCAO -InfoNonlocal::InfoNonlocal() -{ -} -InfoNonlocal::~InfoNonlocal() -{ -} -LCAO_Orbitals::LCAO_Orbitals() -{ -} -LCAO_Orbitals::~LCAO_Orbitals() -{ -} -#endif -Magnetism::Magnetism() -{ - this->tot_magnetization = 0.0; - this->abs_magnetization = 0.0; - this->start_magnetization = nullptr; -} -Magnetism::~Magnetism() -{ - delete[] this->start_magnetization; -} - -/************************************************ - * unit test of sltk_atom_input.cpp - * *********************************************/ - -/** - * - Tested Functions: - * - Constructor: - * - the constructor of Atom_input do almost everything - * - (1) it finds the boundary of atomic coordinates (x_min, x_max, etc.), - * - (2) expands the unitcell according to the boundary condition - * and the searching radius, which is (2 * rcutmax_Phi + 0.01) in - * gamma_only calculation and (2 * (rcutmax_Phi +rcutmax_Beta) + 0.01) - * in multi-k calculation (see atom_arrange::set_sr_NL). And determine - * the number of positive and negative layers in each direction in - * Check_Expand_Condition() - * - (3) records the amount of atoms after expansion (d_amount_expand) and - * sets the expanded lattice grids, including their coordinates and their - * atomic coordinates inside in Expand_Grid() - * - (4) and calculate the number of unitcells in x, y, z directions - * in calculate_cells() - * - Getters: - * - get the values obtained in constructor, for example - * - Clength0() is to get the length in x direction of the expanded big cell in - * unit of lat0 - * - minX() will return x_min if no expanding, or -glayerX_minus in case of - * expansion - * - getCellX() will return the number of unitcells in x direction - * - getCellXLength() will return search radius if no expanding, or 1 in case of - * expansion - * - getRadius() will return the search radius - * - getLatNow() will return lat0 in Bohr - * - getAmount() will return the total number of lattice grids after expansion - */ - -void SetGlobalV() -{ - PARAM.input.test_grid = 0; -} - -class SltkAtomInputTest : public ::testing::Test -{ - protected: - UnitCell* ucell; - UcellTestPrepare utp = UcellTestLib["Si"]; - std::ofstream ofs; - std::ifstream ifs; - bool pbc = true; - double radius = ((8 + 5.01) * 2.0 + 0.01) / 10.2; - int test_atom_in = 0; - std::string output; - void SetUp() - { - SetGlobalV(); - ucell = utp.SetUcellInfo(); - } - void TearDown() - { - delete ucell; - } -}; - -using SltkAtomInputDeathTest = SltkAtomInputTest; - -TEST_F(SltkAtomInputTest, Constructor) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 2; - PARAM.input.test_grid = 1; - Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in); - ofs.close(); - ifs.open("test.out"); - std::string str((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); - EXPECT_THAT(str, testing::HasSubstr("ntype = 1")); - EXPECT_THAT(str, testing::HasSubstr("Amount(atom number) = 2")); - EXPECT_THAT(str, testing::HasSubstr("Periodic_boundary = 1")); - EXPECT_THAT(str, testing::HasSubstr("Searching radius(lat0) = 2.55")); - // EXPECT_THAT(str, testing::HasSubstr("CellLength(unit: lat0) = [ 0.707107, 0.707107, 0.707107 ]")); - EXPECT_THAT(str, testing::HasSubstr("min_tau = [ -0.75, 0, 0 ]")); - EXPECT_THAT(str, testing::HasSubstr("max_tau = [ 0, 0.75, 0.75 ]")); - EXPECT_THAT(str, testing::HasSubstr("glayer+ = [ 6, 6, 6 ]")); - EXPECT_THAT(str, testing::HasSubstr("glayer- = [ 5, 5, 5 ]")); - EXPECT_THAT(str, testing::HasSubstr("CellDim = [ 11, 11, 11 ]")); - ifs.close(); - remove("test.out"); -} - -TEST_F(SltkAtomInputTest, Getters) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 2; - Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in); - EXPECT_TRUE(Atom_inp.getExpandFlag()); - EXPECT_EQ(Atom_inp.getBoundary(), 1); - EXPECT_NEAR(Atom_inp.getRadius(), 2.55196, 1e-5); - EXPECT_DOUBLE_EQ(Atom_inp.minX(), -0.75); - EXPECT_DOUBLE_EQ(Atom_inp.minY(), 0); - EXPECT_DOUBLE_EQ(Atom_inp.minZ(), 0); - EXPECT_EQ(Atom_inp.getCell_nX(), 11); - EXPECT_EQ(Atom_inp.getCell_nY(), 11); - EXPECT_EQ(Atom_inp.getCell_nZ(), 11); - EXPECT_EQ(Atom_inp.getGrid_layerX(), 6); - EXPECT_EQ(Atom_inp.getGrid_layerY(), 6); - EXPECT_EQ(Atom_inp.getGrid_layerZ(), 6); - EXPECT_EQ(Atom_inp.getGrid_layerX_minus(), 5); - EXPECT_EQ(Atom_inp.getGrid_layerY_minus(), 5); - EXPECT_EQ(Atom_inp.getGrid_layerZ_minus(), 5); - ofs.close(); - remove("test.out"); -} - -TEST_F(SltkAtomInputDeathTest, ConstructorWarning1) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 1; - radius = -1; - testing::internal::CaptureStdout(); - EXPECT_EXIT(Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in), - ::testing::ExitedWithCode(1), - ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("search radius < 0,forbidden")); - ofs.close(); - remove("test.out"); -} -/* -TEST_F(SltkAtomInputDeathTest, ConstructorWarning2) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 1; - ucell->atoms[0].taud[1].x = -0.25; - testing::internal::CaptureStdout(); - EXPECT_EXIT(Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in), - ::testing::ExitedWithCode(1), - ""); - // output = testing::internal::GetCapturedStdout(); - // EXPECT_THAT(output, testing::HasSubstr("dminX<0.0")); - ofs.close(); - remove("test.out"); -} - -TEST_F(SltkAtomInputDeathTest, ConstructorWarning3) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 1; - ucell->atoms[0].taud[1].y = -0.25; - testing::internal::CaptureStdout(); - EXPECT_EXIT(Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in), - ::testing::ExitedWithCode(1), - ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("dminY<0.0")); - ofs.close(); - remove("test.out"); -} - -TEST_F(SltkAtomInputDeathTest, ConstructorWarning4) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 1; - ucell->atoms[0].taud[1].z = -0.25; - testing::internal::CaptureStdout(); - EXPECT_EXIT(Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in), - ::testing::ExitedWithCode(1), - ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("dminZ<0.0")); - ofs.close(); - remove("test.out"); -}*/ - -TEST_F(SltkAtomInputTest, ConstructorNoExpand) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 1; - PARAM.input.test_grid = 1; - // this is a bug if radius is too small - // because the expand_flag will be false! - radius = 0; - Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in); - EXPECT_FALSE(Atom_inp.getExpandFlag()); - ofs.close(); - ifs.open("test.out"); - std::string str((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); - EXPECT_THAT(str, testing::HasSubstr("ntype = 1")); - EXPECT_THAT(str, testing::HasSubstr("Amount(atom number) = 2")); - EXPECT_THAT(str, testing::HasSubstr("Periodic_boundary = 1")); - EXPECT_THAT(str, testing::HasSubstr("Searching radius(lat0) = 0")); - // EXPECT_THAT(str, testing::HasSubstr("CellLength(unit: lat0) = [ 0.707107, 0.707107, 0.707107 ]")); - EXPECT_THAT(str, testing::HasSubstr("min_tau = [ -0.75, 0, 0 ]")); - EXPECT_THAT(str, testing::HasSubstr("max_tau = [ 0, 0.75, 0.75 ]")); - EXPECT_THAT(str, testing::HasSubstr("glayer+ = [ 2, 2, 2 ]")); - EXPECT_THAT(str, testing::HasSubstr("glayer- = [ 0, 0, 0 ]")); - ifs.close(); - remove("test.out"); -} - -TEST_F(SltkAtomInputTest, ConstructorSmallSearchRadius) -{ - ofs.open("test.out"); - ucell->check_dtau(); - test_atom_in = 1; - PARAM.input.test_grid = 1; - radius = 0.5; - Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in); - EXPECT_TRUE(Atom_inp.getExpandFlag()); - ofs.close(); - ifs.open("test.out"); - std::string str((std::istreambuf_iterator(ifs)), std::istreambuf_iterator()); - EXPECT_THAT(str, testing::HasSubstr("ntype = 1")); - EXPECT_THAT(str, testing::HasSubstr("Amount(atom number) = 2")); - EXPECT_THAT(str, testing::HasSubstr("Periodic_boundary = 1")); - EXPECT_THAT(str, testing::HasSubstr("Searching radius(lat0) = 0.5")); - // EXPECT_THAT(str, testing::HasSubstr("CellLength(unit: lat0) = [ 0.707107, 0.707107, 0.707107 ]")); - EXPECT_THAT(str, testing::HasSubstr("min_tau = [ -0.75, 0, 0 ]")); - EXPECT_THAT(str, testing::HasSubstr("max_tau = [ 0, 0.75, 0.75 ]")); - EXPECT_THAT(str, testing::HasSubstr("glayer+ = [ 2, 2, 2 ]")); - EXPECT_THAT(str, testing::HasSubstr("glayer- = [ 2, 2, 2 ]")); - EXPECT_THAT(str, testing::HasSubstr("CellDim = [ 4, 4, 4 ]")); - ifs.close(); - remove("test.out"); -} \ No newline at end of file diff --git a/source/module_cell/module_neighbor/test/sltk_grid_test.cpp b/source/module_cell/module_neighbor/test/sltk_grid_test.cpp index af9492e4bc..bdd29428c9 100644 --- a/source/module_cell/module_neighbor/test/sltk_grid_test.cpp +++ b/source/module_cell/module_neighbor/test/sltk_grid_test.cpp @@ -81,10 +81,8 @@ TEST_F(SltkGridTest, Init) ucell->check_dtau(); test_atom_in = 2; PARAM.input.test_grid = 1; - Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in); Grid LatGrid(PARAM.input.test_grid); - LatGrid.init(ofs, *ucell, Atom_inp); - EXPECT_TRUE(LatGrid.init_cell_flag); + LatGrid.init(ofs, *ucell, radius, pbc); EXPECT_EQ(LatGrid.getCellX(), 11); EXPECT_EQ(LatGrid.getCellY(), 11); EXPECT_EQ(LatGrid.getCellZ(), 11); @@ -102,32 +100,21 @@ TEST_F(SltkGridTest, InitSmall) test_atom_in = 2; PARAM.input.test_grid = 1; radius = 0.5; - Atom_input Atom_inp(ofs, *ucell, ucell->nat, ucell->ntype, pbc, radius, test_atom_in); Grid LatGrid(PARAM.input.test_grid); - LatGrid.setMemberVariables(ofs, Atom_inp); - EXPECT_EQ(LatGrid.pbc, Atom_inp.getBoundary()); + LatGrid.init(ofs, *ucell, radius, pbc); + LatGrid.setMemberVariables(ofs); + EXPECT_EQ(LatGrid.pbc, true); EXPECT_TRUE(LatGrid.pbc); - EXPECT_DOUBLE_EQ(LatGrid.sradius2, Atom_inp.getRadius() * Atom_inp.getRadius()); EXPECT_DOUBLE_EQ(LatGrid.sradius2, 0.5 * 0.5); - EXPECT_DOUBLE_EQ(LatGrid.sradius, Atom_inp.getRadius()); EXPECT_DOUBLE_EQ(LatGrid.sradius, 0.5); - // minimal value of x, y, z - EXPECT_DOUBLE_EQ(LatGrid.d_minX, Atom_inp.minX()); - EXPECT_DOUBLE_EQ(LatGrid.d_minY, Atom_inp.minY()); - EXPECT_DOUBLE_EQ(LatGrid.d_minZ, Atom_inp.minZ()); EXPECT_DOUBLE_EQ(LatGrid.true_cell_x, 2); EXPECT_DOUBLE_EQ(LatGrid.true_cell_y, 2); EXPECT_DOUBLE_EQ(LatGrid.true_cell_z, 2); - // number of cells in x, y, z - EXPECT_EQ(LatGrid.cell_nx, Atom_inp.getCell_nX()); - EXPECT_EQ(LatGrid.cell_ny, Atom_inp.getCell_nY()); - EXPECT_EQ(LatGrid.cell_nz, Atom_inp.getCell_nZ()); EXPECT_EQ(LatGrid.cell_nx, 4); EXPECT_EQ(LatGrid.cell_ny, 4); EXPECT_EQ(LatGrid.cell_nz, 4); // init cell flag - EXPECT_TRUE(LatGrid.init_cell_flag); ofs.close(); remove("test.out"); diff --git a/source/module_esolver/esolver_lj.cpp b/source/module_esolver/esolver_lj.cpp index 379aaced17..ccf8d33c30 100644 --- a/source/module_esolver/esolver_lj.cpp +++ b/source/module_esolver/esolver_lj.cpp @@ -88,15 +88,6 @@ void ESolver_LJ::runner(UnitCell& ucell, const int istep) lj_virial(i, j) /= (2.0 * ucell.omega); } } -#ifdef __MPI - atom_arrange::delete_vector( - GlobalV::ofs_running, - PARAM.inp.search_pbc, - grid_neigh, - ucell, - search_radius, - PARAM.inp.test_atom_input); -#endif } double ESolver_LJ::cal_energy() diff --git a/source/module_hamilt_lcao/module_deepks/test/Makefile.Objects b/source/module_hamilt_lcao/module_deepks/test/Makefile.Objects index 9e62c5551a..d962d507a3 100644 --- a/source/module_hamilt_lcao/module_deepks/test/Makefile.Objects +++ b/source/module_hamilt_lcao/module_deepks/test/Makefile.Objects @@ -70,7 +70,6 @@ ORB_gaunt_table.o\ OBJS_NEIGHBOR=sltk_atom_arrange.o\ sltk_atom.o\ -sltk_atom_input.o\ sltk_grid.o\ sltk_grid_driver.o diff --git a/source/module_md/test/CMakeLists.txt b/source/module_md/test/CMakeLists.txt index 16051470fd..2b5e89dc37 100644 --- a/source/module_md/test/CMakeLists.txt +++ b/source/module_md/test/CMakeLists.txt @@ -36,7 +36,6 @@ list(APPEND depend_files ../../module_base/math_integral.cpp ../../module_cell/module_neighbor/sltk_atom_arrange.cpp ../../module_cell/module_neighbor/sltk_atom.cpp - ../../module_cell/module_neighbor/sltk_atom_input.cpp ../../module_cell/module_neighbor/sltk_grid.cpp ../../module_cell/module_neighbor/sltk_grid_driver.cpp ../../module_io/output.cpp From 6e6f1967ef818a687f26cbe670d63777a621e1f4 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 06:49:20 +0000 Subject: [PATCH 02/10] temp remove some case in neighbor --- source/module_cell/module_neighbor/sltk_grid.h | 3 --- 1 file changed, 3 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index 9f93cd40f9..cbac82e5cb 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -111,9 +111,6 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo void Construct_Adjacent_expand_periodic(const int i, const int j, const int k, FAtom& fatom); - void Construct_Adjacent_begin(); - void Construct_Adjacent_nature(const int i, const int j, const int k, FAtom& fatom1); - void Construct_Adjacent_periodic(const int i, const int j, const int k, FAtom& fatom1); void Construct_Adjacent_final(const int i, const int j, const int k, From 34f1ec811de82f4c1be940132befe895b94f738d Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 08:13:49 +0000 Subject: [PATCH 03/10] fix a typo --- source/module_cell/module_neighbor/sltk_grid.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index ce7b254146..8b94a02acb 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -152,8 +152,8 @@ void Grid::setMemberVariables(std::ofstream& ofs_in) } this->true_cell_x = glayerX_minus; - this->true_cell_y = glayerX_minus; - this->true_cell_z = glayerX_minus; + this->true_cell_y = glayerY_minus; + this->true_cell_z = glayerZ_minus; } void Grid::setBoundaryAdjacent(std::ofstream& ofs_in) From ae1ea76850fccf5d92f0e7defc18d9db3ac1fb60 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 08:24:28 +0000 Subject: [PATCH 04/10] remove some variable and format --- .../module_cell/module_neighbor/sltk_atom.cpp | 6 +- .../module_cell/module_neighbor/sltk_atom.h | 30 +--- .../module_neighbor/sltk_atom_arrange.cpp | 132 ++++++++---------- .../module_cell/module_neighbor/sltk_grid.cpp | 20 +-- .../module_neighbor/sltk_grid_driver.cpp | 28 +--- .../module_neighbor/sltk_grid_driver.h | 20 --- .../module_neighbor/test/sltk_atom_test.cpp | 18 +-- 7 files changed, 94 insertions(+), 160 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_atom.cpp b/source/module_cell/module_neighbor/sltk_atom.cpp index 843271cdf1..c6c003f7e0 100644 --- a/source/module_cell/module_neighbor/sltk_atom.cpp +++ b/source/module_cell/module_neighbor/sltk_atom.cpp @@ -4,9 +4,9 @@ /*** Constructors and destructor ***/ FAtom::FAtom() { - d_x = 0.0; - d_y = 0.0; - d_z = 0.0; + x = 0.0; + y = 0.0; + z = 0.0; type = 0; natom = 0; } diff --git a/source/module_cell/module_neighbor/sltk_atom.h b/source/module_cell/module_neighbor/sltk_atom.h index 89789534d1..7c53f346ff 100644 --- a/source/module_cell/module_neighbor/sltk_atom.h +++ b/source/module_cell/module_neighbor/sltk_atom.h @@ -10,10 +10,10 @@ // the type and the index, class FAtom { -private: - double d_x; - double d_y; - double d_z; +public: + double x; + double y; + double z; std::vector adjacent; int type; @@ -22,19 +22,15 @@ class FAtom int cell_x; int cell_y; int cell_z; -public: -//========================================================== -// Default Constructor and deconstructor -//========================================================== FAtom(); FAtom(const double& x_in, const double& y_in, const double& z_in, const int& type_in, const int& natom_in, const int& cell_x_in, const int& cell_y_in, const int& cell_z_in) { - d_x = x_in; - d_y = y_in; - d_z = z_in; + x = x_in; + y = y_in; + z = z_in; type = type_in; natom = natom_in; cell_x = cell_x_in; @@ -52,18 +48,6 @@ class FAtom } const std::vector& getAdjacent() const { return adjacent; } void clearAdjacent() { adjacent.clear(); } -//========================================================== -// MEMBER FUNCTION : -// EXPLAIN : get value -//========================================================== - const double& x() const { return d_x; } - const double& y() const { return d_y; } - const double& z() const { return d_z; } - const int& getType() const { return type;} - const int& getNatom() const { return natom;} - const int& getCellX() const { return cell_x; } - const int& getCellY() const { return cell_y; } - const int& getCellZ() const { return cell_z; } }; #endif diff --git a/source/module_cell/module_neighbor/sltk_atom_arrange.cpp b/source/module_cell/module_neighbor/sltk_atom_arrange.cpp index 68143a5c7e..c9347d5370 100644 --- a/source/module_cell/module_neighbor/sltk_atom_arrange.cpp +++ b/source/module_cell/module_neighbor/sltk_atom_arrange.cpp @@ -1,10 +1,11 @@ #include "sltk_atom_arrange.h" + +#include "module_base/timer.h" #include "module_parameter/parameter.h" #include "sltk_grid.h" #include "sltk_grid_driver.h" -#include "module_base/timer.h" -// update the followig class in near future +// update the followig class in near future #include "module_cell/unitcell.h" atom_arrange::atom_arrange() @@ -15,92 +16,81 @@ atom_arrange::~atom_arrange() { } -double atom_arrange::set_sr_NL( - std::ofstream &ofs_in, - const std::string &output_level, - const double &rcutmax_Phi, - const double &rcutmax_Beta, - const bool gamma_only_local) +double atom_arrange::set_sr_NL(std::ofstream& ofs_in, + const std::string& output_level, + const double& rcutmax_Phi, + const double& rcutmax_Beta, + const bool gamma_only_local) { - ModuleBase::TITLE("atom_arrange","set_sr_NL"); - - if(output_level != "m") //xiaohui add 'output_level', 2015-09-16 - { - ofs_in << "\n\n\n\n"; - ofs_in << " >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>" << std::endl; - ofs_in << " | |" << std::endl; - ofs_in << " | Search adjacent atoms: |" << std::endl; - ofs_in << " | Set the adjacent atoms for each atom and set the periodic boundary |" << std::endl; - ofs_in << " | condition for the atoms on real space FFT grid. For k-dependent |" << std::endl; - ofs_in << " | algorithm, we also need to set the sparse H and S matrix element |" << std::endl; - ofs_in << " | for each atom. |" << std::endl; - ofs_in << " | |" << std::endl; - ofs_in << " <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<" << std::endl; - ofs_in << "\n\n\n\n"; - } - - - //xiaohui add 'output_level' line, 2015-09-16 - if(output_level != "m") { ofs_in << "\n SETUP SEARCHING RADIUS FOR PROGRAM TO SEARCH ADJACENT ATOMS" << std::endl; -} - if(output_level != "m") { ofs_in << std::setprecision(3); -} - if(output_level != "m") { ModuleBase::GlobalFunc::OUT(ofs_in,"longest orb rcut (Bohr)",rcutmax_Phi); -} + ModuleBase::TITLE("atom_arrange", "set_sr_NL"); + + if (output_level != "m") // xiaohui add 'output_level', 2015-09-16 + { + ofs_in << "\n\n\n\n"; + ofs_in << " >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>" << std::endl; + ofs_in << " | |" << std::endl; + ofs_in << " | Search adjacent atoms: |" << std::endl; + ofs_in << " | Set the adjacent atoms for each atom and set the periodic boundary |" << std::endl; + ofs_in << " | condition for the atoms on real space FFT grid. For k-dependent |" << std::endl; + ofs_in << " | algorithm, we also need to set the sparse H and S matrix element |" << std::endl; + ofs_in << " | for each atom. |" << std::endl; + ofs_in << " | |" << std::endl; + ofs_in << " <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<" << std::endl; + ofs_in << "\n\n\n\n"; + ofs_in << "\n SETUP SEARCHING RADIUS FOR PROGRAM TO SEARCH ADJACENT ATOMS" << std::endl; + ofs_in << std::setprecision(3); + ModuleBase::GlobalFunc::OUT(ofs_in, "longest orb rcut (Bohr)", rcutmax_Phi); + ModuleBase::GlobalFunc::OUT(ofs_in, "longest nonlocal projector rcut (Bohr)", rcutmax_Beta); + } -// std::cout << " LONGEST NL PROJ RCUT : " << longest_nl_proj_rcut << std::endl; - if(output_level != "m") { ModuleBase::GlobalFunc::OUT(ofs_in,"longest nonlocal projector rcut (Bohr)", rcutmax_Beta); -} + // check in use_overlap_matrix, + double sr = 0.0; + if (gamma_only_local) + { + sr = 2 * rcutmax_Phi + 0.01; + } + else + { + sr = 2 * (rcutmax_Phi + rcutmax_Beta) + 0.01; // 0.01 is added to make safe. + // sr = 2 * longest_orb_rcut + 0.01; + } - // check in use_overlap_matrix, - double sr = 0.0; - if(gamma_only_local) - { - sr = 2 * rcutmax_Phi + 0.01; - } - else - { - sr = 2 * (rcutmax_Phi +rcutmax_Beta) + 0.01; // 0.01 is added to make safe. - //sr = 2 * longest_orb_rcut + 0.01; - } - - return sr; + return sr; } -void atom_arrange::search( - const bool pbc_flag, - std::ofstream &ofs_in, - Grid_Driver &grid_d, - const UnitCell &ucell, - const double &search_radius_bohr, - const int &test_atom_in, - const bool test_only) +void atom_arrange::search(const bool pbc_flag, + std::ofstream& ofs_in, + Grid_Driver& grid_d, + const UnitCell& ucell, + const double& search_radius_bohr, + const int& test_atom_in, + const bool test_only) { - ModuleBase::TITLE("atom_arrange", "search"); - ModuleBase::timer::tick("atom_arrange","search"); - assert( search_radius_bohr > 0.0 ); + ModuleBase::TITLE("atom_arrange", "search"); + ModuleBase::timer::tick("atom_arrange", "search"); + assert(search_radius_bohr > 0.0); - ModuleBase::GlobalFunc::OUT(ofs_in,"searching radius is (Bohr))", search_radius_bohr); - ModuleBase::GlobalFunc::OUT(ofs_in,"searching radius unit is (Bohr))",ucell.lat0); + ModuleBase::GlobalFunc::OUT(ofs_in, "searching radius is (Bohr))", search_radius_bohr); + ModuleBase::GlobalFunc::OUT(ofs_in, "searching radius unit is (Bohr))", ucell.lat0); - assert(ucell.nat > 0); + assert(ucell.nat > 0); - const double radius_lat0unit = search_radius_bohr / ucell.lat0; + const double radius_lat0unit = search_radius_bohr / ucell.lat0; ModuleBase::timer::tick("atom_arrange", "grid_d.init"); - grid_d.init(ofs_in, ucell, radius_lat0unit, pbc_flag); + grid_d.init(ofs_in, ucell, radius_lat0unit, pbc_flag); ModuleBase::timer::tick("atom_arrange", "grid_d.init"); ModuleBase::timer::tick("atom_arrange", "search"); - // test the adjacent atoms and the box. - if(test_only) - { - std::cout << "radius_lat0unit = " << radius_lat0unit << std::endl; - std::cout << "search_radius_bohr = " << search_radius_bohr << std::endl; + // test the adjacent atoms and the box. + if (test_only) + { + std::cout << "radius_lat0unit = " << radius_lat0unit << std::endl; + std::cout << "search_radius_bohr = " << search_radius_bohr << std::endl; - ofs_in << " " << std::setw(5) << "Type" << std::setw(5) << "Atom" << std::setw(8) << "AdjNum" << std::endl; + ofs_in << " " << std::setw(5) << "Type" << std::setw(5) << "Atom" << std::setw(8) << "AdjNum" << std::endl; std::cout << std::setw(8) << "Labels" << std::setw(15) << "tau.x" << std::setw(15) << "tau.y" << std::setw(15) << "tau.z" << std::setw(8) << "box.x" << std::setw(8) << "box.y" << std::setw(8) << "box.z" << std::endl; diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index 8b94a02acb..b5d8a1ba70 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -199,11 +199,11 @@ void Grid::Build_Hash_Table(const UnitCell& ucell) double z = ucell.atoms[i].tau[j].z + vec1[2] * ix + vec2[2] * iy + vec3[2] * iz; FAtom atom(x, y, z, i, j, ix, iy, iz); int a, b, c; - a = atom.getCellX() + true_cell_x; - b = atom.getCellY() + true_cell_y; - c = atom.getCellZ() + true_cell_z; + a = ix + glayerX_minus; + b = iy + glayerY_minus; + c = iz + glayerZ_minus; - this->Cell[a][b][c].atom_map[atom.getType()][atom.getNatom()] = atom; + this->Cell[a][b][c].atom_map[atom.type][atom.natom] = atom; } } } @@ -260,12 +260,12 @@ void Grid::Construct_Adjacent_final(const int i, const int k2, FAtom& fatom2) { - const double x = fatom1.x(); - const double y = fatom1.y(); - const double z = fatom1.z(); - double x2 = fatom2.x(); - double y2 = fatom2.y(); - double z2 = fatom2.z(); + const double x = fatom1.x; + const double y = fatom1.y; + const double z = fatom1.z; + double x2 = fatom2.x; + double y2 = fatom2.y; + double z2 = fatom2.z; double delta_x = x - x2; double delta_y = y - y2; diff --git a/source/module_cell/module_neighbor/sltk_grid_driver.cpp b/source/module_cell/module_neighbor/sltk_grid_driver.cpp index 239443f12d..9d55a40256 100644 --- a/source/module_cell/module_neighbor/sltk_grid_driver.cpp +++ b/source/module_cell/module_neighbor/sltk_grid_driver.cpp @@ -40,20 +40,15 @@ void Grid_Driver::Find_atom( AdjacentAtomInfo* local_adjs = adjs == nullptr ? &this->adj_info : adjs; local_adjs->clear(); const std::vector & all_atom = Cell[this->true_cell_x][this->true_cell_y][this->true_cell_z].atom_map[ntype][nnumber].getAdjacent(); - //std::cout << "ntype = "<< ntype << " atom size = " << all_atom.size() << std::endl; - - ModuleBase::Vector3 vec1(ucell.latvec.e11, ucell.latvec.e12, ucell.latvec.e13); - ModuleBase::Vector3 vec2(ucell.latvec.e21, ucell.latvec.e22, ucell.latvec.e23); - ModuleBase::Vector3 vec3(ucell.latvec.e31, ucell.latvec.e32, ucell.latvec.e33); for(const FAtom * atom : all_atom) { // std::cout << "atom type = " << atom.getType() << " number = " << atom.getNatom() << " box = " << atom.getCellX() << " " << atom.getCellY() << " " << atom.getCellZ() // << " tau = " << atom.x() << " " << atom.y() << " " << atom.z() << std::endl; - local_adjs->ntype.push_back(atom->getType()); - local_adjs->natom.push_back(atom->getNatom()); - local_adjs->box.push_back(ModuleBase::Vector3(atom->getCellX(), atom->getCellY(), atom->getCellZ())); - local_adjs->adjacent_tau.push_back(ModuleBase::Vector3(atom->x(), atom->y(), atom->z())); + local_adjs->ntype.push_back(atom->type); + local_adjs->natom.push_back(atom->natom); + local_adjs->box.push_back(ModuleBase::Vector3(atom->cell_x, atom->cell_y, atom->cell_z)); + local_adjs->adjacent_tau.push_back(ModuleBase::Vector3(atom->x, atom->y, atom->z)); local_adjs->adj_num++; } local_adjs->ntype.push_back(ntype); @@ -61,25 +56,10 @@ void Grid_Driver::Find_atom( local_adjs->box.push_back(ModuleBase::Vector3(0, 0, 0)); local_adjs->adjacent_tau.push_back(ModuleBase::Vector3(cartesian_pos.x, cartesian_pos.y, cartesian_pos.z)); - ModuleBase::timer::tick("Grid_Driver","Find_atom"); return; } -ModuleBase::Vector3 Grid_Driver::Calculate_adjacent_site(const double x, const double y, const double z, - const double &box11, const double &box12, const double &box13, - const double &box21, const double &box22, const double &box23, - const double &box31, const double &box32, const double &box33, - const short box_x, const short box_y, const short box_z) const -{ - ModuleBase::Vector3 adjacent_site(0, 0, 0); - adjacent_site.x = x + box_x * box11 + box_y * box12 + box_z * box13; - adjacent_site.y = y + box_x * box21 + box_y * box22 + box_z * box23; - adjacent_site.z = z + box_x * box31 + box_y * box32 + box_z * box33; - - return adjacent_site; -} - // filter_adjs delete not adjacent atoms in adjs void filter_adjs(const std::vector& is_adj, AdjacentAtomInfo& adjs) { diff --git a/source/module_cell/module_neighbor/sltk_grid_driver.h b/source/module_cell/module_neighbor/sltk_grid_driver.h index 69e7be13fc..980db70602 100644 --- a/source/module_cell/module_neighbor/sltk_grid_driver.h +++ b/source/module_cell/module_neighbor/sltk_grid_driver.h @@ -103,26 +103,6 @@ class Grid_Driver : public Grid mutable AdjacentAtomInfo adj_info; const int test_deconstructor; // caoyu reconst 2021-05-24 - - //========================================================== - // MEMBER FUNCTIONS : - // NAME : Calculate_adjacent_site - //========================================================== - ModuleBase::Vector3 Calculate_adjacent_site(const double x, - const double y, - const double z, - const double& box11, - const double& box12, - const double& box13, - const double& box21, - const double& box22, - const double& box23, - const double& box31, - const double& box32, - const double& box33, - const short box_x, // three dimensions of the target box - const short box_y, - const short box_z) const; }; namespace GlobalC diff --git a/source/module_cell/module_neighbor/test/sltk_atom_test.cpp b/source/module_cell/module_neighbor/test/sltk_atom_test.cpp index 57e09c9a7c..41e5db2208 100644 --- a/source/module_cell/module_neighbor/test/sltk_atom_test.cpp +++ b/source/module_cell/module_neighbor/test/sltk_atom_test.cpp @@ -32,7 +32,7 @@ TEST_F(SltkAtomTest, AllocateAdjacentSet) test.clearAdjacent(); FAtom test_temp(1.0, 2.0, 3.0, 4, 5, 0, 1, 2); test.addAdjacent(test_temp); - EXPECT_EQ(test.getAdjacent().front()->getType(), 4); + EXPECT_EQ(test.getAdjacent().front()->type, 4); } @@ -40,12 +40,12 @@ TEST_F(SltkAtomTest, SetterGetters) { FAtom test_temp(1.0, 2.0, 3.0, 4, 5, 0, 1, 2); - EXPECT_DOUBLE_EQ(test_temp.x(), 1.0); - EXPECT_DOUBLE_EQ(test_temp.y(), 2.0); - EXPECT_DOUBLE_EQ(test_temp.z(), 3.0); - EXPECT_EQ(test_temp.getType(), 4); - EXPECT_EQ(test_temp.getNatom(), 5); - EXPECT_EQ(test_temp.getCellX(), 0); - EXPECT_EQ(test_temp.getCellY(), 1); - EXPECT_EQ(test_temp.getCellZ(), 2); + EXPECT_DOUBLE_EQ(test_temp.x, 1.0); + EXPECT_DOUBLE_EQ(test_temp.y, 2.0); + EXPECT_DOUBLE_EQ(test_temp.z, 3.0); + EXPECT_EQ(test_temp.type, 4); + EXPECT_EQ(test_temp.natom, 5); + EXPECT_EQ(test_temp.cell_x, 0); + EXPECT_EQ(test_temp.cell_y, 1); + EXPECT_EQ(test_temp.cell_z, 2); } From 5697844f72c0d7974e2123c30e61d1175a99ef90 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 08:32:41 +0000 Subject: [PATCH 05/10] add all adj info --- source/module_cell/module_neighbor/sltk_atom.h | 8 -------- source/module_cell/module_neighbor/sltk_grid.cpp | 10 +++++++++- source/module_cell/module_neighbor/sltk_grid.h | 2 +- .../module_cell/module_neighbor/sltk_grid_driver.cpp | 2 +- .../module_neighbor/test/sltk_atom_test.cpp | 9 --------- 5 files changed, 11 insertions(+), 20 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_atom.h b/source/module_cell/module_neighbor/sltk_atom.h index 7c53f346ff..6acc9cf122 100644 --- a/source/module_cell/module_neighbor/sltk_atom.h +++ b/source/module_cell/module_neighbor/sltk_atom.h @@ -14,7 +14,6 @@ class FAtom double x; double y; double z; - std::vector adjacent; int type; int natom; @@ -39,15 +38,8 @@ class FAtom } ~FAtom() { - adjacent.clear(); } - void addAdjacent(FAtom& atom_in) - { - adjacent.push_back( &atom_in); - } - const std::vector& getAdjacent() const { return adjacent; } - void clearAdjacent() { adjacent.clear(); } }; #endif diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index b5d8a1ba70..7587093a1e 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -53,6 +53,13 @@ void Grid::init(std::ofstream& ofs_in, const UnitCell& ucell, const double radiu } this->Check_Expand_Condition(ucell); this->setMemberVariables(ofs_in); + + all_adj_info.resize(ucell.ntype); + for (int i = 0; i < ucell.ntype; i++) + { + all_adj_info[i].resize(ucell.atoms[i].na); + } + this->Build_Hash_Table(ucell); this->setBoundaryAdjacent(ofs_in); } @@ -154,6 +161,7 @@ void Grid::setMemberVariables(std::ofstream& ofs_in) this->true_cell_x = glayerX_minus; this->true_cell_y = glayerY_minus; this->true_cell_z = glayerZ_minus; + } void Grid::setBoundaryAdjacent(std::ofstream& ofs_in) @@ -275,7 +283,7 @@ void Grid::Construct_Adjacent_final(const int i, if (dr != 0.0 && dr <= this->sradius2) { - fatom1.addAdjacent(fatom2); + all_adj_info[fatom1.type][fatom1.natom].push_back(&fatom2); } } diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index cbac82e5cb..d6c00ded17 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -69,7 +69,7 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo int glayerZ_minus; std::vector>> Cell; // dx , dy ,dz is cell number in each direction,respectly. - + std::vector>> all_adj_info; int getCellX() const { return cell_nx; diff --git a/source/module_cell/module_neighbor/sltk_grid_driver.cpp b/source/module_cell/module_neighbor/sltk_grid_driver.cpp index 9d55a40256..ec678cb6ab 100644 --- a/source/module_cell/module_neighbor/sltk_grid_driver.cpp +++ b/source/module_cell/module_neighbor/sltk_grid_driver.cpp @@ -39,7 +39,7 @@ void Grid_Driver::Find_atom( // store result in member adj_info when parameter adjs is NULL AdjacentAtomInfo* local_adjs = adjs == nullptr ? &this->adj_info : adjs; local_adjs->clear(); - const std::vector & all_atom = Cell[this->true_cell_x][this->true_cell_y][this->true_cell_z].atom_map[ntype][nnumber].getAdjacent(); + const std::vector & all_atom = all_adj_info[ntype][nnumber]; for(const FAtom * atom : all_atom) { diff --git a/source/module_cell/module_neighbor/test/sltk_atom_test.cpp b/source/module_cell/module_neighbor/test/sltk_atom_test.cpp index 41e5db2208..bbebaf7def 100644 --- a/source/module_cell/module_neighbor/test/sltk_atom_test.cpp +++ b/source/module_cell/module_neighbor/test/sltk_atom_test.cpp @@ -27,15 +27,6 @@ class SltkAtomTest : public testing::Test FAtom test; }; -TEST_F(SltkAtomTest, AllocateAdjacentSet) -{ - test.clearAdjacent(); - FAtom test_temp(1.0, 2.0, 3.0, 4, 5, 0, 1, 2); - test.addAdjacent(test_temp); - EXPECT_EQ(test.getAdjacent().front()->type, 4); -} - - TEST_F(SltkAtomTest, SetterGetters) { FAtom test_temp(1.0, 2.0, 3.0, 4, 5, 0, 1, 2); From 1456a8c1023ea3c8424ea9a94ed63466b0a66388 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 08:35:29 +0000 Subject: [PATCH 06/10] remove cellset --- .../module_cell/module_neighbor/sltk_grid.cpp | 20 +++++-------------- .../module_cell/module_neighbor/sltk_grid.h | 9 +-------- 2 files changed, 6 insertions(+), 23 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index 7587093a1e..f6a9eab6ee 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -5,16 +5,6 @@ #include "module_base/memory.h" #include "module_base/timer.h" -//================== -// Class CellSet -//================== -CellSet::CellSet() -{ - in_grid[0] = 0; - in_grid[1] = 0; - in_grid[2] = 0; -} - Grid::Grid(const int& test_grid_in) : test_grid(test_grid_in) { } @@ -180,10 +170,10 @@ void Grid::Build_Hash_Table(const UnitCell& ucell) { for (int k = 0; k < cell_nz; k++) { - Cell[i][j][k].atom_map.resize(ucell.ntype); + Cell[i][j][k].resize(ucell.ntype); for (int it = 0; it < ucell.ntype; ++it) { - Cell[i][j][k].atom_map[it].resize(ucell.atoms[it].na); + Cell[i][j][k][it].resize(ucell.atoms[it].na); } } } @@ -211,7 +201,7 @@ void Grid::Build_Hash_Table(const UnitCell& ucell) b = iy + glayerY_minus; c = iz + glayerZ_minus; - this->Cell[a][b][c].atom_map[atom.type][atom.natom] = atom; + this->Cell[a][b][c][atom.type][atom.natom] = atom; } } } @@ -224,7 +214,7 @@ void Grid::Construct_Adjacent_expand(const int true_i, const int true_j, const i { ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); - for (auto& atom_vector: this->Cell[true_i][true_j][true_k].atom_map) + for (auto& atom_vector: this->Cell[true_i][true_j][true_k]) { for (auto& fatom: atom_vector) { @@ -245,7 +235,7 @@ void Grid::Construct_Adjacent_expand_periodic(const int true_i, const int true_j { for (int k = 0; k < this->cell_nz; k++) { - for (auto& atom_vector: this->Cell[i][j][k].atom_map) + for (auto& atom_vector: this->Cell[i][j][k]) { for (auto& fatom2: atom_vector) { diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index d6c00ded17..7f592b6a79 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -12,13 +12,6 @@ typedef std::vector> AtomMap; -struct CellSet -{ - AtomMap atom_map; - int in_grid[3]; - CellSet(); -}; - //========================================================== // CLASS NAME : // Atom_input : defined elsewhere @@ -68,7 +61,7 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo int glayerZ; int glayerZ_minus; - std::vector>> Cell; // dx , dy ,dz is cell number in each direction,respectly. + std::vector>> Cell; // dx , dy ,dz is cell number in each direction,respectly. std::vector>> all_adj_info; int getCellX() const { From 91f7da2521932b239bcc14f57ab34177487c7a48 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 08:54:26 +0000 Subject: [PATCH 07/10] fix a dont understand bug. and unit test fix --- .../module_cell/module_neighbor/sltk_grid.cpp | 31 +++---------------- .../module_cell/module_neighbor/sltk_grid.h | 14 +-------- .../module_neighbor/test/sltk_grid_test.cpp | 12 +++---- 3 files changed, 12 insertions(+), 45 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index f6a9eab6ee..4996a97ab6 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -163,21 +163,6 @@ void Grid::Build_Hash_Table(const UnitCell& ucell) { ModuleBase::timer::tick("Grid", "Build_Hash_Table"); - // TODO in case expand == false, the following code is over malloc - for (int i = 0; i < cell_nx; i++) - { - for (int j = 0; j < cell_ny; j++) - { - for (int k = 0; k < cell_nz; k++) - { - Cell[i][j][k].resize(ucell.ntype); - for (int it = 0; it < ucell.ntype; ++it) - { - Cell[i][j][k][it].resize(ucell.atoms[it].na); - } - } - } - } ModuleBase::Vector3 vec1(ucell.latvec.e11, ucell.latvec.e12, ucell.latvec.e13); ModuleBase::Vector3 vec2(ucell.latvec.e21, ucell.latvec.e22, ucell.latvec.e23); ModuleBase::Vector3 vec3(ucell.latvec.e31, ucell.latvec.e32, ucell.latvec.e33); @@ -201,7 +186,7 @@ void Grid::Build_Hash_Table(const UnitCell& ucell) b = iy + glayerY_minus; c = iz + glayerZ_minus; - this->Cell[a][b][c][atom.type][atom.natom] = atom; + this->Cell[a][b][c].push_back(atom); } } } @@ -214,12 +199,9 @@ void Grid::Construct_Adjacent_expand(const int true_i, const int true_j, const i { ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); - for (auto& atom_vector: this->Cell[true_i][true_j][true_k]) + for (auto& fatom: this->Cell[true_i][true_j][true_k]) { - for (auto& fatom: atom_vector) - { - Construct_Adjacent_expand_periodic(true_i, true_j, true_k, fatom); - } + Construct_Adjacent_expand_periodic(true_i, true_j, true_k, fatom); } ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); } @@ -235,12 +217,9 @@ void Grid::Construct_Adjacent_expand_periodic(const int true_i, const int true_j { for (int k = 0; k < this->cell_nz; k++) { - for (auto& atom_vector: this->Cell[i][j][k]) + for (auto& fatom2: this->Cell[i][j][k]) { - for (auto& fatom2: atom_vector) - { - Construct_Adjacent_final(true_i, true_j, true_k, fatom, i, j, k, fatom2); - } + Construct_Adjacent_final(true_i, true_j, true_k, fatom, i, j, k, fatom2); } } } diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index 7f592b6a79..e5467656e7 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -10,18 +10,6 @@ #include #include -typedef std::vector> AtomMap; - -//========================================================== -// CLASS NAME : -// Atom_input : defined elsewhere -//========================================================== - -//========================================================== -// CLASS NAME : -// Grid : -//========================================================== - class Grid { public: @@ -61,7 +49,7 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo int glayerZ; int glayerZ_minus; - std::vector>> Cell; // dx , dy ,dz is cell number in each direction,respectly. + std::vector>>> Cell; // dx , dy ,dz is cell number in each direction,respectly. std::vector>> all_adj_info; int getCellX() const { diff --git a/source/module_cell/module_neighbor/test/sltk_grid_test.cpp b/source/module_cell/module_neighbor/test/sltk_grid_test.cpp index bdd29428c9..ad20e36488 100644 --- a/source/module_cell/module_neighbor/test/sltk_grid_test.cpp +++ b/source/module_cell/module_neighbor/test/sltk_grid_test.cpp @@ -108,12 +108,12 @@ TEST_F(SltkGridTest, InitSmall) EXPECT_DOUBLE_EQ(LatGrid.sradius2, 0.5 * 0.5); EXPECT_DOUBLE_EQ(LatGrid.sradius, 0.5); - EXPECT_DOUBLE_EQ(LatGrid.true_cell_x, 2); - EXPECT_DOUBLE_EQ(LatGrid.true_cell_y, 2); - EXPECT_DOUBLE_EQ(LatGrid.true_cell_z, 2); - EXPECT_EQ(LatGrid.cell_nx, 4); - EXPECT_EQ(LatGrid.cell_ny, 4); - EXPECT_EQ(LatGrid.cell_nz, 4); + EXPECT_DOUBLE_EQ(LatGrid.true_cell_x, 1); + EXPECT_DOUBLE_EQ(LatGrid.true_cell_y, 1); + EXPECT_DOUBLE_EQ(LatGrid.true_cell_z, 1); + EXPECT_EQ(LatGrid.cell_nx, 3); + EXPECT_EQ(LatGrid.cell_ny, 3); + EXPECT_EQ(LatGrid.cell_nz, 3); // init cell flag ofs.close(); From e7261dc3746c664f8b0c585a401c670e612ea8d6 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 09:08:05 +0000 Subject: [PATCH 08/10] reconstuct expand condition --- .../module_cell/module_neighbor/sltk_grid.cpp | 65 +++++++------------ .../module_cell/module_neighbor/sltk_grid.h | 12 +--- 2 files changed, 27 insertions(+), 50 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index 4996a97ab6..ef2c74f318 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -51,7 +51,7 @@ void Grid::init(std::ofstream& ofs_in, const UnitCell& ucell, const double radiu } this->Build_Hash_Table(ucell); - this->setBoundaryAdjacent(ofs_in); + this->Construct_Adjacent_expand(ucell); } void Grid::Check_Expand_Condition(const UnitCell& ucell) { @@ -61,42 +61,33 @@ void Grid::Check_Expand_Condition(const UnitCell& ucell) { return; } + double lat03 = ucell.lat0 * ucell.lat0 * ucell.lat0; double a23_1 = ucell.latvec.e22 * ucell.latvec.e33 - ucell.latvec.e23 * ucell.latvec.e32; double a23_2 = ucell.latvec.e21 * ucell.latvec.e33 - ucell.latvec.e23 * ucell.latvec.e31; double a23_3 = ucell.latvec.e21 * ucell.latvec.e32 - ucell.latvec.e22 * ucell.latvec.e31; double a23_norm = sqrt(a23_1 * a23_1 + a23_2 * a23_2 + a23_3 * a23_3); - double extend_v = a23_norm * sradius; - double extend_d1 = extend_v / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; - int extend_d11 = static_cast(extend_d1); - // 2016-09-05, LiuXh - if (extend_d1 - extend_d11 > 0.0) - { - extend_d11 += 1; - } + double extend_d1 = a23_norm * sradius; + extend_d1 = extend_d1 / ucell.omega; + extend_d1 = extend_d1 * lat03; + int extend_d11 = std::ceil(extend_d1); double a31_1 = ucell.latvec.e32 * ucell.latvec.e13 - ucell.latvec.e33 * ucell.latvec.e12; double a31_2 = ucell.latvec.e31 * ucell.latvec.e13 - ucell.latvec.e33 * ucell.latvec.e11; double a31_3 = ucell.latvec.e31 * ucell.latvec.e12 - ucell.latvec.e32 * ucell.latvec.e11; double a31_norm = sqrt(a31_1 * a31_1 + a31_2 * a31_2 + a31_3 * a31_3); - double extend_d2 = a31_norm * sradius / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; - int extend_d22 = static_cast(extend_d2); - // 2016-09-05, LiuXh - if (extend_d2 - extend_d22 > 0.0) - { - extend_d22 += 1; - } + double extend_d2 = a31_norm * sradius; + extend_d2 = extend_d2 / ucell.omega; + extend_d2 = extend_d2 * lat03; + int extend_d22 = std::ceil(extend_d2); double a12_1 = ucell.latvec.e12 * ucell.latvec.e23 - ucell.latvec.e13 * ucell.latvec.e22; double a12_2 = ucell.latvec.e11 * ucell.latvec.e23 - ucell.latvec.e13 * ucell.latvec.e21; double a12_3 = ucell.latvec.e11 * ucell.latvec.e22 - ucell.latvec.e12 * ucell.latvec.e21; double a12_norm = sqrt(a12_1 * a12_1 + a12_2 * a12_2 + a12_3 * a12_3); - double extend_d3 = a12_norm * sradius / ucell.omega * ucell.lat0 * ucell.lat0 * ucell.lat0; - int extend_d33 = static_cast(extend_d3); - // 2016-09-05, LiuXh - if (extend_d3 - extend_d33 > 0.0) - { - extend_d33 += 1; - } + double extend_d3 = a12_norm * sradius; + extend_d3 = extend_d3 / ucell.omega; + extend_d3 = extend_d3* lat03; + int extend_d33 = std::ceil(extend_d3); glayerX = extend_d11 + 1; glayerY = extend_d22 + 1; @@ -154,11 +145,6 @@ void Grid::setMemberVariables(std::ofstream& ofs_in) } -void Grid::setBoundaryAdjacent(std::ofstream& ofs_in) -{ - this->Construct_Adjacent_expand(true_cell_x, true_cell_y, true_cell_z); -} - void Grid::Build_Hash_Table(const UnitCell& ucell) { ModuleBase::timer::tick("Grid", "Build_Hash_Table"); @@ -195,18 +181,22 @@ void Grid::Build_Hash_Table(const UnitCell& ucell) ModuleBase::timer::tick("Grid", "Build_Hash_Table"); } -void Grid::Construct_Adjacent_expand(const int true_i, const int true_j, const int true_k) +void Grid::Construct_Adjacent_expand(const UnitCell & ucell) { ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); - for (auto& fatom: this->Cell[true_i][true_j][true_k]) + for (int i = 0; i < ucell.ntype; i++) { - Construct_Adjacent_expand_periodic(true_i, true_j, true_k, fatom); + for (int j = 0; j < ucell.atoms[i].na; j++) + { + FAtom fatom(ucell.atoms[i].tau[j].x, ucell.atoms[i].tau[j].y, ucell.atoms[i].tau[j].z, i, j, 0, 0, 0); + Construct_Adjacent_expand_periodic(fatom); + } } ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand"); } -void Grid::Construct_Adjacent_expand_periodic(const int true_i, const int true_j, const int true_k, FAtom& fatom) +void Grid::Construct_Adjacent_expand_periodic(const FAtom& fatom) { // if (test_grid)ModuleBase::TITLE(ofs_running, "Grid", "Construct_Adjacent_expand_periodic"); ModuleBase::timer::tick("Grid", "Construct_Adjacent_expand_periodic"); @@ -219,7 +209,7 @@ void Grid::Construct_Adjacent_expand_periodic(const int true_i, const int true_j { for (auto& fatom2: this->Cell[i][j][k]) { - Construct_Adjacent_final(true_i, true_j, true_k, fatom, i, j, k, fatom2); + Construct_Adjacent_final(fatom, fatom2); } } } @@ -228,14 +218,7 @@ void Grid::Construct_Adjacent_expand_periodic(const int true_i, const int true_j } -void Grid::Construct_Adjacent_final(const int i, - const int j, - const int k, - FAtom& fatom1, - const int i2, - const int j2, - const int k2, - FAtom& fatom2) +void Grid::Construct_Adjacent_final(const FAtom& fatom1, FAtom& fatom2) { const double x = fatom1.x; const double y = fatom1.y; diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index e5467656e7..88e84a1adc 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -88,17 +88,11 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo //========================================================== - void Construct_Adjacent_expand(const int i, const int j, const int k); + void Construct_Adjacent_expand(const UnitCell & ucell); - void Construct_Adjacent_expand_periodic(const int i, const int j, const int k, FAtom& fatom); + void Construct_Adjacent_expand_periodic(const FAtom& fatom); - void Construct_Adjacent_final(const int i, - const int j, - const int k, - FAtom& fatom1, - const int i2, - const int j2, - const int k2, + void Construct_Adjacent_final(const FAtom& fatom1, FAtom& fatom2); }; From 3cecc5f89745791f4e3b83908b3ada18ccb5acaf Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 09:16:45 +0000 Subject: [PATCH 09/10] remove true cell index --- .../module_cell/module_neighbor/sltk_grid.cpp | 5 -- .../module_cell/module_neighbor/sltk_grid.h | 48 +++++++++++-------- .../module_neighbor/test/sltk_grid_test.cpp | 6 --- .../hamilt_lcaodft/spar_dh.cpp | 16 ++----- source/module_io/cal_r_overlap_R.cpp | 16 ++----- source/module_io/to_wannier90_lcao.cpp | 21 ++------ 6 files changed, 40 insertions(+), 72 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index ef2c74f318..8df7c1136d 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -138,11 +138,6 @@ void Grid::setMemberVariables(std::ofstream& ofs_in) Cell[i][j].resize(cell_nz); } } - - this->true_cell_x = glayerX_minus; - this->true_cell_y = glayerY_minus; - this->true_cell_z = glayerZ_minus; - } void Grid::Build_Hash_Table(const UnitCell& ucell) diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index 88e84a1adc..3ed90723f3 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -31,9 +31,6 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo int cell_ny; int cell_nz; - int true_cell_x; - int true_cell_y; - int true_cell_z; double x_min; double y_min; double z_min; @@ -51,43 +48,52 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo std::vector>>> Cell; // dx , dy ,dz is cell number in each direction,respectly. std::vector>> all_adj_info; - int getCellX() const + + int getGlayerX() { - return cell_nx; + return glayerX; } - int getCellY() const + int getGlayerY() { - return cell_ny; + return glayerY; } - int getCellZ() const + int getGlayerZ() { - return cell_nz; + return glayerZ; } - int getTrueCellX() const + int getGlayerX_minus() { - return true_cell_x; + return glayerX_minus; } - int getTrueCellY() const + int getGlayerY_minus() { - return true_cell_y; + return glayerY_minus; } - int getTrueCellZ() const + int getGlayerZ_minus() { - return true_cell_z; + return glayerZ_minus; + } + + + int getCellX() const + { + return cell_nx; + } + int getCellY() const + { + return cell_ny; + } + int getCellZ() const + { + return cell_nz; } private: const int test_grid; void setMemberVariables(std::ofstream& ofs_in); - - void setBoundaryAdjacent(std::ofstream& ofs_in); - - //========================================================== void Build_Hash_Table(const UnitCell& ucell); - //========================================================== - void Construct_Adjacent_expand(const UnitCell & ucell); void Construct_Adjacent_expand_periodic(const FAtom& fatom); diff --git a/source/module_cell/module_neighbor/test/sltk_grid_test.cpp b/source/module_cell/module_neighbor/test/sltk_grid_test.cpp index ad20e36488..238ea8b3c3 100644 --- a/source/module_cell/module_neighbor/test/sltk_grid_test.cpp +++ b/source/module_cell/module_neighbor/test/sltk_grid_test.cpp @@ -86,9 +86,6 @@ TEST_F(SltkGridTest, Init) EXPECT_EQ(LatGrid.getCellX(), 11); EXPECT_EQ(LatGrid.getCellY(), 11); EXPECT_EQ(LatGrid.getCellZ(), 11); - EXPECT_EQ(LatGrid.getTrueCellX(), 5); - EXPECT_EQ(LatGrid.getTrueCellY(), 5); - EXPECT_EQ(LatGrid.getTrueCellZ(), 5); ofs.close(); remove("test.out"); } @@ -108,9 +105,6 @@ TEST_F(SltkGridTest, InitSmall) EXPECT_DOUBLE_EQ(LatGrid.sradius2, 0.5 * 0.5); EXPECT_DOUBLE_EQ(LatGrid.sradius, 0.5); - EXPECT_DOUBLE_EQ(LatGrid.true_cell_x, 1); - EXPECT_DOUBLE_EQ(LatGrid.true_cell_y, 1); - EXPECT_DOUBLE_EQ(LatGrid.true_cell_z, 1); EXPECT_EQ(LatGrid.cell_nx, 3); EXPECT_EQ(LatGrid.cell_ny, 3); EXPECT_EQ(LatGrid.cell_nz, 3); diff --git a/source/module_hamilt_lcao/hamilt_lcaodft/spar_dh.cpp b/source/module_hamilt_lcao/hamilt_lcaodft/spar_dh.cpp index 9710ec4b26..58e0bd092c 100644 --- a/source/module_hamilt_lcao/hamilt_lcaodft/spar_dh.cpp +++ b/source/module_hamilt_lcao/hamilt_lcaodft/spar_dh.cpp @@ -68,21 +68,13 @@ void sparse_format::cal_dH(const UnitCell& ucell, void sparse_format::set_R_range(std::set>& all_R_coor, Grid_Driver& grid) { - const int RminX = int(-grid.getTrueCellX()); - const int RminY = int(-grid.getTrueCellY()); - const int RminZ = int(-grid.getTrueCellZ()); - - const int Rx = grid.getCellX(); - const int Ry = grid.getCellY(); - const int Rz = grid.getCellZ(); - - for (int ix = 0; ix < Rx; ix++) + for (int ix = -grid.getGlayerX_minus(); ix < grid.getGlayerX(); ix++) { - for (int iy = 0; iy < Ry; iy++) + for (int iy = -grid.getGlayerY_minus(); iy < grid.getGlayerY(); iy++) { - for (int iz = 0; iz < Rz; iz++) + for (int iz = -grid.getGlayerZ_minus(); iz < grid.getGlayerZ(); iz++) { - Abfs::Vector3_Order temp_R(ix + RminX, iy + RminY, iz + RminZ); + Abfs::Vector3_Order temp_R(ix, iy, iz); all_R_coor.insert(temp_R); } } diff --git a/source/module_io/cal_r_overlap_R.cpp b/source/module_io/cal_r_overlap_R.cpp index b43c181576..80cdf5a477 100644 --- a/source/module_io/cal_r_overlap_R.cpp +++ b/source/module_io/cal_r_overlap_R.cpp @@ -245,23 +245,15 @@ void cal_r_overlap_R::out_rR(const UnitCell& ucell, const int& istep) ModuleBase::timer::tick("cal_r_overlap_R", "out_rR"); int step = istep; - // set R coor range - int R_minX = int(-GlobalC::GridD.getTrueCellX()); - int R_minY = int(-GlobalC::GridD.getTrueCellY()); - int R_minZ = int(-GlobalC::GridD.getTrueCellZ()); - - int R_x = GlobalC::GridD.getCellX(); - int R_y = GlobalC::GridD.getCellY(); - int R_z = GlobalC::GridD.getCellZ(); std::set> all_R_coor; - for (int ix = 0; ix < R_x; ix++) + for (int ix = -GlobalC::GridD.getGlayerX_minus(); ix < GlobalC::GridD.getGlayerX(); ix++) { - for (int iy = 0; iy < R_y; iy++) + for (int iy = -GlobalC::GridD.getGlayerY_minus(); iy < GlobalC::GridD.getGlayerY(); iy++) { - for (int iz = 0; iz < R_z; iz++) + for (int iz = -GlobalC::GridD.getGlayerZ_minus(); iz < GlobalC::GridD.getGlayerZ(); iz++) { - Abfs::Vector3_Order temp_R(ix + R_minX, iy + R_minY, iz + R_minZ); + Abfs::Vector3_Order temp_R(ix, iy, iz); all_R_coor.insert(temp_R); } } diff --git a/source/module_io/to_wannier90_lcao.cpp b/source/module_io/to_wannier90_lcao.cpp index 72651d6991..b3705abd87 100644 --- a/source/module_io/to_wannier90_lcao.cpp +++ b/source/module_io/to_wannier90_lcao.cpp @@ -305,26 +305,15 @@ void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) void toWannier90_LCAO::set_R_coor(const UnitCell& ucell) { - int R_minX = int(-GlobalC::GridD.getTrueCellX()); - int R_minY = int(-GlobalC::GridD.getTrueCellY()); - int R_minZ = int(-GlobalC::GridD.getTrueCellZ()); - - int R_x = GlobalC::GridD.getCellX(); - int R_y = GlobalC::GridD.getCellY(); - int R_z = GlobalC::GridD.getCellZ(); - - int R_num = R_x * R_y * R_z; - R_coor_car.resize(R_num); - int count = 0; - for (int ix = 0; ix < R_x; ix++) + for (int ix = -GlobalC::GridD.getGlayerX_minus(); ix < GlobalC::GridD.getGlayerX(); ix++) { - for (int iy = 0; iy < R_y; iy++) + for (int iy = -GlobalC::GridD.getGlayerY_minus(); iy < GlobalC::GridD.getGlayerY(); iy++) { - for (int iz = 0; iz < R_z; iz++) + for (int iz = -GlobalC::GridD.getGlayerZ_minus(); iz < GlobalC::GridD.getGlayerZ(); iz++) { - ModuleBase::Vector3 tmpR(ix + R_minX, iy + R_minY, iz + R_minZ); - R_coor_car[count] = tmpR * ucell.latvec; + Abfs::Vector3_Order temp_R(ix, iy, iz); + R_coor_car.push_back(temp_R * ucell.latvec); count++; } } From 3407b72ee08442d882984ba1d3568a67edd64452 Mon Sep 17 00:00:00 2001 From: haochong zhang Date: Thu, 5 Dec 2024 09:41:52 +0000 Subject: [PATCH 10/10] fix a memory bug --- source/module_cell/module_neighbor/sltk_grid.cpp | 2 +- source/module_cell/module_neighbor/sltk_grid.h | 2 +- .../module_cell/module_neighbor/sltk_grid_driver.cpp | 12 ++++++------ 3 files changed, 8 insertions(+), 8 deletions(-) diff --git a/source/module_cell/module_neighbor/sltk_grid.cpp b/source/module_cell/module_neighbor/sltk_grid.cpp index 8df7c1136d..512d7c9a38 100644 --- a/source/module_cell/module_neighbor/sltk_grid.cpp +++ b/source/module_cell/module_neighbor/sltk_grid.cpp @@ -230,7 +230,7 @@ void Grid::Construct_Adjacent_final(const FAtom& fatom1, FAtom& fatom2) if (dr != 0.0 && dr <= this->sradius2) { - all_adj_info[fatom1.type][fatom1.natom].push_back(&fatom2); + all_adj_info[fatom1.type][fatom1.natom].push_back(fatom2); } } diff --git a/source/module_cell/module_neighbor/sltk_grid.h b/source/module_cell/module_neighbor/sltk_grid.h index 3ed90723f3..057f2ef52e 100644 --- a/source/module_cell/module_neighbor/sltk_grid.h +++ b/source/module_cell/module_neighbor/sltk_grid.h @@ -47,7 +47,7 @@ void init(std::ofstream& ofs_in, const UnitCell& ucell, const double radius, boo int glayerZ_minus; std::vector>>> Cell; // dx , dy ,dz is cell number in each direction,respectly. - std::vector>> all_adj_info; + std::vector>> all_adj_info; int getGlayerX() { diff --git a/source/module_cell/module_neighbor/sltk_grid_driver.cpp b/source/module_cell/module_neighbor/sltk_grid_driver.cpp index ec678cb6ab..9900517490 100644 --- a/source/module_cell/module_neighbor/sltk_grid_driver.cpp +++ b/source/module_cell/module_neighbor/sltk_grid_driver.cpp @@ -39,16 +39,16 @@ void Grid_Driver::Find_atom( // store result in member adj_info when parameter adjs is NULL AdjacentAtomInfo* local_adjs = adjs == nullptr ? &this->adj_info : adjs; local_adjs->clear(); - const std::vector & all_atom = all_adj_info[ntype][nnumber]; + const std::vector & all_atom = all_adj_info[ntype][nnumber]; - for(const FAtom * atom : all_atom) + for(const FAtom atom : all_atom) { // std::cout << "atom type = " << atom.getType() << " number = " << atom.getNatom() << " box = " << atom.getCellX() << " " << atom.getCellY() << " " << atom.getCellZ() // << " tau = " << atom.x() << " " << atom.y() << " " << atom.z() << std::endl; - local_adjs->ntype.push_back(atom->type); - local_adjs->natom.push_back(atom->natom); - local_adjs->box.push_back(ModuleBase::Vector3(atom->cell_x, atom->cell_y, atom->cell_z)); - local_adjs->adjacent_tau.push_back(ModuleBase::Vector3(atom->x, atom->y, atom->z)); + local_adjs->ntype.push_back(atom.type); + local_adjs->natom.push_back(atom.natom); + local_adjs->box.push_back(ModuleBase::Vector3(atom.cell_x, atom.cell_y, atom.cell_z)); + local_adjs->adjacent_tau.push_back(ModuleBase::Vector3(atom.x, atom.y, atom.z)); local_adjs->adj_num++; } local_adjs->ntype.push_back(ntype);