# 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])
