Loading...
Searching...
No Matches
CollectiveIDPPObjectiveFunction Class Reference

#include <IDPPObjectiveFunction.hpp>

Inheritance diagram for CollectiveIDPPObjectiveFunction:

Public Member Functions

 CollectiveIDPPObjectiveFunction (std::vector< Matter > &pathRef, const Parameters &paramsPassed)
double getEnergy () override
VectorXd getGradient (bool fdstep=false) override
void setPositions (const VectorXd &x) override
VectorXd getPositions () override
int degreesOfFreedom () override
bool isConverged () override
double getConvergence () override
VectorXd difference (const VectorXd &a, const VectorXd &b) override
Public Member Functions inherited from eonc::ObjectiveFunction
 ObjectiveFunction (const Parameters &paramsPassed)
virtual ~ObjectiveFunction ()

Private Member Functions

MatrixXd getDistanceMatrix (const Matter &m)
MatrixXd getIDPPForces (const Matter &m, const MatrixXd &dTarget)

Private Attributes

std::vector< Matter > & path
MatrixXd dInit
MatrixXd dFinal
double lastMaxForce = 100.0

Additional Inherited Members

Protected Attributes inherited from eonc::ObjectiveFunction
const Parametersparams

Detailed Description

Definition at line 79 of file IDPPObjectiveFunction.hpp.

Constructor & Destructor Documentation

◆ CollectiveIDPPObjectiveFunction()

eonc::CollectiveIDPPObjectiveFunction::CollectiveIDPPObjectiveFunction ( std::vector< Matter > & pathRef,
const Parameters & paramsPassed )
inline

Definition at line 81 of file IDPPObjectiveFunction.hpp.

83 : ObjectiveFunction(paramsPassed),
84 path(pathRef) {
85
86 // Initialize distances for endpoints
87 dInit = getDistanceMatrix(path.front());
89 }
ObjectiveFunction(const Parameters &paramsPassed)

Member Function Documentation

◆ degreesOfFreedom()

Implements eonc::ObjectiveFunction.

Definition at line 119 of file IDPPObjectiveFunction.hpp.

119 {
120 return 3 * path[0].numberOfAtoms() * (path.size() - 2);
121 }

◆ difference()

VectorXd eonc::CollectiveIDPPObjectiveFunction::difference ( const VectorXd & a,
const VectorXd & b )
inlineoverridevirtual

Implements eonc::ObjectiveFunction.

Definition at line 130 of file IDPPObjectiveFunction.hpp.

130 {
131 // Simple difference for this purpose, assuming pre-aligned or handling PBC
132 // inside
133 return a - b;
134 }

◆ getConvergence()

◆ getDistanceMatrix()

Definition at line 84 of file IDPPObjectiveFunction.cpp.

84 {
85 int natoms = m.numberOfAtoms();
86 MatrixXd d(natoms, natoms);
87 auto pos = m.getPositions();
88 for (int i = 0; i < natoms; ++i) {
89 for (int j = 0; j < natoms; ++j) {
90 d(i, j) = m.pbc(pos.row(i) - pos.row(j)).norm();
91 }
92 }
93 return d;
94}
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic, eOnStorageOrder > MatrixXd
Definition Eigen.h:33
AtomMatrix pbc(const AtomMatrix &diff) const
Definition Matter.h:163
long int numberOfAtoms() const
Definition Matter.cpp:209
const AtomMatrix & getPositions() const
Definition Matter.cpp:236

◆ getEnergy()

double eonc::CollectiveIDPPObjectiveFunction::getEnergy ( )
inlineoverridevirtual

Implements eonc::ObjectiveFunction.

Definition at line 93 of file IDPPObjectiveFunction.hpp.

93{ return 0.0; }

◆ getGradient()

VectorXd CollectiveIDPPObjectiveFunction::getGradient ( bool fdstep = false)
overridevirtual

Implements eonc::ObjectiveFunction.

Definition at line 124 of file IDPPObjectiveFunction.cpp.

