Source Edit

Differentiation

This module implements various finite difference methods for approximating derivatives.

For multi-valued functions, it is expected that both the inputs and outputs are 1D Tensors.

Here is an example with a scalar multi-variate function:

Example:

import src/numericalnim/differentiate
import arraymancer, std/math
# f(x, y) = sin(x) + 2*sin(y)
proc f(x: Tensor[float]): float =
  sin(x[0]) + 2*sin(x[1])

let x = [1.0, 1.0].toTensor
let gradient = tensorGradient(f, x)

Here's an example with a multi-valued function:

Example:

import src/numericalnim/differentiate
import arraymancer, std/math
# f(x, y) = [sin(x*y), cos(x*y)]
proc f(x: Tensor[float]): Tensor[float] =
  result = newTensor[float](2)
  let arg = x[0] * x[1]
  result[0] = sin(arg)
  result[1] = cos(arg)

let x = [1.0, 1.0].toTensor
let jacobian = tensorJacobian(f, x)

Procs

proc checkGradient[U; T: not Tensor](f: proc (x: Tensor[U]): T; fGrad: proc (x: Tensor[U]): Tensor[T]; x0: Tensor[U]; tol: T): bool: Checks if the provided gradient function fGrad gives the same values as numeric gradient. Source Edit
proc checkGradient[U; T](f: proc (x: Tensor[U]): Tensor[T]; fGrad: proc (x: Tensor[U]): Tensor[T]; x0: Tensor[U]; tol: T): bool: Checks if the provided gradient function fGrad gives the same values as numeric gradient. Source Edit
proc diff1dBackward[U, T](f: proc (x: U): T; x0: U; h: U = U(0.000001)): T: Numerically calculate the derivative of f(x) at x0 using a step size h. Uses backward difference which has accuracy O(h) Source Edit
proc diff1dCentral[U, T](f: proc (x: U): T; x0: U; h: U = U(0.000001)): T: Numerically calculate the derivative of f(x) at x0 using a step size h. Uses central difference which has accuracy O(h^2) Source Edit
proc diff1dForward[U, T](f: proc (x: U): T; x0: U; h: U = U(0.000001)): T: Numerically calculate the derivative of f(x) at x0 using a step size h. Uses forward difference which has accuracy O(h) Source Edit
proc mixedDerivative[U, T](f: proc (x: Tensor[U]): T; x0: var Tensor[U]; indices: (int, int); h: U = U(0.000001)): T: Used internally in tensorHessian. Calculates the mixed derivative d²f/(dx_i dx_j). Modifies x0. Source Edit
proc secondDiff1dBackward[U, T](f: proc (x: U): T; x0: U; h: U = U(0.000001)): T: Numerically calculate the second derivative of f(x) at x0 using a step size h. Source Edit
proc secondDiff1dCentral[U, T](f: proc (x: U): T; x0: U; h: U = U(0.000001)): T: Numerically calculate the second derivative of f(x) at x0 using a step size h. Uses central difference which has accuracy O(h^2) Source Edit
proc secondDiff1dForward[U, T](f: proc (x: U): T; x0: U; h: U = U(0.000001)): T: Numerically calculate the second derivative of f(x) at x0 using a step size h. Source Edit
proc tensorGradient[U, T](f: proc (x: Tensor[U]): Tensor[T]; x0: Tensor[U]; h: U = U(0.000001); fastMode: bool = false): Tensor[T]: Calculates the gradient of multi-valued f(x) w.r.t vector x at x0 using step size h. f(x) is expected to take as input and return a 1D Tensor. Every column is the gradient of one component of f. By default it uses central difference for approximating the derivatives. This requires two function evaluations per derivative. When fastMode is true it will instead use the forward difference which only uses 1 function evaluation per derivative but is less accurate.

Returns:

The gradient as a 2D Tensor with shape (nInputs, nOutputs). This is the transpose of the Jacobian.

Source Edit
proc tensorGradient[U; T: not Tensor](f: proc (x: Tensor[U]): T; x0: Tensor[U]; h: U = U(0.000001); fastMode: bool = false): Tensor[ T]: Calculates the gradient of scalar f(x) w.r.t vector x at x0 using step size h. By default it uses central difference for approximating the derivatives. This requires two function evaluations per derivative. When fastMode is true it will instead use the forward difference which only uses 1 function evaluation per derivative but is less accurate.

Returns:

The gradient as a 1D Tensor.

Source Edit
proc tensorHessian[U; T: not Tensor](f: proc (x: Tensor[U]): T; x0: Tensor[U]; h: U = U(0.000001)): Tensor[T]: Calculates the Hessian of a scalar function f(x) using finite differences at x0. f(x) should accept a 1D Tensor of input values.

Returns:

The Hessian as a 2D Tensor of shape (nInputs, nInputs).

Source Edit
proc tensorJacobian[U, T](f: proc (x: Tensor[U]): Tensor[T]; x0: Tensor[U]; h: U = U(0.000001); fastMode: bool = false): Tensor[T]: Calculates the jacobian of multi-valued f(x) w.r.t vector x at x0 using step size h. f(x) is expected to take as input and return a 1D Tensor. Every row is the gradient of one component of f. By default it uses central difference for approximating the derivatives. This requires two function evaluations per derivative. When fastMode is true it will instead use the forward difference which only uses 1 function evaluation per derivative but is less accurate. Returns:
The Jacobian as a 2D Tensor with shape (nOutputs, nInputs).

Source Edit