集体编程智慧--第五章：优化算法

import math
people=[('Seymour','BOS'),
('Franny','DAL'),
('Zooey','CAK'),
('Walt','MIA'),
('Buddy','ORD'),
('Les','OMA')]
#New York的LaGuardia机场
destination='LGA'
#将数据载入到一个字典中，以起止点为键，以可能的航班详情明细为值
flights={}
#
for line in file('E:/code/schedule.txt'):
origin,dest,depart,arrive,price=line.strip().split(',')
flights.setdefault((origin,dest),[])
#将航班详情添加到航班列表中
flights[(origin,dest)].append((depart,arrive,int(price)))
def getminutes(t):
x=time.strptime(t,'%H:%M')
return x[3]*60+x[4]
def printschedule(r):
for d in range(len(r)/2):
name=people[d][0]
origin=people[d][1]
out=flights[(origin,destination)][r[2*d]]
ret=flights[(destination,origin)][r[2*d+1]]
print '%10s%10s %5s-%5s $%3s %5s-%5s $%3s' % (name,origin,
out[0],out[1],out[2],
ret[0],ret[1],ret[2])

>>> import optimization

>>> s=[1,4,3,2,7,3,6,3,2,4,5,3]

>>> optimization.printschedule(s)
   Seymour       BOS  8:04-10:11 $ 95 12:08-14:05 $142

    Franny       DAL 10:30-14:57 $290  9:49-13:51 $229

     Zooey       CAK 17:08-19:08 $262 10:32-13:16 $139

      Walt       MIA 15:34-18:11 $326 11:08-14:38 $262

     Buddy       ORD  9:42-11:32 $169 12:08-14:47 $231

       Les       OMA 13:37-15:08 $250 11:07-13:24 $171

def schedulecost(sol):
totalprice=0
laterstarrival=0
earliestdep=24*60
for d in range(len(sol)/2):
#得到往程航班和返程航班
origin=people[d][1]
outbound=flights[(origin,destination)][int(sol[2*d])]
returnf=flights[(destination,origin)][int(sol[(2*d)+1])]
#总价格等于所有往程航班和返程航班价格之和
totalprice+=outbound[2]
totalprice+=returnf[2]
#记录最晚达到时间和最早离开时间
if laterstarrival<getminutes(outbound[1]):laterstarrival=getminutes(outbound[1])
if earliestdep>getminutes(returnf[0]):earliestdep=getminutes(returnf[0])
#每个人必须在机场等待直到最后一个人到达为止
#他们也必须在相同时间到达，并等候他们的返程航班
totalwait=0
for d in range(len(sol)/2):
origin=people[d][1]
outbound=flights[(origin,destination)][int(sol[2*d])]
returnf=flights[(destination,origin)][int(sol[2*d+1])]
totalwait+=laterstarrival-getminutes(outbound[1])
totalwait+=getminutes(returnf[0])-earliestdep
#这个题解要求多付一天的汽车租用费用吗？如何是，则费用为50美元
if laterstarrival<earliestdep:totalprice+=50
return totalprice+totalwait

>>> reload(optimization)

<module 'optimization' from 'E:\python\optimization.py'>

>>> s=[1,4,3,2,7,3,6,3,2,4,5,3]

>>> optimization.schedulecost(s)
3242
随机尝试法：

def randomoptimize(domian,costf):
best=999999999
bestr=None
for i in range(1000):
#创建一个随机解
r=[random.randint(domian[i][0],domian[i][1])
for i in range(len(domian))]
#得到成本
cost=costf(r)
#与到目前为止的最优解进行比较
if cost<best:
best=cost
bestr=r
return r>>> reload(optimization)

<module 'optimization' from 'E:\python\optimization.pyc'>

>>> domain=[(0,9)]*(len(optimization.people)*2)

>>> s=optimization.randomoptimize(domain,optimization.schedulecost)

>>> s=[1,4,3,2,7,3,6,3,2,4,5,3]

>>> optimization.schedulecost(s)
3242

