Kieran's Components Logo

Kieran's Components

Explore the world of UI

Numerical ODE Methods

22 Apr 2020




Mathematical Background

Euler’s method is a simple but effective strategy for numerically approximating ODE’s. The method takes advantage of the formal definition of the derivative and Taylor series.

This is called the Forward difference approximation. Changing the value from $x+h$ to $x-h$ creates what’s called the backward difference approximation

This method is easily applied to systems of ODE’s that frequently occur in classical dynamics.

How To Use

  1. Starting with some differential equation algebraically solve for $f’(x)$ so that computing the RHS gives you a numerical value $f’(x)$.
  2. Using a known starting value of $f(x_0) = f_0$, compute the next value using some step size  $h$. $f(x+h) \approx f(x) + h f’(x)$
  3. Repeat successively using the last $f(x+h)$ and the next value for $f(x)$.
  4. Continue until all needed values are obtained.

Important Notes

What This Means For A Programmer

In programming terminology every single expression for a first order derivative must be programmed in as a function. So after reducing your differential equations and solving for the first order derivatives, you must create a function for each expression. That function will look something like this:

func someDerivative(independentValue: Double, otherDependencies: [Double]) -> Double {
    // Your expression wont just be some simple multiplication
   return independentValue*otherDependencies.first!
}

Since the expression for the derivative likely depends not only on the dependent variable but possibly also upon value of the unknown function itself or some other dynamic values.

Using this format for expressing first order derivatives as functions, The amount of code needed to be written for solving larger systems can be significantly reduced. This means that an Euler solver can be written just once instead of for each derivative function.


Simple Euler


Uses simple finite difference with no trial steps Each call to this function performs one step. This implies that the simpleEuler function is meant to be used as part of a larger solve.

func simpleEuler(_ stepSize: Double,
                 _ independentVariable: Double,
                 _ dependentVariables: [Double],
                 functions: [(Double, [Double]) -> Double]) -> [Double] {
    var newValues = [Double]()
    functions.enumerated().forEach { (i, f) in
        newValues.append( dependentVariables[i] + f(independentVariable, dependentVariables)*stepSize)
    }
    return newValues
}

This version of Euler’s method is considered to be the simplest and least effective approximation scheme in terms of convergence. Many other schemes make use of “trial steps”, such as a modified version of Euler’s method or the popular Runge-Kutta 4th Order scheme.


Improved Euler


Makes use of a trial step and then averages the trial and real step values.

func improvedEuler(_ stepSize: Double,
                   _ independentVariable: Double,
                   _ dependentVariables: [Double],
                   functions: [(Double, [Double]) -> Double]) -> [Double] {
    var trialValues = [Double]()
    functions.enumerated().forEach { (i, f) in
        trialValues.append( dependentVariables[i] + f(independentVariable, dependentVariables)*stepSize)
    }
    var newValues = [Double]()
    functions.enumerated().forEach { (i, f) in
        newValues.append(dependentVariables[i] + (f(independentVariable, trialValues) + f(independentVariable, dependentVariables))*stepSize/2)
    }
    return newValues
}


Runge-Kutta 4th Order


By far the most overused albeit useful numerical scheme is the Runge-Kutta 4th Order scheme.

Approximates the solutions to differential equations using the 4th order Runge-Kutta numerical scheme.

func rK4(_ stepSize: Double,
         _ independentVariable: Double,
         _ dependentVariables: [Double],
         functions: [(Double, [Double]) -> Double]) -> [Double] {
    var k1: [Double] = []
    functions.forEach { (f) in
        k1.append(stepSize*f(independentVariable, dependentVariables))
    }
    var k2: [Double] = []
    functions.forEach { (f) in
        let newD = dependentVariables.enumerated().map { (index, nd)  in
            nd + k1[index]/2
        }
        k2.append(stepSize*f(independentVariable + stepSize/2, newD))
    }
    var k3: [Double] = []
    functions.forEach { (f) in
        let newD = dependentVariables.enumerated().map { (index, nd)  in
            nd + k2[index]/2
        }
        k3.append(stepSize*f(independentVariable + stepSize/2, newD))
    }
    var k4: [Double] = []
    functions.forEach { (f) in
        let newD = dependentVariables.enumerated().map { (index, nd)  in
            nd + k3[index]
        }
        k4.append(stepSize*f(independentVariable + stepSize, newD))
    }

    var kTot: [Double] = []
    for index in 0..<dependentVariables.count {
        kTot.append((k1[index] + 2*k2[index] + 2*k3[index] + k4[index])/6)
    }

    var newValues: [Double] = []

    for index in 0..<dependentVariables.count {
        newValues.append(dependentVariables[index] + kTot[index])
    }

    return newValues
}

Wrapping Up

To string all of our schemes together we can make an enumeration to make picking between the easier

enum ODEScheme: Int, CaseIterable, Identifiable, Hashable {
    case rungeKutta
    case euler
    case improvedEuler

    var scheme: (Double, Double, [Double], [(Double, [Double]) -> Double]) -> [Double] {
        switch self {
        case .rungeKutta:
            return rK4(_:_:_:functions:)
        case .euler:
            return simpleEuler(_:_:_:functions:)
        case .improvedEuler:
            return improvedEuler(_:_:_:functions:)
        }
    }
    var name: String {
        switch self {
        case .rungeKutta:
            return "Runge-Kutta"
        case .euler:
            return "Euler"
        case .improvedEuler:
            return "Improved Euler"
        }
    }
    var id: Int { rawValue }
}

Numerical ODE Schemes