Lunar landing game madness

[robve]

So I got one of those "Mathemagician" (APF 1500) calculators from 1980 last week to try and could not stop myself like a lunatic to get perfect Lunar Landing game landings.

The game is fairly simple, like this in HP-71B BASIC:

Code:

100 INTEGER F,V,H,A @ DELAY 0,0
110 F=99 @ V=20 @ H=300
120 DISP F;V;H
130 A=VAL(KEY$)
140 F=F-A @ IF F<0 THEN A=A+F @ F=0
150 H=H-V @ V=V-A+5
160 IF H>0 THEN 120
170 IF V=0 AND H=0 THEN DISP "PERFECT LANDING" @ END
185 REM DISCLAIMER: NOT SURE WHAT MATHEMAGICIAN USES HERE:
180 IF V-H<=5 THEN DISP "CRASH LANDED, REPAIRABLE" @ END
190 DISP "CRASH"

Surely there must be a better way to perform a perfect landing by punching the right sequence of rocket burns to burn as little fuel as possible? It's not that easy (at least not with the Mathemagician version of the game).

I came up with a relatively simple method to find the perfect sequence of burns, which doesn't involve a numerical solver, but does take some time to compute a solution using a brute-force approach. This approach can surely be improved or replaced with a more elegant and faster algorithm

Let's think of this problem as if we're taking off. A landing is a time-reversed takeoff. The best way to preserve fuel is to go full burn when taking off to reach a high velocity. Likewise, when landing, a full burn to the surface from a set altitude and velocity uses less fuel overall. Worst way to take off (or land) is to stall and burn fuel without making any progress.

For example, if the lander is at an altitude (height) h0=100 and velocity v0=20, we set the burn rate to a maximum (say a=9) to land on a lunar surface. Assuming the moon's gravity g=5 (which must be lower than a, of course!), then we have an effective acceleration (time-reversed deceleration) of 4, giving the following velocity and height:

0 (at the surface)
4     4
8     12
12    24
16    40
20    60
24    84
28    112

This consumed 63 units of fuel. If we set a=7 then we consume 70 units of fuel to get to an altitude of 110. Computed in Lua:

g = 5; a = 9; v = 0; h = 0; f = 99
while h >= 0 and h <= h0 do
  f = f - a
  -- out of fuel?
  if f < 0 then a = f + a; f = 0 end
  v = v + a - g
  h = h + v
  print(v, h)
end

In this example we already passed the desired height h0=100 and also our final speed exceeds v0=20. The tricky part is to determine when and where to reduce the burn to reach the exact height h0 with velocity v0.

Let's add an autopilot of sorts that computes a burn rate b that is required to reach the desired velocity at the end of the remaining distance to altitude h0. It may overshoot some, but we want to get the velocity under control. Such burn rate b can be computed with:

d = h0 - h
t = d / v
b = (v0 - v) / t + g

if b is small or negative, then we really have to slow down. In this case we found the point at which the full burn should stop to allow the autopilot to take over control. Let's do this by setting mode flag m=true:

w = 1; m = false
g = 5; a = 9; v = 0; h = 0; f = 99
while h >= 0 and h <= h0 do
  f = f - a
  -- out of fuel?
  if f < 0 then a = f + a; f = 0 end
  v = v + a - g
  h = h + v
  -- remaining time steps t to h0 and burn b to slow to v0
  d = h0 - h
  t = d / v
  b = (v0 - v) / t + g
  -- switch to autopilot?
  if b<=w and nor m then
    m = true
    r = d + v
    s = v - v0
    n = 0
    print("now on autopilot")
  end
  print(v, h)
  if m then
    -- on autopilot
    a = math.max(math.floor(b), 0)
    n = n + 1
  end
end

For h0=300 and v0=20 for the Mathemagician Lunar Lander game we get:

4      4
8      12
12     24
16     40
20     60
24     84
28    112
32    144
now on autopilot
36    180
31    211
27    238
23    261
21    282
19    301

Not bad. But not good enough, because it overshoots h0 by 1 unit and the velocity is one unit off from the target v0. To find the perfect burn sequence, we can use a bit of brute-force to find the solutions by computing the required burns after the full blast stage. For the above conditions, there are nine remaining burns each with five steps for n=6 (n is computed as above) that are perfect:

3 6 0 0 0
4 4 1 0 0
5 2 2 0 0
5 3 0 1 0
6 0 3 0 0
6 1 1 1 0
6 2 0 0 1
7 0 0 2 0
7 0 1 0 1

All of these nine burn sequences of five steps give a perfect landing using the same amount of fuel. For example 3, 6, 0, 0, 0, 9, 9, 9, 9, 9, 9, 9, 9 lands the lunar lander with 9 fuel units remaining (did I mention the Mathemagician version of this game has tight tolerances?).

To brute force compute the solutions we need:

-- the remaining velocities sum to k
k = (n - 1) * g - s

We want to find a burn sequence a[1], a[2], ..., a[n-1] such that

\( \sum_{i=1}^{n-1} a_i = k \)

and

\( v_0 + \sum_{i=1}^{n-1} \sum_{j=1}^i g-a_i = v_0 + \sum_{i=1}^{n-1} (n-i)(g-a_i) = r \)

