Download Solving Harmonic Oscillator Problem using Various Algorithms: A Comparative Analysis - Pro and more Assignments Mathematics in PDF only on Docsity! CMSC/AMSC 660-001 Fall 2004 Homework 4b In this problem we investigate various algorithms for solving our harmonic oscillator problem. Consider System 1 of ordinary differential equations q′ = ωp p′ = −ωq where ω = 1, q(0) = 0 and p(0) = 5. We can add a conservation principle to transform the system to System 2: q′ = ωp − ωqz p′ = −ωq − ωpz z′ = ω(p2 + q2) − 25 where z(0) = 0 Problem 1. (10) Solve System 1 by 4 methods: • ode45 • ode113 • ode15s • ode23t Then solve System 2 using these same 4 methods. In each case, draw a phase plot of the solution on the interval t ∈ [0, 3000], print the time taken by the algorithm, and print the norm of the vector [q(3000), p(3000)], which, for the exact solution, is 5. Answer: See Figures 1 and 2, the posted Matlab files, and the data in the following table, obtained on a Sun 4 Unix machine. System Method Solution norm Time ODE ode45 4.583719 5.19 sec ode113 4.900270 8.07 sec ode15s 11.211087 14.61 sec ode23t 4.997557 37.44 sec DAE ode45 5.000000 116.14 sec ode113 4.996717 46.65 sec ode15s 4.999156 94.11 sec ode23t 5.000007 1273.40 sec Problem 2. (10) • a. Explain what the true solution looks like and what the plots show. • b. What numerical integration methods do each of the four ode solvers use and when should each solver be used? • c. Which method is best for System 1? Which is best for System 2? • d. What are the disadvantages of using System 2? 1 Answer: a. From the class notes we know that true solution is periodic, and the norm of [q(t), p(t)] is 5 for all t. See Figure 3. The true phase plot is just a circle of radius 5 centered at the origin. For System 1 (see Figure 1), the ode15s solution grows, while ode45 is somewhat dissipative; this is also called damping. The solution norm for ode113 oscillates a bit, and that for ode23t is preserved quite nicely, so the plot is nearly a circle. For system 2 (see Figure 2), the solution norm is preserved rather well for each method. b. From http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ode113.html: Name Method Mathwork’s advice ode45 Runge-Kutta of order 4 and 5. Advised for general use. ode113 PECE Adams method. Advised if high accuracy is needed or if f is expensive. ode15s Stiff, NDF method Advised if ode45 fails or system is DAE. ode23t Trapezoidal rule Advised if stiffness is moderate and damping is unacceptable. Note that if you quote from the website, you should indicate that. c. ode113 works well for System 1 – it gives accurate results for a moderate cost, although ode23t gives 2 more significant digits for 4.6 times the cost. All of the DAE methods are expensive, with ode113 the cheapest but least accurate and ode45 giving the best results at roughly 2.5 times the cost. d. System 2 produces very good results but requires slightly more complicated analysis (to design the conservation principle) and programming (to include the conservation principle). It also require significantly more computation time; the penalty ranged from a factor of about 6 for ode113 and ode15s to a factor of 23 for ode45 and 34 for ode23t, because the timestep needs to be small to capture oscillations in z. 2
