|
n2p2 - A neural network potential package
|
Implementation of the Kalman filter method. More...
#include <KalmanFilter.h>


Public Types | |
| enum | KalmanType { KT_STANDARD , KT_FADINGMEMORY } |
| Enumerate different Kalman filter types. More... | |
Public Member Functions | |
| KalmanFilter (std::size_t const sizeState, KalmanType const type) | |
| Kalman filter class constructor. | |
| virtual | ~KalmanFilter () |
| Destructor. | |
| void | setSizeObservation (std::size_t const sizeObservation) |
| Set observation vector size. | |
| void | setState (double *state) |
| Set pointer to current state. | |
| void | setError (double const *const error) |
| Set pointer to current error vector. | |
| void | setError (double const *const error, std::size_t const size) |
| Set pointer to current error vector. | |
| void | setJacobian (double const *const jacobian) |
| Set pointer to current Jacobi matrix. | |
| void | setJacobian (double const *const jacobian, std::size_t const columns) |
| Set pointer to current Jacobi matrix. | |
| void | update () |
| Update error covariance matrix and state vector. | |
| void | update (std::size_t const sizeObservation) |
| Update error covariance matrix and state vector. | |
| void | setParametersStandard (double const epsilon, double const q0, double const qtau, double const qmin, double const eta0, double const etatau, double const etamax) |
| Set parameters for standard Kalman filter. | |
| void | setParametersFadingMemory (double const epsilon, double const q0, double const qtau, double const qmin, double const lambda, double const nu) |
| Set parameters for fading memory Kalman filter. | |
| std::string | status (std::size_t epoch) const |
| Status report. | |
| std::vector< std::string > | statusHeader () const |
| Header for status report file. | |
| std::vector< std::string > | info () const |
| Information about Kalman filter settings. | |
| KalmanType | getType () const |
| Getter for type. | |
| std::size_t | getSizeObservation () const |
| Getter for sizeObservation. | |
| std::size_t | getNumUpdates () const |
| Getter for numUpdates. | |
| double | getEta () const |
| Getter for eta. | |
| double | getEpsilon () const |
| Getter for epsilon. | |
| double | getQ0 () const |
| Getter for q0. | |
| double | getQtau () const |
| Getter for qtau. | |
| double | getQmin () const |
| Getter for qmin. | |
| double | getLambda () const |
| Getter for lambda. | |
| double | getNu () const |
| Getter for nu. | |
| double | getGamma () const |
| Getter for gamma. | |
Public Member Functions inherited from nnp::Updater | |
| virtual void | setupTiming (std::string const &prefix="upd") |
| Activate detailed timing. | |
| virtual void | resetTimingLoop () |
| Start a new timing loop (e.g. | |
| virtual std::map< std::string, Stopwatch > | getTiming () const |
| Return timings gathered in stopwatch map. | |
Private Attributes | |
| KalmanType | type |
| Kalman filter type. | |
| std::size_t | sizeObservation |
| Size of observation (measurement) vector. | |
| std::size_t | numUpdates |
| Total number of updates performed. | |
| double | epsilon |
| Error covariance initialization parameter \(\epsilon\). | |
| double | q |
| Process noise \(q\). | |
| double | q0 |
| Process noise initial value \(q_0\). | |
| double | qtau |
| Process noise exponential decay parameter \(q_{\tau}\). | |
| double | qmin |
| Process noise minimum value \(q_{\text{min}}\). | |
| double | eta |
| Learning rate \(\eta\). | |
| double | eta0 |
| Learning rate initial value \(\eta_0\). | |
| double | etatau |
| Learning rate exponential increase parameter \(\eta_{\tau}\). | |
| double | etamax |
| Learning rate maximum value \(\eta_{\text{max}}\). | |
| double | lambda |
| Forgetting factor for fading memory Kalman filter. | |
| double | nu |
| Parameter for fading memory Kalman filter. | |
| double | gamma |
| Forgetting gain factor gamma for fading memory Kalman filter. | |
| Eigen::Map< Eigen::VectorXd > * | w |
| State vector. | |
| Eigen::Map< Eigen::VectorXd const > * | xi |
| Error vector. | |
| Eigen::Map< Eigen::MatrixXd const > * | H |
| Derivative matrix. | |
| Eigen::MatrixXd | P |
| Error covariance matrix. | |
| Eigen::MatrixXd | K |
| Kalman gain matrix. | |
| Eigen::MatrixXd | X |
| Intermediate result X = P . H. | |
Additional Inherited Members | |
Protected Member Functions inherited from nnp::Updater | |
| Updater (std::size_t const sizeState) | |
| Constructor. | |
Protected Attributes inherited from nnp::Updater | |
| bool | timing |
| Whether detailed timing is enabled. | |
| bool | timingReset |
| Internal loop timer reset switch. | |
| std::size_t | sizeState |
| Number of neural network connections (weights + biases). | |
| std::string | prefix |
| Prefix for timing stopwatches. | |
| std::map< std::string, Stopwatch > | sw |
| Stopwatch map for timing. | |
Implementation of the Kalman filter method.
Definition at line 31 of file KalmanFilter.h.
Enumerate different Kalman filter types.
| Enumerator | |
|---|---|
| KT_STANDARD | Regular Kalman filter. |
| KT_FADINGMEMORY | Kalman filtering with fading memory modification. |
Definition at line 35 of file KalmanFilter.h.
| KalmanFilter::KalmanFilter | ( | std::size_t const | sizeState, |
| KalmanType const | type ) |
Kalman filter class constructor.
| [in] | sizeState | Size of the state vector w. |
| [in] | type | Kalman filter type (see KalmanType). |
Definition at line 27 of file KalmanFilter.cpp.
References epsilon, eta, eta0, etamax, etatau, gamma, H, K, KT_FADINGMEMORY, KT_STANDARD, lambda, nu, numUpdates, P, q, q0, qmin, qtau, sizeObservation, nnp::Updater::sizeState, type, nnp::Updater::Updater(), w, and xi.

|
virtual |
| void KalmanFilter::setSizeObservation | ( | std::size_t const | sizeObservation | ) |
Set observation vector size.
| [in] | size | Size of the observation vector. |
If the size of the observation vector is known in advance it can be set here.
Definition at line 76 of file KalmanFilter.cpp.
References sizeObservation.
Referenced by nnp::Training::update().

|
virtual |
Set pointer to current state.
| [in,out] | state | Pointer to state vector (weights vector), will be changed in-place upon calling update(). |
Implements nnp::Updater.
Definition at line 83 of file KalmanFilter.cpp.
References nnp::Updater::sizeState, and w.
| void KalmanFilter::setError | ( | double const *const | error | ) |
Set pointer to current error vector.
| [in] | error | Pointer to error (difference between reference and neural network potential output). |
Definition at line 90 of file KalmanFilter.cpp.
References setError(), and sizeObservation.
Referenced by setError().


|
virtual |
Set pointer to current error vector.
| [in] | error | Pointer to error (difference between reference and neural network potential output). |
| [in] | size | Number of error vector entries. |
Implements nnp::Updater.
| void KalmanFilter::setJacobian | ( | double const *const | jacobian | ) |
Set pointer to current Jacobi matrix.
| [in] | jacobian | Derivatives of error with respect to weights. |
Definition at line 104 of file KalmanFilter.cpp.
References setJacobian(), and sizeObservation.
Referenced by setJacobian().


|
virtual |
Set pointer to current Jacobi matrix.
| [in] | jacobian | Derivatives of error with respect to weights. |
| [in] | columns | Number of gradients provided. |
Implements nnp::Updater.
|
virtual |
Update error covariance matrix and state vector.
Implements nnp::Updater.
Definition at line 119 of file KalmanFilter.cpp.
References sizeObservation, and update().
Referenced by update().


| void nnp::KalmanFilter::update | ( | std::size_t const | sizeObservation | ) |
| void KalmanFilter::setParametersStandard | ( | double const | epsilon, |
| double const | q0, | ||
| double const | qtau, | ||
| double const | qmin, | ||
| double const | eta0, | ||
| double const | etatau, | ||
| double const | etamax ) |
Set parameters for standard Kalman filter.
| [in] | epsilon | Error covariance initialization parameter \(\epsilon\). |
| [in] | q0 | Process noise initial value \(q_0\). |
| [in] | qtau | Process noise exponential decay parameter \(q_{\tau}\). |
| [in] | qmin | Process noise minimum value \(q_{\text{min}}\). |
| [in] | eta0 | Initial learning rate \(\eta_0\). |
| [in] | etatau | Learning rate exponential increase parameter \(\eta_{\tau}\). |
| [in] | etamax | Learning rate maximum value \(\eta_{\text{max}}\). |
Definition at line 199 of file KalmanFilter.cpp.
References epsilon, eta, eta0, etamax, etatau, P, q, q0, qmin, and qtau.
Referenced by nnp::Training::setupTraining().

| void KalmanFilter::setParametersFadingMemory | ( | double const | epsilon, |
| double const | q0, | ||
| double const | qtau, | ||
| double const | qmin, | ||
| double const | lambda, | ||
| double const | nu ) |
Set parameters for fading memory Kalman filter.
| [in] | epsilon | Error covariance initialization parameter \(\epsilon\). |
| [in] | q0 | Process noise initial value \(q_0\). |
| [in] | qtau | Process noise exponential decay parameter \(q_{\tau}\). |
| [in] | qmin | Process noise minimum value \(q_{\text{min}}\). |
| [in] | lambda | Forgetting factor \(\lambda\). |
| [in] | nu | Fading memory parameter \(\nu\). |
Definition at line 222 of file KalmanFilter.cpp.
References epsilon, gamma, lambda, nu, P, q, q0, qmin, and qtau.
Referenced by nnp::Training::setupTraining().

|
virtual |
Status report.
| [in] | epoch | Current epoch. |
Implements nnp::Updater.
Definition at line 243 of file KalmanFilter.cpp.
References eta, gamma, K, KT_FADINGMEMORY, KT_STANDARD, lambda, numUpdates, P, q, nnp::Updater::sizeState, nnp::strpr(), and type.

|
virtual |
Header for status report file.
Implements nnp::Updater.
Definition at line 268 of file KalmanFilter.cpp.
References nnp::createFileHeader(), KT_FADINGMEMORY, KT_STANDARD, and type.

|
virtual |
Information about Kalman filter settings.
Implements nnp::Updater.
Definition at line 321 of file KalmanFilter.cpp.
References epsilon, eta0, etamax, etatau, KT_FADINGMEMORY, KT_STANDARD, lambda, nu, q0, qmin, qtau, sizeObservation, nnp::Updater::sizeState, nnp::strpr(), and type.

|
inline |
|
inline |
Getter for sizeObservation.
Definition at line 248 of file KalmanFilter.h.
References sizeObservation.
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
| double nnp::KalmanFilter::getGamma | ( | ) | const |
Getter for gamma.
|
private |
Kalman filter type.
Definition at line 196 of file KalmanFilter.h.
Referenced by getType(), info(), KalmanFilter(), status(), and statusHeader().
|
private |
Size of observation (measurement) vector.
Definition at line 198 of file KalmanFilter.h.
Referenced by getSizeObservation(), info(), KalmanFilter(), setError(), setJacobian(), setSizeObservation(), update(), and update().
|
private |
Total number of updates performed.
Definition at line 200 of file KalmanFilter.h.
Referenced by getNumUpdates(), KalmanFilter(), and status().
|
private |
Error covariance initialization parameter \(\epsilon\).
Definition at line 202 of file KalmanFilter.h.
Referenced by getEpsilon(), info(), KalmanFilter(), setParametersFadingMemory(), setParametersStandard(), and update().
|
private |
Process noise \(q\).
Definition at line 204 of file KalmanFilter.h.
Referenced by KalmanFilter(), setParametersFadingMemory(), setParametersStandard(), and status().
|
private |
Process noise initial value \(q_0\).
Definition at line 206 of file KalmanFilter.h.
Referenced by getQ0(), info(), KalmanFilter(), setParametersFadingMemory(), setParametersStandard(), and update().
|
private |
Process noise exponential decay parameter \(q_{\tau}\).
Definition at line 208 of file KalmanFilter.h.
Referenced by getQtau(), info(), KalmanFilter(), setParametersFadingMemory(), setParametersStandard(), and update().
|
private |
Process noise minimum value \(q_{\text{min}}\).
Definition at line 210 of file KalmanFilter.h.
Referenced by getQmin(), info(), KalmanFilter(), setParametersFadingMemory(), setParametersStandard(), and update().
|
private |
Learning rate \(\eta\).
Definition at line 212 of file KalmanFilter.h.
Referenced by getEta(), KalmanFilter(), setParametersStandard(), and status().
|
private |
Learning rate initial value \(\eta_0\).
Definition at line 214 of file KalmanFilter.h.
Referenced by info(), KalmanFilter(), setParametersStandard(), and update().
|
private |
Learning rate exponential increase parameter \(\eta_{\tau}\).
Definition at line 216 of file KalmanFilter.h.
Referenced by info(), KalmanFilter(), setParametersStandard(), and update().
|
private |
Learning rate maximum value \(\eta_{\text{max}}\).
Definition at line 218 of file KalmanFilter.h.
Referenced by info(), KalmanFilter(), setParametersStandard(), and update().
|
private |
Forgetting factor for fading memory Kalman filter.
Definition at line 220 of file KalmanFilter.h.
Referenced by getLambda(), info(), KalmanFilter(), setParametersFadingMemory(), status(), and update().
|
private |
Parameter for fading memory Kalman filter.
Definition at line 222 of file KalmanFilter.h.
Referenced by getNu(), info(), KalmanFilter(), setParametersFadingMemory(), and update().
|
private |
Forgetting gain factor gamma for fading memory Kalman filter.
Definition at line 224 of file KalmanFilter.h.
Referenced by KalmanFilter(), setParametersFadingMemory(), and status().
|
private |
State vector.
Definition at line 226 of file KalmanFilter.h.
Referenced by KalmanFilter(), and setState().
|
private |
|
private |
|
private |
Error covariance matrix.
Definition at line 232 of file KalmanFilter.h.
Referenced by KalmanFilter(), setParametersFadingMemory(), setParametersStandard(), and status().
|
private |
Kalman gain matrix.
Definition at line 234 of file KalmanFilter.h.
Referenced by KalmanFilter(), and status().
|
private |
Intermediate result X = P . H.
Definition at line 236 of file KalmanFilter.h.