# Author: James Collins # Purpose: To allocate shares of cost amongst # a group of friends who benifit by mutual cooperation. from __future__ import division from random import randint from math import factorial from copy import copy # Travelling salesman solver # global variables for tsp qstart = 0 qend = 0 queue = [] bestCost = 0 bestTrack = '' def tsp(nodes, start): global qstart global qend global queue global bestCost global bestTrack n = len(nodes) if n == 0: return [], 0 qstart = 0 qend = n - 1 queue = copy(nodes) #range(n) bestCost = 9999999 # a very large number bestTrack = '' def recur(head, last, CostSoFar, track): global bestCost global bestTrack for j in range( n-head): node = gett() track2 = track + node.initial NewCost = CostSoFar + metric(last, node) if NewCost < bestCost: if head < n-2: recur(head+1, node, NewCost, track2) else: node2 = gett() finalCost = NewCost + metric(node, node2) if finalCost < bestCost: bestCost = finalCost bestTrack = track2 + node2.initial putt(node2) putt(node) def gett(): global qstart x = queue[qstart] qstart += 1 if qstart == n: qstart = 0 return x def putt(x): global queue global qend qend += 1 if qend == n: qend = 0 queue[qend] = x recur(0, start, 0, '') # exhaustive search of all permutations return bestTrack, bestCost # MAIN PROGRAM people = [] class person: def __init__(self, name, x, y): self.name = name self.initial = name[0] self.x = x # house location North/South self.y = y # house location East/West self.bitplace = pow(2, len(people)) def isElementOf(self, n): # If the ith bit in the index of # powerset (base 2)is 1 then it contains the ith person, # or else it dosnt. if self.bitplace & n <> 0: return True return False def metric(a,b): # city grid return abs(a.x - b.x) + abs(a.y - b.y) # random locations are being used for test purposes. # Cost of travel is asssumed to be proportional to distance pub = person('Shakespear', 5, 5) people.append(person('Anne', randint(1, 10), randint(1, 10))) people.append(person('Brian', randint(1, 10), randint(1, 10))) people.append(person('Cindy', randint(1, 10), randint(1, 10))) people.append(person('David', randint(1, 10), randint(1, 10))) people.append(person('Elaine', randint(1, 10), randint(1, 10))) #people.append(person('Francis', randint(1, 10), randint(1, 10))) #people.append(person('Girlee', randint(1, 10), randint(1, 10))) #people.append(person('Him', randint(1, 10), randint(1, 10))) #people.append(person('Indy', randint(1, 10), randint(1, 10))) #people.append(person('Jim', randint(1, 10), randint(1, 10))) #people.append(person('Kathy', randint(1, 10), randint(1, 10))) # On my 3Gz Dual Core this takes 3min, mostly solving tsp lenPeople = len(people) class solution: def __init__(self, coalition): self.coalition = coalition self.size = len(coalition) self.index = len(powerset) self.split = [] if len(coalition) == 1: # individual self.cost = metric(coalition[0], pub) self.trk = coalition[0].initial coalition[0].cost = self.cost else: self.trk, self.cost = tsp(coalition, pub) print 'route: ', self.trk, 'cost>',self.cost powerset = [] # The size of the problem increases exponentialy # w.r.t the number of people. # Fortunately only about four people can fit in a taxi. for i in range(pow(2, lenPeople)): # itterating over the power set newset = [] for j in people: if j.isElementOf(i): newset.append(j) powerset.append(solution(newset)) def ShapleySub(q, p, setindex, personbit): s = q.size t = factorial(s) * factorial(lenPeople - s - 1) v = powerset[setindex + personbit] return t * ( v.cost - q.cost) # Note: The Shapley value method assumes that there is no # incentive for any party to defect from the grand # coalition. This may not be true in this scenario as it # may be more efficient for parties to travel seperately. # In this aproach I calculate the minimum possible cost for the grand coalition # allowing for seperate taxi travel and assign cost portions to each person # as if there were no strategic politicing. # For example if A can travel with B or C, but # A, B and C cannot travel together. Both B and C are advantaged and A pays less # than choosing to travel with either. So even if B ends up travelling home alone # his travel is discounted by the others because of his ability to dispute # the status quoe. # test if union cost > sum cost of partitions for setSize in range(2, lenPeople + 1): for p in powerset: if p.size == setSize: for i in range(1, pow(2, p.size-1)): # itterating all possible splits subset0s = 0 subset1s = 0 for j in range(p.size): k = pow(2, j) if k & i == 0: # note this coresponds to the bit places of the coalition not to the people list. subset0s += p.coalition[j].bitplace # bitplace coresponds to the people list else: subset1s += p.coalition[j].bitplace if powerset[subset0s].cost + powerset[subset1s].cost < p.cost: # subsets gain more by travelling seperately print 'New costing for', p.trk, powerset[subset0s].cost + powerset[subset1s].cost, 'down from', p.cost p.cost = powerset[subset0s].cost + powerset[subset1s].cost p.split = [subset0s, subset1s] print 'Shapley values' for p in people: Shapley = 0.0 personbit = p.bitplace for q in powerset: setindex = q.index if not(p.isElementOf( setindex)): Shapley += ShapleySub(q, p, setindex, personbit) Shapley = Shapley / factorial(lenPeople) diff = p.cost - Shapley if p.cost > 0: perc = int(diff / p.cost * 100) else: perc = 0 print "%s pays %6.2f saves %6.2f saving %d percent." %(p.name, Shapley, diff, perc) # Draws map print print '#######################' for i in range(1, 11): line = '# ' for j in range(1, 11): match = False if i == pub.x and j == pub.y: line += pub.initial match = True for p in people: if p.x == i and p.y == j: line += p.initial match = True if not match: line += '. ' else: line += ' ' line += '#' print line print '#######################' def AllocateTaxis(pset): if len(pset.split) == 0: print pset.trk, pset.cost else: AllocateTaxis(powerset[pset.split[0]]) AllocateTaxis(powerset[pset.split[1]]) print print 'Allocating taxis' AllocateTaxis(powerset[len(powerset)-1])