Portability GHC only? experimental bjorn.buckwalter@gmail.com

Description

Forward Automatic Differentiation via overloading to perform nonstandard interpretation that replaces original numeric type with corresponding generalized dual number type.

Credits:

Authors: Copyright 2008, Barak A. Pearlmutter (barak@cs.nuim.ie) & Jeffrey Mark Siskind (qobi@purdue.edu)

Work started as stripped-down version of higher-order tower code published by Jerzy Karczmarczuk (jerzy.karczmarczuk@info.unicaen.fr) which used a non-standard standard prelude.

Initial perturbation-confusing code is a modified version of http://cdsmith.wordpress.com/2007/11/29/some-playing-with-derivatives/

Notes:

Each invocation of the differentiation function introduces a distinct perturbation, which requires a distinct dual number type. In order to prevent these from being confused, tagging, called branding in the Haskell community, is used. This seems to prevent perturbation confusion, although it would be nice to have an actual proof of this. The technique does require adding invocations of lift at appropriate places when nesting is present.

Synopsis

# Higher-Order Dual Numbers

data Dual tag a Source

The `Dual` type is a concrete representation of a higher-order Dual number, meaning a number augmented with a tower of derivatives. These generalize the Dual numbers of Clifford (1873), which hold only a first derivative. They can be converted to formal power series via division by the sequence of factorials.

Instances

 (Enum a, Num a) => Enum (Dual tag a) (Eq a, Num a) => Eq (Dual tag a) Floating a => Floating (Dual tag a) Fractional a => Fractional (Dual tag a) Num a => Num (Dual tag a) (Ord a, Num a) => Ord (Dual tag a) Real a => Real (Dual tag a) (RealFloat a, RealFrac a) => RealFloat (Dual tag a) RealFrac a => RealFrac (Dual tag a) Show a => Show (Dual tag a)

lift :: Num a => a -> Dual tag aSource

The `lift` function injects a primal number into the domain of dual numbers, with a zero tower. If dual numbers were a monad, `lift` would be `return`.

# First-Order Differentiation Operators

These have two-letter suffices for the arity of the input and output of the passed functions: U for univariate, meaning a number, M for multivariate, meaning a list of numbers.

When the input is multivariate a directional derivative is calculated; this requires an additional "direction" parameter. The multivariate case is treated as a list (on input) and as a functor of arbitrary shape, which includes lists as a special case, on output.

Naming convention:

`diff{U/M}{U/F}`
Derivative-taking operators that return a primal/first-derivative pair, for all combinations of scalar/nonscalar input & output
`diff2{U/M}{U/F}`
Derivative-taking operators that calculate a (primal, first-derivative) pair, for all combinations of scalar/nonscalar input & output

diffUU :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> bSource

The `diffUU` function calculates the first derivative of a scalar-to-scalar function.

diffUF :: (Num a, Num b, Functor f) => (forall tag. Dual tag a -> f (Dual tag b)) -> a -> f bSource

The `diffUF` function calculates the first derivative of scalar-to-nonscalar function.

diffMU :: (Num a, Num b) => (forall tag. [Dual tag a] -> Dual tag b) -> [a] -> [a] -> bSource

The `diffMU` function calculate the product of the Jacobian of a nonscalar-to-scalar function with a given vector. Aka: directional derivative.

diffMF :: (Num a, Num b, Functor f) => (forall tag. [Dual tag a] -> f (Dual tag b)) -> [a] -> [a] -> f bSource

The `diffMF` function calculates the product of the Jacobian of a nonscalar-to-nonscalar function with a given vector. Aka: directional derivative.

diff2UU :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> (b, b)Source

The `diff2UU` function calculates the value and derivative, as a pair, of a scalar-to-scalar function.

diff2UF :: (Num a, Num b, Functor f) => (forall tag. Dual tag a -> f (Dual tag b)) -> a -> (f b, f b)Source

The `diffUF2` function calculates the value and derivative, as a pair, of a scalar-to-nonscalar function.

diff2MU :: (Num a, Num b) => (forall tag. [Dual tag a] -> Dual tag b) -> [a] -> [a] -> (b, b)Source

The `diffMU2` function calculates the value and directional derivative, as a pair, of a nonscalar-to-scalar function.

diff2MF :: (Num a, Num b, Functor f) => (forall tag. [Dual tag a] -> f (Dual tag b)) -> [a] -> [a] -> (f b, f b)Source

The `diffMF2` function calculates the value and directional derivative, as a pair, of a nonscalar-to-nonscalar function.

# Higher-Order Differentiation Operators

Naming convention:

`diffs{U/M}{U/F}`
: Derivative-taking operators that return a list [primal, first-derivative, 2nd-derivative, ...], for all combinations of scalar/nonscalar input & output.

diffsUU :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> [b]Source

The `diffsUU` function calculates a list of derivatives of a scalar-to-scalar function. The 0-th element of the list is the primal value, the 1-st element is the first derivative, etc.

diffsUF :: (Num a, Num b, Functor f, Foldable f) => (forall tag. Dual tag a -> f (Dual tag b)) -> a -> [f b]Source

The `diffsUF` function calculates an infinite list of derivatives of a scalar-to-nonscalar function. The 0-th element of the list is the primal value, the 1-st element is the first derivative, etc.

