Numerical integration of stochastic differential equations (SDE)

## Project description

## Overview

sdeint is a collection of numerical algorithms for integrating Ito and Stratonovich stochastic ordinary differential equations (SODEs). It has simple functions that can be used in a similar way to `scipy.integrate.odeint()` or MATLAB’s `ode45`.

There already exist some python and MATLAB packages providing Euler-Maruyama and Milstein algorithms, and a couple of others. So why am I bothering to make another package?

It is because there has been 25 years of further research with better methods but for some reason I can’t find any open source reference implementations. Not even for those methods published by Kloeden and Platen way back in 1992. So I will aim to gradually add some improved methods here.

This is prototype code in python, so not aiming for speed. Later can always rewrite these with loops in C when speed is needed.

Warning: this is an early pre-release. Wait for version 1.0. Bug reports are very welcome!

## functions

`itoint(f, G, y0, tspan)`for Ito equation dy = f(y,t)dt + G(y,t)dW

`stratint(f, G, y0, tspan)`for Stratonovich equation dy = f(y,t)dt + G(y,t)∘dW

These work with scalar or vector equations. They will choose an algorithm for you. Or you can use a specific algorithm directly:

## specific algorithms:

`itoEuler(f, G, y0, tspan)`: the Euler-Maruyama algorithm for Ito equations.

`stratHeun(f, G, y0, tspan)`: the Stratonovich Heun algorithm for Stratonovich equations.

`itoSRI2(f, G, y0, tspan)`: the Rößler2010 order 1.0 strong Stochastic Runge-Kutta algorithm SRI2 for Ito equations.

`itoSRI2(f, [g1,...,gm], y0, tspan)`: as above, with G matrix given as a separate function for each column (gives speedup for large m or complicated G).

`stratSRS2(f, G, y0, tspan)`: the Rößler2010 order 1.0 strong Stochastic Runge-Kutta algorithm SRS2 for Stratonovich equations.

`stratSRS2(f, [g1,...,gm], y0, tspan)`: as above, with G matrix given as a separate function for each column (gives speedup for large m or complicated G).

`stratKP2iS(f, G, y0, tspan)`: the Kloeden and Platen two-step implicit order 1.0 strong algorithm for Stratonovich equations.

### utility functions:

`deltaW(N, m, h)`: Generate increments of m independent Wiener processes for each of N time intervals of length h.

`Ikpw(dW, h, n=5)`: Approximate repeated Ito integrals.

`Jkpw(dW, h, n=5)`: Approximate repeated Stratonovich integrals.

`Iwik(dW, h, n=5)`: Approximate repeated Ito integrals.

`Jwik(dW, h, n=5)`: Approximate repeated Stratonovich integrals.

## Examples:

`x0 = 0.1`

import numpy as np import sdeint a = 1.0 b = 0.8 tspan = np.linspace(0.0, 5.0, 5001) x0 = 0.1 def f(x, t): return -(a + x*b**2)*(1 - x**2) def g(x, t): return b*(1 - x**2) result = sdeint.itoint(f, g, x0, tspan)

`x = (x1, x2)`,

`dW = (dW1, dW2)`and with initial condition

`x0 = (3.0, 3.0)`

import numpy as np import sdeint A = np.array([[-0.5, -2.0], [ 2.0, -1.0]]) B = np.diag([0.5, 0.5]) # diagonal, so independent driving Wiener processes tspan = np.linspace(0.0, 10.0, 10001) x0 = np.array([3.0, 3.0]) def f(x, t): return A.dot(x) def G(x, t): return B result = sdeint.itoint(f, G, x0, tspan)

## References for these algorithms:

`itoEuler`:

`stratHeun`:

`itoSRI2, stratSRS2`:

`stratKP2iS`:

`Ikpw, Jkpw`:

`Iwik, Jwik`:

## TODO

- Rewrite
`Iwik()`and`Jwik()`so they don’t waste so much memory. - Fix
`stratKP2iS()`. In the unit tests it is currently less accurate than`itoEuler()`and this is likely due to a bug. - Implement the Ito version of the Kloeden and Platen two-step implicit alogrithm.
- Add more strong stochastic Runge-Kutta algorithms. Perhaps starting with Burrage and Burrage (1996)
- Currently prioritizing those algorithms that work for very general d-dimensional systems with arbitrary noise coefficient matrix, and which are derivative free. Eventually will add special case algorithms that give a speed increase for systems with certain symmetries. That is, 1-dimensional systems, systems with scalar noise, diagonal noise or commutative noise, etc. The idea is that
`itoint()`and`stratint()`will detect these situations and dispatch to the most suitable algorithm. - Eventually implement the main loops in C for speed.
- Some time in the dim future, implement support for stochastic delay differential equations (SDDEs).

## See also:

`nsim`: Framework that uses this `sdeint` library to enable massive parallel simulations of SDE systems (using multiple CPUs or a cluster) and provides some tools to analyze the resulting timeseries. https://github.com/mattja/nsim

## Project details

## Download files

Download the file for your platform. If you're not sure which to choose, learn more about installing packages.