124 {
125 int nImgs = path.size() - 2; // Exclude fixed endpoints
126 int natoms = path[0].numberOfAtoms();
127 VectorXd totalGradient(3 * natoms * nImgs);
128 double maxForce = 0.0;
129
130 // 1. Compute Raw IDPP Forces and Tangents
131 std::vector<AtomMatrix> rawForces(path.size());
132 std::vector<AtomMatrix> tangents(path.size());
133
134 // We compute for 1..N (moving images)
135 for (size_t i = 1; i <= nImgs; ++i) {
136 // Interpolate Target
137 double xi = static_cast<double>(i) / (nImgs + 1);
138 MatrixXd dTarget = (1.0 - xi) * dInit + xi * dFinal;
139
140 rawForces[i] = getIDPPForces(path[i], dTarget);
141
142 // Simple Tangent: Next - Prev
143 AtomMatrix nextPos = path[i + 1].getPositions();
144 AtomMatrix prevPos = path[i - 1].getPositions();
145 tangents[i] = path[i].pbc(nextPos - prevPos);
146 tangents[i].normalize(); // Unit tangent
147 }
148
149 // 2. Project Forces and Add Springs (The "NEB" part of IDPP-NEB)
150 double k = params.neb_options.spring.constant;
151
152 for (size_t i = 1; i <= nImgs; ++i) {
153 AtomMatrix f = rawForces[i];
154 AtomMatrix t = tangents[i];
155
156 // Perpendicular Force (IDPP optimization)
157 double f_dot_t = matDot(f, t);
158 AtomMatrix f_perp = f - f_dot_t * t;
159
160 // Spring Force (Spacing optimization)
161 double distNext = path[i].distanceTo(path[i + 1]);
162 double distPrev = path[i].distanceTo(path[i - 1]);
163 AtomMatrix f_spring = k * (distNext - distPrev) * t;
164
165 // Total NEB Force
166 AtomMatrix f_neb = f_perp + f_spring;
167
168 // Store as Gradient (-Force)
169 totalGradient.segment(3 * natoms * (i - 1), 3 * natoms) =
170 VectorXd::Map(f_neb.data(), 3 * natoms) * -1.0;
171
172 // Tracking convergence
173 maxForce = std::max(maxForce, f_neb.template lpNorm<Eigen::Infinity>());
174 }
175
176 lastMaxForce = maxForce;
177 return totalGradient;
178}
double matDot(const AtomMatrix &a, const AtomMatrix &b)
SIMD-optimized dot product for contiguous Eigen matrices.
Definition Eigen.h:50
Eigen::Matrix< double, Eigen::Dynamic, 3, eOnStorageOrder > AtomMatrix
Definition Eigen.h:37
MatrixXd getIDPPForces(const Matter &m, const MatrixXd &dTarget)
const Parameters & params

◆ getIDPPForces()

Definition at line 97 of file IDPPObjectiveFunction.cpp.

98 {
99 int natoms = m.numberOfAtoms();
100 AtomMatrix forces = AtomMatrix::Zero(natoms, 3);
101 auto pos = m.getPositions();
102
103 for (int i = 0; i < natoms; ++i) {
104 for (int j = i + 1; j < natoms; ++j) {
105 Eigen::RowVector3d dr_vec = m.pbc(pos.row(i) - pos.row(j));
106 double r = dr_vec.norm();
107 if (r < 1e-4)
108 r = 1e-4;
109
110 double diff = r - dTarget(i, j);
111 // SOTA Weighting: w = 1/r^4. Gradient logic matches Smidstrup/ASE
112 double r2 = r * r;
113 double r5 = r2 * r2 * r;
114 double dEdr = (diff / r5) * (2.0 * dTarget(i, j) - r);
115
116 Eigen::RowVector3d f = -dEdr * (dr_vec / r); // -dE/dr * r_hat
117 forces.row(i) += f;
118 forces.row(j) -= f;
119 }
120 }
121 return forces;
122}

◆ getPositions()

VectorXd eonc::CollectiveIDPPObjectiveFunction::getPositions ( )
inlineoverridevirtual

Implements eonc::ObjectiveFunction.

Definition at line 107 of file IDPPObjectiveFunction.hpp.

107 {
108 int atoms = path[0].numberOfAtoms();
109 int n_free_images = path.size() - 2;
110 VectorXd pos(3 * atoms * n_free_images);
111
112 for (size_t i = 1; i < path.size() - 1; ++i) {
113 pos.segment(3 * atoms * (i - 1), 3 * atoms) =
114 VectorXd::Map(path[i].getPositions().data(), 3 * atoms);
115 }
116 return pos;
117 }

◆ isConverged()

Implements eonc::ObjectiveFunction.

Definition at line 124 of file IDPPObjectiveFunction.hpp.

124 {
125 return getConvergence() < params.neb_options.initialization.force_tolerance;
126 }

◆ setPositions()

void eonc::CollectiveIDPPObjectiveFunction::setPositions ( const VectorXd & x)
inlineoverridevirtual

Implements eonc::ObjectiveFunction.

Definition at line 98 of file IDPPObjectiveFunction.hpp.

98 {
99 int atoms = path[0].numberOfAtoms();
100 // Skip endpoints (0 and N+1)
101 for (size_t i = 1; i < path.size() - 1; ++i) {
102 path[i].setPositions(AtomMatrix::Map(
103 x.segment(3 * atoms * (i - 1), 3 * atoms).data(), atoms, 3));
104 }
105 }

Member Data Documentation

◆ dFinal

◆ dInit

◆ lastMaxForce

Definition at line 139 of file IDPPObjectiveFunction.hpp.

◆ path

Definition at line 137 of file IDPPObjectiveFunction.hpp.


The documentation for this class was generated from the following files: