Recipe11.2.Applying the RungeKutta Method to SecondOrder Initial Value Problems
Recipe 11.2. Applying the RungeKutta Method to SecondOrder Initial Value ProblemsProblem
You need to
SolutionThis is a standard initial value problem, and you can implement any of a number of standard numerical integration techniques to solve it using Excel and VBA. DiscussionConsider this equation and initial conditions:
Physically, this equation represents the equation of motion for an object moving under some applied
To solve this equation, it helps to rewrite it in the form of two firstorder equations. Let v = ds / dt represent the velocity of the object. Now we have the following:
These equations represent two
I decided to use the RungeKutta method for this example. This method is based on taking more terms in the Taylor series expansion of a function, as I explained in Recipe 11.1. The caveat here is that you need to apply further Taylor series expansions to estimate the higherorder derivatives. The RungeKutta approach
The RungeKutta integration formulas are a little more complicated than the simple formula we used with Euler's basic method. The general equations are as
y
' represents
dy
/
dx
in these equations. You can see that every step requires four intermediate computations prior to predicting the
Example 113 shows the VBA subroutine I wrote implementing these formulas. Example 113. Implementation of RungeKutta method
This is a relatively long subroutine compared to the one shown in Example 112 implementing Euler's method, so I commented this one to help you navigate through it.
As usual, several local variables are declared at the beginning of the subroutine. There are a bunch of them here, so I labeled each one in Example 113 to
After the variables are declared, several of them are
Figure 114. Solution using the RungeKutta method
Values retrieved from the spreadsheet include the timestep size, thrust, mass, drag coefficient, number of iterations, and number of output rows. The block of code in Example 113 containing the
With
statement
{% if main.adsdop %}{% include 'adsenceinline.tpl' %}{% endif %} The next little section of code simply initializes a few variables. This is where the initial values for displacement and velocity specified in the problem are set.
Next, the
For
loop iterates computing estimates for displacement and velocity at each step. Taking a closer look at the code within the
For
loop, you can see how I applied the RungeKutta formulas to compute
v
_{
t
}
_{
+
}
_{
dt
}
. To make things a little more intuitive (at least to me), I broke the computations of the
k
values into a few steps, using the analogy of Newton's second law of motion.
F
in these computations stands for the sum of forces. For this problem, the sum of forces is
I find it easier to think in terms of physical values like force and acceleration whenever possible. Knowing that acceleration is equal to dv / dt , you can readily see that acceleration corresponds to y ' in the RungeKutta formulas for this particular problem. This is because we broke the original equation up and are dealing with the first one involving velocity. After all of the k values are obtained, v _{ t } _{ + } _{ dt } is estimated using the RungeKutta formula discussed earlier. Next, a simple Euler estimate is made for the corresponding displacement, s _{ t } _{ + } _{ dt } . Finally, the results are sent to the active spreadsheet in the same manner I showed you in Recipe 11.1. Figure 114 shows the results, including a plot of velocity and displacement versus time, using a time step of 0.005. You should play around with this time step to see how the accuracy of the results is affected. (If this were a class, I'd assign that exercise for homework.)
I should also mention that I used a button control, as you can see in Figure 114, to call the subroutine when the button is clicked. Refer back to Recipe 9.3 to learn how to add this
