Distance between orbits

[Coelum]

Forum,
I'm addressing the problem of computing the minimum possible distance between two non-interacting bodies on elliptical orbits. From a general point of view, it looks like a minimization problem of a function of two variables, e.g. in the domain [0,2*pi)*[0,2*pi). This problem can be numerically addressed in a standard fashion, e.g. by a conjugate gradient method. But I wonder if an analytical approach exists that can simplify the problem - maybe reducing it to unidimensional - and significantly speed-up the computation.

[Jenab]

Forum,
I'm addressing the problem of computing the minimum possible distance between two non-interacting bodies on elliptical orbits. From a general point of view, it looks like a minimization problem of a function of two variables, e.g. in the domain [0,2*pi)*[0,2*pi). This problem can be numerically addressed in a standard fashion, e.g. by a conjugate gradient method. But I wonder if an analytical approach exists that can simplify the problem - maybe reducing it to unidimensional - and significantly speed-up the computation.
For this problem, you don't need to know where the planet is in its orbit at any particular time, so you can ignore the epoch of mean anomaly and time of perihelion passage. But you do still need the other elements ( a, e, i, L, w ).

You probably know this method already, judging by what you wrote. I used the plain old peck-peck-peck method, except I added an outer loop for homing in on the part of the barnyard where the feed is the thickest. It might speed things up a little.

In the procedure below, for variables having two subscripts, the first subscript will designate which orbit (either 0 or 1) and the second subscript will designate either the beginning (0) or the end (1) of an eccentric anomaly search interval. For variables having only one subscript, the subscript will specify which orbit. All eccentric anomalies, and all other angles, are used only in radians.

a : semimajor axis
e : eccentricity
i : inclination
L : longitude of ascending node
w : argument of perihelion

--------------------------------------------------------------------------------

BEGIN

u00 = 0
u01 = 2*pi
u10 = 0
u11 = 2*pi
drmin = 9.9E+99
drgrandmin = 9.9E+99
count = 0

Repeat

count = count + 1
du0 = (u01 - u00)/100
du1 = (u11 - u10)/100
u0 = u00 - du0

Repeat

u0 = u0 + du0

x0''' = a0 (cos u0 - e0)
y0''' = a0 sin u0 (1 - e0^2 )^0.5

x0'' = x0''' cos w0 - y0''' sin w0
y0'' = x0''' sin w0 + y0''' cos w0

x0' = x0''
y0' = y0'' cos i0
z0' = y0'' sin i0

x0 = x0' cos L0 - y0' sin L0
y0 = x0' sin L0 + y0' cos L0
z0 = z0'

u1 = u10 - du1

Repeat

u1 = u1 + du1

x1''' = a1 (cos u1 - e1)
y1''' = a1 sin u1 (1 - e1^2 )^0.5

x1'' = x1''' cos w1 - y1''' sin w1
y1'' = x1''' sin w1 + y1''' cos w1

x1' = x1''
y1' = y1'' cos i1
z1' = y1'' sin i1

x1 = x1' cos L1 - y1' sin L1
y1 = x1' sin L1 + y1' cos L1
z1 = z1'

dx = x1 - x0
dy = y1 - y0
dz = z1 - z0

dr = (dx^2 + dy^2 + dz^2)^0.5

If dr (is less than) drmin then
begin
drmin = dr
u0,min = u0
u1,min = u1
x0,min = x0
y0,min = y0
z0,min = z0
x1,min = x1
y1,min = y1
z1,min = z1
end

Until { u1 (is greater than) u11 } or { du1 (is less than) 1E-12 }

Until { u0 (is greater than) u01 } or { du0 (is less than) 1E-12 }

u00 = u0,min - du0
u01 = u0,min + du0
u10 = u1,min - du1
u11 = u1,min + du1

Q = abs(drgrandmin - drmin) / drmin

If drmin (is less than) drgrandmin then drgrandmin = drmin

Until { Q (is less than) 1e-12 } or { count (is greater than) 20 }

END.


Jerry Abbott

[Coelum]

Jenab,
thanks for your reply. I'll check carefully your code.
Best regards,

Francesco