Kepler's 2nd. Law

[Ángel Martin]

From a recent conversation with a friend's on the Kepler's laws - his son's subject for a High school paper. The goal was to calculate the value of the swept area between two instants, as determined by the azimuth angles (a1, a2) of the segments linking the focus of the ellipse with the planet at those moments.

Initially I thought the formulas would involve the Elliptic functions, as elliptical sectors were involved - but it appears that's not the case when the coordinates are centered at the focal point, instead of at the center of the ellipse. I found that fact interesting, as it only involves trigonometric functions (not even hyperbolic).

Here's the reference I followed to program it, a good article that describes an ingenious approach - avoiding painful integration steps. It may not be the simplest way to get this done, chime in if you know a better one.

http://www.bado-shanai.net/Platonic%20Dr...Sector.htm

And here's the FOCAL listing for a plain HP-41 - no extensions whatsoever. The result is the area swept between the two positions defined by the angles a1 and a2; a2 > a1. The parameters a,b are the semi-axis of the ellipse, a>b.

Code:

1    LBL "K2+"    
2    RAD    
3    "a^b=?"    
4    PROMPT    
5    X>Y?       b > a ?
6    X<>Y       yes, swap
7    STO 02     b
8    X^2        b^2
9    X<>Y    
10   STO 01     a
11   X^2        a^2
12   -          a^2 - b^2
13   ABS    
14   SQRT    
15   RCL 01     a
16   /    
17   STO 00     e
18   LBL 00    
19   "<)1^<)2=?"    
20   PROMPT    
21   X<Y?       a2 < a1 ?
22   X<>Y       yes, swap
23   STO 04     a2
24   RDN    
25   STO 03     a1    
26   PI       
27   X<=Y?      a1 >= pi ?
28   GTO 01     yes, -> case 2
29   RCL 04     a2
30   X<=Y?      a2 <= pi ?
31   GTO 01     yes, -> case1
32   X<>Y       no, it's case 3
33   XEQ 02     first between [pi, a2]
34   STO O      partial result
35   RCL 03     a1
36   PI         now between [a1, pi]
37   XEQ 02    
38   RCL O    
39   +    
40   RTN    
41   GTO 00    
42   LBL 01    
43   RCL 04     a2
44   RCL 03     a1
45   LBL 02    
46   XEQ 05    
47   STO  M     f1
48   X<>Y       a2
49   XEQ 05    
50   STO  N     f2
51   RCL  M    
52   -    
53   2    
54   /          (f2 - f1) /2
55   ENTER^    
56   SIN    
57   RCL  M     f1
58   RCL  N     f2
59   +    
60   2    
61   /          (f2 + f1) /2
62   COS    
63   *    
64   RCL 00     e
65   *    
66   -    
67   RCL 01     a
68   X^2        a^2
69   *    
70   1    
71   RCL 00    
72   X^2    
73   -          1-e^2
74   SQRT    
75   *        
76   ABS    
77   RTN    
78   GTO 00    
79   LBL 05    
80   COS        cos a
81   RCL 00    
82   X<>Y    
83   +          e + cos a
84   LASTX      cos a
85   RCL 00     e
86   *          e.cos a
87   1    
88   +    
89   /          cos f =  (e + cos a)/(1 + e.cos a)
90   ACOS       f
91   END

(*) the symbol "<)" is for the angle character.

Example: calculate the area swept between a1 = pi/4 and a2 = 3.pi/4, if the ellipse parameters are a= 2, b= 3

Solution: A = 1.989554087

Once the area is obtained, and knowing the period of the orbiting (T), it's straightforward to determine the time taken by the planet to travel between the two positions, with the direct application of Kepler's 2nd. law:

t = T . Area / pi.a.b

Cheers,
ÁM