3.  Numerical integration of ODEs

• read Numerical Recipes, Chapter 17 introduction
• initial-value problems
• integration of Newton's equations of motion
• standard form of an ODE

• some standard integration schemes
  • read Numerical Recipes, Sec. 17.1

  • the Euler method
  • examples:
    • euler1.py: first-order
    • euler2.py: second-order

  • application: the Duffing oscillator
    • starting code: duffing1.py, y0 = 1.3, vy0 = 1.5
    • Exercise 3.1: Run the script and stop the calculation after one orbit. Interpolate the the velocity to determine its value at y = y0. Tune the step dx until the error in vy at y = y0 is less than 0.1%
    • solution: ex3.1.py
    • Exercise 3.2: Repeat Exercise 3.1 using this dx, and plot vy versus y for each of the following initial velocities: vy0 = 0.25, 0.5, 1.0, 1.25, 1.5, 2.0.
    • solution: ex3.2.py
    • Exercise 3.3: Repeat Exercise 3.2 for dx = 1.e-4 and initial velocities vy0 = 1.38, 1.39, 1.40
    • solution: ex3.3.py

  • Taylor series and order of the Euler method

  • problems with Euler
    • low-order, low accuracy
    • unstable: see euler_instability.py
    • implicit Euler: sample script

  • Exercise 3.4: repeat Exercise 3.2 using the implicit Euler method
  • solution: ex3.4.py

  • the Mid-point method
    • example: midpoint2.py

  • Exercise 3.5:
    • modify the program in Exercise 3.2 to use the Mid-point method
    • repeat Exercise 3.2 using Mid-point
    • solution: ex3.5.py

  • Fourth-order Runge-Kutta integrator
    • example: see duffing.py

  • sample code for some Runge-Kutta schemes

  • Exercise 3.6:
    • Use one of these improved integrators (i.e. NOT Euler) to explore how the period of the Duffing oscillator depends on the amplitude of the motion. Use the same parameters as before, choose y0 = 1.3, and compute the time taken to return to the same value of y as a measure of the period. Plot the period as a function of vy0 for vy0 = 0.1, ..., 10, in steps of 0.01.
    • solution: ex3.6.py

• Butcher Tableaux

• Higher-order Runge-Kutta schemes

• The damped, driven, nonlinear oscillator (here's a bare-bones solver as a starting point)

  • Exercise 3.7:
    • Explore how the solution to the problem changes as first damping and then driving are introduced.
      (1) Set alpha = -1, beta = 1, y0 = 1.5, v0 = 0, and modify the script to show two plots: the time sequence (y versus x) at the left and the phase portrait (v versus y) at right.
      (2) Then set delta = 0.3 and 0.6. As expected, the dissipation causes the oscillation to decay on a time scale inversely proportional to delta.
      (3) Now, with delta = 0.3, set omega = 1.2 and gamma = 0.2. Note how the initial oscillation dies away and the system ends up oscillating at the driving frequency omega. From here on, we will discard the decay of the initial conditions by only plotting results after some time. Modify the program to start plotting at x = 200 and continue until x = 500.
      (4) Now successively try gamma = 0.28, 0.29, 0.35, 0.5, 0.65. You should see period doubling and a transition into and out of chaotic motion. How does the behavior you see depend on the choice of integration scheme?

• Homework #2

• some alternative integrators: the second-order predictor-corrector integration scheme
  • time reversibility of the predictor-corrector scheme
  • the Leapfrog and velocity Verlet schemes
• symplectic integrators: the drift-kick-drift and kick-drift-kick schemes

• Exercise 3.8: Apply the drift-kick-drift and kick-drift-kick schemes to the same (undamped, undriven) Duffing system studied in Exercise 3.2, except that now we will switch notation to use t (not x) as the independent variable (time) and pos and vel as position and velocity. Use a time step of dt = 0.05 in each case. Stop when the system has completed one orbit and plot the phase portrait vel versus pos. Starting integrators can be found here.
• solution

• higher-order symplectic schemes
• sample code for some symplectic schemes