diffsMU :: (Num a, Num b) => (forall tag. [Dual tag a] -> Dual tag b) -> [[a]] -> [b]Source

The `diffsMU` function calculates an infinite list of derivatives of a nonscalar-to-scalar function. The 0-th element of the list is the primal value, the 1-st element is the first derivative, etc. The input is a (possibly truncated) list of the primal, first derivative, etc, of the input.

diffsMF :: (Num a, Num b, Functor f, Foldable f) => (forall tag. [Dual tag a] -> f (Dual tag b)) -> [[a]] -> [f b]Source

The `diffsMF` function calculates an infinite list of derivatives of a nonscalar-to-nonscalar function. The 0-th element of the list is the primal value, the 1-st element is the first derivative, etc. The input is a (possibly truncated) list of the primal, first derivative, etc, of the input.

diffs0UU :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> [b]Source

The `diffs0UU` function is like `diffsUU` except the output is zero padded.

diffs0UF :: (Num a, Num b, Functor f, Foldable f) => (forall tag. Dual tag a -> f (Dual tag b)) -> a -> [f b]Source

The `diffs0UF` function is like `diffsUF` except the output is zero padded.

diffs0MU :: (Num a, Num b) => (forall tag. [Dual tag a] -> Dual tag b) -> [[a]] -> [b]Source

The `diffs0MU` function is like `diffsMU` except the output is zero padded.

diffs0MF :: (Num a, Num b, Functor f, Foldable f) => (forall tag. [Dual tag a] -> f (Dual tag b)) -> [[a]] -> [f b]Source

The `diffs0MF` function is like `diffsMF` except the output is zero padded.

# Common access patterns

diff :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> bSource

The `diff` function is a synonym for `diffUU`.

diff2 :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> (b, b)Source

The `diff2` function is a synonym for `diff2UU`.

diffs :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> [b]Source

The `diffs` function is a synonym for `diffsUU`.

diffs0 :: (Num a, Num b) => (forall tag. Dual tag a -> Dual tag b) -> a -> [b]Source

The `diffs0` function is a synonym for `diffs0UU`.

grad :: (Num a, Num b) => (forall tag. [Dual tag a] -> Dual tag b) -> [a] -> [b]Source

The `grad` function calculates the gradient of a nonscalar-to-scalar function, using n invocations of forward AD, where n is the input dimmensionality. NOTE: this is O(n) inefficient as compared to reverse AD.

jacobian :: (Num a, Num b) => (forall tag. [Dual tag a] -> [Dual tag b]) -> [a] -> [[b]]Source

The `jacobian` function calcualtes the Jacobian of a nonscalar-to-nonscalar function, using n invocations of forward AD, where n is the input dimmensionality.

# Optimization Routines

zeroNewton :: Fractional a => (forall tag. Dual tag a -> Dual tag a) -> a -> [a]Source

The `zeroNewton` function finds a zero of a scalar function using Newton's method; its output is a stream of increasingly accurate results. (Modulo the usual caveats.)

TEST CASE: `take 10 \$ zeroNewton (\x->x^2-4) 1 -- converge to 2.0`

TEST CASE :module Data.Complex Numeric.FAD `take 10 \$ zeroNewton ((+1).(^2)) (1 :+ 1) -- converge to (0 :+ 1)`

inverseNewton :: Fractional a => (forall tag. Dual tag a -> Dual tag a) -> a -> a -> [a]Source

The `inverseNewton` function inverts a scalar function using Newton's method; its output is a stream of increasingly accurate results. (Modulo the usual caveats.)

TEST CASE: `take 10 \$ inverseNewton sqrt 1 (sqrt 10) -- converge to 10`

fixedPointNewton :: Fractional a => (forall tag. Dual tag a -> Dual tag a) -> a -> [a]Source

The `fixedPointNewton` function find a fixedpoint of a scalar function using Newton's method; its output is a stream of increasingly accurate results. (Modulo the usual caveats.)

extremumNewton :: Fractional a => (forall tag tag1. Dual tag1 (Dual tag a) -> Dual tag1 (Dual tag a)) -> a -> [a]Source

The `extremumNewton` function finds an extremum of a scalar function using Newton's method; produces a stream of increasingly accurate results. (Modulo the usual caveats.)

argminNaiveGradient :: (Fractional a, Ord a) => (forall tag. [Dual tag a] -> Dual tag a) -> [a] -> [[a]]Source

The `argminNaiveGradient` function performs a multivariate optimization, based on the naive-gradient-descent in the file `stalingrad/examples/flow-tests/pre-saddle-1a.vlad` from the VLAD compiler Stalingrad sources. Its output is a stream of increasingly accurate results. (Modulo the usual caveats.) The gradient is calculated using Forward AD, which is O(n) inefficient as compared to Reverse AD, where n is the input dimensionality.

# Miscellaneous

taylor :: Fractional a => (forall tag. Dual tag a -> Dual tag a) -> a -> a -> [a]Source

The `taylor` function evaluate a Taylor series of the given function around the given point with the given delta. It returns a list of increasingly higher-order approximations.

EXAMPLE: `taylor exp 0 1`