Python程序设计期末作业(共6页).doc

《Python程序设计期末作业(共6页).doc》由会员分享，可在线阅读，更多相关《Python程序设计期末作业(共6页).doc（6页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、精选优质文档-倾情为你奉上球-弹簧物理模型模拟程序一、 编写目的1.模拟一个非线性动力系统，体会系统参数的混沌现象。2.熟悉VPython的使用。二、 物理模型描述模型空间的四周及顶面为刚性壁，底为一有质量的刚性平板，平板由一根轻弹簧支撑。模型空间内有数个刚性球作自由运动。球的初始位置在空间内均匀分布，初始速度满足正态分布。所有碰撞均为完全弹性碰撞。并假设平板只能作Z方向一维运动。三、 设计思想由物理学中的能量守恒和动量守恒推导碰撞前后物体的速度变化。每次迭代先使物体状态按时间小量dt发展，再判断是否有碰撞发生，若有则由前面推导的公式改变相关物体的速度。以适当的速率重复上述步骤，再用visua

2、l库绘制出对应图形，则得到这个物理模型的模拟演示动画。四、 使用方法直接运行程序，则开始模拟并实时绘制。也可修改程序开始处的几个全局变量，改变模型的参数。五、 参考资料VPython官方文档六、 说明在Python2.7下调试通过。本程序中没有一行代码直接来自网络、书籍等处，参考资料仅为VPython官方文档。附截图及源代码。程序运行效果截图from visual import *from random import randomfrom math import sin,cos,sqrt,log,pi#物理模型的参数以全局变量形式定义g=9.81#重力加速度k=10#弹簧胡克系数m1=1#球质

3、量m2=10#平板质量N=10#球数R=1#球半径#物体基类，提供设置、返回位置、速度信息，及迭代一个时间小量dt的方法class object: def setVel(self,v): self.v=v def setPos(self,p): self.p=p def getPos(self): return self.p def getVel(self): return self.v def go(self,dt): self.v+=self.a*dt self.p+=self.v*dt#球，继承于物体类class C_Ball(object): def _init_(self,p,v):

4、 self.p=p self.v=v self.a=vector(0,0,-g)#平板，继承于物体类，有自己的迭代方法，即加速度与位置有关class C_Board(object): def _init_(self): self.p=vector(0,0,0) self.v=vector(0,0,0) self.a=vector(0,0,0) def go(self,dt): self.a=self.p*-k self.v+=self.a*dt self.p+=self.v*dt#生成随机位置向量（模型空间内均匀分布）def randPos(): x=(random()-0.5)*2*(15-

5、R-0.5) y=(random()-0.5)*2*(15-R-0.5) z=random()*(30-R*2-0.5)+(R+0.5) return vector(x,y,z)#生成三维标准正态分布向量，使用Box-Muller公式def randVel(): u1,u2=random(),random() x=sqrt(-2*log(u1)*cos(2*pi*u2) y=sqrt(-2*log(u1)*sin(2*pi*u2) u1,u2=random(),random() z=sqrt(-2*log(u1)*cos(2*pi*u2) return vector(x,y,z)#物理模型，含

6、有一个平板和N个球class C_System: #放置球，并检查冲突 def _init_(self): self.ball= for i in range(N): self.ball.append(C_Ball(randPos(),randVel()*3) j=0 while ji:#若有两个球干涉则重新放置 n=self.balli.getPos()-self.ballj.getPos() if mag(n)=2*R: self.ball-1=C_Ball(randPos(),randVel()*3) j=-1 j+=1 self.board=C_Board() #迭代一个时间小量dt

7、def go(self,dt): for i in range(N): self.balli.go(dt) self.board.go(dt) for i in range(N): #检查与平板的碰撞 if self.balli.getPos().z-R=30-0.5: self.balli.go(-dt) v=self.balli.getVel() n=vector(0,0,1) l=dot(n,v)*2*n self.balli.setVel(v-l) #检查与前后壁的碰撞 if abs(self.balli.getPos().x)+R=15: self.balli.go(-dt) v=s

8、elf.balli.getVel() n=vector(1,0,0) l=dot(n,v)*2*n self.balli.setVel(v-l) #检查与左右壁的碰撞 if abs(self.balli.getPos().y)+R=15: self.balli.go(-dt) v=self.balli.getVel() n=vector(0,1,0) l=dot(n,v)*2*n self.balli.setVel(v-l) #检查两球之间的碰撞 for j in range(i): n=self.balli.getPos()-self.ballj.getPos() if mag(n)=2*R

9、: self.balli.go(-dt) self.ballj.go(-dt) n=norm(n) v1=self.balli.getVel() v2=self.ballj.getVel() l1=dot(n,v1)*2*n l2=dot(n,v2)*2*n self.balli.setVel(v1-l1) self.ballj.setVel(v2-l2) #返回第i个球的位置向量 def getBallPos(self,i): return self.balli.getPos() #返回平板位置向量 def getBoardPos(self): return self.board.getPo

10、s()mySys=C_System()scene=display(title=Animation, x=300,y=100, width=480,height=640, center=(0,0,0), autocenter=False, range=(40,40,60), up=(0,0,1), forward=(-8,-1,-2)B=#球数组for i in range(N): B.append(sphere(pos=mySys.getBallPos(i),radius=R,material=materials.plastic)P=box(pos=mySys.getBoardPos(),le

11、ngth=28.5,height=28.5,width=1,material=materials.wood)#平板H=helix(pos=(0,0,-30),axis=(0,0,1),coils=10,radius=7.5,length=30,material=materials.chrome)#弹簧#生成三面壁和上壁WallL=box(pos=(0,-15,0),length=30,height=1,width=59,material=materials.rough)WallR=box(pos=(0,15,0),length=30,height=1,width=59,material=mat

12、erials.rough)WallF=box(pos=(-15,0,0),length=1,height=30,width=59,material=materials.rough)WallU=box(pos=(0,0,30),length=31,height=31,width=1,material=materials.rough)while True: for i in range(100): mySys.go(4e-4)#迭代100次，每次0.4毫秒，则系统状态更新40毫秒 #根据物理模型重新绘制物体 for i in range(N): Bi.pos=mySys.getBallPos(i) P.pos=mySys.getBoardPos() H.length=P.pos.z+30 rate(50)专心-专注-专业