>>> optimization.printschedule(s)
   Seymour       BOS  8:04-10:11 $ 95 12:08-14:05 $142

    Franny       DAL 10:30-14:57 $290  9:49-13:51 $229

     Zooey       CAK 17:08-19:08 $262 10:32-13:16 $139

      Walt       MIA 15:34-18:11 $326 11:08-14:38 $262

     Buddy       ORD  9:42-11:32 $169 12:08-14:47 $231

       Les       OMA 13:37-15:08 $250 11:07-13:24 $171

爬山法
def hillclimb(domain,costf):
#创建一个随机解：
sol=[random.randint(domain[i][0],domain[i][1])
for i in range(len(domain))]
#主循环
while 1:
#创建相邻的列表
neighbors=[]
for j in range(len(domain)):
#在每个方向上相对于原值偏离一点
if sol[j]>domain[j][0]:
neighbors.append(sol[0:j]+[sol[j]-1]+sol[j+1:])
if sol[j]<domain[j][1]:
neighbors.append(sol[0:j]+[sol[j]+1]+sol[j+1:])
#在相邻解中寻找最优解
current=costf(sol)
best=current
for j in range(len(neighbors)):
cost=costf(neighbors[j])
if cost<best:
best=cost
sol=neighbors[j]
#如果没有更好的解，退出循环
if best==current:
break
return sol>>> reload(optimization)

<module 'optimization' from 'E:\python\optimization.py'>

>>> domain=[(0,9)]*(len(optimization.people)*2)

>>> s=optimization.hillclimb(domain,optimization.schedulecost)

>>> optimization.schedulecost(s)
1821

>>> optimization.printschedule(s)

   Seymour       BOS 20:17-22:22 $102  6:39- 8:09 $ 86

    Franny       DAL  6:12-10:22 $230  9:49-13:51 $229

     Zooey       CAK  8:27-10:45 $139 13:37-15:33 $142

      Walt       MIA  9:15-12:29 $225 15:23-18:49 $150

     Buddy       ORD 14:22-16:32 $126  7:50-10:08 $164

       Les       OMA  9:15-12:03 $ 99 15:07-17:21 $129

        爬山法的缺陷是在局部范围内寻找最小值。解决方法是随机重复爬山法，生成多个随机的初始解，希望有一个解能逼近全局的最小值：
模拟退火算法

def annealingoptimize(domain,costf,T=10000.0,cool=0.98,step=1):
#随机初始化值
vec=[float(random.randint(domain[i][0],domain[i][1]))
for i in range(len(domain))]
while T>0.1:
#选择一个索引值
i=random.randint(0,len(domain)-1)
#选择一个改变索引值的方向
dir=random.randint(-step,step)
#创建一个代表题解的新列表，改变其中一个值
vecb=vec[:]
vecb[i]+=dir
if vecb[i]<domain[i][0]:vecb[i]=domain[i][0]
elif vecb[i]>domain[i][1]:vecb[i]=domain[i][1]
#计算当前成本和新的成本
ea=costf(vec)
eb=costf(vecb)
#它是更好的解吗？或者是趋向最优解的可能的临界解吗？
if (eb<ea or random.random()<pow(math.e,-(eb-ea)/T)):
vec=vecb
#降低温度
T=T*cool
return vec>>> reload(optimization)

<module 'optimization' from 'E:\python\optimization.py'>

>>> domain=[(0,9)]*(len(optimization.people)*2)

>>> s=optimization.annealingoptimize(domain,optimization.schedulecost)

>>> optimization.schedulecost(s)
2094

>>> optimization.printschedule(s)

Traceback (most recent call last):

  File "<pyshell#34>", line 1, in <module>

    optimization.printschedule(s)

  File "E:\python\optimization.py", line 29, in printschedule

    ret=flights[(destination,origin)][r[2*d+1]]

TypeError: list indices must be integers, not float

