File NumericalHessianCalculator.h

Copyright

This code is licensed under the 3-clause BSD license.

Copyright ETH Zurich, Laboratory for Physical Chemistry, Reiher Group.

See LICENSE.txt for details.

namespace Scine

This header file contains functions that allow for common notation for common things that can be done at a different degree of derivatives.

This header contains alias definitions defining which classes to use for the different degrees of derivatives.

namespace Utils
class NumericalHessianCalculator
#include <NumericalHessianCalculator.h>

This class calculates the Hessian and, optionally, the dipole gradient (semi-)numerically.

Useful if analytic second derivatives are not available.

calculateFromEnergyDifferences uses the following, where D = delta / 2: d^2E/dx^2 = 1 / delta^2 * (E(x-2D) - 2 * E(x) + E(x+2D)) d^2E/dxdy = 1 / delta^2 * (E(x+D,y+D) - E(x+D,y-D) - E(x-D,y+D) + E(x-D,y-D))

calculateFromGradientDifferences uses the following, where D = delta / 2: d^2E/dxdy = 1 / (2*delta) * (g_j(i+D,j) - g_j(i-D,j) + g_i(i,j+D) - g_i(i,j-D))

The second formulation is more stable numerically and is used as default.

In order to calculate the dipole gradient, from each displacement the dipole is calculated. dmu/dx_i = 1 / (2*delta) * (mu_x(x_i + delta) - mu_x(x_i - delta)) Only 3N calculations are needed to fill the 3Nx3 dipole gradient matrix, as the dipole comes as a x,y,z vector at each single point calculation.

Public Functions

NumericalHessianCalculator(Core::Calculator &calculator)
Results calculate(double delta = defaultDelta)

Calcualtes the hessian matrix and, if needed, the dipole gradient.

Return

A Results class with an hessian matrix and, if needed, a dipole gradient.

Parameters

  • delta: The step size for the numerical derivative.

HessianMatrix calculateFromEnergyDifferences(double delta)

This method can ONLY calculate the hessian matrix, not the dipole gradient.

Results calculateFromGradientDifferences(double delta)

This method calculated the hessian matrix from the gradient differences and the dipole gradient as dipole difference for each displacement.

void requiredDipoleGradient(bool dipoleGradient)

Sets whether the dipole gradient is also calculated.

Private Functions

double hessianElementSameFromEnergy(int i, const PositionCollection &referencePositions, double delta)
double hessianElementDifferentFromEnergy(int i, int j, const PositionCollection &referencePositions, double delta)
Eigen::VectorXd addGradientContribution(DipoleGradient &dipoleDiff, int i, const Utils::PositionCollection &referencePositions, double delta, std::shared_ptr<Core::Calculator> calculator, std::shared_ptr<State> state)

Private Members

Core::Calculator &calculator_
bool dipoleGradient_ = {false}

Private Static Attributes

constexpr double defaultDelta = 1e-2