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
5 changes: 4 additions & 1 deletion source/source_main/driver_run.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -129,7 +129,10 @@ void Driver::driver_run()
{
ModuleBase::WARNING_QUIT("Driver::driver_run", "cannot recognize the 'calculation' command");
}
p_esolver->after_all_runners(ucell);
if (cal != "scf" && cal != "relax" && cal != "cell-relax" && cal != "nscf")
{
p_esolver->after_all_runners(ucell);
}
}

delete p_esolver;
Expand Down
22 changes: 15 additions & 7 deletions source/source_relax/ions_move_bfgs2.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
#include "source_cell/print_cell.h" // lanshuyue add 2025-06-19
#include "source_cell/update_cell.h"
#include "source_io/module_parameter/parameter.h"
#include "source_relax/relax_criteria.h"

//! initialize H0、H、pos0、force0、force
void Ions_Move_BFGS2::allocate(const int _size)
Expand Down Expand Up @@ -38,15 +39,17 @@ void Ions_Move_BFGS2::reset()
is_initialized = false;
}

bool Ions_Move_BFGS2::relax_step(const ModuleBase::matrix& _force, UnitCell& ucell, std::ofstream& ofs_running)
bool Ions_Move_BFGS2::relax_step(const ModuleBase::matrix& _force,
UnitCell& ucell,
std::ofstream& ofs_running,
const Relax_Criteria& criteria)
{
if (!is_initialized)
{
allocate(ucell.nat);
}
GetPos(ucell, pos);
GetPostaud(ucell, pos_taud);
ucell.ionic_position_updated = true;
assert(_force.nr == force.size() && _force.nc == 3);
for (int i = 0; i < _force.nr; i++)
{
Expand Down Expand Up @@ -79,14 +82,19 @@ bool Ions_Move_BFGS2::relax_step(const ModuleBase::matrix& _force, UnitCell& uce
k += ucell.atoms[i].na;
}

this->CalculateLargestGrad(_force, ucell);
if (this->IsRestrain(criteria.force_thr_ev))
{
return true;
}

this->PrepareStep(force, pos, H, pos0, force0, steplength, dpos, size, ucell);
this->DetermineStep(steplength, dpos, maxstep);
this->UpdatePos(ucell);
this->CalculateLargestGrad(_force, ucell);
bool converged = this->IsRestrain();
ucell.ionic_position_updated = true;
// print out geometry information during bfgs_trad relax
unitcell::print_tau(ucell.atoms, ucell.Coordinate, ucell.ntype, ucell.lat0, ofs_running);
return converged;
return false;
}

void Ions_Move_BFGS2::GetPos(UnitCell& ucell, std::vector<ModuleBase::Vector3<double>>& pos)
Expand Down Expand Up @@ -370,7 +378,7 @@ void Ions_Move_BFGS2::CalculateLargestGrad(const ModuleBase::matrix& _force, Uni
}
}

bool Ions_Move_BFGS2::IsRestrain()
bool Ions_Move_BFGS2::IsRestrain(const double force_thr_ev)
{
return Ions_Move_Basic::largest_grad * ModuleBase::Ry_to_eV / ModuleBase::BOHR_TO_A < PARAM.inp.force_thr_ev;
return Ions_Move_Basic::largest_grad * ModuleBase::Ry_to_eV / ModuleBase::BOHR_TO_A < force_thr_ev;
}
7 changes: 5 additions & 2 deletions source/source_relax/ions_move_bfgs2.h
Original file line number Diff line number Diff line change
Expand Up @@ -11,14 +11,17 @@
#include <tuple>
#include <vector>

struct Relax_Criteria;

class Ions_Move_BFGS2
{
public:
void allocate(const int _size); // initialize parameters
void reset();
bool relax_step(const ModuleBase::matrix& _force,
UnitCell& ucell,
std::ofstream& ofs_running); // a full iteration step
std::ofstream& ofs_running,
const Relax_Criteria& criteria); // a full iteration step

