The Magic Angles of Projectile Motion

Mathematica is used to express the length of a projectile's trajectory in a vacuum as a function of the area it encompasses.  The author shows that independent of kinematics, there is a unique pair of initial angles, (38.82^°) and (78.18^°), resulting in two different parabolas identical in length and encompassing the same areas.

by Haiduke Sarafian

The problem of projectile motion in a vacuum has been discussed in physics texts [e.g. Halliday, Resnick and Walker] and the scope of interest was extended in two articles by  [Sarafian, Tan].  Most recently, with the help of  Mathematica we discovered additional new features associated with the projectile problem.  We update the current status of the projectile problem, and we show that a science project can interactively be completed with Mathematica, avoiding the routine paper and pencil approach.

The traditional topics of  interest of the projectile problem in a vacuum under the influence of gravity are certain geometrical properties of its parabolic trajectory.  Such properties of interest are expressions for the maximum height and the range, the time of flight and ancillary information such as a (45^°) angle to maximize the range.

In this article we first derive the trajectory's length from elementary considerations   We investigate the properties of the length, leading to the discovery of a few new features about the projectile problem.

In the absence of air resistance, a projected point-size projectile moves in a vertical plane.  The parametric, time dependent coordinate of the projectile in a cartesian coordinate system is:  x(t) = (Subscript[v, 0]) cos(θ) t and y(t) = - (1/2)g (t^2)+(Subscript[v, 0]) sin(θ) t,  here (Subscript[v, 0]) is the initial speed,  θ is the initial projection angle above the horizontal and g is the gravity constant.  The position vector of the projectile is  \(TraditionalForm\`\(r\& -> \)\ \((t)\)\  = \ \(x(t)\)\ i\&^ + \(y(
              t)\) j\&^\),   the differential length of the trajectory is  d L = \(TraditionalForm\`\@\(d\ \(\(r\& -> \)(t)\) . d\ \(\(r\& -> \)(t)\)\)\).  Integration of this expression is the length of the trajectory and can be written as L = (∫_0^T) \(TraditionalForm\`\@\(\(\(d\ \(\(r\& -> \)(t)\)\)\/dt\) . \(\(d\ \(\(r\& -> \)(t)\)\)\/dt\)\)\)dt;  T = 2 ((Subscript[v, 0])/g)sin(θ) is the flight time.

A typical parabolic trajectory is shown in Figure 1.  This is a trajectory of a 20 (m/s) projectile thrown at  θ = (30^°) above the horizontal and the gravity constant is 9.8 (m/(s^2)).   The initial values and the position vector are

Off[General::Spell1]

initialVals = {g -> 9.8, 
          (Subscript[v, 0]) -> 20, θ -> (30^°)};

posnvector[
        t_] := {(Subscript[v, 0]) Cos[θ] t, - (1/2) g (t^2) + 
          (Subscript[v, 0]) Sin[θ] t}

plotposnvector = 
        ParametricPlot[
          Evaluate[posnvector[t] /. initialVals], {t, 0, 
            Evaluate[2 ((Subscript[v, 0])/g) Sin[θ] /. initialVals]}, 
          PlotStyle -> {Thickness[0.008]}, 
          AxesLabel -> {"x(m)", "y(m)"}, 
          GridLines -> Automatic];

FIGURE 1: A typical parabolic trajectory in a vacuum.

To evaluate the length of the trajectory, we begin with the integrand of the length equation  \(TraditionalForm\`\@\(\(\(d\ \(\(r\& -> \)(t)\)\)\/dt\) . \(\(d\ \(\(r\& -> \)(t)\)\)\/dt\)\)\).  This is the speed of the projectile and we evaluate the velocity (Overscript[v,  ->]) (t) = (d/dt) (Overscript[r,  ->]) (t) and then the speed.  We had a different approach in [Sarafian].

velocity[t_] := Evaluate[D[posnvector[t], t]]

speed = (Sqrt[velocity[t] . velocity[t]])

In order to more efficiently fulfill our objectives and arrive at an appropriate format for the length equation we replace  gt with gt + (Subscript[v, 0]) sin(θ),  which yields

2 (∫_0^((T/2)))Evaluate[
                    speed /. g t -> 
                        g t  +  (Subscript[v, 0]) Sin[θ]] \[DifferentialD]t /. 
            T -> 2  ((Subscript[v, 0])/g) Sin[θ] // Simplify // 
        PowerExpand // Simplify

The length equation is a product of a kinematic factor (((Subscript[v, 0])^2)/g)and an angular function.  We trace this feature of the length to the separable variable characteristics of the projectile's coordinate.  To study the generic features of the length, we drop the kinematic factor (((Subscript[v, 0])^2)/g), and by applying the logarithm laws we simplify the length equation, which yields

L = (% /. Log[a_ b_] :> Log[a] + Log[b] /. 
            Log[a_] -  Log[b_] :> Log[(a/b)]) (g/((Subscript[v, 0])^2))

The limits of the length for two extreme measures of  θ = (0^°)and  (90^°) are the expected intuitive values.  Recall the maximum height of the projectile is  (Subscript[y, max])= (((Subscript[v, 0])^2)/(2 g))sin(θ).  Twice the value of (Subscript[y, max]) for θ = (90^°) with suppressed kinematic factor (((Subscript[v, 0])^2)/g) is the length of the trajectory and equals unity.

{L /. θ -> 0  °, 
      L /. θ -> N[π]/2}

For further exploration we graph L versus θ.

\(\(plotLength = 
        Plot[Evaluate[L /. θ -> θ\ Degree\ ], {θ, 
            0, 90\ }, PlotRange -> {0, 1.5}, PlotStyle -> {Thickness[0.008]}, 
          AxesLabel -> {\*""\!\(θ\^°\)"", "L"}, 
          GridLines -> Automatic];\)\)

FIGURE 2: Kinematic independent length of a trajectory versus initial projected angle.

Figure 2 shows that independent of the initial speed, the length is at its maximum value if the projectile is thrown at about (55^°).  The limits of the length for θ = {0, (90^°)} corresponding to L = {0,1} are shown in Figure 2.

To evaluate the precise value of the initial projectile angle that maximizes the length we apply the FindMinimum command on the negative length equation, which  yields

maxLength = FindMinimum[Evaluate[- L], {θ, 1. }]

The corresponding angle is (56.4658^°)

In projectile motion the overlooked length-related geometrical quantity is the area encompassed by the trajectory.  One would intuitively expect the longest length to correspond to the largest area; curiously this is not the case.  To explore the reason, we evaluate the area.  First we eliminate t between the coordinate of the projectile.  We then solve the result to express y as a function of x:

yx = 
      Eliminate[{posnvector[
                t][[1]] == 
            x, posnvector[
                t][[2]] == y}, 
        t]

soly = Solve[yx, y]

To evaluate the area we integrate y over x with a limit of 0 to R = (((Subscript[v, 0])^2)/g)sin(2θ).  By suppressing the kinematic factor, as we did earlier, the genuine angle-dependent area becomes

Area = ((∫_0^R)(y /. soly ) \[DifferentialD]x /. 
              R -> (((Subscript[v, 0])^2)/g) Sin[2 θ] // 
            Simplify)  ((((g/((Subscript[v, 0])^2))))^2)

The graph of the area, along with the length is shown in Figure 3.

plotArea = 
        Plot[Evaluate[Area /. θ -> θ Degree], {θ, 0, 
            90}, PlotStyle -> {Dashing[{0.02, 0.02}], Thickness[0.008]}, 
          DisplayFunction -> Identity];

\(\(Show[{plotLength, plotArea}, 
        AxesLabel -> {\*""(θ^°)"", "Length,Area"}, GridLines -> {Range[0, 90, 5], Automatic}, 
        DisplayFunction -> $DisplayFunction];\)\)

FIGURE 3: The length (solid line) and the area encompassed by the length (dashed line) versus the initial projectile angle.

Clearly, the maxima of these two functions occur at different angles.  The angle that maximizes the area is evaluated below.

maxArea = 
      FindMinimum[- 
          Area[[1]], {θ, 
          1. }]

The corresponding angle is (60.0^°).  The angles maximizing the length and the area are not the same; they are off by (3.53^°).  To explore further, we form the algebraic relationship between the length and the area.  In doing so, we eliminate the common angle θ between the length and the area.  Starting with the area equation,  we square the area and replace (sin^2)(θ) with z, forming an algebraic z-dependent symbolic equation and solve the equation for z.

a = (((Area[[1]] /. Cos[q_] :> (Sqrt[1 - Sin[q]^2])))^2) /. {Sin[q_]^2 -> z, 
            Sin[q_]^6 -> (z^3)};

solz = Solve[a == (A^2), z];

The acceptable roots of this equation are solz[[3]] and solz[[4]]; the former is given below.   The corresponding angular domain of validity of these roots are discussed later in the text.

solz[[3]]

To express length L in terms of area A, we manipulate L accordingly.  We suppress the output because of the lengthiness of the expression.  However, we depict the corresponding equation in Figure 4.

L /. Cos[q_] :> (Sqrt[1 - Sin[q]^2]) /. {Sin[q_]^2 -> z, 
            Sin[q_] -> (Sqrt[z])} /. 
        solz[[3]];

Plot[L /. Cos[q_] :> (Sqrt[1 - Sin[q]^2]) /. {Sin[q_]^2 -> z, 
              Sin[q_] -> (Sqrt[z])} /. 
          solz[[3]], {A, (1. /(10^7)), 
          0.216}, PlotRange -> {0. , 1.2}, 
        PlotStyle -> {Thickness[0.008]}, 
        AxesLabel -> {"Area", "Length"}, 
        GridLines -> Automatic];

FIGURE 4:  The length versus the area encompassed by the length for the initial angle (0^°)≤ θ ≤ (60^°).

The maximum range of  A is 0.216.  This corresponds to the maximum value of A(θ = (60^°)); hence solz[[3]] corresponds to 0 ≤ θ ≤ (60^°).   The second real root of solz is solz[[4]] and corresponds to  (60^°)≤ θ ≤ (90^°).  We suppress the analytic output of the length L equation corresponding to solz[[4]];  however, we depict the corresponding L in Figure 5.

L /. Cos[q_] :> (Sqrt[1 - Sin[q]^2]) /. {Sin[q_]^2 -> z, 
            Sin[q_] -> (Sqrt[z])} /. 
        solz[[4]];

Plot[L /. Cos[q_] :> (Sqrt[1 - Sin[q]^2]) /. {Sin[q_]^2 -> z, 
              Sin[q_] -> (Sqrt[z])} /. 
          solz[[4]], {A, (1. /(10^7)), 
          0.216}, PlotRange -> {0. , 1.2}, 
        PlotStyle -> {Thickness[0.008]}, 
        AxesLabel -> {"Area", "Length"}, 
        GridLines -> Automatic];

FIGURE 5:  The length versus the area encompassed by the length for the initial angle (60^°)≤ θ ≤ (90^°).

Therefore, one can overlay Figure 4 on Figure 5 to display length L versus Area for 0 ≤ θ ≤ 
        (90^°).  Alternatively, we plot length L versus Area directly using L(θ)  and Area(θ), which is shown in Figure 6.

ParametricPlot[
        Evaluate[{Area[[1]], 
            L}], {θ, 0. , 89.9 (π/180)}, 
        PlotRange -> {0. , 1.2}, PlotStyle -> {Thickness[0.008]}, 
        AxesLabel -> {"Area", "Length"}, 
        GridLines -> Automatic];

FIGURE 6:  Parametric plot of  L(θ) versus Area(θ)  for  (0^°)≤ θ ≤ (90^°).

According to Figure 6 there is no projectile angle that equates the length to the area. However, there are limited choices of horizontal lines with ordinates greater than unity that will cross the curve at two different points.  The abscissas of the coordinates of the intersecting points are the values of the areas encompassed by the ordinates of the intersecting points.  Using the abscissas, we can identify the values of the corresponding initial projectile angles.  Similarly, but with unlimited choices, there are countless pairs of angles that correspond to the same value of the area with two different values of length; the vertical gridlines have these properties.  For instance, consider the first vertical gridline and solve A = 0.05 for the corresponding angles.

solθ = 
        NSolve[Area[[1]] == 
            0.05, θ];

The acceptable solutions of this equation are the fifth and the sixth roots.  These roots and the corresponding areas and the lengths are:

{{(θ /. 
              solθ[[5]]) (180/π), Area[[1]] /. θ -> θ /. 
          solθ[[5]], L /. θ -> θ /. 
          solθ[[5]]}, {(θ /. 
              solθ[[6]]) (180/π), Area[[1]] /. θ -> θ /. 
          solθ[[6]], L /. θ -> θ /. 
          solθ[[6]]}}

The parabolic trajectories associated with these two angles for projectiles with initial speed of 20 (m/s) are shown in Figure 7.

\(\(Show[{ParametricPlot[
            Evaluate[\(\(posnvector[t] /. 
                    Drop[initialVals, {3}]\) /. θ -> θ\) /. 
                solθ[[5]]], {t, 
              0, Evaluate[\(\(2\ \(v\_0\/g\) Sin[θ] /. 
                      Drop[initialVals, {3}]\) /. θ -> θ\) /. solθ[[5]]]}, 
            PlotStyle -> {Thickness[0.008]}, 
            AxesLabel -> {"x(m)", "y(m)"}, 
            GridLines -> Automatic, DisplayFunction -> Identity], 
          ParametricPlot[
            Evaluate[\(\(posnvector[t] /. 
                    Drop[initialVals, {3}]\) /. θ -> θ\) /. 
                solθ[[6]]], {t, 
              0, Evaluate[\(\(2\ \(v\_0\/g\) Sin[θ] /. 
                      Drop[initialVals, {3}]\) /. θ -> θ\) /. solθ[[6]]]}, 
            PlotStyle -> {Thickness[0.008]}, 
            AxesLabel -> {"x(m)", "y(m)"}, 
            GridLines -> Automatic, DisplayFunction -> Identity]}, 
        DisplayFunction -> $DisplayFunction];\)\)

FIGURE 7: Two typical trajectories with different lengths and the same encompassed areas.

The crossing coordinate of L versus A in Figure 6 has peculiar properties.  At this unique point, there are two special angles at which the corresponding pair of areas and lengths are the same.  To evaluate the value of these magic angles we solve the corresponding A  = 0.128 equation;  0.128 is the abscissa of the crossed point.  As in the previous case this equation has two acceptable solutions, the fifth and the sixth roots.

solmagicθ = 
        NSolve[Area[[1]] == 
            0.128, θ];

magicangles = {(θ /. 
              solmagicθ[[5]]) (180/π), (θ /. 
              solmagicθ[[6]]) (180/π)}

The corresponding areas and lengths of these two angles are

\(TableForm[
      Table[{\((θ /. 
                solmagicθ[[
                  n]])\)\ 180\/π, \(Area /. θ -> θ\) /. 
            solmagicθ[[
              n]], \(L /. θ -> θ\) /. 
            solmagicθ[[n]]}, {n, 
          5, 6}], TableHeadings -> {Automatic, {\*""Magic (Angle^°)"", "Area", "Length"}}]\)

The associated parabolic trajectories are shown in Figure 8.  These parabolas have the same lengths and the same areas.

\(\(Show[{ParametricPlot[
            Evaluate[\(\(posnvector[t] /. 
                    Drop[initialVals, {3}]\) /. θ -> θ\) /. 
                solmagicθ[[5]]], {t, 0, Evaluate[\(\(2\ \(v\_0\/g\) Sin[θ] /. 
                      Drop[initialVals, {3}]\) /. θ -> θ\) /. solmagicθ[[5]]]}, 
            PlotStyle -> {Thickness[0.01]}, 
            AxesLabel -> {"x(m)", "y(m)"}, 
            GridLines -> Automatic, DisplayFunction -> Identity], 
          ParametricPlot[
            Evaluate[\(\(posnvector[t] /. 
                    Drop[initialVals, {3}]\) /. θ -> θ\) /. 
                solmagicθ[[6]]], {t, 0, Evaluate[\(\(2\ \(v\_0\/g\) Sin[θ] /. 
                      Drop[initialVals, {3}]\) /. θ -> θ\) /. solmagicθ[[6]]]}, 
            PlotStyle -> {Thickness[0.01]}, 
            AxesLabel -> {"x(m)", "y(m)"}, 
            GridLines -> Automatic, DisplayFunction -> Identity]}, 
        DisplayFunction -> $DisplayFunction];\)\)

FIGURE 8: A pair of  unique trajectories with the same lengths and the same encompassed areas.

The length can also be written in terms of time of flight.  Since T = 2 ((Subscript[v, 0])/g)sin(θ), we replace sin(θ) in terms of T and because  the kinematic factor can no longer be suppressed, we evaluate the length for (Subscript[v, 0]) = 20 ((m )/s), which yields

length[
        T_] := (L /. {Cos[q_] :> (Sqrt[1 - Sin[q]^2])} /. 
              Sin[q_] -> (g/(2 (Subscript[v, 0]))) T) (((Subscript[v, 0])^2)/g) /. 
        Drop[initialVals, {3}]

The graph of the length versus time of flight is shown in Figure 9.

plotlengthT = 
        Plot[Evaluate[length[T]], {T, 0, 
            2 ((Subscript[v, 0])/g) /. Drop[initialVals, {3}]}, 
          PlotStyle -> {Thickness[0.008]}, 
          AxesLabel -> {"T,sec", "Length,m"}, 
          GridLines -> Automatic, 
          DisplayFunction -> Identity];

plotmagicθ = 
        ListPlot[{{2 ((Subscript[v, 0])/g) Sin[θ Degree] /. 
                  Drop[initialVals, {3}] /. θ -> 
                  magicangles[[1]], Evaluate[length[T]] /. T -> 2 ((Subscript[v, 0])/g) Sin[θ Degree] /. 
                  Drop[initialVals, {3}] /. θ -> 
                  magicangles[[1]]}, {2 ((Subscript[v, 0])/g) Sin[θ Degree] /. 
                  Drop[initialVals, {3}] /. θ -> 
                  magicangles[[2]], Evaluate[length[T]] /. T -> 2 ((Subscript[v, 0])/g) Sin[θ Degree] /. 
                  Drop[initialVals, {3}] /. θ -> 
                  magicangles[[2]]}}, 
          PlotStyle -> {PointSize[0.03]}, 
          DisplayFunction -> Identity];

Show[{plotlengthT, plotmagicθ}, 
        DisplayFunction -> $DisplayFunction];

FIGURE 9: The length versus time of flight for a 20 (m/s) projectile. The coordinates of the dots on the curve correspond to the magic angles.

According to Figure 9 and contrary to our intuitive expectation, the longest length does not correspond to the longest time.  Furthermore, the graph reveals that trajectories with identical lengths differ in time.  The dots in Figure 9 are two such cases corresponding to the pair of magic angles.

In closing, in a right-handed orthogonal  \({e\&^\_θ\),\(e\&^\_T\),\(e\&^\_L\)} coordinate system we construct a special 3D vector  (Overscript[R,  ->]) (θ) = {θ,T(θ),L(θ)}.  The components of this vector are the initial projectile angle θ, time of flight T and the projectile's length L. In this coordinate system (Overscript[R,  ->])  will trace a 3D curve.  In Figure 10 one such curve for (Subscript[v, 0])= 20 (m/s) is shown.  Three shadows of this curve on the three planes of the coordinate system are the discussed topics of interest.

R[θ_] := { θ,  
        2 ((Subscript[v, 0])/g) Sin[θ], (((Subscript[v, 0])^2)/g) L}

\(\(\(ParametricPlot3D[
        Evaluate[{Append[\(R[θ] /. θ -> θ\ Degree\) /. Drop[initialVals, {3}], {Thickness[
                  0.015]}], 
Append[\({\ θ, 
                    2\ \(v\_0\/g\) Sin[θ], 
                    0} /. θ -> θ\ Degree\) /. 
                Drop[initialVals, {3}], {GrayLevel[0.4], Thickness[0.01], 
                Dashing[{0.01, 0.08}]}], 
            Append[\({θ, 0, 
                    v\_0\^2\/g\ L} /. θ -> θ\ Degree\) /. 
                Drop[initialVals, {3}], {Thickness[0.01], 
                Dashing[{0.001, 0.02}]}], 
Append[\({0, 
                    2\ v\_0\/g\ Sin[θ], 
                    v\_0\^2\/g\ L} /. θ -> θ\ Degree\) /. 
                Drop[initialVals, {3}], {Thickness[0.01], 
                Dashing[{0.001, 0.02}]}]}], 
{θ, 0, 
          89.9\ }, 
        AxesLabel -> {"θ,Rad", "T,sec", "Length,m"}, 
        FaceGrids -> All, AspectRatio -> 1, 
        BoxRatios -> {1, 1, 1}, 
        BoxStyle -> {Thickness[0.007]}];\)\(\ \)\)\)

FIGURE 10:  The 3D curve (thick solid line) is the curve traced by vector (Overscript[R,  ->]) .  This curve has three shadows:  its shadow on the front plane is L versus θ (short-dashed line), its shadow on the left plane is L versus T (long-dashed line) and its shadow on the bottom plane is T versus θ (thin solid line).

REFERENCES

Halliday, D., Resnick, R., Walker, J., Fundamental of Physics, 5th ed. (Wiley and Sons Inc.), New York, 1997.
Sarafian, H., On Projectile Motion, The Physics Teacher, 37, 2 (1999) 86.
Tan, A., Extreme Problems in Projectile Motion, Theta 2, 2 (1988) 24.

About the Author

Haiduke Sarafian
The Pennsylvania State University
York,  PA 17403

Electronic Subscriptions

Included in the distribution for each electronic subscription is the file sarafianAngles.nb containing Mathematica code for the material described in this article.