def geneticoptmize(domain,costf,popsize=50,step=1,mutprob=0.2,elite=0.2,maxiter=100):
#变异操作
def mutate(vec):
i=random.randint(0,len(domain)-1)
if random.random()<0.5 and vec[i]>domain[i][0]:
return vec[0:i]+[vec[i]-step]+vec[i+1:]
elif vec[i]<domain[i][1]:
return vec[0:i]+[vec[i]+step]+vec[i+1:]
#交叉操作
def crossover(r1,r2):
i=random.randint(1,len(domain)-2)
return r1[0:i]+r2[i:]
#构造初始种群
pop=[]
for i in range(popsize):
vec=[random.randint(domain[i][0],domain[i][1])
for i in range(len(domain))]
pop.append(vec)
#每一代中有多少胜出者？
topelite=int(elite*popsize)
#主循环
for i in range(maxiter):
scores=[(costf(v),v) for v in pop]
scores.sort()
ranked=[v for (s,v) in scores]
#从纯粹的胜出者开始
pop=ranked[0:topelite]
#添加变异和配对后的胜出者
while len(pop)<popsize:
if random.random()<mutprob:
#变异
c=random.randint(0,topelite)
pop.append(mutate(ranked[c]))
else:
#交叉
c1=random.randint(0,topelite)
c2=random.randint(0,topelite)
pop.append(crossover(ranked[c1],ranked[c2]))
#打印当前最优值
print scores[0][0]
return scores[0][1]

>>> reload(optimization)

<module 'optimization' from 'E:\python\optimization.py'>

>>> domain=[(0,9)]*(len(optimization.people)*2)

>>> s=optimization.annealingoptimize(domain,optimization.schedulecost)

>>> optimization.printschedule(s)

Traceback (most recent call last):

  File "<pyshell#38>", line 1, in <module>

    optimization.printschedule(s)

  File "E:\python\optimization.py", line 29, in printschedule

    ret=flights[(destination,origin)][r[2*d+1]]

TypeError: list indices must be integers, not float

学生宿舍问题
#-*- coding:utf-8 -*-
import random
import math
#代表宿舍，每个宿舍有两个可用的隔阂
dorms=['A','B','C','D','E']
#代表学生及其首选和次选
prefs=[('Toby',('D','C')),('Steve',('A','E')),('Andrea',('B','A')),('Sarah',('A','E')),
('Dave',('B','D')),('Jeff',('C','E')),('Fred',('E','B')),('Suzie',('D','C')),('Laura',('D','C')),
('Neil',('C','B'))]
#[(0,9),(0,8),(0,7),(0,6),...,(0,0)]
domain=[(0,(len(dorms)*2)-i-1) for i in range(0,len(dorms)*2)]
def printsolution(vec):
slots=[]
#为每个宿舍建两个槽
for i in range(len(dorms)):slots+=[i,i]
#遍历每一名学生的安置情况
for i in range(len(vec)):
x=int(vec[i])
#从剩余槽中选择
dorm=dorms[slots[x]]
#输出学生及其被分配的宿舍
print prefs[i][0],dorm
#删除该槽
del slots[x]>>> import dorm

>>> dorm.printsolution([0,0,0,0,0,0,0,0,0,0])
Toby A

Steve A

Andrea B

Sarah B

Dave C

Jeff C

Fred D

Suzie D

Laura E

Neil E

def dormcost(vec):
cost=0
#建立一个槽序列
slots=[0,0,1,1,2,2,3,3,4,4]
#遍历每一名学生
for i in range(len(vec)):
x=int(vec[i])
dorm=dorms[slots[x]]
pref=prefs[i][1]
#首选成本值为0，次选成本值为1
if pref[0]==dorm:cost+=0
elif pref[1]==dorm:cost+=1
else:cost+=3
#不在选择之列则成本值加3
#删除选中的槽
del slots[x]
return cost

>>> reload(dorm)

<module 'dorm' from 'E:\python\dorm.py'>

>>> import optimization

>>> s=optimization.randomoptimize(dorm.domain,dorm.dormcost)

>>> dorm.dormcost(s)
19

>>> optimization.geneticoptimize(dorm.domain,dorm.dormcost)

Traceback (most recent call last):

  File "<pyshell#9>", line 1, in <module>

    optimization.geneticoptimize(dorm.domain,dorm.dormcost)

