PHYS 220 — HW #5 [Solution]


  1. Which of the following ode could give rise to chaos? Why? Note that all of the y/n responses below only allow for the possibility of chaos, they do not guarantee it. Also note that all of them are essentially oscillators of some kind.
    1. y'' = y - ay' + b cos t
      No — lacking a non-linear term in y
    2. d²φ/dt² = sin φ
      No — lacking a driving or dissipative term
    3. x'' = sin x - ax' + b cos t - cx3
      Yes — nonlinear in x, damped & driven
    4. d²y/dt² = cos(y + π/3) - a dy/dt
      No — no driving term.


  2. In an electromagnetic simulation, an rk4 ode integrator produces an estimated local error of 0.8 V at a particular time when the time step is set to Δt = 15 minutes. If you require a precision of at least 0.1 V, what time step should you use for this particular time?
    For this problems, I use my MATLAB function to calculate the new time step.
    new_dt(15,0.8,0.1,5) ≈ 9.89 min (note round down!)



  3. In the problem above, now assume you find that the errors remain roughly constant in time throughout the simulation. If you require the simulation to run for a time of 17 weeks, and can tolerate a global error or of no more than 5 V at the end of the run, what time step (Δt) should you use?
    First, a conversion: 17 wk = 171360 min. Now we can set up a proportion: 171360 min / Δt1 = 5 V / el1. We also know that Δt1 = Δto ( el1 / el0)1/5. Since we know from the results of our last problem that el0 = 0.1 V for Δt0 = 9.89 min, we have two equations and two unknowns.
    The solution requires a little algebra and is given by Δt1 = (5Δto5 / 171360 / el0)1/4 ≈ 2.29 min.


  4. Orbital mechanics problems are governed by a central force law: FG = GMm / r². Write both the figurative construction of vars = [ ... ], and the literal line of MATLAB code corresponding to dvarsdt = [ ... ] for the Cartesian form of this central force law.
    Figurative:

    y = [
        x
        y
        vx
        vy
    ]

    Literal:

    r3 = sqrt(y(1)^2 + y(2)^2)^3;
    dydx = [
        y(3)
        y(4)
        G*M*m*y(1) / r3
        G*M*m*y(2) / r3
    ]

    Note that the trick of x/r3 has the effect of 1/r² in the x-direction; similarly for the y-component.