22 {
23
24 saddle->relax(
false,
true,
false,
"displacementmin");
27
28
29 double eproduct, ereactant, de;
30 eproduct =
product->getPotentialEnergy();
31 ereactant =
reactant->getPotentialEnergy();
32 de = eproduct - ereactant;
33 double kB =
params.constants.kB;
34 double Temperature =
params.main_options.temperature;
35 double arg = -de / (kB * Temperature);
36 double p = std::exp(arg);
38 if (ereactant < eproduct) {
39 if (r > p) {
40 return 1;
41 }
42 }
43
46 neb.path[0]->matter2con("neb_initial_band.con", false))) {
47 QUILL_LOG_WARNING(
log,
"Failed to write neb_initial_band.con");
48 }
49 for (int j = 1; j < neb.numImages; j++) {
50 if (!
eonc::io::io_ok(neb.path[j]->matter2con(
"neb_initial_band",
true))) {
51 QUILL_LOG_WARNING(
log,
"Failed to append neb_initial_band frame");
52 }
53 }
54 neb.compute();
55
56 double Emax = -1e100;
57 int HighestImage = 0;
58
59 for (int i = 1; i < neb.numImages; i++) {
60 double Etest = neb.path[i]->getPotentialEnergy();
61 QUILL_LOG_DEBUG(
log,
"i: {} Etest: {:.1f}", i, Etest);
62 if (Etest > Emax) {
63 Emax = Etest;
64 HighestImage = i;
65 }
66 }
67
68
69 AtomMatrix r_1 = neb.path[HighestImage - 1]->getPositions();
70 AtomMatrix r_2 = neb.path[HighestImage]->getPositions();
71 AtomMatrix r_3 = neb.path[HighestImage + 1]->getPositions();
73 MinModeSaddleSearch dim(neb.path[HighestImage], direction.normalized(),
75 dim.run();
76 *
saddle = *neb.path[HighestImage];
79 return 0;
80}
Eigen::Matrix< double, Eigen::Dynamic, 3, eOnStorageOrder > AtomMatrix
std::shared_ptr< Matter > product
std::shared_ptr< Potential > pot
const Parameters & params
double random(long newSeed=0)
constexpr bool io_ok(IoStatus s) noexcept