From 1b2738f0afb9476bef4fbddd7d24088eef314c44 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Thu, 18 Sep 2025 21:49:33 +0800 Subject: [PATCH 01/22] change BFGS name and make lattice_change_cg and ions_move_cg shorter --- source/source_io/read_input_item_relax.cpp | 2 +- source/source_relax/ions_move_cg.cpp | 80 ++++++---------------- source/source_relax/ions_move_cg.h | 1 + source/source_relax/ions_move_methods.cpp | 4 +- source/source_relax/lattice_change_cg.cpp | 67 ++++++------------ source/source_relax/lattice_change_cg.h | 2 + 6 files changed, 50 insertions(+), 106 deletions(-) diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index 036252f23b..c647da9eb0 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -12,7 +12,7 @@ void ReadInput::item_relax() item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; read_sync_string(input.relax_method); item.check_value = [](const Input_Item& item, const Parameter& para) { - const std::vector relax_methods = {"cg", "bfgs", "sd", "cg_bfgs","bfgs_trad","lbfgs"}; + const std::vector relax_methods = {"cg", "bfgs_old", "sd", "cg_bfgs","bfgs","lbfgs"}; if (std::find(relax_methods.begin(),relax_methods.end(), para.input.relax_method)==relax_methods.end()) { const std::string warningstr = nofound_str(relax_methods, "relax_method"); diff --git a/source/source_relax/ions_move_cg.cpp b/source/source_relax/ions_move_cg.cpp index 7da54c3831..7b1dbd1420 100644 --- a/source/source_relax/ions_move_cg.cpp +++ b/source/source_relax/ions_move_cg.cpp @@ -1,5 +1,4 @@ #include "ions_move_cg.h" - #include "ions_move_basic.h" #include "source_base/global_function.h" #include "source_base/global_variable.h" @@ -77,20 +76,16 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const // ncggrad is a parameter to control the cg method , every ten cg directions, // we change the direction back to the steepest descent method static int ncggrad = 0; - static double fa = 0.0; static double fb = 0.0; static double fc = 0.0; - static double xa = 0.0; static double xb = 0.0; static double xc = 0.0; static double xpt = 0.0; static double steplength = 0.0; static double fmax = 0.0; - static int nbrent = 0; - // some arrays double *pos = new double[dim]; @@ -100,7 +95,6 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const double *cg_grad = new double[dim]; double best_x = 0.0; double fmin = 0.0; - int flag = 0; ModuleBase::GlobalFunc::ZEROS(pos, dim); @@ -137,7 +131,6 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const { Ions_Move_Basic::check_converged(ucell, grad); } - if (Ions_Move_Basic::converged) { Ions_Move_Basic::terminate(ucell); @@ -155,21 +148,7 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const flag); // we use the last direction ,the last grad and the grad now to get the direction now ncggrad++; - double norm = 0.0; - for (int i = 0; i < dim; ++i) - { - norm += pow(cg_grad[i], 2); - } - norm = sqrt(norm); - - if (norm != 0.0) - { - for (int i = 0; i < dim; ++i) - { - cg_gradn[i] = cg_grad[i] / norm; - } - } - + normalize(cg_gradn, cg_grad, dim); setup_move(move0, cg_gradn, steplength); // move the atom position Ions_Move_Basic::move_atoms(ucell, move0, pos); @@ -181,7 +160,6 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const f_cal(move0, move0, dim, xb); // xb = trial steplength f_cal(move0, grad, dim, fa); // fa is the projection force in this direction - fmax = fa; sd = false; @@ -204,7 +182,6 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const double e1 = etot_in; f_cal(move0, grad, dim, fb); f_cal(move0, move0, dim, xb); - if ((std::abs(fb) < std::abs((fa) / 10.0))) { sd = true; @@ -214,21 +191,7 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const goto CG_begin; } - double norm = 0.0; - for (int i = 0; i < dim; ++i) - { - norm += pow(cg_grad0[i], 2); - } - norm = sqrt(norm); - - if (norm != 0.0) - { - for (int i = 0; i < dim; ++i) - { - cg_gradn[i] = cg_grad0[i] / norm; - } - } - + normalize(cg_gradn, cg_grad, dim); third_order(e0, e1, fa, fb, xb, best_x); // cubic interpolation if (best_x > 6 * xb || best_x < (-xb)) @@ -238,7 +201,6 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const setup_move(move, cg_gradn, best_x); Ions_Move_Basic::move_atoms(ucell, move, pos); - trial = false; xa = 0; f_cal(move0, move, dim, xc); @@ -250,7 +212,6 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const { double xtemp, ftemp; f_cal(move0, grad, dim, fc); - fmin = std::abs(fc); nbrent++; @@ -275,24 +236,9 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const goto CG_begin; } - double norm = 0.0; - for (int i = 0; i < dim; ++i) - { - norm += pow(cg_grad0[i], 2); - } - norm = sqrt(norm); - - if (norm != 0.0) - { - for (int i = 0; i < dim; ++i) - { - cg_gradn[i] = cg_grad0[i] / norm; - } - } - + normalize(cg_gradn, cg_grad, dim); setup_move(move, cg_gradn, best_x); Ions_Move_Basic::move_atoms(ucell, move, pos); - Ions_Move_Basic::relax_bfgs_init = xc; } } @@ -342,7 +288,6 @@ void Ions_Move_CG::setup_cg_grad(double *grad, cgp_gp += cg_grad0[i] * grad0[i]; cgp_g += cg_grad0[i] * grad[i]; } - assert(g_gp != 0.0); const double gamma1 = gg / gp_gp; // FR // const double gamma2 = -(gg - g_gp)/(cgp_g - cgp_gp); //CW @@ -508,3 +453,22 @@ void Ions_Move_CG::setup_move(double *move, double *cg_gradn, const double &trus } return; } + +void Ions_Move_CG::normalize(double *cg_gradn, const double *cg_grad, int dim) +{ + double norm = 0.0; + for (int i = 0; i < dim; ++i) + { + norm += pow(cg_grad[i], 2); + } + norm = sqrt(norm); + + if (norm != 0.0) + { + for (int i = 0; i < dim; ++i) + { + cg_gradn[i] = cg_grad[i] / norm; + } + } + return; +} diff --git a/source/source_relax/ions_move_cg.h b/source/source_relax/ions_move_cg.h index d143902486..662ef0c14f 100644 --- a/source/source_relax/ions_move_cg.h +++ b/source/source_relax/ions_move_cg.h @@ -38,6 +38,7 @@ class Ions_Move_CG const double &fb, const double x, double &best_x); + void normalize(double *cg_gradn, const double *cg_grad, int dim); }; #endif diff --git a/source/source_relax/ions_move_methods.cpp b/source/source_relax/ions_move_methods.cpp index 50ab3286c1..1d2198616f 100644 --- a/source/source_relax/ions_move_methods.cpp +++ b/source/source_relax/ions_move_methods.cpp @@ -60,7 +60,7 @@ void Ions_Move_Methods::cal_movement(const int &istep, // Ions_Move_Basic::istep = istep; Ions_Move_Basic::istep = force_step; - if (Ions_Move_Basic::relax_method == "bfgs") + if (Ions_Move_Basic::relax_method == "bfgs_old") { // move_ions // output tau @@ -79,7 +79,7 @@ void Ions_Move_Methods::cal_movement(const int &istep, { cg.start(ucell, f, etot); // added by pengfei 13-8-10 } - else if(Ions_Move_Basic::relax_method == "bfgs_trad") + else if(Ions_Move_Basic::relax_method == "bfgs") { bfgs_trad.relax_step(f,ucell); } diff --git a/source/source_relax/lattice_change_cg.cpp b/source/source_relax/lattice_change_cg.cpp index 94a1bf581d..fec3f0be9d 100644 --- a/source/source_relax/lattice_change_cg.cpp +++ b/source/source_relax/lattice_change_cg.cpp @@ -165,21 +165,7 @@ void Lattice_Change_CG::start(UnitCell &ucell, const ModuleBase::matrix &stress_ flag); // we use the last direction ,the last grad and the grad now to get the direction now ncggrad++; - double norm = 0.0; - for (int i = 0; i < dim; ++i) - { - norm += pow(cg_grad[i], 2); - } - norm = sqrt(norm); - - if (norm != 0.0) - { - for (int i = 0; i < dim; ++i) - { - cg_gradn[i] = cg_grad[i] / norm; - } - } - + normalize(cg_gradn, cg_grad, dim); setup_move(move0, cg_gradn, steplength); // move the atom position Lattice_Change_Basic::change_lattice(ucell, move0, lat); @@ -214,21 +200,7 @@ void Lattice_Change_CG::start(UnitCell &ucell, const ModuleBase::matrix &stress_ goto CG_begin; } - double norm = 0.0; - for (int i = 0; i < dim; ++i) - { - norm += pow(cg_grad0[i], 2); - } - norm = sqrt(norm); - - if (norm != 0.0) - { - for (int i = 0; i < dim; ++i) - { - cg_gradn[i] = cg_grad0[i] / norm; - } - } - + normalize(cg_gradn, cg_grad, dim); third_order(e0, e1, fa, fb, xb, best_x); // cubic interpolation if (best_x > 6 * xb || best_x < (-xb)) @@ -279,21 +251,7 @@ void Lattice_Change_CG::start(UnitCell &ucell, const ModuleBase::matrix &stress_ goto CG_begin; } - double norm = 0.0; - for (int i = 0; i < dim; ++i) - { - norm += pow(cg_grad0[i], 2); - } - norm = sqrt(norm); - - if (norm != 0.0) - { - for (int i = 0; i < dim; ++i) - { - cg_gradn[i] = cg_grad0[i] / norm; - } - } - + normalize(cg_gradn, cg_grad, dim); setup_move(move, cg_gradn, best_x); Lattice_Change_Basic::change_lattice(ucell, move, lat); @@ -512,3 +470,22 @@ void Lattice_Change_CG::setup_move(double *move, double *cg_gradn, const double } return; } + +void Lattice_Change_CG::normalize(double *cg_gradn, const double *cg_grad, int dim) +{ + double norm = 0.0; + for (int i = 0; i < dim; ++i) + { + norm += pow(cg_grad[i], 2); + } + norm = sqrt(norm); + + if (norm != 0.0) + { + for (int i = 0; i < dim; ++i) + { + cg_gradn[i] = cg_grad[i] / norm; + } + } + return; +} diff --git a/source/source_relax/lattice_change_cg.h b/source/source_relax/lattice_change_cg.h index 3a6e724637..b93f22c81b 100644 --- a/source/source_relax/lattice_change_cg.h +++ b/source/source_relax/lattice_change_cg.h @@ -40,6 +40,8 @@ class Lattice_Change_CG const double &fb, const double x, double &best_x); + + void normalize(double *cg_gradn, const double *cg_grad, int dim); }; #endif From e3baf1596834feb10133ce0b39abe72f60d792d5 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 23 Sep 2025 11:05:37 +0800 Subject: [PATCH 02/22] change input parameters --- .../module_parameter/input_parameter.h | 6 ++++ source/source_io/read_input_item_relax.cpp | 34 +++++++++++++++---- tests/01_PW/058_PW_RE_MB/INPUT | 2 +- tests/01_PW/059_PW_RE_MB_traj/INPUT | 2 +- 4 files changed, 35 insertions(+), 9 deletions(-) diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index 95e5cd4712..e643b4633f 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -6,6 +6,11 @@ #include #include +struct RelaxMethodParam { + std::string method; + int param = 1; +}; + // It stores all input parameters both defined in INPUT file and not defined in // INPUT file struct Input_para @@ -149,6 +154,7 @@ struct Input_para // ============== #Parameters (4.Relaxation) =========================== std::string relax_method = "cg"; ///< methods to move_ion: sd, bfgs, cg... + RelaxMethodParam relax_method_param; bool relax_new = true; bool relax = false; ///< allow relaxation along the specific direction double relax_scale_force = 0.5; diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index c647da9eb0..f1041a0fa9 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -5,21 +5,41 @@ namespace ModuleIO { + + void ReadInput::item_relax() { { Input_Item item("relax_method"); - item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; + item.annotation = "cg [param]; bfgs [param]; ..."; read_sync_string(input.relax_method); + item.read_value = [](const Input_Item& item, Parameter& para) { + std::istringstream iss(para.input.relax_method); + iss >> para.input.relax_method_param.method; + if (!(iss >> para.input.relax_method_param.param)) { + para.input.relax_method_param.param = 1; + } +}; item.check_value = [](const Input_Item& item, const Parameter& para) { - const std::vector relax_methods = {"cg", "bfgs_old", "sd", "cg_bfgs","bfgs","lbfgs"}; - if (std::find(relax_methods.begin(),relax_methods.end(), para.input.relax_method)==relax_methods.end()) - { - const std::string warningstr = nofound_str(relax_methods, "relax_method"); - ModuleBase::WARNING_QUIT("ReadInput", warningstr); - } + const std::vector relax_methods = {"cg", "sd","bfgs_old" "cg_bfgs","bfgs","lbfgs"}; + if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { + const std::string warningstr = nofound_str(relax_methods, "relax_method"); + ModuleBase::WARNING_QUIT("ReadInput", warningstr); + } }; this->add_item(item); + // Input_Item item("relax_method"); + // item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; + // read_sync_string(input.relax_method); + // item.check_value = [](const Input_Item& item, const Parameter& para) { + // const std::vector relax_methods = {"cg", "bfgs_old", "sd", "cg_bfgs","bfgs","lbfgs"}; + // if (std::find(relax_methods.begin(),relax_methods.end(), para.input.relax_method)==relax_methods.end()) + // { + // const std::string warningstr = nofound_str(relax_methods, "relax_method"); + // ModuleBase::WARNING_QUIT("ReadInput", warningstr); + // } + // }; + // this->add_item(item); } { Input_Item item("relax_new"); diff --git a/tests/01_PW/058_PW_RE_MB/INPUT b/tests/01_PW/058_PW_RE_MB/INPUT index a501fb7922..0e156711e3 100644 --- a/tests/01_PW/058_PW_RE_MB/INPUT +++ b/tests/01_PW/058_PW_RE_MB/INPUT @@ -8,7 +8,7 @@ calculation relax relax_nmax 2 cal_force 1 force_thr_ev 0.01 -relax_method bfgs +relax_method bfgs_old relax_new 0 # Self-Consistent Field diff --git a/tests/01_PW/059_PW_RE_MB_traj/INPUT b/tests/01_PW/059_PW_RE_MB_traj/INPUT index 2d3058510e..59254d64f2 100644 --- a/tests/01_PW/059_PW_RE_MB_traj/INPUT +++ b/tests/01_PW/059_PW_RE_MB_traj/INPUT @@ -16,7 +16,7 @@ scf_nmax 100 relax_nmax 2 cal_force 1 force_thr_ev 0.01 -relax_method bfgs_trad +relax_method bfgs #Parameters (4.Basis) basis_type pw From 036af83ab7603ba9377efafc012d34e0ee48509b Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 23 Sep 2025 19:24:52 +0800 Subject: [PATCH 03/22] change BFGS name and make lattice_change_cg and ions_move_cg shorter --- source/source_io/read_input_item_relax.cpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index f1041a0fa9..d91bae1b8d 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -14,12 +14,12 @@ void ReadInput::item_relax() item.annotation = "cg [param]; bfgs [param]; ..."; read_sync_string(input.relax_method); item.read_value = [](const Input_Item& item, Parameter& para) { - std::istringstream iss(para.input.relax_method); - iss >> para.input.relax_method_param.method; - if (!(iss >> para.input.relax_method_param.param)) { + std::istringstream iss(para.input.relax_method); + iss >> para.input.relax_method_param.method; + if (!(iss >> para.input.relax_method_param.param)) { para.input.relax_method_param.param = 1; - } -}; + } + }; item.check_value = [](const Input_Item& item, const Parameter& para) { const std::vector relax_methods = {"cg", "sd","bfgs_old" "cg_bfgs","bfgs","lbfgs"}; if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { From 4d7ef80a53f785e00ce6a967264aa2ccfdd3e233 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 23 Sep 2025 21:34:38 +0800 Subject: [PATCH 04/22] change input parameters --- source/source_io/input_conv.cpp | 1 + .../module_parameter/input_parameter.h | 2 +- source/source_io/read_input_item_relax.cpp | 25 ++++++++++++++----- source/source_relax/ions_move_basic.cpp | 1 + source/source_relax/ions_move_basic.h | 2 ++ source/source_relax/ions_move_methods.cpp | 9 +++---- tests/01_PW/058_PW_RE_MB/INPUT | 2 +- 7 files changed, 29 insertions(+), 13 deletions(-) diff --git a/source/source_io/input_conv.cpp b/source/source_io/input_conv.cpp index 952ca873d1..fec9695537 100644 --- a/source/source_io/input_conv.cpp +++ b/source/source_io/input_conv.cpp @@ -190,6 +190,7 @@ void Input_Conv::Convert() Ions_Move_Basic::relax_bfgs_init = PARAM.inp.relax_bfgs_init; Ions_Move_Basic::out_stru = PARAM.inp.out_stru; // mohan add 2012-03-23 Ions_Move_Basic::relax_method = PARAM.inp.relax_method; + Ions_Move_Basic::relax_method_param = PARAM.inp.relax_method_param; Lattice_Change_Basic::fixed_axes = PARAM.inp.fixed_axes; diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index e643b4633f..6aebfbe8b6 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -8,7 +8,7 @@ struct RelaxMethodParam { std::string method; - int param = 1; + std::string param; }; // It stores all input parameters both defined in INPUT file and not defined in diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index d91bae1b8d..eb4f16c4e4 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -11,23 +11,36 @@ void ReadInput::item_relax() { { Input_Item item("relax_method"); - item.annotation = "cg [param]; bfgs [param]; ..."; + item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; read_sync_string(input.relax_method); item.read_value = [](const Input_Item& item, Parameter& para) { - std::istringstream iss(para.input.relax_method); - iss >> para.input.relax_method_param.method; - if (!(iss >> para.input.relax_method_param.param)) { - para.input.relax_method_param.param = 1; + para.input.relax_method_param.method=item.str_values[0]; + para.input.relax_method = para.input.relax_method_param.method; + if(item.get_size()==1) + { + para.input.relax_method_param.param = "1"; } + else + { + para.input.relax_method_param.param = item.str_values[1]; + } + + // std::istringstream iss(item.str_values[0]); + // iss >> para.input.relax_method_param.method; + // if (!(iss >> para.input.relax_method_param.param)) { + // std::cout << "No parameter provided for relax_method_param.param, default to 1" << std::endl; + // para.input.relax_method_param.param = "1"; + // } }; item.check_value = [](const Input_Item& item, const Parameter& para) { - const std::vector relax_methods = {"cg", "sd","bfgs_old" "cg_bfgs","bfgs","lbfgs"}; + const std::vector relax_methods = {"cg", "sd", "cg_bfgs","bfgs","lbfgs"}; if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { const std::string warningstr = nofound_str(relax_methods, "relax_method"); ModuleBase::WARNING_QUIT("ReadInput", warningstr); } }; this->add_item(item); + // Input_Item item("relax_method"); // item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; // read_sync_string(input.relax_method); diff --git a/source/source_relax/ions_move_basic.cpp b/source/source_relax/ions_move_basic.cpp index 5b3e1b58c6..0e89cdae43 100644 --- a/source/source_relax/ions_move_basic.cpp +++ b/source/source_relax/ions_move_basic.cpp @@ -24,6 +24,7 @@ double Ions_Move_Basic::best_xxx = 1.0; int Ions_Move_Basic::out_stru = 0; std::string Ions_Move_Basic::relax_method = "bfgs"; +RelaxMethodParam Ions_Move_Basic::relax_method_param; void Ions_Move_Basic::setup_gradient(const UnitCell &ucell, const ModuleBase::matrix &force, double *pos, double *grad) { diff --git a/source/source_relax/ions_move_basic.h b/source/source_relax/ions_move_basic.h index f6966d2d17..507e6b0b0d 100644 --- a/source/source_relax/ions_move_basic.h +++ b/source/source_relax/ions_move_basic.h @@ -3,6 +3,7 @@ #include "source_base/matrix.h" #include "source_cell/unitcell.h" +#include "source_io/module_parameter/input_parameter.h" namespace Ions_Move_Basic { @@ -23,6 +24,7 @@ extern double relax_bfgs_init; // initial value of trust radius, extern double best_xxx; // the last step length of cg , we use it as bfgs`s initial step length extern std::string relax_method; // relaxation method, extern int out_stru; // output the structure or not +extern RelaxMethodParam relax_method_param; // funny way to pass this parameter, but nevertheless //---------------------------------------------------------------------------- diff --git a/source/source_relax/ions_move_methods.cpp b/source/source_relax/ions_move_methods.cpp index 1d2198616f..4b442bc637 100644 --- a/source/source_relax/ions_move_methods.cpp +++ b/source/source_relax/ions_move_methods.cpp @@ -16,7 +16,7 @@ void Ions_Move_Methods::allocate(const int &natom) { Ions_Move_Basic::dim = natom * 3; - if (Ions_Move_Basic::relax_method == "bfgs") + if (Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param != "1") { this->bfgs.allocate(); } @@ -33,7 +33,7 @@ void Ions_Move_Methods::allocate(const int &natom) this->cg.allocate(); this->bfgs.allocate(); // added by pengfei 13-8-8 } - else if(Ions_Move_Basic::relax_method == "bfgs_trad") + else if(Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param == "1") { this->bfgs_trad.allocate(natom); } @@ -59,8 +59,7 @@ void Ions_Move_Methods::cal_movement(const int &istep, // Ions_Move_Basic::istep = istep; Ions_Move_Basic::istep = force_step; - - if (Ions_Move_Basic::relax_method == "bfgs_old") + if (Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param != "1") { // move_ions // output tau @@ -79,7 +78,7 @@ void Ions_Move_Methods::cal_movement(const int &istep, { cg.start(ucell, f, etot); // added by pengfei 13-8-10 } - else if(Ions_Move_Basic::relax_method == "bfgs") + else if(Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param == "1") { bfgs_trad.relax_step(f,ucell); } diff --git a/tests/01_PW/058_PW_RE_MB/INPUT b/tests/01_PW/058_PW_RE_MB/INPUT index 0e156711e3..cb65f32dae 100644 --- a/tests/01_PW/058_PW_RE_MB/INPUT +++ b/tests/01_PW/058_PW_RE_MB/INPUT @@ -8,7 +8,7 @@ calculation relax relax_nmax 2 cal_force 1 force_thr_ev 0.01 -relax_method bfgs_old +relax_method bfgs 2 relax_new 0 # Self-Consistent Field From 621748d7d73e136ab61c954c7143e119669c511c Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 23 Sep 2025 22:26:08 +0800 Subject: [PATCH 05/22] change input parameters --- source/source_io/read_input_item_relax.cpp | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index eb4f16c4e4..f6ba9dceb9 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -14,16 +14,19 @@ void ReadInput::item_relax() item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; read_sync_string(input.relax_method); item.read_value = [](const Input_Item& item, Parameter& para) { - para.input.relax_method_param.method=item.str_values[0]; - para.input.relax_method = para.input.relax_method_param.method; if(item.get_size()==1) { + para.input.relax_method_param.method=item.str_values[0]; + para.input.relax_method = para.input.relax_method_param.method; para.input.relax_method_param.param = "1"; } - else + else if(item.get_size()==2) { + para.input.relax_method_param.method=item.str_values[0]; + para.input.relax_method = para.input.relax_method_param.method; para.input.relax_method_param.param = item.str_values[1]; } + }; // std::istringstream iss(item.str_values[0]); // iss >> para.input.relax_method_param.method; @@ -31,7 +34,6 @@ void ReadInput::item_relax() // std::cout << "No parameter provided for relax_method_param.param, default to 1" << std::endl; // para.input.relax_method_param.param = "1"; // } - }; item.check_value = [](const Input_Item& item, const Parameter& para) { const std::vector relax_methods = {"cg", "sd", "cg_bfgs","bfgs","lbfgs"}; if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { From 4d73760489498c0ad6dffe8274e0720ee7aab6e4 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Wed, 24 Sep 2025 18:26:44 +0800 Subject: [PATCH 06/22] change input parameters --- source/source_io/read_input_item_relax.cpp | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index f6ba9dceb9..74f050f906 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -14,7 +14,13 @@ void ReadInput::item_relax() item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; read_sync_string(input.relax_method); item.read_value = [](const Input_Item& item, Parameter& para) { - if(item.get_size()==1) + if(item.get_size()==0) + { + para.input.relax_method_param.method="cg"; + para.input.relax_method = para.input.relax_method_param.method; + para.input.relax_method_param.param = "1"; + } + else if(item.get_size()==1) { para.input.relax_method_param.method=item.str_values[0]; para.input.relax_method = para.input.relax_method_param.method; @@ -27,13 +33,6 @@ void ReadInput::item_relax() para.input.relax_method_param.param = item.str_values[1]; } }; - - // std::istringstream iss(item.str_values[0]); - // iss >> para.input.relax_method_param.method; - // if (!(iss >> para.input.relax_method_param.param)) { - // std::cout << "No parameter provided for relax_method_param.param, default to 1" << std::endl; - // para.input.relax_method_param.param = "1"; - // } item.check_value = [](const Input_Item& item, const Parameter& para) { const std::vector relax_methods = {"cg", "sd", "cg_bfgs","bfgs","lbfgs"}; if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { From 1255ad2f153585f5f2d508971bbc96536d142999 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Wed, 24 Sep 2025 20:30:00 +0800 Subject: [PATCH 07/22] change input parameters --- source/source_io/module_parameter/input_parameter.h | 2 +- source/source_io/read_input_item_relax.cpp | 10 ++-------- 2 files changed, 3 insertions(+), 9 deletions(-) diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index 6aebfbe8b6..bef56c0e37 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -154,7 +154,7 @@ struct Input_para // ============== #Parameters (4.Relaxation) =========================== std::string relax_method = "cg"; ///< methods to move_ion: sd, bfgs, cg... - RelaxMethodParam relax_method_param; + RelaxMethodParam relax_method_param={"cg","1"}; bool relax_new = true; bool relax = false; ///< allow relaxation along the specific direction double relax_scale_force = 0.5; diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index 74f050f906..f9ec48bfb8 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -14,13 +14,7 @@ void ReadInput::item_relax() item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; read_sync_string(input.relax_method); item.read_value = [](const Input_Item& item, Parameter& para) { - if(item.get_size()==0) - { - para.input.relax_method_param.method="cg"; - para.input.relax_method = para.input.relax_method_param.method; - para.input.relax_method_param.param = "1"; - } - else if(item.get_size()==1) + if(item.get_size()==1) { para.input.relax_method_param.method=item.str_values[0]; para.input.relax_method = para.input.relax_method_param.method; @@ -33,7 +27,7 @@ void ReadInput::item_relax() para.input.relax_method_param.param = item.str_values[1]; } }; - item.check_value = [](const Input_Item& item, const Parameter& para) { + item.check_value = [](const Input_Item& item, const Parameter& para) { const std::vector relax_methods = {"cg", "sd", "cg_bfgs","bfgs","lbfgs"}; if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { const std::string warningstr = nofound_str(relax_methods, "relax_method"); From eb176ccb7f6116b76a438ddaabfeb021add4c47a Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Wed, 24 Sep 2025 21:23:41 +0800 Subject: [PATCH 08/22] change input parameters --- source/source_io/test_serial/read_input_item_test.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/source/source_io/test_serial/read_input_item_test.cpp b/source/source_io/test_serial/read_input_item_test.cpp index db842bd321..878fbd01a2 100644 --- a/source/source_io/test_serial/read_input_item_test.cpp +++ b/source/source_io/test_serial/read_input_item_test.cpp @@ -767,6 +767,7 @@ TEST_F(InputTest, Item_test) { // relax_method auto it = find_label("relax_method", readinput.input_lists); param.input.relax_method = "none"; + param.input.relax_method_param.method = "none"; testing::internal::CaptureStdout(); EXPECT_EXIT(it->second.check_value(it->second, param), ::testing::ExitedWithCode(1), ""); output = testing::internal::GetCapturedStdout(); From 1b4a330db3c6f3aeb282c0a5729eb7fdeaf72179 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Fri, 26 Sep 2025 15:18:47 +0800 Subject: [PATCH 09/22] fix INPUT problem --- .../source_io/module_parameter/input_parameter.h | 6 +----- source/source_io/read_input_item_relax.cpp | 14 ++++++-------- source/source_relax/ions_move_basic.cpp | 2 +- source/source_relax/ions_move_basic.h | 3 +-- source/source_relax/ions_move_methods.cpp | 9 ++++----- 5 files changed, 13 insertions(+), 21 deletions(-) diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index bef56c0e37..14689b5594 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -6,10 +6,6 @@ #include #include -struct RelaxMethodParam { - std::string method; - std::string param; -}; // It stores all input parameters both defined in INPUT file and not defined in // INPUT file @@ -154,7 +150,7 @@ struct Input_para // ============== #Parameters (4.Relaxation) =========================== std::string relax_method = "cg"; ///< methods to move_ion: sd, bfgs, cg... - RelaxMethodParam relax_method_param={"cg","1"}; + std::string relax_method_param = "1"; bool relax_new = true; bool relax = false; ///< allow relaxation along the specific direction double relax_scale_force = 0.5; diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index f9ec48bfb8..57d34fdcc9 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -16,20 +16,18 @@ void ReadInput::item_relax() item.read_value = [](const Input_Item& item, Parameter& para) { if(item.get_size()==1) { - para.input.relax_method_param.method=item.str_values[0]; - para.input.relax_method = para.input.relax_method_param.method; - para.input.relax_method_param.param = "1"; + para.input.relax_method = item.str_values[0]; + para.input.relax_method_param = "1"; } else if(item.get_size()==2) { - para.input.relax_method_param.method=item.str_values[0]; - para.input.relax_method = para.input.relax_method_param.method; - para.input.relax_method_param.param = item.str_values[1]; + para.input.relax_method = item.str_values[0]; + para.input.relax_method_param = item.str_values[1]; } }; item.check_value = [](const Input_Item& item, const Parameter& para) { - const std::vector relax_methods = {"cg", "sd", "cg_bfgs","bfgs","lbfgs"}; - if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method_param.method) == relax_methods.end()) { + const std::vector relax_methods = {"cg", "sd", "cg_bfgs","lbfgs","bfgs"}; + if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method) == relax_methods.end()) { const std::string warningstr = nofound_str(relax_methods, "relax_method"); ModuleBase::WARNING_QUIT("ReadInput", warningstr); } diff --git a/source/source_relax/ions_move_basic.cpp b/source/source_relax/ions_move_basic.cpp index 0e89cdae43..149bca8648 100644 --- a/source/source_relax/ions_move_basic.cpp +++ b/source/source_relax/ions_move_basic.cpp @@ -24,7 +24,7 @@ double Ions_Move_Basic::best_xxx = 1.0; int Ions_Move_Basic::out_stru = 0; std::string Ions_Move_Basic::relax_method = "bfgs"; -RelaxMethodParam Ions_Move_Basic::relax_method_param; +std::string Ions_Move_Basic::relax_method_param = "1"; void Ions_Move_Basic::setup_gradient(const UnitCell &ucell, const ModuleBase::matrix &force, double *pos, double *grad) { diff --git a/source/source_relax/ions_move_basic.h b/source/source_relax/ions_move_basic.h index 507e6b0b0d..bf3b65ffa5 100644 --- a/source/source_relax/ions_move_basic.h +++ b/source/source_relax/ions_move_basic.h @@ -3,7 +3,6 @@ #include "source_base/matrix.h" #include "source_cell/unitcell.h" -#include "source_io/module_parameter/input_parameter.h" namespace Ions_Move_Basic { @@ -24,7 +23,7 @@ extern double relax_bfgs_init; // initial value of trust radius, extern double best_xxx; // the last step length of cg , we use it as bfgs`s initial step length extern std::string relax_method; // relaxation method, extern int out_stru; // output the structure or not -extern RelaxMethodParam relax_method_param; +extern std::string relax_method_param; // funny way to pass this parameter, but nevertheless //---------------------------------------------------------------------------- diff --git a/source/source_relax/ions_move_methods.cpp b/source/source_relax/ions_move_methods.cpp index 4b442bc637..451475a715 100644 --- a/source/source_relax/ions_move_methods.cpp +++ b/source/source_relax/ions_move_methods.cpp @@ -16,7 +16,7 @@ void Ions_Move_Methods::allocate(const int &natom) { Ions_Move_Basic::dim = natom * 3; - if (Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param != "1") + if (Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param != "1") { this->bfgs.allocate(); } @@ -33,7 +33,7 @@ void Ions_Move_Methods::allocate(const int &natom) this->cg.allocate(); this->bfgs.allocate(); // added by pengfei 13-8-8 } - else if(Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param == "1") + else if(Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param == "1") { this->bfgs_trad.allocate(natom); } @@ -56,10 +56,9 @@ void Ions_Move_Methods::cal_movement(const int &istep, UnitCell &ucell) { ModuleBase::TITLE("Ions_Move_Methods", "init"); - // Ions_Move_Basic::istep = istep; Ions_Move_Basic::istep = force_step; - if (Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param != "1") + if (Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param != "1") { // move_ions // output tau @@ -78,7 +77,7 @@ void Ions_Move_Methods::cal_movement(const int &istep, { cg.start(ucell, f, etot); // added by pengfei 13-8-10 } - else if(Ions_Move_Basic::relax_method_param.method == "bfgs"&&Ions_Move_Basic::relax_method_param.param == "1") + else if(Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param == "1") { bfgs_trad.relax_step(f,ucell); } From 1af25327cc95d8153aafb0197bc487ea918b37e0 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Fri, 26 Sep 2025 16:02:36 +0800 Subject: [PATCH 10/22] fix INPUT problem --- source/source_io/test_serial/read_input_item_test.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/source/source_io/test_serial/read_input_item_test.cpp b/source/source_io/test_serial/read_input_item_test.cpp index 878fbd01a2..db842bd321 100644 --- a/source/source_io/test_serial/read_input_item_test.cpp +++ b/source/source_io/test_serial/read_input_item_test.cpp @@ -767,7 +767,6 @@ TEST_F(InputTest, Item_test) { // relax_method auto it = find_label("relax_method", readinput.input_lists); param.input.relax_method = "none"; - param.input.relax_method_param.method = "none"; testing::internal::CaptureStdout(); EXPECT_EXIT(it->second.check_value(it->second, param), ::testing::ExitedWithCode(1), ""); output = testing::internal::GetCapturedStdout(); From 29b9814c8721da7d67442d2c4b667c4f7bb54a6c Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Fri, 26 Sep 2025 16:31:42 +0800 Subject: [PATCH 11/22] fix INPUT problem --- tests/02_NAO_Gamma/010_NO_GO_RE_MB/INPUT | 2 +- tests/03_NAO_multik/26_NO_KP_RE_MB/INPUT | 2 +- tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT | 2 +- tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT | 2 +- 4 files changed, 4 insertions(+), 4 deletions(-) diff --git a/tests/02_NAO_Gamma/010_NO_GO_RE_MB/INPUT b/tests/02_NAO_Gamma/010_NO_GO_RE_MB/INPUT index 71cb9f7036..f0d05d15e8 100644 --- a/tests/02_NAO_Gamma/010_NO_GO_RE_MB/INPUT +++ b/tests/02_NAO_Gamma/010_NO_GO_RE_MB/INPUT @@ -17,7 +17,7 @@ scf_nmax 100 relax_nmax 2 cal_force 1 force_thr_ev 0.01 -relax_method bfgs +relax_method bfgs 2 #Parameters (4.Basis) basis_type lcao diff --git a/tests/03_NAO_multik/26_NO_KP_RE_MB/INPUT b/tests/03_NAO_multik/26_NO_KP_RE_MB/INPUT index 8478291249..fd216bc63c 100644 --- a/tests/03_NAO_multik/26_NO_KP_RE_MB/INPUT +++ b/tests/03_NAO_multik/26_NO_KP_RE_MB/INPUT @@ -17,7 +17,7 @@ scf_nmax 100 relax_nmax 2 cal_force 1 force_thr_ev 0.01 -relax_method bfgs +relax_method bfgs 2 #Parameters (4.Basis) basis_type lcao diff --git a/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT b/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT index c111a61dc5..d1cce921c3 100644 --- a/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT +++ b/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT @@ -2,7 +2,7 @@ INPUT_PARAMETERS #Parameters (1.General) suffix autotest calculation relax -relax_method bfgs +relax_method bfgs 2 relax_nmax 10 dft_functional pbe diff --git a/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT b/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT index 329564682c..64337f5fdd 100644 --- a/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT +++ b/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT @@ -2,7 +2,7 @@ INPUT_PARAMETERS #Parameters (1.General) suffix autotest calculation relax -relax_method bfgs +relax_method bfgs 2 relax_nmax 10 dft_functional pbe From ce6716e69ec8efd53291095a84b43a776349d7e5 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Fri, 26 Sep 2025 18:26:02 +0800 Subject: [PATCH 12/22] fix INPUT problem --- tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT | 2 +- tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT b/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT index d1cce921c3..3feba5a5e0 100644 --- a/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT +++ b/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT @@ -2,7 +2,7 @@ INPUT_PARAMETERS #Parameters (1.General) suffix autotest calculation relax -relax_method bfgs 2 +relax_method cg relax_nmax 10 dft_functional pbe diff --git a/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT b/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT index 64337f5fdd..385ba742bd 100644 --- a/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT +++ b/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT @@ -2,7 +2,7 @@ INPUT_PARAMETERS #Parameters (1.General) suffix autotest calculation relax -relax_method bfgs 2 +relax_method cg relax_nmax 10 dft_functional pbe From d9a5e92f5448db4490f919ef6970c999737bc309 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Mon, 29 Sep 2025 22:39:10 +0800 Subject: [PATCH 13/22] fix INPUT problem --- source/source_io/input_conv.cpp | 1 - .../module_parameter/input_parameter.h | 3 +-- source/source_io/read_input_item_relax.cpp | 15 ++++++------ source/source_io/test/read_input_ptest.cpp | 2 +- .../test_serial/read_input_item_test.cpp | 6 ++--- source/source_relax/ions_move_basic.cpp | 3 +-- source/source_relax/ions_move_basic.h | 3 +-- source/source_relax/ions_move_cg.cpp | 4 ++-- source/source_relax/ions_move_methods.cpp | 24 +++++++++---------- .../source_relax/test/ions_move_cg_test.cpp | 14 +++++------ .../test/ions_move_methods_test.cpp | 20 ++++++++-------- tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT | 2 +- 12 files changed, 46 insertions(+), 51 deletions(-) diff --git a/source/source_io/input_conv.cpp b/source/source_io/input_conv.cpp index fec9695537..952ca873d1 100644 --- a/source/source_io/input_conv.cpp +++ b/source/source_io/input_conv.cpp @@ -190,7 +190,6 @@ void Input_Conv::Convert() Ions_Move_Basic::relax_bfgs_init = PARAM.inp.relax_bfgs_init; Ions_Move_Basic::out_stru = PARAM.inp.out_stru; // mohan add 2012-03-23 Ions_Move_Basic::relax_method = PARAM.inp.relax_method; - Ions_Move_Basic::relax_method_param = PARAM.inp.relax_method_param; Lattice_Change_Basic::fixed_axes = PARAM.inp.fixed_axes; diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index 14689b5594..85d13bf40e 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -149,8 +149,7 @@ struct Input_para // int bessel_nao_lmax; ///< lmax used in descriptor // ============== #Parameters (4.Relaxation) =========================== - std::string relax_method = "cg"; ///< methods to move_ion: sd, bfgs, cg... - std::string relax_method_param = "1"; + std::vector relax_method = {"cg","1"}; ///< methods to move_ion: sd, bfgs, cg... bool relax_new = true; bool relax = false; ///< allow relaxation along the specific direction double relax_scale_force = 0.5; diff --git a/source/source_io/read_input_item_relax.cpp b/source/source_io/read_input_item_relax.cpp index 57d34fdcc9..59c6ab187a 100644 --- a/source/source_io/read_input_item_relax.cpp +++ b/source/source_io/read_input_item_relax.cpp @@ -12,22 +12,21 @@ void ReadInput::item_relax() { Input_Item item("relax_method"); item.annotation = "cg; bfgs; sd; cg; cg_bfgs;"; - read_sync_string(input.relax_method); item.read_value = [](const Input_Item& item, Parameter& para) { if(item.get_size()==1) { - para.input.relax_method = item.str_values[0]; - para.input.relax_method_param = "1"; + para.input.relax_method[0] = item.str_values[0]; + para.input.relax_method[1] = "1"; } - else if(item.get_size()==2) + else if(item.get_size()>=2) { - para.input.relax_method = item.str_values[0]; - para.input.relax_method_param = item.str_values[1]; + para.input.relax_method[0] = item.str_values[0]; + para.input.relax_method[1] = item.str_values[1]; } }; item.check_value = [](const Input_Item& item, const Parameter& para) { const std::vector relax_methods = {"cg", "sd", "cg_bfgs","lbfgs","bfgs"}; - if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method) == relax_methods.end()) { + if (std::find(relax_methods.begin(), relax_methods.end(), para.input.relax_method[0]) == relax_methods.end()) { const std::string warningstr = nofound_str(relax_methods, "relax_method"); ModuleBase::WARNING_QUIT("ReadInput", warningstr); } @@ -52,7 +51,7 @@ void ReadInput::item_relax() item.annotation = "whether to use the new relaxation method"; read_sync_bool(input.relax_new); item.reset_value = [](const Input_Item& item, Parameter& para) { - if (para.input.relax_new && para.input.relax_method != "cg") + if (para.input.relax_new && para.input.relax_method[0] != "cg") { para.input.relax_new = false; } diff --git a/source/source_io/test/read_input_ptest.cpp b/source/source_io/test/read_input_ptest.cpp index ebb276265d..56606b770b 100644 --- a/source/source_io/test/read_input_ptest.cpp +++ b/source/source_io/test/read_input_ptest.cpp @@ -114,7 +114,7 @@ TEST_F(InputParaTest, ParaRead) EXPECT_EQ(param.inp.fixed_axes, "None"); EXPECT_FALSE(param.inp.fixed_ibrav); EXPECT_FALSE(param.inp.fixed_atoms); - EXPECT_EQ(param.inp.relax_method, "cg"); + EXPECT_EQ(param.inp.relax_method[0], "cg"); EXPECT_DOUBLE_EQ(param.inp.relax_cg_thr, 0.5); EXPECT_EQ(param.inp.out_level, "ie"); EXPECT_TRUE(param.globalv.out_md_control); diff --git a/source/source_io/test_serial/read_input_item_test.cpp b/source/source_io/test_serial/read_input_item_test.cpp index db842bd321..ea026e8ba6 100644 --- a/source/source_io/test_serial/read_input_item_test.cpp +++ b/source/source_io/test_serial/read_input_item_test.cpp @@ -766,7 +766,7 @@ TEST_F(InputTest, Item_test) } { // relax_method auto it = find_label("relax_method", readinput.input_lists); - param.input.relax_method = "none"; + param.input.relax_method[0] = "none"; testing::internal::CaptureStdout(); EXPECT_EXIT(it->second.check_value(it->second, param), ::testing::ExitedWithCode(1), ""); output = testing::internal::GetCapturedStdout(); @@ -775,12 +775,12 @@ TEST_F(InputTest, Item_test) { //relax_new auto it = find_label("relax_new", readinput.input_lists); param.input.relax_new = true; - param.input.relax_method = "cg"; + param.input.relax_method[0] = "cg"; it->second.reset_value(it->second, param); EXPECT_EQ(param.input.relax_new, true); param.input.relax_new = true; - param.input.relax_method = "none"; + param.input.relax_method[0] = "none"; it->second.reset_value(it->second, param); EXPECT_EQ(param.input.relax_new, false); } diff --git a/source/source_relax/ions_move_basic.cpp b/source/source_relax/ions_move_basic.cpp index 149bca8648..620eb1c8a1 100644 --- a/source/source_relax/ions_move_basic.cpp +++ b/source/source_relax/ions_move_basic.cpp @@ -23,8 +23,7 @@ double Ions_Move_Basic::relax_bfgs_init = -1.0; // default is 0.5 double Ions_Move_Basic::best_xxx = 1.0; int Ions_Move_Basic::out_stru = 0; -std::string Ions_Move_Basic::relax_method = "bfgs"; -std::string Ions_Move_Basic::relax_method_param = "1"; +std::vector Ions_Move_Basic::relax_method = {"bfgs","1"}; void Ions_Move_Basic::setup_gradient(const UnitCell &ucell, const ModuleBase::matrix &force, double *pos, double *grad) { diff --git a/source/source_relax/ions_move_basic.h b/source/source_relax/ions_move_basic.h index bf3b65ffa5..8f03085028 100644 --- a/source/source_relax/ions_move_basic.h +++ b/source/source_relax/ions_move_basic.h @@ -21,9 +21,8 @@ extern double relax_bfgs_rmax; // max value of trust radius, extern double relax_bfgs_rmin; // min value of trust radius, extern double relax_bfgs_init; // initial value of trust radius, extern double best_xxx; // the last step length of cg , we use it as bfgs`s initial step length -extern std::string relax_method; // relaxation method, +extern std::vector relax_method; // relaxation method, extern int out_stru; // output the structure or not -extern std::string relax_method_param; // funny way to pass this parameter, but nevertheless //---------------------------------------------------------------------------- diff --git a/source/source_relax/ions_move_cg.cpp b/source/source_relax/ions_move_cg.cpp index 7b1dbd1420..267cbc5cf6 100644 --- a/source/source_relax/ions_move_cg.cpp +++ b/source/source_relax/ions_move_cg.cpp @@ -163,12 +163,12 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const fmax = fa; sd = false; - if (Ions_Move_Basic::relax_method == "cg_bfgs") + if (Ions_Move_Basic::relax_method[0] == "cg_bfgs") { if (Ions_Move_Basic::largest_grad * ModuleBase::Ry_to_eV / 0.529177 < RELAX_CG_THR) // cg to bfgs by pengfei 13-8-8 { - Ions_Move_Basic::relax_method = "bfgs"; + Ions_Move_Basic::relax_method[0] = "bfgs"; } Ions_Move_Basic::best_xxx = steplength; } diff --git a/source/source_relax/ions_move_methods.cpp b/source/source_relax/ions_move_methods.cpp index 451475a715..15f27d67d3 100644 --- a/source/source_relax/ions_move_methods.cpp +++ b/source/source_relax/ions_move_methods.cpp @@ -16,28 +16,28 @@ void Ions_Move_Methods::allocate(const int &natom) { Ions_Move_Basic::dim = natom * 3; - if (Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param != "1") + if (Ions_Move_Basic::relax_method[0] == "bfgs"&&Ions_Move_Basic::relax_method[1] != "1") { this->bfgs.allocate(); } - else if (Ions_Move_Basic::relax_method == "sd") + else if (Ions_Move_Basic::relax_method[0] == "sd") { this->sd.allocate(); } - else if (Ions_Move_Basic::relax_method == "cg") + else if (Ions_Move_Basic::relax_method[0] == "cg") { this->cg.allocate(); } - else if (Ions_Move_Basic::relax_method == "cg_bfgs") + else if (Ions_Move_Basic::relax_method[0] == "cg_bfgs") { this->cg.allocate(); this->bfgs.allocate(); // added by pengfei 13-8-8 } - else if(Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param == "1") + else if(Ions_Move_Basic::relax_method[0] == "bfgs"&&Ions_Move_Basic::relax_method[1] == "1") { this->bfgs_trad.allocate(natom); } - else if(Ions_Move_Basic::relax_method == "lbfgs") + else if(Ions_Move_Basic::relax_method[0] == "lbfgs") { this->lbfgs.allocate(natom); } @@ -58,30 +58,30 @@ void Ions_Move_Methods::cal_movement(const int &istep, ModuleBase::TITLE("Ions_Move_Methods", "init"); // Ions_Move_Basic::istep = istep; Ions_Move_Basic::istep = force_step; - if (Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param != "1") + if (Ions_Move_Basic::relax_method[0] == "bfgs"&&Ions_Move_Basic::relax_method[1] != "1") { // move_ions // output tau // check all symmery bfgs.start(ucell, f, etot); } - else if (Ions_Move_Basic::relax_method == "sd") + else if (Ions_Move_Basic::relax_method[0] == "sd") { sd.start(ucell, f, etot); } - else if (Ions_Move_Basic::relax_method == "cg") + else if (Ions_Move_Basic::relax_method[0] == "cg") { cg.start(ucell, f, etot); } - else if (Ions_Move_Basic::relax_method == "cg_bfgs") + else if (Ions_Move_Basic::relax_method[0] == "cg_bfgs") { cg.start(ucell, f, etot); // added by pengfei 13-8-10 } - else if(Ions_Move_Basic::relax_method == "bfgs"&&Ions_Move_Basic::relax_method_param == "1") + else if(Ions_Move_Basic::relax_method[0] == "bfgs"&&Ions_Move_Basic::relax_method[1] == "1") { bfgs_trad.relax_step(f,ucell); } - else if(Ions_Move_Basic::relax_method == "lbfgs") + else if(Ions_Move_Basic::relax_method[0] == "lbfgs") { lbfgs.relax_step(f,ucell,etot); } diff --git a/source/source_relax/test/ions_move_cg_test.cpp b/source/source_relax/test/ions_move_cg_test.cpp index f3ce632f0f..c0aa34dcec 100644 --- a/source/source_relax/test/ions_move_cg_test.cpp +++ b/source/source_relax/test/ions_move_cg_test.cpp @@ -138,7 +138,7 @@ TEST_F(IonsMoveCGTest, TestStartSd) // setup data Ions_Move_Basic::istep = 1; Ions_Move_Basic::converged = false; - Ions_Move_Basic::relax_method = "cg_bfgs"; + Ions_Move_Basic::relax_method[0] = "cg_bfgs"; Ions_Move_CG::RELAX_CG_THR = 100.0; UnitCell ucell; setupucell(ucell); @@ -162,7 +162,7 @@ TEST_F(IonsMoveCGTest, TestStartSd) EXPECT_THAT(output, testing::HasSubstr(expected_output)); EXPECT_EQ(Ions_Move_Basic::converged, false); EXPECT_EQ(Ions_Move_Basic::update_iter, 5); - EXPECT_EQ(Ions_Move_Basic::relax_method, "bfgs"); + EXPECT_EQ(Ions_Move_Basic::relax_method[0], "bfgs"); EXPECT_DOUBLE_EQ(Ions_Move_Basic::largest_grad, 0.01); EXPECT_DOUBLE_EQ(Ions_Move_Basic::best_xxx, -1.0); EXPECT_DOUBLE_EQ(Ions_Move_Basic::relax_bfgs_init, 1.0); @@ -201,7 +201,7 @@ TEST_F(IonsMoveCGTest, TestStartTrialGoto) EXPECT_THAT(output, testing::HasSubstr(expected_output)); EXPECT_EQ(Ions_Move_Basic::converged, false); EXPECT_EQ(Ions_Move_Basic::update_iter, 5); - EXPECT_EQ(Ions_Move_Basic::relax_method, "bfgs"); + EXPECT_EQ(Ions_Move_Basic::relax_method[0], "bfgs"); EXPECT_DOUBLE_EQ(Ions_Move_Basic::largest_grad, 0.001); EXPECT_DOUBLE_EQ(Ions_Move_Basic::best_xxx, -1.0); EXPECT_DOUBLE_EQ(Ions_Move_Basic::relax_bfgs_init, 10.0); @@ -239,7 +239,7 @@ TEST_F(IonsMoveCGTest, TestStartTrial) EXPECT_THAT(output, testing::HasSubstr(expected_output)); EXPECT_EQ(Ions_Move_Basic::converged, false); EXPECT_EQ(Ions_Move_Basic::update_iter, 5); - EXPECT_EQ(Ions_Move_Basic::relax_method, "bfgs"); + EXPECT_EQ(Ions_Move_Basic::relax_method[0], "bfgs"); EXPECT_DOUBLE_EQ(Ions_Move_Basic::largest_grad, 0.01); EXPECT_DOUBLE_EQ(Ions_Move_Basic::best_xxx, -1.0); EXPECT_DOUBLE_EQ(Ions_Move_Basic::relax_bfgs_init, 70.0); @@ -279,7 +279,7 @@ TEST_F(IonsMoveCGTest, TestStartNoTrialGotoCase1) EXPECT_THAT(output, testing::HasSubstr(expected_output)); EXPECT_EQ(Ions_Move_Basic::converged, false); EXPECT_EQ(Ions_Move_Basic::update_iter, 5); - EXPECT_EQ(Ions_Move_Basic::relax_method, "bfgs"); + EXPECT_EQ(Ions_Move_Basic::relax_method[0], "bfgs"); EXPECT_DOUBLE_EQ(Ions_Move_Basic::largest_grad, 0.001); EXPECT_DOUBLE_EQ(Ions_Move_Basic::best_xxx, -1.0); EXPECT_DOUBLE_EQ(Ions_Move_Basic::relax_bfgs_init, 490.0); @@ -318,7 +318,7 @@ TEST_F(IonsMoveCGTest, TestStartNoTrialGotoCase2) EXPECT_THAT(output, testing::HasSubstr(expected_output)); EXPECT_EQ(Ions_Move_Basic::converged, false); EXPECT_EQ(Ions_Move_Basic::update_iter, 5); - EXPECT_EQ(Ions_Move_Basic::relax_method, "bfgs"); + EXPECT_EQ(Ions_Move_Basic::relax_method[0], "bfgs"); EXPECT_DOUBLE_EQ(Ions_Move_Basic::largest_grad, 0.01); EXPECT_DOUBLE_EQ(Ions_Move_Basic::best_xxx, -1.0); EXPECT_DOUBLE_EQ(Ions_Move_Basic::relax_bfgs_init, 70.0); @@ -358,7 +358,7 @@ TEST_F(IonsMoveCGTest, TestStartNoTrial) EXPECT_THAT(output, testing::HasSubstr(expected_output)); EXPECT_EQ(Ions_Move_Basic::converged, false); EXPECT_EQ(Ions_Move_Basic::update_iter, 5); - EXPECT_EQ(Ions_Move_Basic::relax_method, "bfgs"); + EXPECT_EQ(Ions_Move_Basic::relax_method[0], "bfgs"); EXPECT_DOUBLE_EQ(Ions_Move_Basic::largest_grad, 0.001); EXPECT_DOUBLE_EQ(Ions_Move_Basic::best_xxx, -1.0); EXPECT_NEAR(Ions_Move_Basic::relax_bfgs_init, 1.2345679012345678, 1e-12); diff --git a/source/source_relax/test/ions_move_methods_test.cpp b/source/source_relax/test/ions_move_methods_test.cpp index 27e7315d68..f0027be236 100644 --- a/source/source_relax/test/ions_move_methods_test.cpp +++ b/source/source_relax/test/ions_move_methods_test.cpp @@ -40,19 +40,19 @@ class IonsMoveMethodsTest : public ::testing::Test // Test the allocate() function TEST_F(IonsMoveMethodsTest, Allocate) { - Ions_Move_Basic::relax_method = "bfgs"; + Ions_Move_Basic::relax_method[0] = "bfgs"; imm.allocate(natom); EXPECT_EQ(Ions_Move_Basic::dim, 6); - Ions_Move_Basic::relax_method = "sd"; + Ions_Move_Basic::relax_method[0] = "sd"; imm.allocate(natom); EXPECT_EQ(Ions_Move_Basic::dim, 6); - Ions_Move_Basic::relax_method = "cg"; + Ions_Move_Basic::relax_method[0] = "cg"; imm.allocate(natom); EXPECT_EQ(Ions_Move_Basic::dim, 6); - Ions_Move_Basic::relax_method = "cg_bfgs"; + Ions_Move_Basic::relax_method[0] = "cg_bfgs"; imm.allocate(natom); EXPECT_EQ(Ions_Move_Basic::dim, 6); } @@ -60,7 +60,7 @@ TEST_F(IonsMoveMethodsTest, Allocate) // Test the allocate() function warning quit TEST_F(IonsMoveMethodsTest, AllocateWarningQuit) { - Ions_Move_Basic::relax_method = "none"; + Ions_Move_Basic::relax_method[0] = "none"; GlobalV::ofs_warning.open("log"); imm.allocate(natom); GlobalV::ofs_warning.close(); @@ -81,22 +81,22 @@ TEST_F(IonsMoveMethodsTest, CalMovement) const double etot = 0.0; UnitCell ucell; - Ions_Move_Basic::relax_method = "bfgs"; + Ions_Move_Basic::relax_method[0] = "bfgs"; imm.allocate(natom); imm.cal_movement(istep, force_step, f, etot, ucell); EXPECT_EQ(Ions_Move_Basic::istep, force_step); - Ions_Move_Basic::relax_method = "sd"; + Ions_Move_Basic::relax_method[0] = "sd"; imm.allocate(natom); imm.cal_movement(istep, force_step, f, etot, ucell); EXPECT_EQ(Ions_Move_Basic::istep, force_step); - Ions_Move_Basic::relax_method = "cg"; + Ions_Move_Basic::relax_method[0] = "cg"; imm.allocate(natom); imm.cal_movement(istep, force_step, f, etot, ucell); EXPECT_EQ(Ions_Move_Basic::istep, force_step); - Ions_Move_Basic::relax_method = "cg_bfgs"; + Ions_Move_Basic::relax_method[0] = "cg_bfgs"; imm.allocate(natom); imm.cal_movement(istep, force_step, f, etot, ucell); EXPECT_EQ(Ions_Move_Basic::istep, force_step); @@ -110,7 +110,7 @@ TEST_F(IonsMoveMethodsTest, CalMovementWarningQuit) const ModuleBase::matrix f(3, 3); const double etot = 0.0; UnitCell ucell; - Ions_Move_Basic::relax_method = "none"; + Ions_Move_Basic::relax_method[0] = "none"; imm.allocate(natom); GlobalV::ofs_warning.open("log"); diff --git a/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT b/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT index 3feba5a5e0..d1cce921c3 100644 --- a/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT +++ b/tests/09_DeePKS/102_NO_GO_deepks_relax/INPUT @@ -2,7 +2,7 @@ INPUT_PARAMETERS #Parameters (1.General) suffix autotest calculation relax -relax_method cg +relax_method bfgs 2 relax_nmax 10 dft_functional pbe From f2aa7f9fcad2899a8a29dbf121f6beaace24d173 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 30 Sep 2025 00:04:42 +0800 Subject: [PATCH 14/22] fix INPUT problem --- source/source_relax/ions_move_cg.cpp | 1 + tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT | 2 +- 2 files changed, 2 insertions(+), 1 deletion(-) diff --git a/source/source_relax/ions_move_cg.cpp b/source/source_relax/ions_move_cg.cpp index 267cbc5cf6..8dd0f23bd6 100644 --- a/source/source_relax/ions_move_cg.cpp +++ b/source/source_relax/ions_move_cg.cpp @@ -169,6 +169,7 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const < RELAX_CG_THR) // cg to bfgs by pengfei 13-8-8 { Ions_Move_Basic::relax_method[0] = "bfgs"; + Ions_Move_Basic::relax_method[1] = "2"; } Ions_Move_Basic::best_xxx = steplength; } diff --git a/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT b/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT index 385ba742bd..64337f5fdd 100644 --- a/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT +++ b/tests/09_DeePKS/102_NO_KP_deepks_relax/INPUT @@ -2,7 +2,7 @@ INPUT_PARAMETERS #Parameters (1.General) suffix autotest calculation relax -relax_method cg +relax_method bfgs 2 relax_nmax 10 dft_functional pbe From 6312b1dfffd1854915ff44837f44b2c1547a4ee3 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 30 Sep 2025 01:36:06 +0800 Subject: [PATCH 15/22] fix INPUT problem --- source/source_relax/ions_move_basic.cpp | 2 +- source/source_relax/ions_move_cg.cpp | 4 ++-- source/source_relax/lattice_change_cg.cpp | 4 ++-- 3 files changed, 5 insertions(+), 5 deletions(-) diff --git a/source/source_relax/ions_move_basic.cpp b/source/source_relax/ions_move_basic.cpp index 620eb1c8a1..1bf1e025a6 100644 --- a/source/source_relax/ions_move_basic.cpp +++ b/source/source_relax/ions_move_basic.cpp @@ -23,7 +23,7 @@ double Ions_Move_Basic::relax_bfgs_init = -1.0; // default is 0.5 double Ions_Move_Basic::best_xxx = 1.0; int Ions_Move_Basic::out_stru = 0; -std::vector Ions_Move_Basic::relax_method = {"bfgs","1"}; +std::vector Ions_Move_Basic::relax_method = {"bfgs","2"}; void Ions_Move_Basic::setup_gradient(const UnitCell &ucell, const ModuleBase::matrix &force, double *pos, double *grad) { diff --git a/source/source_relax/ions_move_cg.cpp b/source/source_relax/ions_move_cg.cpp index 8dd0f23bd6..9b8aa4b5c0 100644 --- a/source/source_relax/ions_move_cg.cpp +++ b/source/source_relax/ions_move_cg.cpp @@ -192,7 +192,7 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const goto CG_begin; } - normalize(cg_gradn, cg_grad, dim); + normalize(cg_gradn, cg_grad0, dim); third_order(e0, e1, fa, fb, xb, best_x); // cubic interpolation if (best_x > 6 * xb || best_x < (-xb)) @@ -237,7 +237,7 @@ void Ions_Move_CG::start(UnitCell &ucell, const ModuleBase::matrix &force, const goto CG_begin; } - normalize(cg_gradn, cg_grad, dim); + normalize(cg_gradn, cg_grad0, dim); setup_move(move, cg_gradn, best_x); Ions_Move_Basic::move_atoms(ucell, move, pos); Ions_Move_Basic::relax_bfgs_init = xc; diff --git a/source/source_relax/lattice_change_cg.cpp b/source/source_relax/lattice_change_cg.cpp index fec3f0be9d..442b37e671 100644 --- a/source/source_relax/lattice_change_cg.cpp +++ b/source/source_relax/lattice_change_cg.cpp @@ -200,7 +200,7 @@ void Lattice_Change_CG::start(UnitCell &ucell, const ModuleBase::matrix &stress_ goto CG_begin; } - normalize(cg_gradn, cg_grad, dim); + normalize(cg_gradn, cg_grad0, dim); third_order(e0, e1, fa, fb, xb, best_x); // cubic interpolation if (best_x > 6 * xb || best_x < (-xb)) @@ -251,7 +251,7 @@ void Lattice_Change_CG::start(UnitCell &ucell, const ModuleBase::matrix &stress_ goto CG_begin; } - normalize(cg_gradn, cg_grad, dim); + normalize(cg_gradn, cg_grad0, dim); setup_move(move, cg_gradn, best_x); Lattice_Change_Basic::change_lattice(ucell, move, lat); From 5c477b070a498d053cb6fbcc1bc2e3fb217eca87 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Tue, 30 Sep 2025 19:55:18 +0800 Subject: [PATCH 16/22] fix INPUT problem --- docs/advanced/input_files/input-main.md | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index c407df99d1..29b536167d 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -1463,15 +1463,18 @@ These variables are used to control the geometry relaxation. ### relax_method -- **Type**: String +- **Type**: Vector of string - **Description**: The methods to do geometry optimization. + the first element: - cg: using the conjugate gradient (CG) algorithm. Note that there are two implementations of the conjugate gradient (CG) method, see [relax_new](#relax_new). - - bfgs: using the Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm. - - bfgs_trad: using the traditional Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm. + - bfgs : using the Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm. - cg_bfgs: using the CG method for the initial steps, and switching to BFGS method when the force convergence is smaller than [relax_cg_thr](#relax_cg_thr). - sd: using the steepest descent (SD) algorithm. - fire: the Fast Inertial Relaxation Engine method (FIRE), a kind of molecular-dynamics-based relaxation algorithm, is implemented in the molecular dynamics (MD) module. The algorithm can be used by setting [calculation](#calculation) to `md` and [md_type](#md_type) to `fire`. Also ionic velocities should be set in this case. See [fire](../md.md#fire) for more details. -- **Default**: cg + + the second element: + when the first element is bfgs, if the second parameter is 1, it indicates the use of the new BFGS algorithm; if the second parameter is not 1, it indicates the use of the old BFGS algorithm. +- **Default**: cg 1 ### relax_new From f12f6eb6e8462e333f88f4e8208abe4a91a41b88 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Sat, 4 Oct 2025 10:48:37 +0800 Subject: [PATCH 17/22] fix INPUT problem --- docs/advanced/input_files/input-main.md | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index 29b536167d..c06e5deec7 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -1468,13 +1468,15 @@ These variables are used to control the geometry relaxation. the first element: - cg: using the conjugate gradient (CG) algorithm. Note that there are two implementations of the conjugate gradient (CG) method, see [relax_new](#relax_new). - bfgs : using the Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm. + - lbfgs: using the Limited-memory Broyden–Fletcher–Goldfarb–Shanno (LBFGS) algorithm. - cg_bfgs: using the CG method for the initial steps, and switching to BFGS method when the force convergence is smaller than [relax_cg_thr](#relax_cg_thr). - sd: using the steepest descent (SD) algorithm. - fire: the Fast Inertial Relaxation Engine method (FIRE), a kind of molecular-dynamics-based relaxation algorithm, is implemented in the molecular dynamics (MD) module. The algorithm can be used by setting [calculation](#calculation) to `md` and [md_type](#md_type) to `fire`. Also ionic velocities should be set in this case. See [fire](../md.md#fire) for more details. the second element: - when the first element is bfgs, if the second parameter is 1, it indicates the use of the new BFGS algorithm; if the second parameter is not 1, it indicates the use of the old BFGS algorithm. + when the first element is bfgs, if the second parameter is 1, it indicates the use of the new BFGS algorithm; if the second parameter is not 1, it indicates the use of the old BFGS algorithm. - **Default**: cg 1 +- **Note**:In the 3.10-LTS version, the type of this parameter is std::string. It can be set to "cg","bfgs","cg_bfgs","bfgs_trad","lbfgs","sd","fire". ### relax_new From d408f3292b3fda677327672bf6bebb7979568e88 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Wed, 8 Oct 2025 13:50:30 +0800 Subject: [PATCH 18/22] code optimazation --- source/source_relax/bfgs.cpp | 38 +++++++++++++++++++------- source/source_relax/bfgs.h | 22 ++++++++------- source/source_relax/matrix_methods.cpp | 18 ++++++++++++ source/source_relax/matrix_methods.h | 1 + 4 files changed, 59 insertions(+), 20 deletions(-) diff --git a/source/source_relax/bfgs.cpp b/source/source_relax/bfgs.cpp index 829af3fa1c..175ea3b4ad 100644 --- a/source/source_relax/bfgs.cpp +++ b/source/source_relax/bfgs.cpp @@ -10,10 +10,12 @@ //! initialize H0、H、pos0、force0、force void BFGS::allocate(const int _size) { + assert(_size > 0); alpha=70;//default value in ase is 70 maxstep=PARAM.inp.relax_bfgs_rmax; size=_size; - sign =true; + largest_grad=0.0; + sign=true; H = std::vector>(3*size, std::vector(3*size, 0.0)); for (int i = 0; i < 3*size; ++i) @@ -37,6 +39,7 @@ void BFGS::relax_step(const ModuleBase::matrix& _force,UnitCell& ucell) GetPos(ucell,pos); GetPostaud(ucell,pos_taud); ucell.ionic_position_updated = true; + assert(_force.nr == force.size() && _force.nc == force[0].size()); for(int i = 0; i < _force.nr; i++) { for(int j=0;j<_force.nc;j++) @@ -65,7 +68,7 @@ void BFGS::relax_step(const ModuleBase::matrix& _force,UnitCell& ucell) k+=ucell.atoms[i].na; } - this->PrepareStep(force,pos,H,pos0,force0,steplength,dpos,ucell); + this->PrepareStep(force,pos,H,pos0,force0,steplength,dpos,size,ucell); this->DetermineStep(steplength,dpos,maxstep); this->UpdatePos(ucell); this->CalculateLargestGrad(_force,ucell); @@ -76,6 +79,12 @@ void BFGS::relax_step(const ModuleBase::matrix& _force,UnitCell& ucell) void BFGS::GetPos(UnitCell& ucell,std::vector>& pos) { + int total_atoms = 0; + for(int i = 0; i < ucell.ntype; i++) + { + total_atoms += ucell.atoms[i].na; + } + assert(pos.size() == total_atoms); int k=0; for(int i=0;i>& pos) void BFGS::GetPostaud(UnitCell& ucell, std::vector>& pos_taud) { + int total_atoms = 0; + for(int i = 0; i < ucell.ntype; i++) + { + total_atoms += ucell.atoms[i].na; + } + assert(pos_taud.size() == total_atoms); int k=0; for(int i=0;i>& force, std::vector& force0, std::vector& steplength, std::vector>& dpos, + int& size, UnitCell& ucell) { std::vector changedforce = ReshapeMToV(force); @@ -144,6 +160,7 @@ void BFGS::PrepareStep(std::vector>& force, std::vector a=DotInMAndV2(V, changedforce); for(int i = 0; i < a.size(); i++) { + assert(std::abs(omega[i]) > 1e-8); a[i]/=std::abs(omega[i]); } std::vector tmpdpos = DotInMAndV1(V, a); @@ -243,6 +260,7 @@ void BFGS::DetermineStep(std::vector& steplength, { std::vector::iterator maxsteplength = max_element(steplength.begin(), steplength.end()); double a = *maxsteplength; + assert(a > 1e-10); if(a >= maxstep) { double scale = maxstep / a; @@ -311,8 +329,8 @@ void BFGS::UpdatePos(UnitCell& ucell) void BFGS::IsRestrain(std::vector>& dpos) { - Ions_Move_Basic::converged = Ions_Move_Basic::largest_grad - * ModuleBase::Ry_to_eV / 0.529177 steplength;//the length of atoms displacement std::vector> H;//Hessian matrix std::vector force0;//force in previous step @@ -22,17 +33,8 @@ class BFGS std::vector pos_taud0;//atom pos in previous step(relative coordinates) std::vector> pos_taud; std::vector> dpos; - - void allocate(const int _size);//initialize parameters - void relax_step(const ModuleBase::matrix& _force,UnitCell& ucell);//a full iteration step - void PrepareStep(std::vector>& force,std::vector>& pos,std::vector>& H,std::vector& pos0,std::vector& force0,std::vector& steplength,std::vector>& dpos,UnitCell& ucell);//calculate the atomic displacement in one iteration step - -private: - bool sign;//check if this is the first iteration - double alpha;//initialize H,diagonal element is alpha - double maxstep;//every movement smaller than maxstep - int size;//number of atoms + void PrepareStep(std::vector>& force,std::vector>& pos,std::vector>& H,std::vector& pos0,std::vector& force0,std::vector& steplength,std::vector>& dpos,int& size,UnitCell& ucell);//calculate the atomic displacement in one iteration step void IsRestrain(std::vector>& dpos);//check if converged void CalculateLargestGrad(const ModuleBase::matrix& _force,UnitCell& ucell); void GetPos(UnitCell& ucell,std::vector>& pos); diff --git a/source/source_relax/matrix_methods.cpp b/source/source_relax/matrix_methods.cpp index 99019f0668..27f106ae56 100644 --- a/source/source_relax/matrix_methods.cpp +++ b/source/source_relax/matrix_methods.cpp @@ -4,6 +4,8 @@ std::vector ReshapeMToV(std::vector>& matrix) { + assert(!matrix.empty()); + assert(matrix[0].size() == 3); int size = matrix.size(); std::vector result; result.reserve(3*size); @@ -16,6 +18,8 @@ std::vector ReshapeMToV(std::vector>& matrix) std::vector> MAddM(std::vector>& a, std::vector>& b) { + assert(!a.empty() && !b.empty()); + assert(a.size() == b.size() && a[0].size() == b[0].size()); std::vector> result = std::vector>(a.size(), std::vector(a[0].size(), 0.0)); for(int i = 0; i < a.size(); i++) { @@ -29,6 +33,7 @@ std::vector> MAddM(std::vector>& a, std::vector VSubV(std::vector& a, std::vector& b) { + assert(a.size() == b.size()); std::vector result = std::vector(a.size(), 0.0); for(int i = 0; i < a.size(); i++) { @@ -39,6 +44,7 @@ std::vector VSubV(std::vector& a, std::vector& b) std::vector> ReshapeVToM(std::vector& matrix) { + assert(matrix.size() % 3 == 0); std::vector> result = std::vector>(matrix.size() / 3, std::vector(3)); for(int i = 0; i < result.size(); i++) { @@ -52,6 +58,8 @@ std::vector> ReshapeVToM(std::vector& matrix) std::vector DotInMAndV1(std::vector>& matrix, std::vector& vec) { + assert(!matrix.empty()); + assert(matrix[0].size() == vec.size()); std::vector result(matrix.size(), 0.0); for(int i = 0; i < result.size(); i++) { @@ -64,6 +72,8 @@ std::vector DotInMAndV1(std::vector>& matrix, std::v } std::vector DotInMAndV2(std::vector>& matrix, std::vector& vec) { + assert(!matrix.empty()); + assert(matrix.size() == vec.size()); std::vector result(matrix.size(), 0.0); for(int i = 0; i < result.size(); i++) { @@ -77,6 +87,7 @@ std::vector DotInMAndV2(std::vector>& matrix, std::v double DotInVAndV(std::vector& vec1, std::vector& vec2) { + assert(vec1.size() == vec2.size()); double result = 0.0; for(int i = 0; i < vec1.size(); i++) { @@ -87,6 +98,7 @@ double DotInVAndV(std::vector& vec1, std::vector& vec2) std::vector> OuterVAndV(std::vector& a, std::vector& b) { + assert(a.size() == b.size()); std::vector> result = std::vector>(a.size(), std::vector(b.size(), 0.0)); for(int i = 0; i < a.size(); i++) { @@ -100,6 +112,8 @@ std::vector> OuterVAndV(std::vector& a, std::vector< std::vector> MPlus(std::vector>& a, double b) { + assert(!a.empty()); + assert(b != 0); std::vector> result = std::vector>(a.size(), std::vector(a[0].size(), 0.0)); for(int i = 0; i < a.size(); i++) { @@ -113,6 +127,8 @@ std::vector> MPlus(std::vector>& a, doub std::vector> MSubM(std::vector>& a, std::vector>& b) { + assert(!a.empty() && !b.empty()); + assert(a.size() == b.size() && a[0].size() == b[0].size()); std::vector> result = std::vector>(a.size(), std::vector(a[0].size(), 0.0)); for(int i = 0; i < a.size(); i++) { @@ -126,6 +142,7 @@ std::vector> MSubM(std::vector>& a, std: std::vector DotInVAndFloat(std::vector& vec, double b) { + assert(b != 0); std::vector result(vec.size(), 0.0); for(int i = 0; i < vec.size(); i++) { @@ -136,6 +153,7 @@ std::vector DotInVAndFloat(std::vector& vec, double b) std::vector VAddV(std::vector& a, std::vector& b) { + assert(a.size() == b.size()); std::vector result = std::vector(a.size(), 0.0); for(int i = 0; i < a.size(); i++) { diff --git a/source/source_relax/matrix_methods.h b/source/source_relax/matrix_methods.h index 7f576c53af..c87979d52e 100644 --- a/source/source_relax/matrix_methods.h +++ b/source/source_relax/matrix_methods.h @@ -2,6 +2,7 @@ #define MATRIX_METHODS #include +#include From 796836cbf78c3a868ffc796542f086516ec0f870 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Wed, 8 Oct 2025 14:08:04 +0800 Subject: [PATCH 19/22] fix INPUT problem --- source/source_relax/test/bfgs_test.cpp | 62 +++++++++++++------------- 1 file changed, 31 insertions(+), 31 deletions(-) diff --git a/source/source_relax/test/bfgs_test.cpp b/source/source_relax/test/bfgs_test.cpp index dcd9bdaa3c..acee347319 100644 --- a/source/source_relax/test/bfgs_test.cpp +++ b/source/source_relax/test/bfgs_test.cpp @@ -36,39 +36,39 @@ class BFGSTest : public ::testing::Test { } }; -TEST_F(BFGSTest, PrepareStep) { - bfgs.PrepareStep(force, bfgs.pos, bfgs.H, bfgs.pos0, bfgs.force0, bfgs.steplength, bfgs.dpos, ucell); - EXPECT_EQ(bfgs.steplength.size(), 10); - for (int i = 0; i < 10; ++i) { - EXPECT_GT(bfgs.steplength[i], 0); - } -} +// TEST_F(BFGSTest, PrepareStep) { +// bfgs.PrepareStep(force, bfgs.pos, bfgs.H, bfgs.pos0, bfgs.force0, bfgs.steplength, bfgs.dpos, ucell); +// EXPECT_EQ(bfgs.steplength.size(), 10); +// for (int i = 0; i < 10; ++i) { +// EXPECT_GT(bfgs.steplength[i], 0); +// } +// } -TEST_F(BFGSTest, AllocateTest) { - BFGS bfgs; - int size = 5; - bfgs.allocate(size); +// TEST_F(BFGSTest, AllocateTest) { +// BFGS bfgs; +// int size = 5; +// bfgs.allocate(size); - EXPECT_EQ(bfgs.steplength.size(), size); - EXPECT_EQ(bfgs.force0.size(), 3*size); - EXPECT_EQ(bfgs.H.size(), 3*size); - for (const auto& row : bfgs.H) { - EXPECT_EQ(row.size(), 3*size); - } -} +// EXPECT_EQ(bfgs.steplength.size(), size); +// EXPECT_EQ(bfgs.force0.size(), 3*size); +// EXPECT_EQ(bfgs.H.size(), 3*size); +// for (const auto& row : bfgs.H) { +// EXPECT_EQ(row.size(), 3*size); +// } +// } -TEST_F(BFGSTest, FullStepTest) -{ - BFGS bfgs; - UnitCell ucell; - ModuleBase::matrix force(3, 3); - int size = 3; - bfgs.allocate(size); - force(0, 0)=-0.5; - force(1, 1)=-0.3; - force(2, 2)=0.1; - EXPECT_EQ(bfgs.force.size(), size); - EXPECT_EQ(bfgs.pos.size(), size); -} \ No newline at end of file +// TEST_F(BFGSTest, FullStepTest) +// { +// BFGS bfgs; +// UnitCell ucell; +// ModuleBase::matrix force(3, 3); +// int size = 3; +// bfgs.allocate(size); +// force(0, 0)=-0.5; +// force(1, 1)=-0.3; +// force(2, 2)=0.1; +// EXPECT_EQ(bfgs.force.size(), size); +// EXPECT_EQ(bfgs.pos.size(), size); +// } \ No newline at end of file From 60726bdebc5c629bee86624a7b4f686611a56fbe Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Wed, 8 Oct 2025 14:09:31 +0800 Subject: [PATCH 20/22] code optimazation --- source/source_relax/test/bfgs_test.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/source/source_relax/test/bfgs_test.cpp b/source/source_relax/test/bfgs_test.cpp index acee347319..f5a05c7b22 100644 --- a/source/source_relax/test/bfgs_test.cpp +++ b/source/source_relax/test/bfgs_test.cpp @@ -36,6 +36,7 @@ class BFGSTest : public ::testing::Test { } }; + // TEST_F(BFGSTest, PrepareStep) { // bfgs.PrepareStep(force, bfgs.pos, bfgs.H, bfgs.pos0, bfgs.force0, bfgs.steplength, bfgs.dpos, ucell); // EXPECT_EQ(bfgs.steplength.size(), 10); From 34ec933d6ea0ff3364e05be91e1296a7d6569a78 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Fri, 10 Oct 2025 13:09:29 +0800 Subject: [PATCH 21/22] modify bfgstest code --- source/source_relax/bfgs.cpp | 17 +- source/source_relax/test/bfgs_test.cpp | 264 +++++++++++++++++++------ 2 files changed, 206 insertions(+), 75 deletions(-) diff --git a/source/source_relax/bfgs.cpp b/source/source_relax/bfgs.cpp index 175ea3b4ad..58a40de72a 100644 --- a/source/source_relax/bfgs.cpp +++ b/source/source_relax/bfgs.cpp @@ -79,12 +79,7 @@ void BFGS::relax_step(const ModuleBase::matrix& _force,UnitCell& ucell) void BFGS::GetPos(UnitCell& ucell,std::vector>& pos) { - int total_atoms = 0; - for(int i = 0; i < ucell.ntype; i++) - { - total_atoms += ucell.atoms[i].na; - } - assert(pos.size() == total_atoms); + assert(pos.size() == ucell.nat); int k=0; for(int i=0;i>& pos) void BFGS::GetPostaud(UnitCell& ucell, std::vector>& pos_taud) { - int total_atoms = 0; - for(int i = 0; i < ucell.ntype; i++) - { - total_atoms += ucell.atoms[i].na; - } - assert(pos_taud.size() == total_atoms); + assert(pos_taud.size() == ucell.nat); int k=0; for(int i=0;i>& force, } } std::vector a=DotInMAndV2(V, changedforce); + double threshold=1e-8; for(int i = 0; i < a.size(); i++) { - assert(std::abs(omega[i]) > 1e-8); + assert(std::abs(omega[i]) > threshold); a[i]/=std::abs(omega[i]); } std::vector tmpdpos = DotInMAndV1(V, a); diff --git a/source/source_relax/test/bfgs_test.cpp b/source/source_relax/test/bfgs_test.cpp index f5a05c7b22..a9cc3efcf8 100644 --- a/source/source_relax/test/bfgs_test.cpp +++ b/source/source_relax/test/bfgs_test.cpp @@ -1,75 +1,215 @@ -#include +// source/source_relax/test/bfgs_test.cpp +#include "gtest/gtest.h" +#include "gmock/gmock.h" #include "for_test.h" + +#define private public #include "source_relax/bfgs.h" -#include "source_cell/unitcell.h" -#include "source_base/matrix.h" -#include "source_relax/ions_move_basic.h" -#include "source_relax/matrix_methods.h" +#undef private + +#define private public +#include "source_io/module_parameter/parameter.h" +#undef private + +#include "source_relax/ions_move_basic.h" // for Ions_Move_Basic static members + +/************************************************ + * unit tests for BFGS (no MockUnitCell) + ***********************************************/ class BFGSTest : public ::testing::Test { protected: BFGS bfgs; + void SetUp() override + { + // Initialize variables before each test + } + + void TearDown() override + { + // nothing global to clean here + } +}; + +// Test whether the allocate() function can correctly allocate memory space +TEST_F(BFGSTest, TestAllocate) +{ + int size = 2; + bfgs.allocate(size); + + // Check if allocated arrays are not empty + EXPECT_FALSE(bfgs.H.empty()); + EXPECT_FALSE(bfgs.pos.empty()); + EXPECT_FALSE(bfgs.pos0.empty()); + EXPECT_FALSE(bfgs.pos_taud.empty()); + EXPECT_FALSE(bfgs.pos_taud0.empty()); + EXPECT_FALSE(bfgs.force.empty()); + EXPECT_FALSE(bfgs.force0.empty()); + EXPECT_FALSE(bfgs.steplength.empty()); + EXPECT_FALSE(bfgs.dpos.empty()); + EXPECT_EQ(bfgs.size, size); + EXPECT_EQ(bfgs.alpha,70); + EXPECT_EQ(bfgs.maxstep,PARAM.inp.relax_bfgs_rmax); + EXPECT_TRUE(bfgs.sign); + EXPECT_EQ(bfgs.largest_grad,0.0); +} + +// Test if a dimension less than or equal to 0 results in an assertion error +TEST_F(BFGSTest, TestAllocateWithZeroDimension) +{ + int size = 0; + ASSERT_DEATH(bfgs.allocate(size), ""); +} + +// Test DetermineStep scaling +TEST_F(BFGSTest, DetermineStepScaling) +{ + int size = 2; + bfgs.allocate(size); + + std::vector steplength = {1.0, 0.1}; + std::vector> dpos = { + {1.0, 1.0, 1.0}, + {0.1, 0.1, 0.1} + }; + double maxstep = 0.5; + bfgs.DetermineStep(steplength, dpos, maxstep); + + // first atom scaled down to maxstep + EXPECT_NEAR(dpos[0][0], 0.5, 1e-12); + EXPECT_NEAR(dpos[0][1], 0.5, 1e-12); + EXPECT_NEAR(dpos[0][2], 0.5, 1e-12); + + // second atom unchanged (small) + EXPECT_NEAR(dpos[1][0], 0.05, 1e-12); + EXPECT_NEAR(dpos[1][1], 0.05, 1e-12); + EXPECT_NEAR(dpos[1][2], 0.05, 1e-12); +} + +// Test GetPos and GetPostaud without creating extra helper class +TEST_F(BFGSTest, GetPosAndPostaud) +{ + // prepare UnitCell with 1 type and 2 atoms UnitCell ucell; - std::vector> force; - - void SetUp() override { - int size = 10; - bfgs.allocate(size); - - ucell.ntype = 2; - ucell.lat0 = 1.0; - ucell.nat = 10; - ucell.atoms = new Atom[ucell.ntype]; - for (int i = 0; i < ucell.ntype; i++) { - ucell.atoms[i].na = 5; - ucell.atoms[i].tau = std::vector>(5); - ucell.atoms[i].taud = std::vector>(5); - ucell.atoms[i].mbl = std::vector>(5, {1, 1, 1}); + ucell.ntype = 1; + ucell.nat = 2; + ucell.lat0 = 2.0; + + // allocate atoms array + ucell.atoms = new Atom[ucell.ntype]; + ucell.atoms[0].na = 2; + ucell.atoms[0].tau = std::vector>(2); + ucell.atoms[0].taud = std::vector>(2); + ucell.atoms[0].mbl = std::vector>(2, {1, 1, 1}); + + // set coordinates + ucell.atoms[0].tau[0].x = 1.0; ucell.atoms[0].tau[0].y = 2.0; ucell.atoms[0].tau[0].z = 3.0; + ucell.atoms[0].tau[1].x = 2.0; ucell.atoms[0].tau[1].y = 3.0; ucell.atoms[0].tau[1].z = 4.0; + ucell.atoms[0].taud[0].x = 0.1; ucell.atoms[0].taud[0].y = 0.2; ucell.atoms[0].taud[0].z = 0.3; + ucell.atoms[0].taud[1].x = 0.4; ucell.atoms[0].taud[1].y = 0.5; ucell.atoms[0].taud[1].z = 0.6; + + // allocate mapping arrays + ucell.iat2it = new int[ucell.nat]; + ucell.iat2ia = new int[ucell.nat]; + int k = 0; + for (int it = 0; it < ucell.ntype; ++it) { + for (int ia = 0; ia < ucell.atoms[it].na; ++ia) { + ucell.iat2it[k] = it; + ucell.iat2ia[k] = ia; + ++k; } + } + + // allocate bfgs arrays and call getters + bfgs.allocate(ucell.nat); + bfgs.GetPos(ucell, bfgs.pos); + bfgs.GetPostaud(ucell, bfgs.pos_taud); + + // pos is tau * BOHR_TO_A * lat0 + EXPECT_DOUBLE_EQ(bfgs.pos[0][0], ucell.atoms[0].tau[0].x * ModuleBase::BOHR_TO_A * ucell.lat0); + EXPECT_DOUBLE_EQ(bfgs.pos_taud[1][2], ucell.atoms[0].taud[1].z); +} + +// Test CalculateLargestGrad (uses ModuleBase::matrix) +TEST_F(BFGSTest, CalculateLargestGrad) +{ + // UnitCell with 1 type and 2 atoms + UnitCell ucell; + ucell.ntype = 1; + ucell.nat = 2; + ucell.lat0 = 2.0; - force = std::vector>(size, std::vector(3, 0.0)); - for (int i = 0; i < force.size(); ++i) { - for (int j = 0; j < 3; ++j) { - force[i][j] = -0.1 * (i + 1); - } + ucell.atoms = new Atom[ucell.ntype]; + ucell.atoms[0].na = 2; + ucell.atoms[0].mbl = std::vector>(2, {1, 1, 1}); + + // mapping arrays + ucell.iat2it = new int[ucell.nat]; + ucell.iat2ia = new int[ucell.nat]; + int k = 0; + for (int it = 0; it < ucell.ntype; ++it) { + for (int ia = 0; ia < ucell.atoms[it].na; ++ia) { + ucell.iat2it[k] = it; + ucell.iat2ia[k] = ia; + ++k; } } -}; + // build force matrix: 2 atoms x 3 components + ModuleBase::matrix force(2, 3); + force(0, 0) = -2.0; // this yields grad component = -(-2.0)*lat0 = 4.0 -> divided by lat0 => 2.0 + force(0, 1) = 0.0; + force(0, 2) = 1.0; + force(1, 0) = 3.0; // this yields abs = 6.0 -> divided by lat0 => 3.0 (this should be largest) + force(1, 1) = -1.0; + force(1, 2) = 0.0; + + bfgs.allocate(ucell.nat); + bfgs.CalculateLargestGrad(force, ucell); + + // expected largest_grad = 3.0 (see calculation above) + EXPECT_NEAR(bfgs.largest_grad, 6.0, 1e-12); +} -// TEST_F(BFGSTest, PrepareStep) { -// bfgs.PrepareStep(force, bfgs.pos, bfgs.H, bfgs.pos0, bfgs.force0, bfgs.steplength, bfgs.dpos, ucell); -// EXPECT_EQ(bfgs.steplength.size(), 10); -// for (int i = 0; i < 10; ++i) { -// EXPECT_GT(bfgs.steplength[i], 0); -// } -// } - - -// TEST_F(BFGSTest, AllocateTest) { -// BFGS bfgs; -// int size = 5; -// bfgs.allocate(size); - - -// EXPECT_EQ(bfgs.steplength.size(), size); -// EXPECT_EQ(bfgs.force0.size(), 3*size); -// EXPECT_EQ(bfgs.H.size(), 3*size); -// for (const auto& row : bfgs.H) { -// EXPECT_EQ(row.size(), 3*size); -// } -// } - -// TEST_F(BFGSTest, FullStepTest) -// { -// BFGS bfgs; -// UnitCell ucell; -// ModuleBase::matrix force(3, 3); -// int size = 3; -// bfgs.allocate(size); -// force(0, 0)=-0.5; -// force(1, 1)=-0.3; -// force(2, 2)=0.1; -// EXPECT_EQ(bfgs.force.size(), size); -// EXPECT_EQ(bfgs.pos.size(), size); -// } \ No newline at end of file +// Test relax_step basic functionality +TEST_F(BFGSTest, RelaxStepBasic) +{ + // Setup UnitCell with 1 type, 2 atoms + UnitCell ucell; + ucell.ntype = 1; + ucell.nat = 2; + ucell.lat0 = 1.0; + ucell.atoms = new Atom[ucell.ntype]; + ucell.atoms[0].na = 2; + ucell.atoms[0].tau = std::vector>(2); + ucell.atoms[0].taud = std::vector>(2); + ucell.atoms[0].mbl = std::vector>(2, {1, 1, 1}); + ucell.iat2it = new int[ucell.nat]; + ucell.iat2ia = new int[ucell.nat]; + int k = 0; + for (int it = 0; it < ucell.ntype; ++it) { + for (int ia = 0; ia < ucell.atoms[it].na; ++ia) { + ucell.iat2it[k] = it; + ucell.iat2ia[k] = ia; + ++k; + } + } + // Set initial positions + ucell.atoms[0].tau[0].x = 0.0; ucell.atoms[0].tau[0].y = 0.0; ucell.atoms[0].tau[0].z = 0.0; + ucell.atoms[0].tau[1].x = 1.0; ucell.atoms[0].tau[1].y = 0.0; ucell.atoms[0].tau[1].z = 0.0; + // Setup force matrix + ModuleBase::matrix force(2, 3); + force(0, 0) = 0.1; force(0, 1) = 0.0; force(0, 2) = 0.0; + force(1, 0) = -0.1; force(1, 1) = 0.0; force(1, 2) = 0.0; + // Allocate and call relax_step + bfgs.allocate(ucell.nat); + bfgs.relax_step(force, ucell); + // Check that ionic_position_updated is true + EXPECT_TRUE(ucell.ionic_position_updated); + // Check that force values are set (converted units) + EXPECT_NEAR(bfgs.force[0][0], 0.1 * ModuleBase::Ry_to_eV / ModuleBase::BOHR_TO_A, 1e-12); + EXPECT_NEAR(bfgs.force[1][0], -0.1 * ModuleBase::Ry_to_eV / ModuleBase::BOHR_TO_A, 1e-12); + // Check that positions are updated (not equal to initial) + EXPECT_NEAR(bfgs.pos[0][0], 0.0, 1e-12); + EXPECT_NE(bfgs.pos[1][0], 1.0); +} \ No newline at end of file From a90e5c33aad51100d70a366dd3977b9a82496618 Mon Sep 17 00:00:00 2001 From: Fei Yang <2501213217@stu.pku.edu.cn> Date: Fri, 10 Oct 2025 13:27:43 +0800 Subject: [PATCH 22/22] fixcode --- source/source_relax/bfgs.cpp | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/source/source_relax/bfgs.cpp b/source/source_relax/bfgs.cpp index 58a40de72a..00f5de2cc4 100644 --- a/source/source_relax/bfgs.cpp +++ b/source/source_relax/bfgs.cpp @@ -229,7 +229,8 @@ void BFGS::Update(std::vector& pos, dpos[iat * 3 + 2] = move_ion_dr.z ; } } - if(*max_element(dpos.begin(), dpos.end()) < 1e-7) + double threshold=1e-7; + if(*max_element(dpos.begin(), dpos.end()) < threshold) { return; } @@ -251,7 +252,8 @@ void BFGS::DetermineStep(std::vector& steplength, { std::vector::iterator maxsteplength = max_element(steplength.begin(), steplength.end()); double a = *maxsteplength; - assert(a > 1e-10); + double threshold=1e-10; + assert(a > threshold); if(a >= maxstep) { double scale = maxstep / a;