Skip to content
Open
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
14 changes: 10 additions & 4 deletions source/source_relax/bfgs_basic.cpp
Original file line number Diff line number Diff line change
@@ -1,9 +1,11 @@
#include "bfgs_basic.h"
#include <algorithm>

#include "ions_move_basic.h"
#include "source_base/global_function.h"
#include "source_base/global_variable.h"

#include <algorithm>

using namespace Ions_Move_Basic;

double BFGS_Basic::relax_bfgs_w1 = -1.0; // default is 0.01
Expand Down Expand Up @@ -31,7 +33,7 @@ void BFGS_Basic::allocate_basic(void)
return;
}

void BFGS_Basic::update_inverse_hessian(const double &lat0, std::ofstream& ofs)
void BFGS_Basic::update_inverse_hessian(const double& lat0, std::ofstream& ofs)
{
// ModuleBase::TITLE("Ions_Move_BFGS","update_inverse_hessian");
assert(dim > 0);
Expand Down Expand Up @@ -157,7 +159,11 @@ void BFGS_Basic::save_bfgs(void)
// a new bfgs step is done
// we have already done well in the previous direction
// we should get a new direction in this case
void BFGS_Basic::new_step(const double &lat0, int& update_iter, std::ofstream& ofs, std::vector<double>& etot_info, const int test_relax_method)
void BFGS_Basic::new_step(const double& lat0,
int& update_iter,
std::ofstream& ofs,
std::vector<double>& etot_info,
const int test_relax_method)
{
ModuleBase::TITLE("BFGS_Basic", "new_step");

Expand Down Expand Up @@ -341,7 +347,7 @@ void BFGS_Basic::compute_trust_radius(std::ofstream& ofs, std::vector<double>& e
return;
}

double BFGS_Basic::check_move(const double &lat0, const double &pos, const double &pos_p)
double BFGS_Basic::check_move(const double& lat0, const double& pos, const double& pos_p)
{
// this must be careful.
// unit is ucell.lat0.
Expand Down
17 changes: 11 additions & 6 deletions source/source_relax/bfgs_basic.h
Original file line number Diff line number Diff line change
@@ -1,9 +1,10 @@
#ifndef BFGS_BASIC
#define BFGS_BASIC

#include "source_base/matrix.h"

#include <fstream>
#include <iostream>
#include "source_base/matrix.h"
#include <vector>

// references
Expand Down Expand Up @@ -129,7 +130,11 @@ class BFGS_Basic

protected:
void allocate_basic(void);
void new_step(const double& lat0, int& update_iter, std::ofstream& ofs, std::vector<double>& etot_info, const int test_relax_method);
void new_step(const double& lat0,
int& update_iter,
std::ofstream& ofs,
std::vector<double>& etot_info,
const int test_relax_method);
void reset_hessian(void);
void save_bfgs(void);

Expand All @@ -146,15 +151,15 @@ class BFGS_Basic
static double relax_bfgs_w2; // fixed: parameters for Wolfe conditions.

protected:
bool save_flag=false;
bool tr_min_hit=false; //.TRUE. if the trust_radius has already been set
// to the minimum value at the previous step
bool save_flag = false;
bool tr_min_hit = false; //.TRUE. if the trust_radius has already been set
// to the minimum value at the previous step

// mohan add 2010-07-27
double check_move(const double& lat0, const double& pos, const double& pos_p);

private:
bool wolfe_flag=false;
bool wolfe_flag = false;
ModuleBase::matrix inv_hess;

int bfgs_ndim;
Expand Down
35 changes: 21 additions & 14 deletions source/source_relax/cg_base.cpp
Original file line number Diff line number Diff line change
@@ -1,10 +1,16 @@
#include "cg_base.h"
#include <cmath>

#include "source_base/global_function.h"

void CG_Base::setup_cg_grad(const int dim, double *grad, const double *grad0,
double *cg_grad, const double *cg_grad0,
const int &ncggrad, int &flag)
#include <cmath>

void CG_Base::setup_cg_grad(const int dim,
double* grad,
const double* grad0,
double* cg_grad,
const double* cg_grad0,
const int& ncggrad,
int& flag)
{
ModuleBase::TITLE("CG_Base", "setup_cg_grad");
double gamma = 0.0;
Expand Down Expand Up @@ -32,7 +38,7 @@ void CG_Base::setup_cg_grad(const int dim, double *grad, const double *grad0,
cgp_g += cg_grad0[i] * grad[i];
}

const double gamma1 = gg / gp_gp; // FR
const double gamma1 = gg / gp_gp; // FR
const double gamma2 = (gg - g_gp) / gp_gp; // PRP

if (gamma1 < 0.5)
Expand All @@ -51,7 +57,7 @@ void CG_Base::setup_cg_grad(const int dim, double *grad, const double *grad0,
}
}

void CG_Base::setup_move(const int dim, double *move, double *cg_gradn, const double &trust_radius)
void CG_Base::setup_move(const int dim, double* move, double* cg_gradn, const double& trust_radius)
{
ModuleBase::TITLE("CG_Base", "setup_move");
for (int i = 0; i < dim; ++i)
Expand All @@ -60,9 +66,7 @@ void CG_Base::setup_move(const int dim, double *move, double *cg_gradn, const do
}
}

void CG_Base::Brent(double &fa, double &fb, double &fc,
double &xa, double &xb, double &xc,
double &best_x, double &xpt)
void CG_Base::Brent(double& fa, double& fb, double& fc, double& xa, double& xb, double& xc, double& best_x, double& xpt)
{
ModuleBase::TITLE("CG_Base", "Brent");
double dmove = 0.0;
Expand Down Expand Up @@ -129,7 +133,7 @@ void CG_Base::Brent(double &fa, double &fb, double &fc,
xc = dmove;
}

void CG_Base::f_cal(const int dim, const double *g0, const double *g1, double &f_value)
void CG_Base::f_cal(const int dim, const double* g0, const double* g1, double& f_value)
{
ModuleBase::TITLE("CG_Base", "f_cal");
double hv0 = 0.0;
Expand All @@ -142,9 +146,12 @@ void CG_Base::f_cal(const int dim, const double *g0, const double *g1, double &f
f_value = hel / sqrt(hv0);
}

void CG_Base::third_order(const double &e0, const double &e1,
const double &fa, const double &fb,
const double x, double &best_x)
void CG_Base::third_order(const double& e0,
const double& e1,
const double& fa,
const double& fb,
const double x,
double& best_x)
{
ModuleBase::TITLE("CG_Base", "third_order");
double k3 = 0.0;
Expand Down Expand Up @@ -193,7 +200,7 @@ void CG_Base::third_order(const double &e0, const double &e1,
best_x = dmove;
}

void CG_Base::normalize(const int dim, double *cg_gradn, const double *cg_grad)
void CG_Base::normalize(const int dim, double* cg_gradn, const double* cg_grad)
{
ModuleBase::TITLE("CG_Base", "normalize");
double norm = 0.0;
Expand Down
39 changes: 22 additions & 17 deletions source/source_relax/cg_base.h
Original file line number Diff line number Diff line change
Expand Up @@ -4,13 +4,14 @@
/**
* @class CG_Base
* @brief Base class for conjugate gradient optimization algorithms.
*
*
* This class provides common CG utility functions that can be shared between
* different optimization implementations (ions movement and lattice change).
* The functions are designed to be independent of specific data sources.
*/
class CG_Base {
public:
class CG_Base
{
public:
/**
* @brief Setup conjugate gradient direction.
* @param dim Dimension of the optimization problem
Expand All @@ -21,9 +22,13 @@ class CG_Base {
* @param ncggrad CG iteration counter
* @param flag Restart flag
*/
void setup_cg_grad(const int dim, double *grad, const double *grad0,
double *cg_grad, const double *cg_grad0,
const int &ncggrad, int &flag);
void setup_cg_grad(const int dim,
double* grad,
const double* grad0,
double* cg_grad,
const double* cg_grad0,
const int& ncggrad,
int& flag);

/**
* @brief Brent's method for one-dimensional minimization.
Expand All @@ -36,9 +41,7 @@ class CG_Base {
* @param best_x Output best step length
* @param xpt Previous best point
*/
void Brent(double &fa, double &fb, double &fc,
double &xa, double &xb, double &xc,
double &best_x, double &xpt);
void Brent(double& fa, double& fb, double& fc, double& xa, double& xb, double& xc, double& best_x, double& xpt);

/**
* @brief Third-order polynomial interpolation for line search.
Expand All @@ -49,9 +52,12 @@ class CG_Base {
* @param x Current step length
* @param best_x Output optimal step length
*/
void third_order(const double &e0, const double &e1,
const double &fa, const double &fb,
const double x, double &best_x);
void third_order(const double& e0,
const double& e1,
const double& fa,
const double& fb,
const double x,
double& best_x);

/**
* @brief Calculate projection of gradient onto search direction.
Expand All @@ -60,7 +66,7 @@ class CG_Base {
* @param g1 Second vector (usually gradient)
* @param f_value Output projection value
*/
void f_cal(const int dim, const double *g0, const double *g1, double &f_value);
void f_cal(const int dim, const double* g0, const double* g1, double& f_value);

/**
* @brief Calculate movement vector from CG direction.
Expand All @@ -69,17 +75,16 @@ class CG_Base {
* @param cg_gradn Normalized CG direction
* @param trust_radius Step size limit
*/
void setup_move(const int dim, double *move, double *cg_gradn, const double &trust_radius);

protected:
void setup_move(const int dim, double* move, double* cg_gradn, const double& trust_radius);

protected:
/**
* @brief Normalize vector.
* @param dim Dimension of vector
* @param cg_gradn Output normalized vector
* @param cg_grad Input vector
*/
void normalize(const int dim, double *cg_gradn, const double *cg_grad);
void normalize(const int dim, double* cg_gradn, const double* cg_grad);
};

#endif
Loading
Loading