Projectile Motion: Finding the Optimal Launch Angle - Report

Pro jectile Motion: Finding the Optimal Launch Angle 1 Abstract If we want to throw a projectile as far as possible, at what angle should it be launched? This paper focuses on how the answer to this question changes depending on the situation. We look at launching projectiles onto di↵erently shaped hills, as well as how varying initial velocity and height a↵ect the launch angle. Finally, we add air resistance to the projectile problem and compare two di↵erent models: air resistance proportional to the projectile’s velocity and air resistance proportional to velocity squared. 2 Contents 1 Introduction 4 2 The projectile problem 4 3 Equations of motion: no air resistance 5 4 The optimal launch angle 6 4.1 The distance function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 4.2 The enveloping parabola . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 4.2.1 Derivation of the enveloping parabola: height maximization . . . . . 8 4.2.2 Derivation of the enveloping parabola: expanding circles . . . . . . . 9 4.3 The projectile problem solution . . . . . . . . . . . . . . . . . . . . . . . . . 12 5 Varying the impact function 12 5.1 The linear impact function . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 5.2 The parabolic impact function 5.3 The semicircular impact function . . . . . . . . . . . . . . . . . . . . . . . . 15 5.4 The sinusoidal impact function . . . . . . . . . . . . . . . . . . . . . . . . . 17 5.4.1 5.5 . . . . . . . . . . . . . . . . . . . . . . . . . 14 Newton’s method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 Varying initial conditions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 6 Air resistance 24 6.1 Equations of motion: linear air resistance . . . . . . . . . . . . . . . . . . . 25 6.2 The level ground impact function . . . . . . . . . . . . . . . . . . . . . . . . 27 6.3 The parabolic impact function 6.4 Quadratic air resistance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 6.4.1 6.5 . . . . . . . . . . . . . . . . . . . . . . . . . 29 The numerical solution . . . . . . . . . . . . . . . . . . . . . . . . . . 31 The drag coefficient . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 7 Conclusion 36 Index 38 3 1 Introduction In this paper, we examine how to find the optimal launch angle, which is the angle at which a projectile is launched that maximizes its horizontal distance traveled. We find that this angle depends on numerous factors, including the projectile’s initial velocity, the e↵ects of air resistance, and the surface upon which the projectile lands. This paper addresses these relationships in three parts: finding a general solution for the optimal launch angle, exploring di↵erent landing surfaces, and adding the e↵ects of air resistance. In the first four sections, we derive a solution to the projectile problem by considering the projectile’s equations of motion. We introduce the importance of the enveloping parabola and derive its equation in two ways. In the fifth section we examine specific examples of the projectile landing on di↵erently shaped hills, specifically linear, parabolic, semicircular and sinusoidal. We also explore the dependence of the optimal launch angle on the projectile’s initial height and velocity. Finally, in the sixth section we add air resistance to the problem, examining both the linear and the quadratic model. We compare these two models and aim to understand why certain models work better than others in a given situation. 2 The projectile problem We define the projectile problem as follows: a projectile is launched from a tower of height h, with initial velocity v, and at an angle ✓ measured with respect to the horizontal. We aim to find ✓m , the launch angle that maximizes horizontal distance. The projectile is launched onto a hill, which is defined by the function , called the impact function. The impact function varies depending on the shape of the surface we want to explore, and in this paper we will look at -functions that are linear, parabolic, semicircular and sinusoidal. Figure 1 shows a typical setup for the projectile problem. 4 Vertical distance p(x) v (x) h Horizontal distance Figure 1: The projectile problem. 3 Equations of motion: no air resistance We first consider the situation of a projectile launched from a tower of height h onto some impact function , ignoring the e↵ect of air resistance. In order to solve for ✓m , we need to find equations for motion in the x- and y-directions. We define ✓ to be the angle above the horizontal at which the projectile is launched. The projectile is launched with an initial velocity v, which has magnitude v, and when broken up into x- and y-components, gives us the initial conditions x(0) = 0; y(0) = h; x0 (0) = v cos ✓; y 0 (0) = v sin ✓. Without air resistance, acceleration in the x-direction is zero, while in the y-direction it is solely due to gravity, where g = 9.8 m/s2 . Thus we can solve the second-order di↵erential equations to find our two motion equations. In each step, we integrate both sides of the 5 equation with respect to t and apply our initial conditions. In the x-direction, we have x00 (t) = 0; x0 (t) = v cos ✓; x(t) = vt cos ✓. (1) The motion in the y-direction is described by y 00 (t) = g; y 0 (t) = y(t) = gt + v sin ✓; 1 2 gt + vt sin ✓ + h. 2 (2) We now have a set of parametric equations for the motion of the projectile as a function of t, but to maximize the projectile’s horizontal distance, we want to find a path function, p, that defines the projectile’s height as a function of horizontal distance, x. Solving for t in (1) and substituting into (2) yields t= x , v cos ✓ and therefore x ⌘ 1 ⇣ x ⌘2 g v cos ✓ 2 v cos ✓ gx2 = h + x tan ✓ sec2 ✓. 2v 2 p(x) = h + v sin ✓ ⇣ (3) We now have one equation that describes the motion of the projectile, which is useful in finding the launch angle that maximizes x. 4 The optimal launch angle Next we will explore the process of finding ✓m , the projectile’s optimal launch an- gle. This angle will vary depending on the height and velocity at which the projectile is 6 launched, as well as what type of surface the projectile lands on. The impact function, is a function of horizontal distance, x. We find that for certain closed-form solution for ✓m is possible, while for other , functions, obtaining a functions, we must instead turn to approximations of ✓m . 4.1 The distance function For each value of ✓, there is a value of x where the projectile hits the hill, which occurs when p(x) = (x). We call this x-value d(✓) since it varies depending on the launch angle. By maximizing d(✓), we can find the angle that maximizes the projectile’s horizontal distance. We note that for every xi within the projectile’s horizontal range, we can find a launch angle corresponding to at least one projectile that hits that xi . This fact, along with the Implicit Function Theorem (see [11]), lets us assume d(✓) is a di↵erentiable function. With this condition, we can use implicit di↵erentiation on our distance function, p, to maximize d(✓) [3]. We have (d(✓)) = p(d(✓)) = h + d(✓) tan ✓ d(✓)2 g sec2 ✓. 2v 2 Di↵erentiation gives 0 (d(✓))d0 (✓) = d(✓) sec2 ✓ + d0 (✓) tan ✓ g d(✓)2 sec ✓(sec ✓ tan ✓) + d(✓)d0 (✓) sec2 ✓ . v2 Since d0 (✓m ) = 0, we find that g d(✓m )2 sec2 ✓m tan ✓m v2 0 = d(✓m ) sec2 ✓m d(✓m ) = v2 cot ✓m , g which tells us that the maximum horizontal distance the projectile travels is dependent upon initial velocity and gravity. We now have a value for horizontal distance in terms of the launch angle ✓. This value is independent of the impact function , but the angle that maximizes x will vary for di↵erent -functions, a connection we explore in Section 5. 7 4.2 The enveloping parabola An enveloping parabola is a path that encloses and intersects all possible projectile paths. For each value of ✓ in [ ⇡, ⇡], the enveloping parabola intersects the projectile path at exactly one point, and at this point the two functions share a common tangent line (see [1]). Figure 2 shows the enveloping parabola that intersects each possible path at exactly one point, where a possible path corresponds to a unique launch angle in the range [ ⇡, ⇡]. We derive the equation for the enveloping parabola in two ways, and show that both methods yield the same answer. (x) Figure 2: The enveloping parabola intersects each possible projectile path at one point. 4.2.1 Derivation of the enveloping parabola: height maximization We first derive the enveloping parabola by maximizing the height of the projectile for a given horizontal distance x, which will give us the path that encloses all possible paths. In Section 3, we derived the path of the projectile for a given launch angle ✓ to be gx2 (1 + tan2 ✓). 2v 2 y = h + x tan ✓ To simplify this equation, we let u = tan ✓, gx2 (1 + u2 ). 2v 2 y = h + ux 8 (4) Now, given any x within the projectile’s horizontal range, we can maximize y as a function of u (see [1]). Taking the derivative of (4) with respect to u and setting y 0 = 0, we solve for u: y0 = x 0=x u= gx2 (2u2 ) 2v 2 gx2 u v2 v2 . gx (5) This is the u-value at which y is maximized given a fixed x. Substituting (5) into (4), we have y =h+ v2 x gx ⇣ v 2 ⌘2 ⌘ gx2 ⇣ v2 1 + = h + 2v 2 gx 2g gx2 . 2v 2 (6) This is the equation for the enveloping parabola, which, as shown in Figure 3, is the function that encloses and intersects all possible projectile paths. (x) Figure 3: The enveloping parabola (black) encloses all possible parabolic paths for a constant initial height and velocity. 4.2.2 Derivation of the enveloping parabola: expanding circles Next we consider a derivation of the enveloping parabola that involves expanding circles. Instead of examining individual projectile paths, we focus instead on circles that are composed of every possible projectile’s position at a time t. These circles have a radius 9 of vt, where v is the initial velocity, since each possible projectile travels this distance in time t. For each value of ✓, there is a unique projectile path that corresponds to a point on this circle. Now that we have the radius of the circles, we need to add in the component of gravity acting in the negative y-direction. The y-position of the center of each circle, yc , is given by yc = y0 + v 0 t 1 2 gt . 2 (7) In our case, the center of each circle is not initially moving, so v0 = 0, and each circle is centered at (0, h), so y0 = h. Thus the center of each circle will be moving in the negative y-direction at a rate of h 1 2 2 gt (see [1]). Adding in this shift for the y-values of each circle, we find the equation of each expanding circle at a given time t is x2 + (y 1 h + gt2 )2 = v 2 t2 . 2 (8) Figure 4 shows some of these circles for varying values of t. Figure 4: Each circle is composed of every possible projectile’s position at a time t, and the enveloping parabola encloses these circles. Note that for small values of t, the circles appear to be centered at the origin, but for larger values of t, the center falls with gravity. Again, these circles represent the positions of all possible projectiles at each given t-value. By enclosing these circles, we are enclosing all projectile paths. We expect each circle to intersect the enveloping parabola at a unique 10 time t, since each point on the enveloping parabola corresponds to exactly one point from one projectile path. The equation for each expanding circle is a quadratic for the variable t2 , so solving this equation for t2 and taking the positive square root, we can find either two, one or zero real solutions for t. Points that lie within the enveloping parabola correspond to two real solutions, since, as shown in Figure 4, each point is crossed by two circles with two di↵erent t-values. Points that lie on the enveloping parabola correspond to a single solution for t, since each point is unique to one projectile path. Since the quadratic equation gives us two solutions when the discriminant is nonzero, we will find a single solution for t when the discriminant of (8) is equal to zero. We have x2 + (y 1 h)gt2 + g 2 t4 = v 2 t2 4 h)2 + (y 1 2 4 g t + ((y 4 h)g v 2 )t2 + x2 + (y ((y h)g v 2 )2 h)2 = 0. The discriminant is ((y h)g)2 2(y 4 h)gv 2 + v 4 ⇣1⌘ 4 (x2 + (y g 2 x2 g 2 (y 2ygv 2 + 2hgv 2 + v 4 h)2 ) = 0 h)2 = g 2 x2 = and solving for y gives 1 2hgv 2 v 4 + g 2 x2 2gv 2 v2 gx2 =h+ . 2g 2v 2 y= (9) Note that (9) is the same equation for the enveloping parabola that we found in (6) using the height maximization method. 11 4.3 The projectile problem solution The enveloping parabola has important significance for solving the projectile problem. Let the enveloping parabola be defined by the function (x) = h + Since v2 2g g 2 x . 2v 2 (10) encloses all possible projectile paths, maximizing the horizontal distance, x, for is equivalent to maximizing all projectile paths. Thus we want to find the point at which the enveloping parabola intersects the impact function , and then find the ✓-value that corresponds to this point on the enveloping parabola. This ✓-value will be the optimal launch angle, ✓m . Recall that maximizing horizontal distance yielded the equation d(✓m ) = v2 cot ✓m , g (11) so the solution to the projectile problem requires first finding the x-value where (c) = (c), then solving for ✓m in (11). 5 Varying the impact function We now apply the method derived in the last section to examples of specific impact functions. Specifically, we explore linear, parabolic, semicircular, and sinusoidal impact functions. We also analyze how changing initial conditions, like h and v, a↵ect the optimal launch angle. 5.1 The linear impact function Consider the linear impact function of the form (x) = mx, where m is a positive constant. Recall that finding the optimal launch angle involves first maximizing horizontal distance by solving for the x-value at which the enveloping parabola and the impact 12 function intersect. We want to find this value c such that (c) = (c). So we have v2 2g g 2 c 2v 2 ⇣ g v2 ⌘ 0 = 2 c2 + mc h+ . 2v 2g mc = h + The solutions to this quadratic equation are c= mv 2 v 2 ± g g r m2 + 2gh + 1. v2 We know the point c = d(✓m ) corresponds to the intersection between the enveloping parabola and the projectile path with launch angle ✓m . This intersection point has the value c = d(✓m ) = g/v 2 cot ✓m , so we have g cot ✓m = v2 and therefore ✓m mv 2 v 2 ± g g r m2 + 2gh + 1, v2 ⇣ mv 4 v 4 r ⌘ 2gh 2+ = arccot ± m + 1 . g2 g2 v2 (12) We find that the angle that lies in the first quadrant indeed maximizes horizontal distance of the projectile by graphing this path alongside the enveloping parabola. (x) p(x) (x) = mx Figure 5: The optimal projectile path and the enveloping parabola plotted with initial height 10 m, initial velocity 20 m/s and impact function (x) = 0.4x. 13 Figure 5 shows that the two paths intersect at one point, and this point is in fact the intersection of three functions – the enveloping parabola, the impact function and the optimal projectile path. Since all paths meet at one point, we know we have found the launch angle for which the projectile’s horizontal path distance reaches a maximum. 5.2 The parabolic impact function Next we consider an impact function of the form (x) = ax2 , where a is a positive constant. As we did in the previous example, we solve the equation (c) = (c) for c. We have ac2 = h + 2 c = c= h+ a s v2 2g g 2 c 2v 2 v2 2g g 2v 2 = 2hgv 2 + v 4 2agv 2 + g 2 2hgv 2 + v 4 . 2agv 2 + g 2 Recall c = v 2 /g cot ✓m , so we find that the optimal initial angle, ✓m , is ✓m g = arccot 2 v s = arccot s 2hgv 2 + v 4 2agv 2 + g 2 ! 2hg 2 + gv 2 . 2av 4 + gv 2 ! (13) We now have a solution for optimal launch angle as a function of the projectile’s initial height and velocity. Figure 6 shows the path of the projectile using this angle. Note that the projectile path, the enveloping parabola, and the impact function all intersect at a single point, which implies we have found the maximized path. 14 (x) p(x) (x)= ax2 Figure 6: The projectile path p intersects the enveloping parabola tion at a single point. 5.3 and the impact func- The semicircular impact function The next example we consider is an impact function that is a semicircle. Let us define p the impact function as (x) = r2 x2 , which is a semicircle that is concave down and centered at the origin. In this example, we still define h to be the vertical distance from the origin to the launch point. The closed-form solution exists but is more complicated than previous examples, as we must square both sides of the equation (c) = (c) in order to isolate c. We have p r2 r2 v2 g 2 c 2g 2v 2 ⇣ v2 g 2 ⌘2 c2 = h + c 2g 2v 2 ⇣ 1 gh ⌘ g2 0 = 4 c4 + c 2 + h2 4v 2 v2 c2 = h + r2 + v2h v4 + 2. g 4g The quadratic formula, once simplified, gives us 4hv 2 c2 = g and therefore c= 4hv 2 g v4 ±2 g2 s v 4 r2 g2 s v4 v 4 r2 ± 2 g2 g2 15 2hv 6 , g3 2hv 6 g3 !1/2 . (14) The semicircular impact function has a property we have not encountered before – there are four possible values for c, two of which are negative x-values, which we can disregard for our situation, but we are left with two intersections of the parabolic path with the impact function. Figure 7 shows an example of this situation, where if we choose the second point of intersection, and its corresponding projectile path, the projectile does not actually reach the second point and our solution would not be the maximized path. Therefore we must choose the smaller of the two positive c-values to avoid this issue. Figure 7: For a semicircle, we run into the issue of two intersections between the enveloping parabola and the impact function. Furthermore, we can find the radius of the semicircle for which the projectile hits at exactly one point. In this case we want there to be only two solutions for c: a positive and a negative value. Hence we must solve for r when the discriminant in (14) is zero. We find that g 2 r2 v4 2gh =0 v2 s 2hv 2 r= . g (15) As shown in Figure 8, the semicircle with this radius intersects the enveloping parabola at exactly one point. For any radius greater than r in (15), the enveloping parabola will intersect the semicircle twice, and for any r less than this value, the enveloping parabola will not hit the impact function. 16 Figure 8: A semicircle with the radius in (15) results in one intersection between the enveloping parabola and the impact function. 5.4 The sinusoidal impact function Next we consider an example of an impact function that does not have a closedform solution: (x) = sin x. We are still able to approximate ✓m by solving for the intersection between the enveloping parabola and the impact function. For a sinusoidal impact function, we have (c) = (c) sin c = h + v2 2g g 2 c . 2v 2 Note that this equation does not have a closed-form solution for c, and this we must turn to a di↵erent method to approximate the optimal launch angle. 5.4.1 Newton’s method Newton’s method is used to approximate the roots of functions. In this case, we apply Newton’s method to the function f (x) = (x) (x), since the zero of f corresponds to the intersection between the enveloping parabola and the impact function (see [12]). Newton’s method requires the following steps: 1. Make an initial x-value guess, x1 , as to where f (x) = 0; 2. Find the tangent line at the point (x1 , f (x1 )); 17 3. Solve for the x-value, x2 , at which this tangent line intersects the x-axis; 4. Repeat steps 2 and 3 using the values x2 , x3 , x4 , . . . , xn . Carrying out step 2, we find that the tangent line at the point (x1 , f (x1 )) takes the form y f (x1 ) = f 0 (x1 )(x x1 ). As stated in step 3, we next find the root of this tangent line, which is done by setting y = 0. Letting x = x2 when y = 0, we have f (x1 ) = f 0 (x1 )(x2 x2 = x1 x1 ) f (x1 ) . f 0 (x1 ) As we repeat steps 2 and 3, the x-values approach a more exact approximation of the root of f . We can generalize our result above to the nth iteration, and thus we have xn+1 = xn Recall that f (x) = (x) f (xn ) . f 0 (xn ) (x), so we have f (x) = sin x (h + f 0 (x) = cos x + g x. v2 v2 2g g 2 x ) 2v 2 Therefore, our general equation for the nth approximation of the root of f is xn+1 = xn g v2 2 2g + 2v 2 xn . cos xn + vg2 xn sin xn h Figure 9 shows Newton’s method in graph form, with initial guess x1 and the process to find x2 . 18 (x 1, f(x 1 )) Enveloping parabola ( ) End goal: approximate this x-value x2 x 1 (guess) Impact function ( ) f(x) = (x) - (x) Tangent line to (x 1 ,f(x 1 )) Figure 9: The first iteration of Newton’s method, where we aim to approximate the zero of f , where and intersect. There are some constraints on our initial guess for where the root occurs. We have to be careful when we make this guess, for if we choose a point where the derivative of has a di↵erent sign than at the zero, our method of using tangent lines will not lead to x-values that approach our zero. Instead the values might oscillate or diverge infinitely and thus Newton’s method will fail. To avoid this issue, we must isolate our x-value guess ⇡ ⇡ 3⇡ to the ⇡-length interval [ 3⇡ 2 + 2⇡n, 2 + 2⇡n] or [ 2 + 2⇡n, 2 + 2⇡n] as shown in Figure 10, where we know the derivative will not change sign. 19 xmin xmax Figure 10: We must narrow down the range for the initial guess. To find the interval specific to the intersection we want to approximate, we first isolate the window in which the enveloping parabola ranges between 1 and 1. We know the intersection of the enveloping parabola and the sine curve occurs in this window since sin x has a maximum of 1 and a minimum of 1. Solving for c when = 1 and = 1 gives us a narrowed range in which the intersection occurs. We then adjust this interval to the ⇡ nearest [ 3⇡ 2 + 2⇡n, 2 + 2⇡n] interval, since as we previously found, this interval guarantees we can apply Newton’s method. Once we have found this interval, we can make our initial guess. (x) p(x) (x) Figure 11: The maximized projectile path onto a sinusoidal impact function. Now that we have our guess, we can apply Newton’s method to approximate the intersection. From there we can find the angle ✓m at which the projectile must be launched 20 to hit this intersection. Figure 11 shows this optimal path, where the launch angle is the output from five Newton’s method iterations. One problem that we come across when finding the intersection is the possibility of multiple intersections between the enveloping parabola and the impact function. Recall that we saw this issue in the semicircular impact function example, where choosing any intersection but the first one results in a ”wrong” solution for ✓m . Figure 12 shows the impact function (x) = 6 sin x, where the sinusoidal nature of the impact function allows for multiple intersections and therefore we must be careful when we narrow in on the range to apply Newton’s method. (x) (x) Figure 12: In this case, the enveloping parabola intersects the impact function three times. Physically, it only makes sense to choose the first intersection, since the projectile will hit the hill at this point and does not reach the other intersection points. For this reason, when deciding where to make our initial guess for Newton’s method, we look to guess in the interval that contains the first intersection. 5.5 Varying initial conditions We now explore the e↵ect of changing the initial height and velocity of the projectile. Since the maximized projectile path is dependent on both height and initial velocity, we expect the initial angle for maximum distance to vary as well. Looking at the parabolic impact function example, recall the optimal angle is given by ✓m = arccot s 21 2hg 2 + gv 2 2av 4 + gv 2 ! . (16) We see that height h is in the numerator of the cotangent argument. Recalling the graph of cot ✓, we know that on the interval (0, ⇡), cot ✓ decreases as ✓ increases. Thus we would expect ✓m to decrease as initial height increases. 35 a = 0.1 a = 0.05 a = 0.2 a = 0.5 1.5 30 Initial launch angle 25 y 20 15 10 5 0 1.0 0.5 (x) = 0.1x 2 0 2 4 6 8 10 12 0.0 14 0 5 10 15 20 25 30 Initial height x Figure 14: Optimal launch angle ✓m varies inversely with initial height h, and shows the same trends for di↵ering values of a in the impact function (x) = ax2 . Figure 13: Maximized projectile paths for varying heights between 0 and 30 m. Figure 13 shows a plot of maximized projectile paths onto a parabolic impact function for varying initial heights. For each initial h in Figure 13, the projectile travels the farthest horizontal distance before intersecting the impact function. Figure 14 shows a plot of each path’s initial angle ✓m that maximizes horizontal distance for a given h, where the di↵erent sets of points correspond to di↵erent a-values (recall a is the coefficient for the parabolic impact function). As we expected, ✓m decreases as h increases. Additionally a steeper parabolic impact function corresponds to an overall larger optimal launch angle, as indicated in Figure 14, where di↵erent colored points correspond to di↵erent coefficients for the parabolic impact function. 22 a = .05 a = .01 a = .001 15 Launch angle 0.8 10 0.6 0.4 5 0.2 (x) = 0.1x 2 0 0 5 10 15 20 25 0.0 30 0 Figure 15: Maximized projectile paths for varying initial velocities between 1 and 15 m/s. 2 4 6 8 10 12 14 Initial velocity Figure 16: Optimal launch angle ✓m varies directly with initial velocity v, and shows the same trends for di↵ering values of a in the impact function (x) = ax2 . Varying the initial velocity leads to opposite results – an increase in v leads to an increase in ✓m . This again is expected by the equation for the optimal launch angle, ✓m = arccot s 2hg 2 + gv 2 2av 4 + gv 2 ! , as we see there is a greater degree of v in the denominator. This means that for increasing v-values, the argument of cotangent decreases and thus ✓m increases. Figure 15 shows the optimal paths for varying velocities and a parabolic impact function. Additionally, changing the coefficient a of the parabola results in a more linear relationship for lower a-values and a more parabolic relationship between v and ✓m for larger a-values. This is shown in Figure 16. For di↵erent impact functions, the relationship between initial height and optimal angle varies, but for our three examples, the two quantities are always inversely related. Figure 17 shows the results of changing initial height for di↵erent impact functions. 23 1.2 Linear Parabolic Sinusoidal 1.0 Optimal launch angle 0.8 0.6 0.4 0.2 0.0 0 5 10 15 20 Initial height Figure 17: We see a decrease in optimal launch angle for increasing initial heights. 6 Air resistance In order to make the projectile problem more realistic, we consider the e↵ects of air resistance. Air resistance is a force, called the drag force, that acts in the direction opposite an object’s motion. Air resistance takes the form of a Taylor series where the terms are powers of the projectile’s velocity. Overall, two of the terms tend to be much larger than all the other terms: the v and the v 2 terms. For this reason, we consider two main types of drag force: linear and quadratic. These forces take the form flin = bv, fquad = cv 2 , where b and c are constants with b = D and c = D2 , where D is the diameter of the object and and are constants that depend on the nature of the medium (see [10]). An object subject to air resistance is modeled by a linear combination of both linear and quadratic drag forces, but in many cases, one term is much larger than the other, indicating the situation can be modeled as either linear or quadratic. In order to understand which drag force will be most realistic for the problem at hand, we estimate the ratio fquad /flin , which tells us if one model of air resistance can be neglected. For our projectile problem, let the projectile be similar to a baseball. Suppose the projectile moves with velocity 20 24 m/s and diameter 8 cm, and assume the projectile is traveling in air at standard pressure = 0.25 Ns2 /m4 , and temperature. At these values, = 1.6 ⇥ 10 4 Ns/m2 (see [10]). We find that fquad ⇡ 2500, flin which indicates that our projectile’s motion is better modeled by quadratic drag force. An example of an object that experiences more linear drag force would be a small drop of oil. In this case, let the droplet of oil have a diameter of 15 mm and move at a rate of 5⇥10 5 m/s. We find that the ratio in this case is fquad ⇡ 7.8 ⇥ 10 flin 6 , which indeed indicates that a linear model would be correct. When solving the projectile problem, we first consider linear air resistance in order to show how air resistance is included in projectile motion. We choose this model because mathematically it is simpler and provides insight into modeling air resistance. Quadratic air resistance is a more realistic model for our situation, but we will see that the equations of motion in this case are not solvable in a closed form using di↵erential equations techniques. We therefore turn to computer approximations to model this motion, but first let us further explore the linear model. 6.1 Equations of motion: linear air resistance Recall that when we explored the case without air resistance, we solved second order di↵erential equations for the x- and y-motion of the projectile as functions of time, t. The method here is the same, except we must consider a force other than gravity acting on the projectile. Drag force acts in the direction opposing motion, so both the x- and yequations will have a component of air resistance. Our second order di↵erential equations will be the x- and y-components of acceleration, which we find from Newton’s 2nd law. Bolded values indicate a vector, meaning the variable has both a magnitude and direction. 25 We define the drag force to be FD and the gravitational force is Fg . We have ma = F = Fg + FD = mgŷ b(x̂ + ŷ), and letting k = b/m, we can separate the above equation into x- and y-equations. We have x00 (t) = kx0 (t), y 00 (t) = g ky 0 (t). Next we will solve the above di↵erential equations using the initial conditions: x(0) = 0; y(0) = h; x0 (0) = v cos ✓; y 0 (0) = v sin ✓; where v is the initial velocity of the projectile. Using separation of variables to solve the x-equation, we obtain x00 (t) = kx0 (t), x0 (t) = Ce x(t) = kt = v cos ✓e v cos ✓ e k kt kt +C = v cos ✓ ⇣ 1 k Similarly, for the y-equation, we have y 00 (t) = g ky 0 (t) dy 0 = ( g ky 0 )dt dy 0 = dt g ky 0 1 ln(g + ky 0 ) = k t+C g + ky 0 = Ce 26 kt . e kt ⌘ . We solve for C using the initial condition y 0 (0) = v sin ✓, g k(v sin ✓) = C. Therefore, g ky 0 = ( g y0 = kv sin ✓)e kt ) 1 ( g + (g + kv sin ✓)e k kt ). Using the initial condition y(0) = h, we integrate once more with respect to t to find the motion equation for y: ⇣g ⌘ g t+ + v sin ✓e kt + C k k 1 1 y = h + (v sin ✓ gt) + 2 (g e k k y= kt (g + kv sin ✓)). Now that we have the x- and y-equations of motion, we want to eliminate the variable t, by solving the x-equation for t and substituting this t-equation into the y-equation. We have ⌘ v cos ✓ ⇣ 1 e kt k ⌘ 1 ⇣ xk t= ln +1 , k v cos ✓ x(t) = which gives us y =h+ ⇣ xk ⌘⌘ 1⇣ 1⇣ v sin ✓ + g ln +1 + 2 g k v cos ✓ k (g + kv sin ✓) ⇣ ⌘⌘ xk +1 . v cos ✓ (17) We now have a function for the projectile’s path in terms of x. 6.2 The level ground impact function Previously, before we included air resistance, our next step would have been to find the enveloping parabola of all possible projectile paths. However, since the projectile paths are 27 no longer parabolic, this method does not work. Instead we explore a di↵erent method; we maximize x at the intersection of the path, represented by (17), and the impact function, . To illustrate this method, consider the case of the horizontal impact function y = 0. We aim to maximize x when y = 0. In other words, we want to find how far the projectile travels horizontally before it hits the ground. The algebra is a bit dense, but the method is this: we simplify our y-equation, set y = 0, and then take the derivative with respect to ✓, since we want to maximize x as a function of the launch angle. Finally we let x0 = 0 and solve for x. We have ⇣ xk ⌘ ⇣ xk ⌘⌘ 1⇣ g ln + 1 + g + kv sin ✓ (g + kv sin ✓) + 1 k2 v cos ✓ v cos ✓ ⇣ v cos ✓ xk ⌘ ⌘ 1⇣ gxk = h + 2 g ln + + xk 2 tan ✓ k v cos ✓ v cos ✓ ⇣ ⇣ ⇣ ⌘ ⌘⌘ d d 1 v cos ✓ xk gxk (0) = h + 2 g ln + + xk 2 tan ✓ d✓ d✓ k v cos ✓ v cos ✓ ⇣ v cos ✓ ⌘⇣ v cos ✓( v sin ✓) + v sin ✓(v cos ✓ xk) ⌘ gxkv sin ✓ =g + 2 + xk 2 sec2 ✓ v cos ✓ xk v 2 cos2 ✓ v cos2 ✓ 0=h+ = vgxk sin ✓ g sin ✓ + kv + . v cos ✓(v cos ✓ xk) cos ✓ Solving for x, we find x = d(✓m ) = v 2 cos ✓m . g sin ✓m + kv (18) This equation for horizontal distance in terms of ✓m represents the farthest distance a projectile will travel when launched onto the line y = 0. Kantrowitz and Neumann present a similar method for this solution involving the Lambert W function, which deals with solving equations with both exponential and linear parts [9]. We found (18) using a di↵erent method, and it is helpful to have found the same solution. Note that our solution is not complete however, since our end goal is to find the launch angle ✓m that maximizes the distance found in (18). However, solving for ✓m in (18) requires introducing the Lambert W function or approximating a solution using a computer, which could be explored in future work but we do not address it here. 28 6.3 The parabolic impact function Let us consider a parabolic impact function of the form (x) = ax2 , with the goal of finding the equation for maximized x in terms of ✓. We begin by setting our impact function equal to the path equation: ax2 = h + ⇣ xk ⌘⌘ 1⇣ 1⇣ v sin ✓ + g ln +1 + 2 g k v cos ✓ k (g + kv sin ✓) ⇣ ⌘⌘ xk +1 . v cos ✓ (19) We di↵erentiate with respect to ✓, set x0 = 0, and solve for x. The presence of an x in every term allows for convenient cancellation that leaves just a single x term in the denominator of the first term on the right side of the equation: ⇣ v cos ✓ xk ⌘ ⌘⌘ d d⇣ 1⇣ gxk (ax2 ) = h + 2 g ln + + xk 2 tan ✓ d✓ d✓ k v cos ✓ v cos ✓ 2axk 2 = 2ak = gxvk sin ✓ gxkv sin ✓ + xk 2 v 2 sin2 ✓ 2 + xk + v 2 cos2 ✓ xkv cos ✓ v 2 cos2 ✓ g sin ✓ g sin ✓ + xkv sin2 ✓ + k + . v cos2 ✓ xk cos ✓ v cos2 ✓ Solving for x gives xmax = v 2 cos ✓m (1 2a cos2 ✓m ) . kv 2akv cos2 ✓m + g sin ✓m (20) Again, this is the equation for maximum distance a projectile will travel when launched at its optimal angle, and we would need di↵erent techniques to find ✓m . 6.4 Quadratic air resistance We now turn to the much more complicated situation of a projectile that is subject to quadratic air resistance. This drag force, FD , has a magnitude of cv 2 , where v is the projectile’s speed and c is a constant. As we did in the case of linear air resistance, we solve equations of motion for the projectile in the x- and y-directions. Beginning with Newton’s 2nd Law, we have F = ma, 29 where again a bold letter indicates a vector quantity. This is important as we next write the left hand side as a sum of the forces acting on the projectile. In this case, we have gravity, acting in the downward ( ŷ) direction, and drag force, acting in the opposite direction of velocity (-v̂). So we have cv 2 v̂ = ma. mgŷ (21) From the definition of a unit vector, we know v̂ = v , v and thus (21) becomes mgŷ cvv = ma. For our situation, v has both a horizontal and a vertical component, so we can write v = vx x̂ + vy ŷ, where vx is the x-component and vy is the y-component of velocity. This gives us mgŷ cv(vx x̂ + vy ŷ) = ma. We also know v, the speed of the projectile, is the magnitude of the velocity vector and q thus we have v = vx2 + vy2 . Our equation becomes mgŷ c q vx2 + vy2 (vx x̂ + vy ŷ) = ma, as shown in [10]. We can break this equation into x- and y-equations of motion, max = mx00 (t) = may = my 00 (t) = q vx2 + vy2 , q mg cvy vx2 + vy2 . cvx 30 (22) (23) Recall that at this point in the linear air resistance case, we solved the two second order di↵erential equations for x(t) and y(t), where our set of equations were x00 (t) = kx0 (t), y 00 (t) = g ky 0 (t). (24) (25) Equations (24) and (25) are linear and we are able to solve them by the separation of variables. However, equations (22) and (23) are not linear; the vx2 and vy2 terms increase the complexity of the di↵erential equations. We also run into issues because they are coupled equations, meaning that in order to solve for vx we need to know vy , or vice versa. The extra v term in (21) results in the coupled terms, where as in the linear case we did not have this problem. Physically, this means that the motion of the projectile cannot be interpreted in the x- and y-directions separately. The linear case allowed us to make assumptions about the projectile’s motion that simplified the solutions, but here the xand y-positions of the projectile depend on both the vertical and horizontal components of velocity and hence cannot be separated. These coupled equations do not have analytic solutions when we consider a projectile that has both x- and y-components of motion (see [9]). There are of course special cases that we are able to solve. These include when a projectile has strictly vertical motion, as in a rock dropped straight o↵ of a cli↵, or strictly horizontal motion, as in a bullet shot from a gun, where gravity’s e↵ects are negligible in the initial stage of motion. In these cases, there is only one component of velocity and only one equation to solve, with which we can use separation of variables techniques. However, the problem we have been looking at has two directions of motion, and thus we cannot solve these equations analytically, but we can approximate their solutions with a computer algebra system. 6.4.1 The numerical solution We will not be able to find closed-form solutions to (24) and (25), so we turn to Mathematica to model the motion of the projectile subject to quadratic air resistance. The NDSolve function in Mathematica finds solutions to di↵erential equations by approximat31 ing the solutions on small subintervals, and is described further in [6]. The output of NDSolve is numerical on an interval, rather than a closed-form function. We can then plot these numerical solutions for x and y and acquire data points at certain values of t. Fitting a curve to these data points will give us the approximation for the projectile’s motion. Below we work through an example. Let us look at a projectile launched from a height of 10 m at an angle of ✓ = 50 with an initial velocity of 10 m/s. This means vx (0) = 10 cos(50 ) and vy (0) = 10 sin(50 ). We will use the constant k = c/m as our drag coefficient, and in this case, let k = 0.008 as this is the constant for a projectile resembling a baseball. Using Mathematica, we can plot the x- and y-solutions to the set of di↵erential equations in this example. We plot x and y as a function of time, t, but cannot plot them parametrically since they do not have an analytic form. We can extract data points from the plot of x and y at the same values t and plot these (x, y) points in Figure 18. 14 12 10 8 6 4 2 0 0 5 10 15 20 Figure 18: Plot of x- and y-positions of the projectile at increasing t-values. We are then able to fit a curve to these data points and compare the path a projectile would take depending on if it is modeled by no air resistance or quadratic air resistance. 32 14 No air resistance Quadratic air resistance 12 10 8 6 4 2 0 0 5 10 15 20 Figure 19: Two models of projectile motion under varying drag forces. Figure 19 shows the three paths plotted all with an initial height of h = 10 m and initial velocity v = 10 m/s. Figure 19 is important because although we cannot plot an exact function for the projectile’s motion under quadratic air resistance, we can still approximate its behavior and compare it to other models. In this example, the projectile takes a much di↵erent path when drag force is proportional to velocity squared, and its horizontal motion is less than the other model. Hence choosing the right type of air resistance for a certain problem is important, and we must consider the physical properties of the projectile, its speed and the medium it travels through when deciding which option is best for a given problem. Due to the difficulty and approximate nature of our numerical solution for quadratic air resistance, finding the optimal launch angle would involve a guess and check method, where we narrow in on an angle by plotting and fitting points from di↵erent launch angles. We will not pursue this solution, as it focuses less on the math behind projectile motion and more on the power of computer systems to find a solution. Figure 20 shows three paths of a projectile under quadratic air resistance launched at theta values of 70 , 50 , and 30 . 33 6 4 = 70° = 50° 2 0 -2 = 30° -4 -6 0 2 4 6 8 10 Figure 20: Three possible projectile paths subject to quadratic air resistance and launched at di↵erent angles. 6.5 The drag coefficient One important reason as to why the linear model does not work well for our situation lies in calculating the drag coefficient k. For linear and quadratic models, this coefficient is defined as klin = where and D , m kquad = D2 , m (26) are coefficients that depend on the medium through which the projectile launches, D is the diameter of the object, and m is its mass. So far we have assumed that our projectiles are similar to baseballs, golfballs, cannon balls. These objects move at relatively large speeds (⇡ 10 30 m/s) and have diameters in the range of a few centimeters to 0.5 m. In contrast, linear drag models apply to objects with small speeds (on the order of 10 3 or less) and are small in size (< 1 mm in diameter). For this reason, linear drag coefficients we calculate for our types of projectiles are small, and thus this model does not end up a↵ecting the path of the projectile by a noticeable factor. Let us consider an example of linear versus quadratic coefficients. Suppose our projectile is a baseball with diameter 10 cm and has a mass of 100 g. In this case, we assume the projectile is launched in air at standard temperature and pressure, and thus 34 = 0.25 Ns2 /m4 , = 1.6 ⇥ 10 4 Ns/m4 (see [10]). Plugging in our known values to (26), we have 1.6 ⇥ 10 4 · 0.1 = 0.00016 0.1 0.25 · 0.12 = = 0.025. 0.1 klin = kquad These coefficients di↵er by a factor of 100, but how do they a↵ect the motion of the projectile? Recall that klin is applied to the air resistance case where drag is proportional to the velocity of the projectile, while kquad is for drag proportional to the square of velocity. Figure 21 shows three projectile paths: no air resistance, linear air resistance and quadratic air resistance. Note that the linear air resistance path, with its very small drag coefficient, is almost exactly the same path as the no air resistance path. This occurs in this example because the linear model does not realistically apply to the objects we want to consider. If we were looking at tiny objects moving at very low speeds, we could apply this model. No air resistance Linear air resistance Quadratic air resistance Figure 21: Three models of a projectile launched at an angle of 30 with initial velocity 15 m/s and height 10 m. Overall, klin -values tend to fall in the range of 10 4 and smaller, where kquad -values tend to be within the range 0.001 2. Of course, drag coefficients are not exact, universally accepted values, since they are calculated experimentally, and many conditions factor into these values, such as size of object, mass of object, temperature and pressure of medium. 35 7 Conclusion We have thoroughly explored the projectile problem, examining the relationships be- tween optimal launch angle and other variables. We found that in the no air resistance case, solving for the optimal launch angle required finding the intersection between the enveloping parabola and the impact function. For certain impact functions, we were able to find closed-form solutions for this launch angle, but for other impact functions such as the sinusoidal one, we turned to approximation techniques that were still helpful in understanding the relationship between launch angle and other factors. Finally, we added air resistance to the problem, and found that linear air resistance did not adequately model the projectile problem. Instead, we used quadratic drag to understand how air resistance a↵ects the path of the projectile. This problem has many applications, and understanding the connections between each aspect allows us to accurately model the projectile and to find a solution for the optimal launch angle. In the future, the projectile problem could be further explored in many directions. We focused on fairly simple impact functions, but it would be interesting to see how the relationships between optimal launch angle and changing initial conditions vary with more complicated impact functions. Also, we were not able to find a complete solution for the projectile subject to linear air resistance, and learning more about the Lambert W function would be helpful in finding the optimal launch angle for this situation. Finally, we could look at more unusual projectiles. This paper focused on projectiles that resembled baseballs, but it would be interesting to see how the projectile problem changes for other types of projectiles. 36 [ n. pag. [ [ 37 Index air resistance multiple intersections, 16, 21 comparing models, 24, 25, 33, 35 parabolic, 14, 15, 22 drag coefficient, 34, 35 parabolic, with air resistance, 29 linear, 25 semicircular, 15–17 quadratic, 29, 31 sinusoidal, 17 initial conditions quadratic approximation, 31, 33 changing height, 22 di↵erential equations, 5, 26, 27 changing velocity, 23 coupled equations, 29, 31 Lambert W function, 27, 36 distance function, 7, 8, 12, 27 Newton’s method, 18, 21 enveloping parabola, 8, 12 initial guess, 19, 20 derivation, 9–12 equations of motion optimal launch angle, 4, 6, 27 linear air resistance, 25–27 general solution, 12 no air resistance, 5 linear, 13 quadratic air resistance, 29 parabola, 15 expanding circles, 10, 11 sinusoidal, 21 impact function, 4, 7 path function, 9 level ground, with air resistance, 27 linear air resistance, 27 linear, 13, 14 no air resistance, 6 38