Central Force Motion: Kepler’s Laws

Lecture L15 - Central Force Motion: Kepler’s Laws When the only force acting on a particle is always directed to­ wards a fixed point, the motion is called central force motion. This type of motion is particularly relevant when studying the orbital movement of planets and satellites. The laws which gov­ ern this motion were first postulated by Kepler and deduced from observation. In this lecture, we will see that these laws are a con­ sequence of Newton’s second law. An understanding of central force motion is necessary for the design of satellites and space vehicles. Kepler’s Problem We consider the motion of a particle of mass m, in an inertial reference frame, under the influence of a force, F , directed towards the origin. We will be particularly interested in the case when the force is inversely proportional to the square of the distance between the particle and the origin, such as the gravitational force. In this case, F =− µ mer , r2 where µ is the gravitational parameter, r is the modulus of the position vector, r, and er = r/r. It can be shown that, in general, Kepler’s problem is equivalent to the two-body problem, in which two masses, M and m, move solely due to the influence of their mutual gravitational attraction. This equivalence is obvious when M � m, since, in this case, the center of mass of the system can be taken to be at M . 1 However, even in the more general case when the two masses are of similar size, we shall show that the problem can be reduced to a ”Kepler” problem. Although most problems in celestial mechanics involve more than two bodies, many problems of practical interest can be accurately solved by just looking at two bodies at a time. When more than two bodies are involved, the problem is considerably more complicated, and, in this case, no general solutions are known. The two body problem was studied by Kepler (1571-1630) who lived before Newton was born. His interest was in describing the motion of planets around the sun. He postulated the following laws: 1.- The orbits of the planets are ellipses with the Sun at one focus 2.- The line joining a planet to the Sun sweeps out equal areas in equal intervals of time 3.- The square of the period of a planet is proportional to the cube of the major axis of its elliptical orbit In this lecture, we will start from Newton’s laws and verify that the above three laws can indeed be derived from Newtonian mechanics. Equivalence between the two-body problem and Kepler’s problem Here we consider the problem of two isolated bodies of masses M and m which interact though gravitational attraction. Let r M and r m denote the position vectors of the two bodies relative to a fixed origin O. Since the only force acting on the bodies is the force of mutual gravitational attraction, the motion is governed by Newton’s law with an equal and opposite force acting on each body. Mm er , r2 Mm = −G 2 er , r M r̈ M = G mr̈ m (1) (2) where r = |r|, er = r/r, and G is the gravitational constant. The position of the center of gravity, G, of the two bodies will be rG = M r M + mr m . M +m 2 (3) Since the two bodies are isolated, we will have, from momentum conservation, that ṙ G =constant, and r̈ G = 0. Therefore, the position of the center of gravity, at all times, can be found trivially from the initial conditions. If the position vector of m as observed by M , r = r m − r M , is known, then the position vectors of M and m could be computed as rM = rG − m r, M +m rm = rG + M r. M +m (4) Therefore, since we know the position of the center of mass rG for all time, we shall show that the problem of determining r M and r m is equivalent to that of determining r, the vector distance between them. The governing equations for rm and rM are given in equation (1) and (2). Subtracting these two expressions, we obtain, r̈ = r̈ m − r̈ M = −G M +m er , r2 (5) or, Mm Mm r̈ = −G 2 er . M +m r (6) The above expression shows that the motion of m relative to M is in fact a Kepler problem in which the force is given by −GM mer /r2 (this is indeed the real force), but the mass of the orbiting body (m in this case), has been replaced by the reduced mass, M m/(M + m). Note that when M >> m, the reduced mass becomes m. However, the above expression is general and applies to general masses M and m. Alternatively, the above expression can be written as mr̈ = −G (M + m)m er , r2 (7) which is again a Kepler problem for an orbiting body of mass m, in which the gravitational parameter µ is given by µ = G(M + m). Example Solution to the Two Body Problem There are two approaches to the solution of the two-body problem. One is a direct numerical attack on equations (1) and (2); the other is to use the analytic solution of the Kepler problem, equation(7), and having found r(t), to use the equation for the position of the center of mass, r G (t) and equation (4) to determine r m (t) and r M (t). The position of the center of mass is determined by the initial conditions (position and velocity) of the bodies. Consider the motion of two bodies as shown in a). The masses of the two bodies are M = 4 and m = 1; for convenience G was set equal to 10. The initial conditions (vector components) are given as r m = (1, 0), ṙ m = (2, 3) and r M = (−2, 0) and ṙ M = (−2, 0). The motion of the two bodies with time is shown in a). From the boundary conditions, we obtain the position of the center of mass with time as r G = (−7/5, 0) + (−6/5, 3/5)t; this position with time is shown in b). The bodies ”orbit” about the instantaneous position of the center of mass. 3 The solution to the ”Kepler” problem for these bodies is shown in c); the solution to the ”Kepler” orbital problem gives the instantaneous position of the relative position of the two bodies, r(t) = r m − r m . The Kepler problem has its origin as the center of mass, which also is the focus of the elliptical orbit. To recover the orbits of the two bodies, we use equation (4). The two orbits are shown in d). These are also the solutions that would be obtained by a direct numerical solution of the two-body problem with boundary conditions chosen to place the center of mass at the origin. The origin serves as the focus for each elliptical orbit. This example shows the importance of formulating the velocity and position boundary conditions so that the center of mass remains fixed at the origin. If this is done, the bodies will orbit about the center of mass, producing the simplest solution to the two-body problem. Equations of Motion The equation of motion (F = ma), is − µm er = mr̈. r2 Since the only force in the system is directed towards point O, the angular momentum of m with respect to the origin will be constant. Therefore, the position and velocity vectors, r and ṙ, will be in a plane 4 orthogonal to the angular momentum vector, and, as a consequence, the motion will be planar. Using cylindrical coordinates, with ez being parallel to the angular momentum vector, we have, − µ er = (r̈ − rθ̇2 )er + (rθ¨ + 2ṙθ̇)eθ . r2 Now, we consider the radial and circumferential components of this equation separately. Circumferential component We have, 0 = rθ¨ + 2ṙθ˙ . Using the following identity, 1 r � � d 2 ˙ (r θ̇) = rθ¨ + 2ṙθ, dt the above equation implies that r2 θ̇ = h ≡ constant. (8) We note that the constant of integration, h, that will be determined by the initial conditions, is precisely the magnitude of the specific angular momentum vector, i.e. h = |r × v|. In a time dt, the area, dA, swept by r will be dA = r rdθ/2. Therefore, dA 1 h = r2 θ̇ = , dt 2 2 which proves Kepler’s second law:The line joining a planet to the Sun sweeps out equal areas in equal intervals of time. Radial component The radial component of the equation of motion reads, − d Since −r2 dt �1� r µ = r̈ − rθ̇2 . r2 = r, ˙ and r2 = h/θ̇ from equation 8, we can write � � � � h d 1 d 1 ṙ = − = −h . dθ r θ̇ dt r Differentiating with respect to time, r̈ = −h d2 dθ2 � � � � 1 h2 d2 1 θ̇ = − 2 2 . r r dθ r 5 (9) Inserting this expression into equation 9, and using equation 8, we obtain the following differential equation for 1/r as a function of θ. d2 dθ2 � � 1 1 µ + = 2. r r h This is a linear second order ordinary differential equation which has a general solution of the form, 1 µ = 2 (1 + e cos(θ + ψ)) , r h where e and ψ are two constants of integration. If we choose θ to be zero when r is minimum, then e will be positive, and ψ = 0. The equation describing the trajectory will be r= h2 /µ . 1 + e cos θ (10) We shall see below that this is the equation of a conic section in polar coordinates. Conic Sections Conic sections are planar curves that are defined as follows: given a line, or directrix, and a point, or focus O, a conic section is the locus of points, P , such that the ratio of the distance between the point and the focus, P O, to the distance between the point and the directrix, P A, is a constant e. That is, e = P O/P A. Since P O = r and P A = p/e − r cos θ, we have r= p . 1 + e cos θ (11) Here, p is the parameter of the conic and is equal to r when θ = ±90o . The constant e ≥ 0 is called the eccentricity, and, depending on its value, the conic surface will be either an open or closed curve. In particular, we have that when e=0 the curve is a circle e<1 the curve is an ellipse e=1 the curve is a parabola e>1 the curve is a hyperbola. 6 Comparing equation(11) which deals solely with the property of a conic section, and equation(10) which provides the solution of the motion of a point mass in a gravitational field, we can identify the properties of the conic section orbits in terms of the physical parameters of the Kepler problem. In particular, we see that the trajectory of a mass under the influence of a central force will be a conic curve with parameter p = h2 /µ. (12) When e < 1, the trajectory is an ellipse, thus proving Kepler’s first law:The orbits of the planets are ellipses with the Sun at one focus. The point in the trajectory which is closest to the focus is called the periapsis and is denoted by π. For elliptical orbits, the point in the trajectory which is farthest away from the focus is called the apoapsis and is denoted by α. When considering orbits around the earth, these points are called the perigee and apogee, whereas for orbits around the sun, these points are called the perihelion and aphelion, respectively. Elliptical Trajectories If a is the semi-major axis of the ellipse, then 2a = rπ + rα . (13) Using equation 11 to evaluate rπ (θ = 0) and rα (θ = π), we obtain a = p/(1 − e2 ). (14) Thus from the geometric properties of an ellipse, rπ = p = a(1 − e), 1+e rα = p = a(1 + e). 1−e Also, the distance between O and the center of the ellipse will be a − rπ = a e. 7 (15) Other geometric properties of the ellipse are that the distance between point D and the directrix will be equal to DO/e, which in turn will be equal to the sum of the distance between the focus and the center of the ellipse, plus the distance between the focus and the directrix. That is, DO/e = ae + p/e. Therefore, DO = a e2 + p = a. Hence, using Pythagoras’ theorem, b2 + (a e)2 = a2 , the semi-minor axis of the ellipse √ will be b = a 1 − e2 . The area of the ellipse is given by A = πab. (16) dA/dt = h/2 (17) A = hτ /2, (18) Also, since is a constant, we have where τ is the period of the orbit. Equating these two expressions and expressing h in terms of the semi-major axis as h2 = µp = µa(1 − e2 ), we have � µ= 2π τ �2 a3 , (19) (20) which proves Kepler’s third law:The square of the period of a planet is proportional to the cube of the major axis of its elliptical orbit. This can be rewritten to obtain the time of flight or period of the orbit. 2π τ = √ a3/2 µ (21) Time of Flight (TOF) in Elliptical Trajectories We have found r(θ), the prediction of the shape of the orbit. However, this solution gives us no direct information about the time behavior of the motions, such as θ(t). In many situations we will need to determine the time of flight between two arbitrary points along the ellipse. In order to do that, we use Kepler’s second law, which states that the motion of the planet sweeps out area at a constant rate. Consider the orbital motion from point 0, to point 1. We would like to determine the time taken t1 . If the motion continues, returning to point 2, the total time taken will be τ . We define the time to reach point 1 as T1 and the time to reach point 2 as t2 . The total time taken is the t1 + t2 = τ , where τ is the total period of the orbit. 8 From Kepler’s second law, equal areas are swept out in equal times. Thus the time taken to reach point 1 is given by t1 /A1 = t2 /A2 = τ /πab (22) where πab is the total area of the ellipse, πab = A1 + A2 . We now construct a more detailed analysis to determine the area Ap swept out by the orbit to a point P. Referring to the figure, we see that the time required to travel between the point π, the periapsis, and an arbitrary point P is proportional to the curved area denoted by AP (AP is the sector defined by O, π, P ). More specifically, since the total period of the orbit is τ and the total area of an ellipse is πab, tP , the time required to travel from π to P , equals the fraction that the area AP represents of the total area of the ellipse. tP = AP τ. πab (23) To find the area AP we construct a circle of radius a with origin at the center of the ellipse. We identify a point P � on the circle to be in a vertical line with the point of interest P on the ellipse intersecting the point O” on the axis. The various geometric quantities of the elliptical orbit have standard definitions: the position angle θ is often called the true anomaly. The radial line of the circle for the origin O� to P � and the major axis of the ellipse 9 major axis define an angle u, which is referred to as the eccentric anomaly. In addition, we define a third anomaly, the mean anomaly M of the point P , as MP = 2πtP . τ (24) Here, tP is the time of flight from the periapsis to the point P . Thus, if we want to determine the time of flight between two points 1 and 2 on the ellipse, we can use equation (24) and write T OF = t2 − t1 = τ A2 − A1 (M2 − M1 ) = τ , 2π πab where A2 − A1 is the area swept out between points 1 and 2. The mean anomaly for point P can also be written as MP = 2π ∗ AP , AT (25) where AP is the area swept out up to the point P . When the area swept out equals the total area of the ellipse AT , the time t equals the period τ and the mean anomaly Mπ = 2 ∗ π. (The subscript π denotes the return to the periapsis π.) Thus the mean anomaly can be thought of as the fraction of the total angle 2π that would be swept out in a time τ by an object reaching point P . The focus is on time not on actual spatial angle. All is needed now is an expression for the mean anomaly M as a function of the orbit parameters. We start by obtaining a relation between θ and u. From simple trigonometry, we have that a cos u − r cos θ = ae (26) or noting that r = a(1 − e2 )/(1 + e cos θ), cos u = e + cos θ , 1 + e cos θ → cos θ = cos u − e . 1 − e cos u (27) We now develop relationships between the various areas indicated on the figure, with the goal to find the formula for the area AP , the area swept out by the point r as it travels from the periapsis π to the point P . The area A1 is the wedge in the circle occupied by the angle u. A1 = a2 u; the area of the large triangle formed by the angle u within the circle is A2 = a2 cos(u)sin(u)/2. Therefore, the area of the large curved segment from O”, P � , π is A1 − A2 = (1/2) × a2 × (u − Cos(u)Sin(u)). (28) The base of the small triangle of area A4 , O, O�� , P , is r × cos(θ) = a × cos(u) − ae by equation (26). The height of the small triangle is b × sin(u). Therefore, the area of the small triangle is A4 = (1/2) × (a ∗ cos(u) − e) × b × sin(u). This area plus the curved segment O”, P, π is the total area swept by the point P . The final step in identifying the area segment swept out between point π and P is to identify the curved segment from O”, P, π, which is then added to the triangle section A4 to form the complete swept area. The curved vertical segment formed from removing the large imbedded triangle A2 from the arc segment of 10 the circle A1 –call it A3 – is geometrically similar to the curved segment formed by removing the small triangle from the area of the swept segment of the ellipse. Since the vertical height of the ellipse is b, and the vertical height of the circle is a, the area of the desired curved segment can by obtained from that of the corresponding segment of the circle by multiplying by b/a. Specifically, AP − A4 = (b/a) × (A1 − A2 ). (29) Therefore, the final result for the area swept out by the point r moving from point π to point P is AP = b/a × (A1 − A2 ) + A4 (30) And the mean anomaly for the point P is MP = 2π ∗ AP 2π × (b/a × (A1 − A2 )) + A4 = AT AT (31) Thus, combining equations (24), (28) and (30), we obtain the mean anomaly for the point P , called Kepler’s equation (It took a Kepler to work this out.) u − e sin u = MP = 2πtP . τ where u is the eccentric anomaly for the point P , defined in the figure. This equation is very easy to use if we want to know the time tP at which the satellite is at position θ. The only thing required, in this case, is the calculation of the eccentric anomaly u using equation (27). On the other hand, if we need to find the position θ of the satellite at a given time t, then, we need to solve Kepler’s equation which is non-linear using an iterative numerical algorithm such as Newton’s method. ADDITIONAL READING J.L. Meriam and L.G. Kraige, Engineering Mechanics, DYNAMICS, 5th Edition 3/13 (except energy analysis) 11