Solving the Reve's puzzle

I'm trying to solve the Reve's puzzle (a variant of the Tower of Hanoi puzzle with four pegs). This is what I've come up with so far:

def solve(disk, source, temppeg1, temppeg2, destination):
    if disk == 1:
        print((source, destination))
    elif disk == 2:
        print((source, temppeg1))
        print((source, destination))
        print((temppeg1, destination))
    else:
        solve(disk - 2, source, temppeg2, destination, temppeg1)
        print((source, temppeg2))
        print((source, destination))
        print((temppeg2, destination))
        solve(disk - 2,  temppeg1, source, temppeg2, destination)

>>> solve(16, 'a', 'b', 'c', 'd')
('a', 'c')
... 763 lines deleted ...
('a', 'd')

While this works fine, I was wondering if there's a better solution (more optimal) than what I have so far. Currently with what I have, I can solve 16 disk in 765 moves.

if you are just trying to find the minimum number of moves and not necessarily a solution you can use the Frame–Stewart algorithm that you linked to earlier

this builds up a solution to the number of moves to achieve a solution.

def FrameStewart(ndisks,npegs):
    if ndisks ==0: #zero disks require zero moves
        return 0
    if  ndisks == 1 and npegs > 1: #if there is only 1 disk it will only take one move
        return 1
    if npegs == 3:#3 pegs is well defined optimal solution of 2^n-1
        return 2**ndisks - 1
    if npegs >= 3 and ndisks > 0:
        potential_solutions = (2*FrameStewart(kdisks,npegs) + FrameStewart(ndisks-kdisks,npegs-1) for kdisks in range(1,ndisks))
        return min(potential_solutions) #the best solution
    #all other cases where there is no solution (namely one peg, or 2 pegs and more than 1 disk)
    return float("inf") 

print FrameStewart(16,4) #prints 161

this tells us that the optimal solution for 16 disks and 4 pegs is 161 moves, note that it does not tell us what those moves are

if you actually need the moves you will have to heavily modify this solution.

start by solving the towersofhanoi with 3 pegs as that is the traditional layout and has well defined algorithms to solve

def towers3(ndisks,start=1,target=3,peg_set=set([1,2,3])):
   if ndisks == 0 or start == target: #if there are no disks, or no move to make
      return [] #no moves
   my_move = "move(%s,%s)"%(start,target) 
   if ndisks == 1: #trivial case if there is only one disk, just move it
      return [my_move]
   helper_peg = peg_set.difference([start,target]).pop()
   moves_to_my_move = towers3(ndisks-1,start,helper_peg)
   moves_after_my_move = towers3(ndisks-1,helper_peg,target)
   return moves_to_my_move + [my_move] + moves_after_my_move

you can easily verify that this is returning optimal solutions by knowing that the optimal solution to the towers of hanoi with 3 pegs is 2ndisks - 1

all thats left is to change FrameStewart to also return moves

def FrameStewartSolution(ndisks,start=1,end=4,pegs=set([1,2,3,4])):
    if ndisks ==0 or start == end: #zero disks require zero moves
        return []
    if  ndisks == 1 and len(pegs) > 1: #if there is only 1 disk it will only take one move
        return ["move(%s,%s)"%(start,end)]  
    if len(pegs) == 3:#3 pegs is well defined optimal solution of 2^n-1
        return towers3(ndisks,start,end,pegs)
    if len(pegs) >= 3 and ndisks > 0:
        best_solution = float("inf")
        best_score = float("inf")
        for kdisks in range(1,ndisks):
            helper_pegs = list(pegs.difference([start,end]))
            LHSMoves = FrameStewartSolution(kdisks,start,helper_pegs[0],pegs)
            pegs_for_my_moves = pegs.difference([helper_pegs[0]]) # cant use the peg our LHS stack is sitting on
            MyMoves = FrameStewartSolution(ndisks-kdisks,start,end,pegs_for_my_moves) #misleading variable name but meh 
            RHSMoves = FrameStewartSolution(kdisks,helper_pegs[0],end,pegs)#move the intermediat stack to 
            if any(move is None for move in [LHSMoves,MyMoves,RHSMoves]):continue #bad path :(
            move_list = LHSMoves + MyMoves + RHSMoves
            if(len(move_list) < best_score):
                best_solution = move_list
                best_score = len(move_list)
        if best_score < float("inf"):       
            return best_solution
    #all other cases where there is no solution (namely one peg, or 2 pegs and more than 1 disk)
    return None

note that this is going to be much slower than the version that does not need to find the actual solution (this being codereview maybe some folks have suggestions to make it run faster)
Timings from this experiment

towers3(16)  # 0.09 secs
FrameStewart(16) #0.04 secs
FrameStewartSolution(16) #67.04 secs!!!

really slow as you can see

you can speed it up alot by memoizing it

import json

def fsMemoizer(f): #just a junky quick memoizer
    cx = {}
    def f2(*args):
        try:
            key= json.dumps(args)
        except:
            key =json.dumps(args[:-1] + (sorted(list(args[-1])),))
        if key not in cx:
            cx[key] = f(*args)
        return cx.get(key)
    return f2
@fsMemoizer
def FrameStewartSolution(ndisks,start=1,end=4,pegs=set([1,2,3,4])):
    ...

after memoization the time to calculate became much faster (less than a second)

StackExchange Code Review Q#42524, answer score: 6