Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 9 additions & 4 deletions docs/advanced/input_files/input-main.md
Original file line number Diff line number Diff line change
Expand Up @@ -1487,15 +1487,20 @@ 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.
- 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.
- **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
- **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

Expand Down
3 changes: 2 additions & 1 deletion source/source_io/module_parameter/input_parameter.h
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
#include <string>
#include <vector>


// It stores all input parameters both defined in INPUT file and not defined in
// INPUT file
struct Input_para
Expand Down Expand Up @@ -150,7 +151,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::vector<std::string> 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;
Expand Down
43 changes: 34 additions & 9 deletions source/source_io/read_input_item_relax.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -5,28 +5,53 @@

namespace ModuleIO
{


void ReadInput::item_relax()
{
{
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<std::string> relax_methods = {"cg", "bfgs", "sd", "cg_bfgs","bfgs_trad","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);
}
item.read_value = [](const Input_Item& item, Parameter& para) {
if(item.get_size()==1)
{
para.input.relax_method[0] = item.str_values[0];
para.input.relax_method[1] = "1";
}
else if(item.get_size()>=2)
{
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<std::string> relax_methods = {"cg", "sd", "cg_bfgs","lbfgs","bfgs"};
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);
}
};
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<std::string> 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");
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;
}
Expand Down
2 changes: 1 addition & 1 deletion source/source_io/test/read_input_ptest.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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);
Expand Down
6 changes: 3 additions & 3 deletions source/source_io/test_serial/read_input_item_test.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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();
Expand All @@ -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);
}
Expand Down
2 changes: 1 addition & 1 deletion source/source_relax/ions_move_basic.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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::string Ions_Move_Basic::relax_method = "bfgs";
std::vector<std::string> Ions_Move_Basic::relax_method = {"bfgs","2"};

void Ions_Move_Basic::setup_gradient(const UnitCell &ucell, const ModuleBase::matrix &force, double *pos, double *grad)
{
Expand Down
2 changes: 1 addition & 1 deletion source/source_relax/ions_move_basic.h
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ 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<std::string> relax_method; // relaxation method,
extern int out_stru; // output the structure or not
// funny way to pass this parameter, but nevertheless

Expand Down
85 changes: 25 additions & 60 deletions source/source_relax/ions_move_cg.cpp
Original file line number Diff line number Diff line change
@@ -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"
Expand Down Expand Up @@ -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];
Expand All @@ -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);
Expand Down Expand Up @@ -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);
Expand All @@ -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);

Expand All @@ -181,16 +160,16 @@ 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;

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::relax_method[1] = "2";
}
Ions_Move_Basic::best_xxx = steplength;
}
Expand All @@ -204,7 +183,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;
Expand All @@ -214,21 +192,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_grad0, dim);
third_order(e0, e1, fa, fb, xb, best_x); // cubic interpolation

if (best_x > 6 * xb || best_x < (-xb))
Expand All @@ -238,7 +202,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);
Expand All @@ -250,7 +213,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++;

Expand All @@ -275,24 +237,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_grad0, dim);
setup_move(move, cg_gradn, best_x);
Ions_Move_Basic::move_atoms(ucell, move, pos);

Ions_Move_Basic::relax_bfgs_init = xc;
}
}
Expand Down Expand Up @@ -342,7 +289,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
Expand Down Expand Up @@ -508,3 +454,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;
}
1 change: 1 addition & 0 deletions source/source_relax/ions_move_cg.h
Original file line number Diff line number Diff line change
Expand Up @@ -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
Loading
Loading