""" Algorithms to solve the rod-cutting problem from CLRS 15.1. Author: Michael Goldwasser 1) cut_rod exponential-time, top-down recursive approach 2) memoized_cut_rod O(n^2)-time top-down dynamic programming approach 3) bottom_up_cut_rod: O(n^2)-time bottom-up dynamic programming approach 4) extended_bottom_up_cut_rod: O(n^2)-time bottom-up dynamic programming approach with additional trace of solution 5) report_cut_rod_solution: O(n)-time function for reporting actual solution given extended trace 6) reconstruct_solution; O(n^2)-time algorithm to reconstruct solution using only table of revenues from original dynamic programming implementation. Use -h flag for documentation on usage. """ #----------------------------------------------------------------------- def cut_rod(p, n): """ Exponential-time, top-down recursive approach. p is assumed to be an array of size n+1 such that p[j] defines price for length j (with p[0] = 0). """ if n == 0: return 0 else: q = float('-inf') # Consider a first cut of length i, for i from 1 to n inclusive for i in range(1, n+1): q = max(q, p[i] + cut_rod(p, n-i)) return q #----------------------------------------------------------------------- def memoized_cut_rod(p, n, r=None): """ Optional parameter r is an existing table of known revenues, such that r[j] is maximum revenue for a rod of length j. """ # create empty table, if none given if r is None: r = [None] * (n+1) # use known answer, if one exists if r[n] is not None: return r[n] # otherwise use original recursive formulation if n == 0: q = 0 else: q = float('-inf') # Consider a first cut of length i, for i from 1 to n inclusive for i in range(1, n+1): q = max(q, p[i] + memoized_cut_rod(p, n-i, r)) # recur on length n-i # memoize answer in table before returning r[n] = q return q #----------------------------------------------------------------------- def bottom_up_cut_rod(p, n, r=None): """ Optional parameter r is an existing table of known revenues, such that r[i] is maximum revenue for a rod of length i. """ if r is None: r = [0] * (n+1) # will be list of maximum revenues. r[0] = 0 for j in range(1, n+1): q = float('-inf') for i in range(1, j+1): q = max(q, p[i] + r[j-i]) r[j] = q return r[n] #----------------------------------------------------------------------- def extended_bottom_up_cut_rod(p, n): """ Returns lists r and s, where r[i] is maximum revenue for rod of length i and s[i] is size of smallest piece used in maximizing revenue for such rod. """ r = [0] * (n+1) # will be list of maximum revenues. r[0] = 0 s = [0] * (n+1) # s[i] is size of smallest piece in optimal solution for i for j in range(1, n+1): q = float('-inf') for i in range(1, j+1): if p[i] + r[j-i] > q: # found a better solution q = p[i] + r[j-i] s[j] = i r[j] = q return (r, s) # return complete tables #----------------------------------------------------------------------- def report_cut_rod_solution(n, s): """Return optimal list of rod sizes given table 's' from extended algorithm.""" pieces = [] while n > 0: cut = s[n] pieces.append(cut) n = n - cut # continue on what remains return pieces #----------------------------------------------------------------------- def reconstruct_solution(p, n, r): """Build optimal list of rod sizes given only revenue table form original dynamic program.""" pieces = [] while n > 0: q = float('-inf') for i in range(1, n+1): if p[i] + r[n-i] > q: # this is a better solution q = p[i] + r[n-i] cut = i pieces.append(cut) n = n - cut return pieces #----------------------------------------------------------------------- def wrap_extended(p, n): r,s = extended_bottom_up_cut_rod(p,n) pieces = report_cut_rod_solution(n,s) return r,pieces #----------------------------------------------------------------------- def wrap_reconstruct(p, n): r = [None] * (n+1) memoized_cut_rod(p,n,r) pieces = reconstruct_solution(p, n, r) return r,pieces #----------------------------------------------------------------------- if __name__ == '__main__': from optparse import OptionParser, OptionGroup import sys import time import random distributions = { 'gauss' : (lambda:random.gauss(options.density, 0.5)), 'pareto' : (lambda:options.density*random.paretovariate(10)), } algorithms = { 'recursion' : (cut_rod, None), 'memoize' : (memoized_cut_rod, wrap_reconstruct), 'bottomup' : (bottom_up_cut_rod, wrap_extended), } def quit(message=None): if message: print(message) sys.exit(1) parser = OptionParser(usage='usage: %prog [options]') parser.add_option('-n', dest='n', metavar='LENGTH', type='int', default=10, help='length of original rod [default: %default]') parser.add_option('-P', dest='echo', default=False, action='store_true', help='echo price list to console [default: %default]') group = OptionGroup(parser, "Algorithm Options", 'Available algorithms: ' + ', '.join(sorted(algorithms.keys()))) group.add_option('-a', dest='algorithm', default='bottomup', choices=sorted(algorithms.keys()), help='algorithm choice [default: %default]') group.add_option('-V', dest='verbose', default=False, action='store_true', help='generate verbose solution [default: %default]') parser.add_option_group(group) group = OptionGroup(parser, "Pricing Options") group.add_option('-f', dest='file', default=None, help='read (whitespace separated) price list from file') group.add_option('-d', dest='density', type='int', default=10, help='average price per unit [default: %default]') group.add_option('-s', dest='seed', type='int', default=None, help='random seed for generating prices [default: %default]') group.add_option('-r', dest='distribution', default='gauss', choices=sorted(distributions.keys()), help='random distribution function for price densities [default: %default]') parser.add_option_group(group) group = OptionGroup(parser, "Available Distributions", ', '.join(sorted(distributions.keys()))) parser.add_option_group(group) (options,args) = parser.parse_args() if options.n <= 0: quit('n must be positive') alg = algorithms[options.algorithm][options.verbose] if alg is None: print('Verbose solutinos not implemented for algorithm {0}'.format(options.algorithm)) sys.exit(1) if options.file is not None: try: raw = open(options.file,'r').read() except IOError: print('Unable to open file {0}'.format(options.file)) sys.exit(1) try: prices = [0] + [int(token) for token in raw.split()] except ValueError: print('Prices must be integers.') sys.exit(1) if len(prices) < 1 + options.n: pad = 1 + options.n - len(prices) prices.extend([0] * pad) # longer pieces worth 0 else: if options.seed is None: options.seed = random.randrange(1000000) print('Using seed: {0}'.format(options.seed)) random.seed(options.seed) prices = [int(k*distributions[options.distribution]()) for k in range(1 + options.n)] if options.echo: print('Prices: {0}'.format(repr(prices))) result = alg(prices, options.n) if not options.verbose: print('optimal revenue is {0}'.format(result)) else: rev,pieces = result print('optimal revenue is {0}'.format(rev[options.n])) print('cut into pieces with size: {0}'.format(', '.join(str(p) for p in pieces)))