SaddleSearchJob.cpp

Go to the documentation of this file.

/*
** This file is part of eOn.
**
** SPDX-License-Identifier: BSD-3-Clause
**
** Copyright (c) 2010--present, eOn Development Team
** All rights reserved.
**
** Repo:
** https://github.com/TheochemUI/eOn
*/
#include "eon/SaddleSearchJob.h"
#ifdef WITH_ARTN
#include "eon/ARTnSaddleSearch.h"
#endif
#include "eon/EonLogger.h"
#include "eon/EpiCenters.h"
#include "eon/HelperFunctions.h"
#include "eon/Potential.h"
 
#include <filesystem>
#include <format>
#include <fstream>
#include <stdexcept>
#include <string>
 
namespace eonc {
 
std::vector<std::string> SaddleSearchJob::run() {
  std::string reactantFilename = eonc::helpers::getRelevantFile("pos.con");
  std::string displacementFilename("displacement.con");
  std::string modeFilename("direction.dat");
 
  if (params.main_options().checkpoint) {
    if (std::filesystem::exists("displacement_cp.con") &&
        std::filesystem::exists("mode_cp.dat")) {
      displacementFilename = "displacement_cp.con";
      modeFilename = "mode_cp.dat";
      QUILL_LOG_DEBUG(log, "Resuming from checkpoint");
    } else {
      QUILL_LOG_DEBUG(log, "No checkpoint files found");
    }
  }
 
  initial = std::make_shared<Matter>(pot, params);
  displacement = std::make_shared<Matter>(pot, params);
  saddle = std::make_shared<Matter>(pot, params);
 
  if (!eonc::io::io_ok(initial->con2matter(reactantFilename))) {
    EONC_LOG_CRITICAL("Failed to load {}", reactantFilename);
    throw std::runtime_error("failed to load " + reactantFilename);
  }
 
  const bool standaloneARTn = params.saddle_search_options().method == "artn";
 
  AtomMatrix mode = AtomMatrix::Zero(initial->numberOfAtoms(), 3);
  if (!standaloneARTn && params.saddle_search_options().displace_type ==
                             eonc::EpiCenters::DISP_LOAD) {
    // Load displacement.con, or synthesize from pos.con + direction.dat (#79).
    if (!eonc::helpers::loadOrSynthesizeDisplacement(
            *saddle, *initial, displacementFilename, modeFilename,
            params.saddle_search_options().displace_magnitude)) {
      EONC_LOG_CRITICAL("Failed to load {} (and no usable {})",
                        displacementFilename, modeFilename);
      throw std::runtime_error("missing displacement.con and direction.dat");
    }
    if (std::filesystem::exists(modeFilename)) {
      mode = eonc::helpers::loadMode(modeFilename, initial->numberOfAtoms());
    }
  } else if (!standaloneARTn && eonc::helpers::applyClientDisplacement(
                                    *saddle, *initial, params, &mode)) {
    // listed_atoms / random / last_atom / least_coordinated / not_fcc_hcp
  } else {
    *saddle = *initial;
  }
 
  if (standaloneARTn && std::filesystem::exists(modeFilename)) {
    mode = eonc::helpers::loadMode(modeFilename, initial->numberOfAtoms());
  }
 
  (void)runPrepared(mode);
  return returnFiles;
}
 
std::shared_ptr<Matter>
SaddleSearchJob::runFromMatter(std::shared_ptr<Matter> seed) {
  if (!seed) {
    throw std::runtime_error("SaddleSearchJob::runFromMatter: null Matter");
  }
  initial = seed;
  initial->setPotential(pot);
  displacement = std::make_shared<Matter>(pot, params);
  saddle = std::make_shared<Matter>(pot, params);
  AtomMatrix mode = AtomMatrix::Zero(initial->numberOfAtoms(), 3);
  if (!eonc::helpers::applyClientDisplacement(*saddle, *initial, params,
                                              &mode)) {
    *saddle = *initial;
  }
  *displacement = *saddle;
  return runPrepared(mode);
}
 
std::shared_ptr<Matter> SaddleSearchJob::runPrepared(const AtomMatrix &mode) {
  const bool useStandaloneARTn =
      params.saddle_search_options().method == "artn";
  const bool useARTnAsMinMode =
      params.saddle_search_options().method == "min_mode" &&
      params.saddle_search_options().minmode_method == "artn";
 
#ifdef WITH_ARTN
  if (useStandaloneARTn || useARTnAsMinMode) {
    saddleSearch =
        std::make_unique<ARTnSaddleSearch>(saddle, pot, mode, params);
  } else
#endif
  {
#ifndef WITH_ARTN
    if (useStandaloneARTn) {
      throw std::runtime_error(
          "saddle_search.method=artn requires a build with ARTn support");
    }
    if (useARTnAsMinMode) {
      throw std::runtime_error(
          "saddle_search.minmode_method=artn requires a build with ARTn "
          "support");
    }
#endif
    saddleSearch = std::make_unique<MinModeSaddleSearch>(
        saddle, mode, initial->getPotentialEnergy(), params, pot);
  }
 
  int status = doSaddleSearch();
  printEndState(status);
  saveData(status);
 
  return saddle;
}
 
int SaddleSearchJob::doSaddleSearch() {
  long status;
  int f1{0};
  f1 = this->pot->forceCallCounter;
  try {
    status = saddleSearch->run();
  } catch (int e) {
    if (e == 100) {
      status = MinModeSaddleSearch::STATUS_POTENTIAL_FAILED;
    } else {
      QUILL_LOG_ERROR(log, "unknown exception: {}", e);
      throw e;
    }
  } catch (const std::exception &e) {
    QUILL_LOG_ERROR(log, "Saddle search potential failed: {}", e.what());
    status = MinModeSaddleSearch::STATUS_POTENTIAL_FAILED;
  }
 
  if (params.saddle_search_options().method == "min_mode" &&
      params.saddle_search_options().minmode_method ==
          LowestEigenmode::MINMODE_GPRDIMER) {
    fCallsSaddle = saddleSearch->getForceCalls();
  } else if (params.saddle_search_options().method == "artn") {
    fCallsSaddle = saddleSearch->getForceCalls();
  } else {
    fCallsSaddle += this->pot->forceCallCounter - f1;
  }
 
  return status;
}
 
void SaddleSearchJob::saveData(int status) {
  std::string resultsFilename("results.dat");
  returnFiles.push_back(resultsFilename);
 
  std::ofstream out(resultsFilename, std::ios::binary);
  if (out) {
    out << std::format("{} termination_reason\n", status);
    out << std::format("{} termination_reason_text\n",
                       saddleSearch->describeStatus(status));
    out << "saddle_search job_type\n";
    out << std::format("{} random_seed\n", params.main_options().randomSeed);
    out << std::format(
        "{} potential_type\n",
        magic_enum::enum_name<PotType>(params.potential_options().potential));
    out << std::format("{} total_force_calls\n",
                       this->pot->forceCallCounter.load());
    out << std::format("{} force_calls_saddle\n", fCallsSaddle);
    out << std::format("{} iterations\n", saddleSearch->getIterationCount());
    if (status != MinModeSaddleSearch::STATUS_POTENTIAL_FAILED) {
      out << std::format("{:f} potential_energy_saddle\n",
                         saddle->getPotentialEnergy());
      out << std::format("{:f} final_eigenvalue\n",
                         saddleSearch->getEigenvalue());
    }
    out << std::format("{:f} potential_energy_reactant\n",
                       initial->getPotentialEnergy());
  }
 
  std::string modeFilename("mode.dat");
  returnFiles.push_back(modeFilename);
  eonc::helpers::saveMode(modeFilename, saddle, saddleSearch->getEigenvector());
 
  std::string saddleFilename("saddle.con");
  returnFiles.push_back(saddleFilename);
  if (!eonc::io::io_ok(saddle->matter2con(saddleFilename))) {
    QUILL_LOG_ERROR(log, "Failed to write {}", saddleFilename);
  }
}
 
void SaddleSearchJob::printEndState(int status) {
  auto msg = saddleSearch->describeStatus(status);
  if (status == MinModeSaddleSearch::STATUS_GOOD) {
    QUILL_LOG_DEBUG(log, "[Saddle Search] {}", msg);
  } else {
    QUILL_LOG_WARNING(log, "[Saddle Search] {}", msg);
  }
}
 
} // namespace eonc