用python做数值计算，spyder 上的调试分形算法

smtmobly · 发表于 2010-11-2 20:53

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。

您需要登录才可以下载或查看，没有账号？我要加入

x

再发一个分形的代码，这是《用python做数值计算》上的例子，是迭代法计算描绘分形图形的。先看一下效果。分形.png

这是spyder上界面，并且窗口都是浮动的，可以任意调节，这个界面是用PyQt写的，pythonxy上有完整的设计。关于pyqt界面设计，也有很多介绍。但是我个人觉得如果是自己做一些简单的使用用tkinter就比较不错。先不说那么多，看看代码吧。

# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as pl
import time
# 蕨类植物叶子的迭代函数和其概率值
eq1 = np.array([[0,0,0],[0,0.16,0]])
p1 = 0.01
eq2 = np.array([[0.2,-0.26,0],[0.23,0.22,1.6]])
p2 = 0.07
eq3 = np.array([[-0.15, 0.28, 0],[0.26,0.24,0.44]])
p3 = 0.07
eq4 = np.array([[0.85, 0.04, 0],[-0.04, 0.85, 1.6]])
p4 = 0.85
def ifs(p, eq, init, n):
"""
进行函数迭代
p: 每个函数的选择概率列表
eq: 迭代函数列表
init: 迭代初始点
n: 迭代次数
返回值：每次迭代所得的X坐标数组， Y坐标数组，计算所用的函数下标
"""
# 迭代向量的初始化
pos = np.ones(3, dtype=np.float)
pos[:2] = init
# 通过函数概率，计算函数的选择序列
p = np.add.accumulate(p)
rands = np.random.rand(n)
select = np.ones(n, dtype=np.int)*(n-1)
for i, x in enumerate(p[::-1]):
select[rands<x] = len(p)-i-1
# 结果的初始化
result = np.zeros((n,2), dtype=np.float)
c = np.zeros(n, dtype=np.float)
for i in xrange(n):
eqidx = select[i] # 所选的函数下标
tmp = np.dot(eq[eqidx], pos) # 进行迭代
pos[:2] = tmp # 更新迭代向量
# 保存结果
result[i] = tmp
c[i] = eqidx
return result[:,0], result[:, 1], c
start = time.clock()
x, y, c = ifs([p1,p2,p3,p4],[eq1,eq2,eq3,eq4], [0,0], 100000)
print time.clock() - start
pl.figure(figsize=(6,6))
pl.subplot(121)
pl.scatter(x, y, s=1, c="g", marker="s", linewidths=0)
pl.axis("equal")
pl.axis("off")
pl.subplot(122)
pl.scatter(x, y, s=1,c = c, marker="s", linewidths=0)
pl.axis("equal")
pl.axis("off")
pl.subplots_adjust(left=0,right=1,bottom=0,top=1,wspace=0,hspace=0)
pl.gcf().patch.set_facecolor("white")
pl.show()

复制代码

使用pythonxy，你会发现，绘图时，你喜欢用什么都可以，甚至可以使用gnuplot这个在tex界大名鼎鼎的绘图包。

smtmobly · 发表于 2010-11-2 22:00

对于matlab如http://forum.vibunion.com/thread-96985-1-1.html
帖子上的问题，我们可以使用如下code得到：

# -*- coding: utf-8 -*-
from numpy import sin,pi
import numpy as np
import pylab as pl
x=np.arange(0,2*pi,0.1)
pl.figure(figsize=(8,6))
for i in range(1,6):
pl.subplot(510+i)
pl.plot(x, i*sin(x+i*pi/2),'.')
pl.legend()
pl.show()

复制代码

Rainyboy · 发表于 2010-11-2 22:26

回复 smtmobly 的帖子

import numpy as np
import pylab as pl

复制代码

为什么要这么import呢，难道还有和pylab有相同接口的其他绘图库类，所以这么做方便重构？

smtmobly · 发表于 2010-11-2 22:31

python的结构是松散的，模块与模块间的联系通过import来连接，除了内建的几个基本的模块外不会加载其他模块，在使用的时候必须加载进去才能使用，比如math中有sin，cos等，numpy里也有同样名称的函数，他们的计算方式运行速度等都不一样。
松散结构可以让程序不会因为加载了不使用的模块而影响速度

Rainyboy · 发表于 2010-11-2 22:44

回复 smtmobly 的帖子

我认为松散的结构是为了使程序不会因为加载了不同的模块而作过多的改写吧？

smtmobly · 发表于 2010-11-2 22:50

当然如果你加载过多的不使用的模块，在改写的时候会麻烦些。改写上，python有更好的优势，因为函数，类，等特别的结构都可以作为变量赋值。
所以只要一句赋值语句就可以改写。
最需要注意的一一个关于list的，比如：x=[1,2,3]
y=x
y[2]=1
那么x也会变成[1,2,1]

Rainyboy · 发表于 2010-12-6 18:40

回复 6 # smtmobly 的帖子

最需要注意的一一个关于list的，比如：x=[1,2,3]
y=x
y[2]=1
那么x也会变成[1,2,1]

现在对你说的这个有体会了，要想不发生这种情况，就要:
y = np.range(1,3)
y[:] = x[:]
呵呵，是吧？

wqsong · 发表于 2010-12-6 19:21

本帖最后由 wqsong 于 2010-12-6 19:47 编辑

呵呵，上次写了一个用到矩阵运算的程序。初始化一个二维数组也发现了一个问题，这里贴出来，一起看看。两种不同的初始化方式：

1、
>>> arr = [[0]*4]*4
>>> arr[0][2] = 6
>>> arr
[[0, 0, 6, 0], [0, 0, 6, 0], [0, 0, 6, 0], [0, 0, 6, 0]]
2、
>>> arr = [[0]*4 for i in range(4)]
>>> arr[0][2] = 6
>>> arr
[[0, 0, 6, 0], [0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0]]

复制代码

第一种相当于C++中常说的位拷贝，第二种就是逻辑拷贝。
这块挺诡异的，无章可循，目前是死记的。。。

Rainyboy · 发表于 2010-12-6 19:28

回复 8 # wqsong 的帖子

我也觉得PYthon序列的寻址或者叫切片吧，反正英文就是Indexing挺盘根错节的……
同死记……

账号		自动登录	找回密码
密码			我要加入

[Python] 用python做数值计算，spyder 上的调试分形算法

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。

评分

本帖被以下淘专辑推荐:

评分