THE UNIVERSITY OF MICHIGAN COLLEGE OF ENGINEERING Department of Aeronautical and Astronautical Engineering Technical Report TRAJECTORY CALCULATIONS FOR SOUNDING ROCKETS By Robert M Ho.we UMRI Project 04003 under contract with: DEPARTMENT OF THE ARMY BALLISTIC RESEARCH LABORATORIES ABERDEEN PROVING GROUND, MARYLAND administered by: THE UNIVERSITY OF MICHIGAN RESEARCH INSTITUTE ANN ARBOR July 1959

TRAJECTORY CALCULATIONS FOR SOUNDING ROCKETS SUMMARY This report describes the derivation of the equations used to calculate sounding-rocket peak altitude, range, and sensitivity to velocity error, flight path elevation error, and yaw angle at burnout. Curves which allow rapid computation of these quantities are shown for rockets with equivalent straight-up altitudes of 0-2000 miles and burnout flight path angles of 60-90 degrees, 1. Basic Free-Fall Equations We will use as an inertial coordinate frame of reference a system with origin fixed at the center of the earth. The axes of this system will be non-rotating, so that the earth rotates with respect to the system. We will neglect perturbations due to the non-spherical earth corrections, so that any free-fall trajectory will be in a plane determined by its velocity vector at any point and the center of the earth, In general, we will use polar coordinates r, 0 as shown in Fig, 1. 1 to locate the rocket in this plane. Y v Since we will find it convenient rO h ~to use energy methods for many of the calculations, we will use Lagrange's equation to set up the usual equations vr of motion for a point mass m (the rocket mass = m) in a central force 0 \ 9 field. Thus the Lagrangian function C I — enter —- L-PW L = T-V, where T = kinetic energy and Center of Earth V = potential energy*. If qi is the ith general coordinate, then the equations of motion can be written as d di L 0 dt - =O_ (1, 1) dt 1 For our polar-coordinate case there are two general coordinates, ql = r, q2 = 8. Hence Kinetic Energy T= T = m (2 + r22) (1.2).Slater and Frank, Mechanics, McGraw-Hill, 1947, pp 72-74

