Using the AlgLib ODE (Runge-Kutta) Solver with Excel

Or to be more correct, the AlgLib Cash-Karp Solver, Cash-Karp being a refinement of the Runge-Kutta method of solving ordinary differential equations.

Solution of differential equations is an iterative process requiring the repeated application of the solver routine followed by evaluation of the differential equation using an external routine.  This is performed in AlgLib by a process known as reverse communication, whereby an initialisation routine creates a “State” object which is used to transfer data between the solver routine and the external equation evaluation routine (which is in this case a VBA function).  The steps required to wrap this process in an Excel User Defined Function (UDF) are:

  1. Read the name of the ODE evaluation routine and the required input data from the spreadsheet.
  2. Convert data in spreadsheet ranges to variant arrays.
  3. Convert variant arrays (which will always be base 1, 2 dimension arrays) into the form required by the AlgLib routines, which will always be base 0, and may be 1 or 2 dimensional.
  4. Evaluate any additional values required by the AlgLib routines
  5. Set up any required output arrays
  6. Call the appropriate AlgLib routine to initialise the State object
  7. Loop through the AlgLib solver routine and the external ODE evaluation routine until the solver routine returns “False”, which indicates that it has finished.  The VBA function named in “FuncName” is called using the Application.Run method.
  8. Extract the desired results from the State object using the appropriate AlgLib routine
  9. If necessary convert the AlgLib output arrays to a form suitable for Excel.
  10. Assign the output array to the function return value.
  11. Import the necessary AlgLib routines; in this case the only required modules are ODESolver.bas and ap.bas

This process is illustrated in the code for the UDF ODE() shown below (may take some time to display correctly). 

 This code and all the related AlgLib and other routines are included in the file ODESolver.zip.

Function ODE(FuncName As String, Initial As Variant, XA As Variant, _
CoeffA As Variant, Optional Eps As Double = 0.000001, _
Optional Step As Double = 0, Optional MaxIts As Long = 100)
Dim M As Long, N As Long, State As ODESolverState, _
YA() As Double, XA2() As Double, YA2() As Double, i As Long, _
Rtn As Boolean, Rep As ODESolverReport, NC As Long

' The VBA function arguments are the name of the VBA
' ODE evaluation routine, followed by the required input data.

' Convert data in spreadsheet ranges to variant arrays
Initial = GetArray(Initial)
XA = GetArray(XA)
CoeffA = GetArray(CoeffA)

' Convert variant arrays (which will always be base
' 1, 2 dimension arrays) into the form required by the AlgLib
' routines, which will always be base 0, and may be 1 or 2 dimensional.
Rtn = VarAtoDouble1D_0(Initial, YA, N, NC)
If N = 1 And NC > 1 Then N = NC
Rtn = VarAtoDouble1D_0(XA, XA2, M, NC)

' Evaluate any additional values required by the AlgLib routines
MaxIts = MaxIts * M

'Set up any required output arrays
ReDim YA2(0 To M - 1, 0 To N - 1)

'Call the appropriate AlgLib routine to initialise the State object
Call ODESolverRKCK(YA(), N, XA2, M, Eps, Step, State)

' Loop through the AlgLib solver routine and the external ODE
' evaluation routine until the solver routine returns "False",
' which indicates that it has finished.
' The VBA function named in "FuncName" is called using
' the Application.Run method.
Rtn = True
i = 0
Do While Rtn = True And i < MaxIts
Rtn = ODESolverIteration(State)
State.DY = Application.Run(FuncName, State.X, State.Y, CoeffA)
i = i + 1
Loop

 ' Extract the desired results from the State
 ' object using the appropriate AlgLib routine
 Call ODESolverResults(State, M, XA2, YA2, Rep)

' If necessary convert the AlgLib output array(s) to
' a form suitable for Excel. In this case YA2 is
' a 2D base 0 array, which may be assigned to the
' function return value without further processing.

' Assign the output array to the function return value
 ODE = YA2

End Function

The screen shots below show examples of the application of the ODE() function (more details of the final example, buckling of a strut, will be given in the next post):

ODE input and output details, click for full size view

Example 1; First Order kinetic process

Example 1 Output

Example 2; Process dependent on X and Y

Example 2 results

Example 3; damped vibration system

Example 3 output

Example 4; Undamped SHM, ball through the centre of the Earth

Example 4 results

Example 5; deflection of a cantilever strut and buckling load

Example 5 results

This entry was posted in AlgLib, Beam Bending, Differential Equations, Excel, Maths, Newton, UDFs, VBA and tagged , , , , , , , , . Bookmark the permalink.

5 Responses to Using the AlgLib ODE (Runge-Kutta) Solver with Excel

  1. Pingback: Daily Download 21: Assorted Solvers | Newton Excel Bach, not (just) an Excel Blog

  2. miro ilias says:

    Dear Doug,

    thanks for this post, which for single ODE replaces the XNumbers ODE solver.

    The question is whether your solver can be applied for system of differential equations…

    Yours, Miro

    Like

  3. Pingback: Rabbits, Foxes, and Lorenz Attractors | Newton Excel Bach, not (just) an Excel Blog

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s