AttributeError: 'module' object has no attribute 'geneticoptimize'

>>> optimization.geneticoptmize(dorm.domain,dorm.dormcost)
4

4

Traceback (most recent call last):

  File "<pyshell#10>", line 1, in <module>

    optimization.geneticoptmize(dorm.domain,dorm.dormcost)

  File "E:\python\optimization.py", line 145, in geneticoptmize

    scores=[(costf(v),v) for v in pop]

  File "E:\python\dorm.py", line 30, in dormcost

    for i in range(len(vec)):

TypeError: object of type 'NoneType' has no len()

>>> dorm.printsolution(s)
Toby C

Steve E

Andrea D

Sarah E

Dave C

Jeff B

Fred B

Suzie A

Laura D

Neil A

网络可视化：质点弹簧算法，各结点彼此向对方施以推力并试图分离，而结点间的连接则试图将关联结点彼此拉近。这样一来网络会逐渐呈现出未关联的结点被推离，而关联的结点则被彼此拉近，却又不会靠的很拢。但无法避免交叉线
#-*- coding:utf-8 -*-
import math
people=['Charlie','Augustus','Veruca','Violet','Mike','Joe','Willy','Miranda']
links=[('Augustus','Willy'),('Mike','Joe'),('Miranda','Mike'),('Violet','Augustus'),
('Miranda','Willy'),('Charlie','Mike'),('Veruca','Joe'),('Miranda','Augustus'),
('Willy','Augustus'),('Joe','Charlie'),('Veruca','Augustus'),('Miranda','Joe')]
def crosscount(v):
#将数字序列转换成一个person:(x,y)的字典
loc=dict([(people[i],(v[i*2],v[i*2+1])) for i in range(0,len(people))])
total=0
#遍历每一对连线
for i in range(len(links)):
for j in range(i+1,len(links)):
#获取坐标位置
(x1,y1),(x2,y2)=loc[links[i][0]],loc[links[i][1]]
(x3,y3),(x4,y4)=loc[links[j][0]],loc[links[j][1]]
den=(y4-y3)*(x2-x1)-(x4-x3)*(y2-y1)
#如果两线平行，则den==0
if den==0:continue
#否则，ua与ub就是两条交叉线的分数值
ua=((x4-x3)*(y1-y3)-(y4-y3)*(x1-x3))/den
ub=((x2-x1)*(y1-y3)-(y2-y1)*(x1-x3))/den
#如果两条线的分数值介于0和1之间，则两线彼此交叉
if ua>0 and ua<1 and ub>0 and ub<1:
total+=1
return total
>>> import socialnetwork

>>> import optimization
>>> domain=[(10,370)]*(len(people)*2)

>>> sol=optimization.geneticoptmize(domain,crosscount)

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

>>> print sol
[10, 193, 150, 137, 10, 284, 220, 37, 165, 56, 16, 64, 72, 236, 332, 309]
from PIL import Image,ImageDraw
def drawnetwork(sol):
#建立image对象
img=Image.new('RGB',(400,400),(255,255,255))
draw=ImageDraw.Draw(img)
#建立标示位置信息的字典
pos=dict([(people[i],(sol[i*2],sol[i*2+1])) for i in range(0,len(people))])
#绘制连线
for (a,b) in links:
draw.line((pos[a],pos[b]),fill=(255,0,0))
#绘制代表人的结点
for n,p in pos.items():
draw.text(p,n,(0,0,0))
img.show(
9e6a
)

>>> reload(socialnetwork)

<module 'socialnetwork' from 'E:\python\socialnetwork.py'>

>>> import optimization

>>> sol=optimization.randomoptimize(socialnetwork.domain,socialnetwork.crosscount)

>>> socialnetwork.crosscount(sol)

4.638490077108166

>>> sol

[149, 288, 155, 290, 307, 320, 226, 256, 330, 219, 316, 205, 87, 173, 172, 97]

>>> print socialnetwork.drawnetwork(sol)

None
这是为什么呢？