/// @brief whether allocate() has run and the optimiser holds live state.
/// reset() clears it; ions_move_methods_test checks that.
Expand Down Expand Up @@ -57,7 +60,7 @@ class Ions_Move_BFGS2
std::vector<ModuleBase::Vector3<double>>& dpos,
int& size,
UnitCell& ucell); // calculate the atomic displacement in one iteration step
bool IsRestrain(); // check if converged
bool IsRestrain(const double force_thr_ev); // check if converged
void CalculateLargestGrad(const ModuleBase::matrix& _force, UnitCell& ucell);
void GetPos(UnitCell& ucell, std::vector<ModuleBase::Vector3<double>>& pos);
void GetPostaud(UnitCell& ucell, std::vector<ModuleBase::Vector3<double>>& pos_taud);
Expand Down
12 changes: 7 additions & 5 deletions source/source_relax/ions_move_lbfgs.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -36,8 +36,6 @@ bool Ions_Move_LBFGS::relax_step(const ModuleBase::matrix _force,
{
get_pos(ucell, pos);
get_pos_taud(ucell, pos_taud);
// solver=p_esolver;
ucell.ionic_position_updated = true;
for (int i = 0; i < _force.nr; i++)
{
for (int j = 0; j < _force.nc; j++)
Expand Down Expand Up @@ -65,14 +63,18 @@ bool Ions_Move_LBFGS::relax_step(const ModuleBase::matrix _force,
}
k += ucell.atoms[i].na;
}
this->calculate_largest_grad(_force, ucell);
if (this->is_restrain(criteria))
{
return true;
}
this->prepare_step(force, pos, H, pos0, force0, dpos, ucell, etot);
this->determine_step(steplength, dpos, maxstep);
this->update_pos(ucell);
this->calculate_largest_grad(_force, ucell);
bool converged = this->is_restrain(criteria);
ucell.ionic_position_updated = true;
// mohan add 2025-06-22
unitcell::print_tau(ucell.atoms, ucell.Coordinate, ucell.ntype, ucell.lat0, ofs_running);
return converged;
return false;
}

void Ions_Move_LBFGS::get_pos(UnitCell& ucell, std::vector<ModuleBase::Vector3<double>>& pos)
Expand Down
2 changes: 1 addition & 1 deletion source/source_relax/ions_move_methods.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -83,7 +83,7 @@ void Ions_Move_Methods::cal_movement(const int& istep,
}
else if (relax_method == "bfgs_trad")
{
converged_ = bfgs_trad.relax_step(f, ucell, ofs);
converged_ = bfgs_trad.relax_step(f, ucell, ofs, criteria);
}
else if (relax_method == "lbfgs")
{
Expand Down
104 changes: 74 additions & 30 deletions source/source_relax/relax_driver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,7 @@ void Relax_Driver::relax_driver(ModuleESolver::ESolver* p_esolver,
{
Socket_Driver socket_driver;
socket_driver.socket_driver(p_esolver, ucell, inp, ofs_running);
p_esolver->after_all_runners(ucell);
ModuleBase::timer::end("Relax_Driver", "relax_driver");
return;
}
Expand All @@ -51,36 +52,73 @@ void Relax_Driver::relax_driver(ModuleESolver::ESolver* p_esolver,
// After relax_step() proposes a new geometry, it is not evaluated until
// the next esolve() call; if we exit the loop early, force/stress are stale.
bool geometry_evaluated = false;
RelaxStopReason stop_reason = RelaxStopReason::Finished;
bool after_all_runners_called = false;

if (inp.relax_nmax == 0 && (inp.calculation == "relax" || inp.calculation == "cell-relax"))
{
stop_reason = RelaxStopReason::DryRun;
}

while (steps[0] < inp.relax_nmax)
{
this->iter_info(steps, inp);
this->esolve(steps[0], p_esolver, ucell, inp, force, stress, etot);
geometry_evaluated = true;
this->stru_out(steps[0], ucell, inp, etot, stress, force);
bool converged = this->relax_step(steps, p_esolver, ucell, inp, force, stress, etot, ofs_running);
if (!converged)
this->json_out(p_esolver, ucell, inp, force, stress);

if (ModuleIO::read_exit_file(my_rank_, "EXIT", ofs_running))
{
geometry_evaluated = false;
stop_reason = RelaxStopReason::UserExit;
break;
}
this->json_out(p_esolver, ucell, inp, force, stress);

const bool laststep =
(inp.calculation == "relax" || inp.calculation == "cell-relax")
&& steps[0] + 1 == inp.relax_nmax;
if (laststep)
{
p_esolver->after_all_runners(ucell);
after_all_runners_called = true;
}

const bool converged = this->relax_step(steps, p_esolver, ucell, inp, force, stress, etot, ofs_running);

// Check stop conditions
if (ucell.ionic_position_updated || ucell.cell_parameter_updated)
{
geometry_evaluated = false;
Comment on lines +89 to 91

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For the default relax_method=cg, these flags remain set after an earlier Relax::move_cell_ions() call; the later converged return does not clear them. That evaluated frame is then marked unevaluated, so final structure output loses its computed energy, stress, and forces. Reset the flags for each optimizer step after the solver consumes prior changes, and cover a move-then-converge sequence.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The reset now clears both flags before the converged return. Non-blocking: ClearsPreviousGeometryUpdate uses calculation="relax", so cell_parameter_updated stays false throughout and that assertion would pass even without its reset. Add a cell-relax move followed by convergence, asserting the cell flag is true after the move and false after convergence.

}

if (converged)
{
// Relaxation converged, exit loop immediately
if (inp.calculation == "relax" || inp.calculation == "cell-relax")
{
stop_reason = RelaxStopReason::Converged;
}
else
{
stop_reason = RelaxStopReason::Finished;
}
break;
Comment thread
19hello marked this conversation as resolved.
}
else if (ModuleIO::read_exit_file(my_rank_, "EXIT", ofs_running))

if (steps[0] + 1 >= inp.relax_nmax)
{
// EXIT file detected, exit loop
stop_reason = RelaxStopReason::MaxSteps;
break;
}

++steps[0];
}

this->final_out(steps[0], ucell, inp, etot, stress, force, geometry_evaluated, ofs_running);
this->final_out(steps[0], stop_reason, ucell, inp, etot, stress, force, geometry_evaluated, ofs_running);

if (!after_all_runners_called)
{
p_esolver->after_all_runners(ucell);
}

ModuleBase::timer::end("Relax_Driver", "relax_driver");
return;
Expand Down Expand Up @@ -260,6 +298,7 @@ void Relax_Driver::json_out(ModuleESolver::ESolver* p_esolver,
}

void Relax_Driver::final_out(const int istep,
const RelaxStopReason stop_reason,
UnitCell& ucell,
const Input_para& inp,
const double etot,
Expand All @@ -268,6 +307,12 @@ void Relax_Driver::final_out(const int istep,
const bool geometry_evaluated,
std::ofstream& ofs_running)
{
int ionic_steps = istep + 1;
if (stop_reason == RelaxStopReason::DryRun)
{
ionic_steps = 0;
}

// Structure final output is effective for scf/nscf/relax/cell-relax;
// relax-specific screen messages remain guarded below.
const bool is_relax = (inp.calculation == "relax" || inp.calculation == "cell-relax");
Expand All @@ -294,37 +339,37 @@ void Relax_Driver::final_out(const int istep,
write_force);
}

// relax_nmax == 0 is a valid dry-run mode: no relaxation step was ever
// taken, so it must not be reported as a failed relaxation.
if (istep == inp.relax_nmax && inp.relax_nmax > 0)
if (stop_reason == RelaxStopReason::MaxSteps)
{
if (is_relax)
{
print_not_converged_summary(inp, ofs_running);
print_not_converged_summary(ionic_steps, inp, ofs_running);

std::cout << "\n ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
std::cout << " Geometry relaxation stops here due to reaching the maximum " << std::endl;
std::cout << " relaxation steps. More steps are needed to converge the results " << std::endl;
std::cout << " ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
}
std::cout << "\n ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
std::cout << " Geometry relaxation stops here due to reaching the maximum " << std::endl;
std::cout << " relaxation steps. More steps are needed to converge the results " << std::endl;
std::cout << " ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
}
else if (inp.relax_nmax > 0)
else if (stop_reason == RelaxStopReason::Converged)
{
if (is_relax)
{
std::cout << "\n ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
std::cout << " Geometry relaxation thresholds are reached within " << istep << " steps." << std::endl;
std::cout << " ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
}
std::cout << "\n ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
std::cout << " Geometry relaxation thresholds are reached within " << ionic_steps << " steps." << std::endl;
std::cout << " ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
}

if (is_relax && inp.relax_nmax == 0)
else if (stop_reason == RelaxStopReason::UserExit)
{
ofs_running << " Geometry relaxation stopped by EXIT after "
<< ionic_steps << " ionic step(s)." << std::endl;
std::cout << " Geometry relaxation stopped by EXIT after "
<< ionic_steps << " ionic step(s)." << std::endl;
}
else if (stop_reason == RelaxStopReason::DryRun)
{
print_dry_run_message();
}
}

void Relax_Driver::print_not_converged_summary(const Input_para& inp, std::ofstream& ofs_running) const
void Relax_Driver::print_not_converged_summary(const int ionic_steps,
const Input_para& inp,
std::ofstream& ofs_running) const
{
// Unified not-converged summary for both relaxation paths so ASE and
// users can read the method, the actual number of steps taken, and
Expand All @@ -333,7 +378,6 @@ void Relax_Driver::print_not_converged_summary(const Input_para& inp, std::ofstr
= inp.uses_simultaneous_relaxation() ? rl.get_max_force_history() : rl_old.get_max_force_history();
const std::vector<double>& stress_hist
= inp.uses_simultaneous_relaxation() ? rl.get_max_stress_history() : rl_old.get_max_stress_history();
const int ionic_steps = static_cast<int>(force_hist.size());
const std::string method = inp.relax_method.empty() ? "unknown" : inp.relax_method;
ofs_running << " Relaxation method: " << method << std::endl;
ofs_running << " Relaxation stopped after " << ionic_steps << " ionic step(s) (relax_nmax = " << inp.relax_nmax
Expand Down
16 changes: 15 additions & 1 deletion source/source_relax/relax_driver.h
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,15 @@ class Relax_Driver
std::ofstream& ofs_running);

private:
enum class RelaxStopReason
{
Converged,
MaxSteps,
UserExit,
DryRun,
Finished
};

/// New relaxation optimizer (Relax class)
Relax rl;
/// Old relaxation optimizer (IonCellOptimizer class)
Expand Down Expand Up @@ -148,6 +157,7 @@ class Relax_Driver
* @brief Output final results after relaxation.
*
* @param istep Final iteration step.
* @param stop_reason Reason why the calculation stopped.
* @param ucell Reference to the unit cell.
* @param inp Input parameters for the calculation.
* @param etot Total energy of the final step.
Expand All @@ -158,6 +168,7 @@ class Relax_Driver
* omitted / marked N/A.
*/
void final_out(const int istep,
const RelaxStopReason stop_reason,
UnitCell& ucell,
const Input_para& inp,
const double etot,
Expand All @@ -173,10 +184,13 @@ class Relax_Driver
* per-step largest force / stress history, so ASE and users can diagnose
* why the relaxation stopped at relax_nmax.
*
* @param ionic_steps Number of evaluated ionic steps.
* @param inp Input parameters (selects the active optimizer and relax_nmax).
* @param ofs_running Output stream for the running log.
*/
void print_not_converged_summary(const Input_para& inp, std::ofstream& ofs_running) const;
void print_not_converged_summary(const int ionic_steps,
const Input_para& inp,
std::ofstream& ofs_running) const;

/**
* @brief Print the dry-run success message to screen.
Expand Down
Loading
Loading