This lab will probably take us two lab periods. The assignment is therefore due at the beginning of lecture on Friday, 23 Feb. Your solution should include all codes and plots, as well as written responses where appropriate. You may use your own ode solver or the generic ODE solver & codes I have posted off of the course outline.

In this lab we'll investigate projectile motion much like you studied in PHYS 211; however, we'll also consider the effects of air resistance. For these problems the Euler method is sufficient to obtain precise solutions; however, you may use a higher order method if you so choose. You may also chose to investigate the differences between methods in arriving at your solutions (e.g., compare Runge-Kutta 4 to Euler).

Throughout this lab we'll consider a few simplifying assumptions. All of our projectiles are fired from a planar surface, in the presence of a uniform gravitational field. We'll also assume that the atmosphere is motionless, isobaric and isothermal. If you're having no trouble with MATLAB programming, you are encouraged to relax some of those assumptions (see 8-10). I expect that everyone should, using my codes if necessary, be able to work through 1-7 by the due date.

From PHYS 211 you should recall that, neglecting air resistance, projectiles launched from/to a plane have identical ranges (r) for complimentary angles of launch. That is, r_θ = r_90°-θ, where r(θ) = v_o² sin 2θ / g. This result isn't too practical for actual weapons systems — mortars typically cannot be depressed to less than 65° and howitzers (modern cannons) cannot be elevated to more than 35° — so our "gun" will have to be somewhat of a hypothetical weapon. We'll use the projectile equations of motion we developed in class for this lab. Keep in mind the points we discussed in class.

Our gun will be a customized 105 mm howitzer firing a high explosive (HE) round. We'll assume the following characteristics for this gun:

• muzzle velocity: 310 m/s (fixed and just subsonic)
• blast radius: 30 m
• elevation angle range: 5° – 60° (this is the only unrealistic characteristic, it's actually -5° – 35°)
• projectile mass: 15 kg
• drag coefficient: 6×10^-4 kg/m

Given these characteristics, complete the following tasks:

1. select two angles and their complements, within the elevation limits of our gun;
2. verify graphically that the ranges for the pairs are identical in vacuum;
3. verify numerically that these results are within half a blast radius of the theoretical range in vacuum;
   
   For §§1-7, I make use of artyShell_1.m as my derivative function. To fire the cannon without air resistance, I simply set b = 0 and re-run the code each time for each angle of fire. Saving the intermediate results at each angle allows me to plot them later. Below you'll see numerical results that correspond to the plot of the four trajectories for §§1-3. I used the Euler method for each. Since that method is relatively imprecise, I had to use a timestep of 0.05 s — using rk4, it could have been as large as 0.5 s, executed in the same amount of CPU time and given an order of magnitude better precision! The results below also make use of a code that performs a linear interpolation to find the impact point.

      theta     comp         ranges (km)     diff (blast r)
      35.0000   55.0000    9.2180    9.2142    0.1269
      15.0000   75.0000    4.9130    4.9021    0.3651
   

   
   
   
4. allow for air resistance and graphically compare the results with the same angles as above;
5. verify numerically that these results are within 1/10-th of a blast radius of one another.
   
   This result is as expected, the rounds fall short and our trajectories are less parabolic on the downward-falling side.
   
   
   
6. allow for air resistance and find the range for θ = 35°;
7. find another angle φ > θ such that r_θ = r_φ when air resistance is a factor;
   
   For §§5-7 I use the rk2 or mid-point method for improved precision with a larger timestep. It turns out that for these problems, rk2 delivers the precision of rk4 at half the computational cost. The table below shows how the search algorithm cuts the steps in θ by 1/2 and changes direction in the search each time I cross past the desired "target", or where the round landed for a 35° shot. Note that the effective complement is shallower than 55° without the effect of air and with Δt = 0.04 s we greatly exceed the required r_blast/10 criterion. A complete code, including this search algorithm (as I detailed in lab last time) is included here; it is exactly what I used to produce the output below.
   

     theta_try   range     target
      37.0000    7.4562    7.3398
      38.0000    7.5019    7.3398
      39.0000    7.5393    7.3398
      40.0000    7.5684    7.3398
      41.0000    7.5893    7.3398
      42.0000    7.6019    7.3398
      43.0000    7.6063    7.3398
      44.0000    7.6026    7.3398
      45.0000    7.5907    7.3398
      46.0000    7.5708    7.3398
      47.0000    7.5427    7.3398
      48.0000    7.5067    7.3398
      49.0000    7.4627    7.3398
      50.0000    7.4108    7.3398
      51.0000    7.3510    7.3398
      52.0000    7.2834    7.3398
      51.5000    7.3182    7.3398
      51.0000    7.3510    7.3398
      51.2500    7.3348    7.3398
      51.1250    7.3430    7.3398
      51.1875    7.3389    7.3398
   
      angle     range(km)
      35.0000    7.3398
    
      comp      range      diff (blast r)
      51.1875    7.3389    0.0295
   

   
   
   

Those students who are feeling comfortable with MATLAB and desire a greater challenge are encouraged to do the following problems. We're going to consider two changes to our modeling of the physics in the problem. The first is the fact that our assumption F_g = -mg was only an approximation to the force due to gravity. We know that this force changes with altitude, so that F_g = -GM_em/r² is a better model. We further assumed that F_drag = -bv², which is a good approximation of the drag at the firing position of the gun. However, since the shell reaches considerable altitude — the atmosphere is not isobaric, nor isothermal — a better approximation would be F^*_drag = ρ/ρ_o F_drag (y = 0), where ρ(y) = ρ_o ( 1 - ay/T_o)^α. All the other assumptions remain intact. For the following problems you should provide code(s), graphical and or numerical evidence that you've verified your results.

8. A practice target is situated 6 km from the firing point. Find the firing angle(s) that cause a shell to "hit" the target (impact within 1/10-th of a blast radius) using the new model. Assume α = 7/2, a = 6.5×10^-3 K/m, and T_o = 300 K.
   
   For §§8-9, I make use of artyShell_2.m as my derivative function, which handles the extra parameters required. You'll note that it takes scalar or vector input for g; this allows me to easily test the difference between the effects of the assumptions tested in §§9. I still use the midpoint method to fire these rounds. Below are the effective complements for variable air resistance that "hit" the target (that is, to better than 1/100^th of a blast radius!). You'll see how the same search algorithm I used in §§5-7 is applied here and how it gets the round within a 20 centimeters of the desired target.
   

    theta_try    range     target
      15.0000    4.3511    6.0000
      20.0000    5.4239    6.0000
      25.0000    6.2978    6.0000
      22.5000    5.8859    6.0000
      23.7500    6.0981    6.0000
      23.1250    5.9936    6.0000
      23.4375    6.0463    6.0000
      23.2812    6.0200    6.0000
      23.1250    5.9936    6.0000
      23.2031    6.0068    6.0000
      23.1641    6.0002    6.0000
      23.1250    5.9936    6.0000
      23.1445    5.9969    6.0000
      23.1641    6.0002    6.0000
   
    theta_try    range     target
      33.1641    7.2889    6.0002
      38.1641    7.6248    6.0002
      43.1641    7.7526    6.0002
      48.1641    7.6718    6.0002
      53.1641    7.3838    6.0002
      58.1641    6.8929    6.0002
      63.1641    6.2070    6.0002
      68.1641    5.3388    6.0002
      65.6641    5.7947    6.0002
      63.1641    6.2070    6.0002
      64.4141    6.0064    6.0002
      65.6641    5.7947    6.0002
      65.0391    5.9019    6.0002
      64.4141    6.0064    6.0002
      64.7266    5.9545    6.0002
      64.5703    5.9806    6.0002
      64.4141    6.0064    6.0002
      64.4922    5.9935    6.0002
      64.4531    6.0000    6.0002
      64.4141    6.0064    6.0002
      64.4336    6.0032    6.0002
      64.4531    6.0000    6.0002
   

   
   
   
9. Considering the problem above, compare the effects of each of the new model components with the model used in §§4-7. We can do this by variously setting a = 0 (isobaric/isothermal atmosphere) and r = r_e (uniform gravitational field). Which assumption, a = 0 or r = r_e, is physically more sound? Justify your response.
   
   To answer §§9, we've got to test the following combinations:
   • variable density & gravitational field
   • constant density & variable gravitational field
   • variable density & constant gravitational field
   • constant density & gravitational field
   I do this and get the results shown in the table & plot below. Notice that in either, it's the assumption of variable air density that has the biggest effect (in terms of θ_o) when we "turn it off". By switching off variable gravity, we see little effect. So even though our high trajectory reaches more than 3 km altitude, it's still the change in the density of air that makes the biggest impact on how we aim the gun. Note that although this is also true for our low trajectory, the differences are an order-of-magnitude reduced – i.e., neither assumption is particularly important for the low angles of fire.

      theta     theta_c    dTh       dTh_c
      23.1641   64.4531         0         0
      23.3203   63.3594    0.1562   -1.0938
      23.2031   64.4141    0.0391   -0.0391
      23.3594   63.2812    0.1953   -1.1719
   

   
   
   
   
10. Automate your approach(es) to finding solutions to any of the problems 1-9.