PHYS 220 — HW #4 [Solution]


The equations of motion for a golf ball in a realistic atmosphere are given by (be careful of signs!):

The drag force is the same "quadratic" drag force we encountered in Lab #3 and lecture. S and ω represent a coefficient used to estimate the Magnus force on a spinning ball and the rotation rate of that ball.

  1. Write out the figurative representation of the vector containing the variables in the problem as vars = [ ... ].
    Figurative:
    y = [
    x
    y
    vx
    vy
    ]


  2. Using valid MATLAB code, write the single line of code for the derivative function (if you defined an auxiliary variable for ease of writing out that line, provide it as well) dvarsdt = [ ... ].
    MATLAB:
    vx = y(3);
    vy = y(4);
    v = sqrt(vx^2 + vy^2);
    dydt = [
    vx
    vy
    -(b*v*vx + S*w*vy)/m
    -g - (b*v*vy - S*w*vx)/m % NOTE SIGN CHANGE!
    ];


  3. List the parameters you need to pass into your derivative function.
    You'll need five parameters this time:
    • m → mass of golf ball
    • g → acceleration due to gravity
    • b → quadratic drag coefficient
    • S → Magnus coefficient
    • w → spin rate of ball



The following problem concerns error analysis/control. You're modelling a non-linear oscillator that has a characteristic amplitude of 45°. You're currently running your simulation with a time step Δto = 0.05 s. For the problem at hand you desire a fractional error of 0.5% of the amplitude scale.
For these problems, I wrote a short MATLAB function to calculate the new time step. Below I show the calls to it.

  1. At some point in your simulation you estimate a local truncation error of -0.09°. Calculate Δt1 required for the following methods:
    1. Euler method
      new_dt(0.05,-0.09,45*0.005,2) = 0.079
    2. Mid-point method (rk2) new_dt(0.05,-0.09,45*0.005,3) = 0.067
    3. Runge-Kutta method (rk4) new_dt(0.05,-0.09,45*0.005,5) = 0.060
  2. Further into the simulation you estimate a local truncation error of 0.9°. Calculate Δt1 required for the following methods:
    1. Euler method
      new_dt(0.05,0.9,45*0.005,2) = 0.025
    2. Mid-point method (rk2)
      new_dt(0.05,0.9,45*0.005,3) = 0.031
    3. Runge-Kutta method (rk4)
      new_dt(0.05,0.9,45*0.005,5) = 0.037
  3. Comment on the results found in 4 & 5, above.
    A few things to take notice of here. In each case I intentionally round down to a value lower than the 5-digit MATLAB result, ensuring that I still meet the error criteria with a smaller time step than is required. In §§4 you see that our time step is unnecessarily small for what we need out of the problem and that we can make it larger. Also, you can get away with larger time steps when switching to Runge-Kutta 4 from Euler, etc. In §§5 we have the opposite problem: too large a time step, so we see it made smaller for each method. Again, rk4 lets us get away with the largest time steps due to its improved precision.