Errors in programming of Runge Kutta method in VBA

[brunette1216]

Hi, I am a Chemical Engineering Graduate student taking a class in programming. I was assigned a problem where a valve is closed quickly in a reservoir piping system during an emergency. I was asked to write a three step(rkck, rkqs, and odeint) program for the fifth order Runge Kutta method with adaptive step size which solves the governing differential equations (dydx(0) and dydx(1)).

I was given an algorithm handout for each of the subroutines required and followed them exactly, However, I keep running into issues when I run the macro at the end that incorporates three subroutines together.

The current issue is that yscal() in the program is zero. From the handout, yscal() is the vector of values against which the errors are scaled, but I am not sure what that means or how to make that work in VBA. I don't actually calculate yscal() and I feel like I should but I don't know how. I attached my excel file with the VBA coding and macros.

Can someone please help???

Thank you.

10852

[Bob Phillips]

I can't say I have followed all of the code, and I have no idea what the Runge Kutta method is, but I did see that you are loading yscal in procedure odeint, using this code


For I = 0 To nvar - 1
yscal(I) = Abs(y(I)) + Abs(h * dydx(I)) + Tiny
Next I


You say yscal is zero, where are you seeing that?

[Kenneth Hobs]

I had to change the Derivs1 to Derivs since you called it so many times.

You need an If() statement or an On Error routine to handle division by zero errors. Either skip that step or set the value in the If().

Use F8 to step through the code to see the computation of yscal as xld noted. The error first occurs in the sub rkqs when you run the main sub RungeKutta.

Add the Compile button to the VBE toolbar and compile your code. Your will see that in the Sub odeint, y() is not defined (dimmed). There are other similar problems such as dydx(). Option Explicit is your friend and shows these problems at Compile time. These Compile errors need to be resolved.

e.g.

Sub odeint(ystart() As Double, nvar As Integer, x1 As Double, x2 As Double, eps As Double, h1 As Double, hmin As Double, NOk As Integer, NBad As Integer)
Dim y() As Double
ReDim y(0 To nvar - 1) As Double