and Potential Energy = V = - Mm (1, 3) where V is the gravitational constant and M is the mass of the earth. Eq, (1.3)can be written in more useful form by noting that at the earth's surface (r = r0) V Mm mg0 = (1. 4) r where go is the acceleration of gravity at the earth's surface, not including the term roC due to the angular rate of earth rotation W. In subsequent calculations we will 2 assume go = 32. 2 ft/sec In terms of Eq. (1.4), Eq. (1.3) becomes 2 V= -mg0 r (1, 5) Thus L becomes 1 2 2.2r 22(. 2 Or L =T - V= 2m( +r22) + mg0 r (1.6) Substituting Eq. (1. 6) into (1. 1) we have 2 Ir- rO2 +go = 0 (1.7) r and d mr21 = 0 (1,8) dt Eq. (1. 7) states that the acceleration due to the trajectory balances the. acceleration due to gravitational attraction. Eq. (1, 8) is simply the equation for constant angular momentum, since after integration 26 mr = p = angular momentum = constant. The time variable can be eliminated from Eqs.(1, 7) and (1.8) to yield an equation describing the trajectory shape (i. e., involving r and 0 only) by defining a new variable u given by U=- (1.9) Then. = dr du _ 1 du = p du rdu dt - "^2 -G y; = m d9 (1, 10) 2

Similarly P) =1 d2u( (1. 11) ^ ^.2 2 r dO From Eqs. (1. 10) and(1. 11) Eq. (1. 7) becomes 2 2 2 du +u= mgrO (1.12) d02 2p But the solution to this equation can be written immediately as 2 2 1 - m g0r0( - - = u = C cos (9 - 90) + (1 13) r 0 + p where 0 and C are arbitrary constants which depend on the initial conditions of the trajectory, or, in this case, specified values of r and 0 at some given time and the angular momentum p. Solving Eq. (1. 13) for r, we have 2 p 2 2 m g0r 0 r = 2m- (1. 14) 1 + C 2 2 cos (0-0) 2 2 m g0r0 The equation for an ellipse of eccentricity C has the form 2 a(l-~2) r (+ cs (1. 15) where a is the semimajor axis and where the origin of polar coordinates is at the right hand foci of the ellipse (see Fig. 1. 2). This Fig. 1. 2 Ellipse 3

can be seen from the law of cosines for the triangle formed by r and r" and from the equation r + rq = 2a - constant, which defines an ellipseo r2 2 2 r92 r + (2a ) + 4r(aE) cos O Next we consider the determination of 6 and a in terms of initial conditions for the trajectory. The solution to Eq, (1. 12) can also be written as 2 2 1 m g0r0 1 u A cos 9 + B sin 0 + (1. 16) r 2 p where A and B are constants determined by initial conditions. Differentiating Eq. (1. 16) with respect to time - L r - A sin 0 + B cos 0 (1.17) At r rl, let ir =1, 9 - 1 Q =0 1 Then Eqs. (1. 16) and (1. 17) can be written as 2 A cos 0 + B sin 0 - 0 (1. 18) - 1 1 r r 0 1 0 1 and' ~1 -A sin 8 + B cos 98 - (1. 19) 1i 1 at 2..1 2 r r!!1 Solving for A and B, we have 2 g r 2 2 1 r and 1 g0r0 B f B sin - + cos (1.21) 1.. 1 2 1 r91 0 r01 Note that r019 - v0 where V0 is the component of total velocity v parallel to the earth's surface at r rl' 2 Comparison of Eqso (1. 14) and(l. 15) shows that ~ = C - 2 But from f->-, n gar0 Eqs. (1. 14) and (1. 16) it is clear that C.=/A2+ B2, Thus 4~~~~~

----- --- 2 4r 2 E = /A2 + B2 m g0r g0r vent escape, C 41, and from Eq. (1. 22) this implies that or / —------------------— r 32 2 V2 (r + (1 < 2 r; 2 (1.22) circula r satellit where orbit, = 0v For the components of rocket velhyperbole and the rocket escapes from the earth's gravitationalfield. Thus to prevent escape, C /. 1, and from Eq. (1. 22) this implies that ocity at r = r1 parallel and perpendicular, respectively, to the earth's surface, *2 2 2 2 2 Hence (rl01) + 1 = v v = v1, where vl = total rocket velocity. At the earth's surface r1 = rO and vp < 2rpg0 to prevent escape. Hence we define the escape velocity vE as vE = 2r gi t1.24) = /2 x 3950 x 32. 2 x 5280 = 36, 700 ft/sec In terms of the escape velocity vE Eq. (1. 22) for e can be rewritten as fl /"BA2\2 r r/ l VO 2 C/(2 1( ) )2-[1 e I rl 2 ] (1. 25) E2 V2 where {1 = (0 at r = r1, and where (O is a dimensionless radius equal to the ratio of actual radius r to earth's radius r0. Thus e r (1.28.) If we know the radius /1 and the velocity components vg and r, then Eq. (1. 23) allows us to calculate e. To calculate the semimajor axis a we note by comparing Eqs. (1, 14) and (l. 15) that 5

a(1-E2 P2 2 2 r0 2 2 1E) or \ 2 2-1 l^ 0o 201 (1 r- (1 27) 1 2 Now that we have developed the basic equations describing trajectories of a rocket in free-flight in the gravitational field of the spherical earth, we turn to consideration of some specific trajectory problems. 2. Vertical Launch Free Flight Trajectory In the case of a vertical launch the problem is complicated by the rotation of the launching platform, namely, the earth. Thus the velocity vector at burnout will be the sum of the rotating-earth component and the propulsion-unit contribution. To simplify the initial analysis, let us assume that this sum is such that at finalstage burnout, i. e., at the beginning of free-flight, the total velocity vector is perpendicular to the earth's surface. Since this will give maximum altitude, small deviations from this perpendicular direction will only cause second order effects in total altitude, Under these conditions a calculation based on energy conservation is the easiest approach. Assume r = r1 at burnout and r2 at the apex of the trajectory. Then the change in potential energy V1 - V2 must equal the change in kinetic energy -T, or from Eq. (1. 5) 2 0'0 1 2 mg0rO 2 - 2- r4 v12 (2.1) where vl is the total velocity at burnout. Replacing go0 yy vE2 /2 in accordance with Eq. (1. 25) and solving for v1 we have 1; VE - (2. 2) where ~1 = r1/r0 and {~2 = r2/r0. If we let h1 equal the altitude at burnout and h2 equal the peak altitude, then 1 = 1 h, / 2 1+ 2. Since hl/r0 <~1 2ne /a 2is 1/ 1 0 1- z /- ad- 1 ) bh 0 r0 r0 under normal circumstances, 1/P1, 1 - h /r, and Eq. (2. 2) becomes 6

v 1v h <r (2. 3 1 VE r0 1+ r 0 (2 3) Eq. (2. 3) can be solved exactly for peak altitude h2. Thus v 2 h +r 1 1 0 2 h 2 2E-, 2 < 1 (2.4) h v 0 rO 2 VE In Fig. 2. 1 the peak altitude h2 in miles for various burnout velocities v1 is plotted for hi = 0. For h1 > 0, the maximum altitude is obtained approximately by adding hn / 1 - (v)1 to the ordinate shown, providing r- < 1. 1E I r 3. Calculation of Peak Altitude for Non-Vertical Launch Next we consider a calculation of the peak altitude reached when the velocity vector at burnout is not vertical but makes an angle / with the plane parallel to the earth's surface directly below the burnout point. It is assumed that the velocity vector is measured with respect to the inertial (non-rotating) frame of reference with the origin at the center of the earth, as in Section 1. Let vl be the velocity at burnout v2 =vo k-2 e2 Max. Alt, v r1 Burnout r2 r/ Fig. 3. 1 Center of \ Earth 1 7

1400 -iiiI 1,300 PEAK ALTITUDE VERSUS VELOCITY FOR A STRAIGHT-UP VACUUM TRAJECTORY, SEA LEVEL LAUNCH 1200 - r0o( V )2 VE 1100 h2 = —. t- 900 800 1800 Z 700 1700 co c: w 800 1800 0 500. / 5 100 u 400 - 1400 w a- 300 1300 200 - 1200 o00 2000 4000 6000 8000 10,000 12,000 14,000 16,000 18,000 20,000 VELOCITY IN FEET/SEC AT SEA LEVEL Figure 3.1

with components v and v parallel ta and perpendicular to, respectively, the 91 r1 earth's surface, as shown in Fig. 3. 1. Let r1 be the polar radius of the rocket at this point. At the apex of the trajectory the velocity v2 = v2, and r- r2. The most direct method to calculate r2 and hence h2, the maximum altitude, is to use equations for constant energy and constant angular momentum. From Eq. (1. 5) we have for constant total energy 2 2 1 2 r0 1 2 r0 — m v -mg0 r- = -mv2 - mg0 r2 (3.1) Replacing g0r0 by VE2/2, where vE is the escape velocity, it follows that 2 2 V1 2 1,1 ) ) — p-+ =0 (3, 2) where 1 = r1/rOs (2 = r2/r0, r = radius of earth. Since the angular momentum is constant, we have mrlv1 cos Y = mr2v2 (3.3) Eliminating v2 in Eq. (3. 1) by means of (3. 3), we obtain 22 tE \2 + W 2 (02 i[l (v ( + ) cos2 = 0 (3.4) Solving for p 2 which represents maximum P, we have 2 V21:_ 2 2v1 2 L~F:11 W @2 L' ( + 1-4 (vv )2e12 COSV[G (3.5) The above expression reduces to'2 = for V = 0 and to Eq. (2, 2) of the previous section for / = 90~ For the case where 4( 2_ cos_ 1 V v which occurs when vE (41 or cos21, Eq. (3. 5) can be written approximately in terms of h2 = r0 ( - 1). Thus 9

h -h sin (3.6) 2 2= /2 where h2 is the maximum altitude given in Eq. (2. 4) for the case of a V= 1712 vertical launching. Thus the effect of burnout angle V different from 90 degrees is,approximately, to reduce the altitude by the factor sin V. The same result is obtained by considering the firing as one straight up with a velocity equal to the vertical component v1 sin V of the velocity v1. The correction factor which must be applied to Eq. (3. 6) for various altitudes and launch angles is shown in Fig. 3, 2. The peak altitude reached is always somewhat higher than that given by Eq. (3. 6) due to the influence of centrifugal acceleration. 4. Calculation of Sounding Rocket Range The calculation of range is of interest not only because it is necessary for the prediction of impact point but also because it provides the equations necessary to calculate the sensitivity of range to burn-out velocity and flight-path angle. The equation for the trajectory of free-flight after burnout is given byEq. (1. 16), which is repeated here for convenience 2 1 _ g0r0 A cos 0 + B sin 0= - - (4. 1) rl 0. where A and B are given in terms of the burnout variables, i. e., r1, 01, 1l' 0l' in Eqs. (1. 20) and (1. 21). Let us arbitrarily assume that the impact occurs at 0 = 0. Also r = r0 at impact. For these values of r and 0 we have from Eqs. (4. 1) and (1.20) 2 ~. ~ 2 - 42 CCos e + I s in 8 =- g0r-o (4. 2) 4 2: 2I 4) 2 r1 1 rl 1 Replacing g0r0 by VE /2, rll1 by v0, and t by vr, we have where l = rl/r0. Finally, vr = v1 sin i, v0 = v1 cos V. Thus Eq. (4.3) becomes 10

1.28 1.26 1.24 CORRECTION FACTOR FOR CALCULATING MAXIMUM ___ SOUNDING-ROCKET ALTITUDE 1.22 Max. Alt. hy =9oo [sin2 ], K 1.20 y = 45~^ 1.18 1.16 1.14 Kf 1.12 1.10 1.08 1.06 1.04 1.02 1.00 0 200 400 600 800 1000 1200 1400 1600 1800 2000 MAXIMUM ALTITUDE IN STATUATE MILES, y=90~ Figure 3.2

1 2\2 ~1 V 2 cos2/ - cos 01 + VE in2V sin9 (Vli cos2- 2~/ (4.4) VE 2 For a given velocity v1, flight-path angle %/, and radius ~1 at burnout, the angle 01 given by Eq. (4. 4) represents the burnout 0 coordinate necessary to give impact at 0 = 0. The range of the rocket is, then simply -r001. Eq. (4. 4) has the general form C1 cos 80 + C2 sin = D (4.5) But this can be rewritten as C + C2 cos ( + ) D (4.6) where C 7= - -tan- 2- (4.7) Hence D C c-1 D1 2 1 =cos 2 2 + tan 2 (4.8) 1 /C 2+C2 C1 By analogy between (4. 4) and (4. 5) 1 V cos -E cos 1 V1 - 2V/ 1 v" 2 For the special case where (1 = 1, i. e., where the burnout is assumed to take place at sealevel, the two angles on the right side of Eq. (4. 9) are equal and 12

v \2 v, 1 sin 2 0 =2 tan1 2 IE, (4. 10) Vl 2 1 -v cos - 2r VE 2 Since the burnout altitude for most sounding rockets is low (perhaps 10-30 miles), Eq. (4. 10) will give a reasonably accurate estimate of sounding-rocket range. For more accurate computation Eq, (4. 9) can always be used, Eq. (4, 10) should, however, be quite adequate for establishing the dependence of range on v1 and /. Note that we have also ignored the reentry trajectory here. If a missile with velocity v1 and flight path angle V/ were launched in a constantgravity field of acceleration go, then the total time of flight would be (2v1sin V )/go and the horizontal distance traveled would be [(2v1 sinY)/g [v1 cosY = (v12 sin 2 V )/go. This makes a convenient norm to compare the actual distance traveled along the earth's surface as computed from Eq. (4. 10). The required correction factor K is shown in Fig. 4.1 as a function of rocket altitude for vertical burnout for various burnout angles V. Thus 2 v sin 2 i Range - K (4. 11) g0 This allows a quick calculation of range for a wide variety of burnout conditions. In Fig. 4. 2 are plots of actual range R for various V versus peak altitude for straight-up launch. 5. Effect of Errors in Burnout Velocity and Angle on Impact Point 5. 1 Effect of Velocity Error In Eq. (4. 9) we derived the launch polar angle 01 for an impact at 0 = 0~ with velocity v1 and flight-path angle vat launch, Eq. (4. 10) is an excellent approximate representation of 01 ( it is exact for the case where launch, or more correctly, burnout occurs at P 1 = 1, the earth's surface), Eq. (4. 10) can be rewritten as -1 sin2 2 0 - -2 tan si (5. 1) 1V 2 -2 cos / 13

CORRECTION FACTOR FOR CALCULATING SOUNDING ROCKET RANGE FOR VARIOUS____ BURNOUT FLIGHT-PATH ELEVATION ANGLES Range 0= Z y K 1.20 / /50~ 1.10 1.00 200 400 600 800 1000 1200 1400 1600 1800 2000 MAXIMUM ALTITUDE IN STATUATE MILES FOR y=90~ Figure 4.1

1400 RANGE OF SOUNDING ROCKETS FOR VARIOUS BURNOUT FLIGHT-PATH ANGLES __ 1200 --— _ _1000 70~ U) w -J w / 5 / (c H~~ Z 600 ____ 400~~~~~~~~~~~~~~~~~~~~~~~~~~~~10 w 0 z 400__ _ _ 85 200 0 200 400 600 800 1000 1200 1400 1600 1800 2000 MAXIMUM ALTITUDE IN STATUATE MILES FORy =90~ Figure 4.2

.26 EFFECT OF BURNOUT VELOCITY ERROR.24 --- ON RANGE 8.22 d Y.=60^^^^' c).20 C) w x 0 I - 81 4 T I- -^^ LL 75 O 1.F w cn r 3~~~~~. w _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ U.1 w.0 2 0 0 00 6080 00 20 140 100 10 H ~Z 10 _______80 00 2014010010 20 MAXIMUM ALTITUDE IN STATUATE MILES FOR =900 Figure 5.1

Taking the partial derivative with respect to vi, we have vE 3 A3 v ) (5c 2) 1 l.1E - _2 cos20 + sin 2 / v1 The change in range R due to changes in vl is simply r 0 a 1/ vi. Thus the error in range A R is ae a R -rO av _ v (5. 3) where A vl is the change in burnout velocity from the norm value. vertical sounding rockets. Thus if (vE/Vl 2 2 cos >/, and (vE/vl) 4 sin2 2', 91 4 Posin 2' / v 1 O - I 1, 4) v vE } For a 1000 mile sounding rocket with a burnout angle of 75 the approximate expression (5. 4) gives r0 a 1/ v1 = -0. 0968 mi/ft/sec compared with the exact value of -0. 099 from Eq. (5. 2). The exact value corresponds to 16. 4 miles error in range for a 1% error in velocity vl. In Fig. 5. 1 are shown plots of r Q01/) v for burnout angles between 60~ and 900. 5. 2 Effect of Elevation Angle Error Next we calculate the effect of errors in flight-path elevation angle P/ at burnout. Taking the partial derivative of Eq. (5. 2) with respect to, we have siii VE 2_2 1 vE 2 2Again a fair approximation is obtained for sounding rockets when (vE/v) 2 cos 0117 (VE/Vl)4>> sin22 V', and|(vE/vl)2 cos 2V ^>sin22 / 17

I: 0 ar, 120|.. w EFFECT OF FLIGHT-PATH ELEVATION ANGLE z ERROR ON RANGE 0 > 100C —_J go o, I: I 80 o -------- 80~_ _ _0 _7 z 20 II 0w 60 co tL4'9 40 z 20 200 400 600 800 1000 1200 1400 1600 1800 2000 MAXIMUM ALTITUDE IN STATUATE MILES FOR y =900 Figure 5.2

1 -4 ( v1i cos 2V (5. 6) 1 E For a 1000 mile sounding rocket with V = 75~, 8 91/ = 0. 706 from the approximate Eq. (5. 6) and 0. 763 from the exact expression (5. 5) the error in range A R is given by AR = r 1 Av^ (5.7) where AV is the error in / at burnout. For the rocket example earlier this turns out to be 48. 7 miles/degree. In Fig. 5. 2 are shown plots of r0 90la/8/for burnout angles between 60~ and 90~ 5. 3 Effect of Yaw Angle Error Let us define a yaw angle ) as measured in the plane through the desired velocity vector at burnout but normal to the desired plane of motion, as shown in Fig. 5.3. This defines a new plane of free-flight Sk^~i ~ \ ~motion which makes an angle jr = )/ cos Vector Vector at Vctual Vell\ I \ Vectore atl, with the original plane, assuming, ~, 1. / \ Burnout Hence as long as the range R is not so /\. \1^ ~long that spherical-surface effects predominate, the error in impact position //~! \ \ ~d at right angles to the plane of the mo/, \ Projection of tion is given by' Vel. Vecotrs on "^' \ Horizontal Plane d= Cos VI (5.8) [~rt^ ~ For the 1000 mile rocket with V = 75~, Fig. 5. 3 R 820 miles and d = 55, 3 miles/deg. 6. Effect of the Earth's Rotation on Sounding Rocket Trajectories 6. 1 Calculation of the Burnout Flight-Path Angle and Velocity in Inertial Coordinates All of the trajectory computations up to now have been made with respect to a non-rotating coordinate system with its origin at the center of the earth. Actually, 19

the rocket is launched from the surface of the earth, which is moving eastward at a velocity v = J r0 cos' = 1515 cosetft/sec, where SJ is the angular velocity of rotation of the earth, r0 is the earth radius, and 0' is the latitude of the launch point. All of the rocket performance calculations are made with respect to the moving-earth's surface at the launch point. To the burnout velocity v1' calculated with respect to this system we must add vectorally the eastward velocity component ve to obtain the burnout velocity v1 with respect to the inertial coordinate system. The simplest case is the one where we fire the rocket directly east or directly west, i. e., where vl' and v lie in the same plane, as shown in Fig. 6. 1. ivf~ Ve Here the horizontal velocity component -...... lx v1 cos / ~ v and the vertical vy l....... lx = ev component vly v'sinv', where 9' is v 1' / / is the burnout angle as calculated with v1 respect to the rotating earth and V is the y / burnout angle as calculated with respect to the inertial coordinates. For an eastVa/ | ward launch we have +ve in the formula v X v for v1x; for a westward launch we have.../- SLY -ve. The total velocity vlin the inertial Fig. 6.1 x 1 coordinate system is given by ly 1 ---— 1 — mately ~ e ~~1' e V1 [1 l~ cos I + - - (6. 2) This approximation formula should be quite accurate for near vertical launches of rockets to several hundred miles or higher. From Fig. 6. 1 it is apparent that J = tan 1 y = tan 1 cosy17'v (6. 3) If 2Pe -y ( (1, then we can write the approximate formula 20

v'I / y?+ e sin / (6.4) vi As a specific example, consider a nominal 1000 mile sounding rocket fired eastward at a latitude of 40~N. Let v1' = 16, 500 ft/sec and / - 75~ Then v 1160 ft/sec, v1 = 16840 ft/sec, and Y= 71. 1~ Thus the effect of launching the rocket east has been to raise the burnout velocity by 340 ft/sec and lower the flight path angle at burnout by 3. 9. If we had neglected this effect the peak altitude for the original 16500 ft/sec velocity and 750~ burnout angle would have yielded a calculated peak altitude of 951 miles. With the earthgs rotational velocity taken into account the peak altitude is 970 miles. For a westward launch with vl = 16500 ft/sec ands/' = 75~ as before, vl = 16240 ft/sec and V = 78. 9~ The correct peak altitude is then 935 miles. Thus an eastward launch yields a slight advantage in peak altitude. For azimuth angles at launch which are approximately eastward or westward, the above results will apply with reasonable accuracy. For other azimuth angles the analysis is more complicated but is equally straightforward. 21