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:
- Read the name of the ODE evaluation routine and the required input data from the spreadsheet.
- Convert data in spreadsheet ranges to variant arrays.
- 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.
- Evaluate any additional values required by the AlgLib routines
- Set up any required output arrays
- Call the appropriate AlgLib routine to initialise the State object
- 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.
- Extract the desired results from the State object using the appropriate AlgLib routine
- If necessary convert the AlgLib output arrays to a form suitable for Excel.
- Assign the output array to the function return value.
- 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):