where \( r \) is the remaining distance to h0 plus one more step (the current velocity v) as set by the autopilot r = d + v.

The Lua code to produce perfect lunar landing sequences:

Code:

-- compute perfect solutions to a lunar lander game
-- assumes burn rate is between 0 and 9

v0 = 20   -- velocity (lunar lander initial velocity)
h0 = 300  -- height (lunar lander initial height)

g = 5     -- gravitational constant 0 < g < 9

a = 9     -- current engine burn (acceleration), we land with max burn
v = 0     -- current velocity
h = 0     -- current height
f = 99    -- current fuel level

m = false -- autopilot when true
w = 1     -- weight factor to go autopilot: higher is smoother ride but burns more

print("==> v0 = " .. v0 .. " m/s")
print("==> h0 = " .. h0 .. " m")
print("fuel burn speed height      t      b")

while h >= 0 and h <= h0 do
  f = f - a
  -- out of fuel?
  if f < 0 then a = f + a; f = 0 end
  v = v + a - g
  h = h + v
  -- remaining time steps t to h0 and burn b to slow to v0
  local d = h0 - h              -- distance remaining
  local t = d / v               -- time steps remaining
  local b = (v0 - v) / t + g    -- autopilot engine burn
  -- switch to autopilot?
  if b <= w and not m then
    m = true
    r = d + v
    s = v - v0
    n = 0
    print("==> r = " .. r .. " m remaining to go on autopilot")
  end
  print(string.format("%4d %4d %5d %6d %6.1f %6.1f", f, a, v, h, t, b))
  if m then
    -- on autopilot
    a = math.max(math.floor(b), 0)
    n = n + 1
  end
end

if f > 0 then

  -- actual fuel left when starting at v0 instead of v
  f = f + v - v0
  print("==> " .. f .. " units of fuel left")

  print("==> " .. n - 1 .. " initial engine burn steps, which can be:")
  -- the remaining velocities sum to k
  k = (n - 1) * g - s
  -- engine burn (accelerations) for the n-1 steps to compute
  a = {}
  j = k
  for i = n - 1, 1, -1 do
    a[i] = j % 10
    j = j // 10
  end

  p = 0

  repeat
    -- sum burns (accelerations)
    local s = 0
    for i = 1, n - 2 do s = s + a[i] end
    a[n - 1] = k - s
    if a[n - 1] >= 0 and a[n - 1] <= 9 then
      -- compute height difference (quadratic)
      local q = n * v0
      for i = 1, n - 1 do q = q + (n - i) * (g - a[i]) end
      -- print solution
      if q == r then
        print(unpack(a))
        p = p + 1
      end
    end
    -- next burn sequence
    local f = true
    for i = n - 2, 1, -1 do
      a[i] = a[i] + 1
      if a[i] <= 9 then f = false; break end
      a[i] = 0
    end
    -- ran out?
  until f

  if p == 0 then
    print("*** BUMMER, NO SOLUTION ***")
  end

else

  print("*** OUT OF FUEL! ***")

end

A version in a generic BASIC:

Code:

100 CLEAR : V0=20 : H0=300
110 G=5
120 A=9 : V=0 : H=0 : F=99
130 M=0 : W=1
140 PRINT "V0=";V0
150 PRINT "H0=";H0

160 F=F-A : IF F<0 THEN LET A=F+A : F=0
170 V=V+A-G : H=H+V
180 D=H0-H : T=D/V : B=(V0-V)/T+G
190 IF B<=W AND M=0 THEN LET M=1 : R=D+V : S=V-V0 : N=0 : PRINT R;" TO GO"
200 PRINT F;A;V;H;T;B
210 IF M THEN LET N=N+1 : A=INT(B) : IF A<0 THEN LET A=0
220 IF H>=0 AND H<=H0 THEN 160

230 IF F=0 THEN PRINT "OUT OF FUEL!" : END
240 DIM A(N-1) : K=(N-1)*G-S : J=K
250 FOR I=N-1 TO 1 STEP -1 : A(I)=J-10*INT(J/10) : J=INT(J/10) : NEXT I
260 P=0

270 S=0 : FOR I=1 TO N-2 : S=S+A(I) : NEXT I : A(N-1)=K-S
280 IF A(N-1)<0 OR A(N-1)>9 THEN 310
290 Q=N*V0 : FOR I=1 TO N-1 : Q=Q+(N-I)*(G-A(I)) : NEXT I
300 IF Q=R THEN P=P+1 : FOR J=1 TO N-1 : PRINT A(J);" "; : NEXT J : PRINT
310 F=1
320 FOR I=N-2 TO 1 STEP -1
330 A(I)=A(I)+1 : IF A(I)<=9 THEN LET F=0 : I=1 : GOTO 350
340 A(I)=0
350 NEXT I
360 IF F=0 THEN 270

370 IF P=0 THEN PRINT "NO SOLUTION"
380 END

I wonder what other methods or solvers are applicable. Because it is a discrete problem with discrete solutions, rather than a physics-based continuous model, it smells like a bin-packing problem (that's NP-hard), but on the other hand, it's not the same thing.

- Rob

EDIT: moved a loop to compute the next burn sequence, so that the first burn sequence is not omitted when it could be a solution. Also added game code.