加速Python中的计算(模拟磁场中的粒子)
作者：互联网
A用Python编写的程序速度有问题.该程序是“模拟磁场中的铁磁颗粒”,更具体地说是磁惰性液体.该程序可以工作,但与用C编写的相同程序相比非常缓慢,但我用Python编写了一个项目来研究.
总的来说,程序代码基于循环,有很多浮点计算.存在随机数量的粒子(产生随机位置),其在磁场的影响下彼此相互作用.
这是最初的立场：
http://i.stack.imgur.com/T15Bb.jpg
决赛：
http://i.stack.imgur.com/0nU5D.jpg
主循环(在SymMain.py中,具有k变量)迭代是时间步长,计算当时粒子在其中的坐标和作用在其上的力(吸引力和小排斥力).为了加快速度,我想使用并行处理来同时计算迭代次数,但这是不可能的,因为一次迭代中的计算取决于前一次迭代中的计算.
我不知道Python比C慢.例如,计算一次性步骤中529个粒子的模拟(在我的计算机上)：
C~0.5s
Python~50s
那么,如何加快程序？请帮忙.
This code is only calculate, not drawing particles as in the pictures. The number of simulate partcile is on SymConst.py, this is nrH*nrL.
SymMain.py
#coding:windows-1250 from os import system from SymCalc import * from SymParticle import * if __name__ == '__main__': App = SymCalc() App.MainLoop()SymParticle.py
#coding:windows-1250 from random import randint from math import * from SymConst import * from SymParticle import * class SymCalc(object): def __init__(self): # declaration lists containing the properties of the particles ParticleList = [] ParticleListTemp = [] t = 0.0 # the initial values of particle for x in range(0, nParticle): ParticleList.append(Particle(x+1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 8e-6)) ParticleListTemp.append(Particle(x+1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 8e-6)) # generating random coordinates X, Y for x in range(0, nParticle): self.Rand(ParticleList[x], x) # time steps for k in range(0, k0): print "Time step = {0}".format(k) # calculation of forces for i in range(0, nParticle): for j in range(0, i-1): self.ParticleCalculate(ParticleList[i], ParticleList[j], 1) for j in range(i+1, nParticle): self.ParticleCalculate(ParticleList[i], ParticleList[j], 1) # display data if(k%Constant == 0): for i in range(0, nParticle): self.Print(k, i, ParticleList[i], t) # changing the position of the particle for i in range(0, nParticle): self.ChangeOfPosition (ParticleList[i], dt) # reset forces for j in range(0, nParticle): self.ParticleCalculate(ParticleList[i], ParticleList[j], 0) self.ParticleCalculate(ParticleList[i], ParticleListTemp[j], 0) # next time step t += dt # random coordinates of the particles def Rand(self, part, lp): l = lp%nrL # współrzędna X pola part.x = (l-nrL/2)*(L0/nrL) + ((randint(0,32767)%100)-50)*(L0/nrL-2*part.r)/99 # X h = (lp+1-(lp)%nrL)/nrL # Y part.y = (h-nrH/2)*(H0/nrH) + ((randint(0,32767)%100)-50)*(H0/nrH-2*part.r)/99 # współrzędna Y cząstki # calculating function of the force acting on the particles def ParticleCalculate(self, part, part_tmp, p): # auxiliary variables # r_4, dx, dy, rr, sum, sx, sy, chi_eff, tmp, frepx, frepy, f0 md = [] if(p == 0): part.frx = 0 part.fry = 0 part.fx = 0 part.fy = 0 if(p == 1): # versor coordinates connecting the geometrical center of the particle dx = part.x - part_tmp.x dy = part.y - part_tmp.y # the distance between the geometric means of the particle rr = sqrt(dx*dx + dy*dy) if(rr < 0.85*(part.r + part_tmp.r)): print "ERROR: Invalid distance between the particles! Simulation aborted..." if(rr >= 10*part.r): # magnetic dipoles chi_eff = (3.*(MI_P - 1.))/((MI_P - 1) + 3.) md.append((4.*MI_0*pi*part.r*part.r*part.r*chi_eff*M_H0)/3.0) md.append((4.*MI_0*pi*part_tmp.r*part_tmp.r*part_tmp.r*chi_eff*M_H0)/3.0) tmp = pow(rr,7) # first member sum = (5.*(md[0]*part.nx*dx + md[0]*part.ny*dy) * (md[1]*part_tmp.nx*dx + md[1]*part_tmp.ny*dy)) / tmp sx = sum*dx sy = sum*dy tmp = tmp / (rr*rr) # second member sum = (md[0]*part.nx*md[1]*part_tmp.nx + md[0]*part.ny*md[1]*part_tmp.ny) / tmp sx -= sum*dx sy -= sum*dy # third member sx -= (md[0]*(md[1]*part_tmp.nx*dx + md[1]*part_tmp.ny*dy)*part.nx + md[1]*(md[0]*part.nx*dx + md[0]*part.ny*dy)*part_tmp.nx)/tmp sy -= (md[0]*(md[1]*part_tmp.nx*dx + md[1]*part_tmp.ny*dy)*part.ny + md[1]*(md[0]*part.nx*dx + md[0]*part.ny*dy)*part_tmp.ny)/tmp # finally tmp = (-3./(4*pi*MI_0)) sx *= tmp sy *= tmp part.fx += sx part.fy += sy # short-range repulsive force tmp = pow(rr,15) r_4 = pow((part.r+part.r),4) f0 = 3.*fabs(md[0])*fabs(md[1])/(4.*pi*MI_0*r_4) frepx = pow((part.r+part.r),15)*dx/(tmp*rr)*f0 frepy = pow((part.r+part.r),15)*dy/(tmp*rr)*f0 part.frx += frepx; part.fry += frepy; # change the position of the particle def ChangeOfPosition(self, part, dt): part.ddx = 0 part.ddy = 0 part.ddx = 1/(6*pi*part.r*eta)*(part.fx+part.frx)*dt part.ddy = 1/(6*pi*part.r*eta)*(part.fy+part.fry)*dt # particles new position value part.x += part.ddx part.y += part.ddy if(part.x < -L0/2): part.x += L0 elif(part.x > L0/2): part.x -= L0 elif(part.y < -H0/2): part.y += H0 elif(part.y > H0/2): part.y -= H0 # display data def Print(self, k, i, part, t): print "-"*50 print "\nParticle {0}".format(i+1) print "The resultant magnetostatic force fx = {0}".format(part.fx - part.frx) print "The resultant magnetostatic force fy = {0}\n".format(part.fy - part.fry) if(i == nParticle-1): print "\n\t\t---t={0}[s]---".format(t) print "-"*50SymParticle.py
#coding:windows-1250 class Particle(object): # generating a particle properties def __init__(self, num, x, y, fx, fy, frx, fry, nx, ny, ddx, ddy, r): self.num = num self.x = x self.y = y self.fx = fx self.fy = fy self.frx = frx self.fry = fry self.nx = nx self.ny = ny self.ddx = ddx self.ddy = ddy self.r = rSymConst.py
#coding:windows-1250 ### Constant M_H0 = 3e4 MI_0 = 12.56e-7 MI_P = 2000 eta = 0.1 k0 = 10001 dt = 1e-6 # time step H0 = 5.95e-4 L0 = 5.95e-4 nrH = 4 nrL = 4 Constant = 100 nParticle = nrH*nrL # number of particle
解决方法:
TL; DR – 使用PyPy！
要研究这个问题,首先使用Python cProfile模块来分析您的代码并找出所有时间的使用位置(将print语句注释掉)：
$python -m cProfile -s time SymMain.py
30264977 function calls in 34.912 seconds
Ordered by: internal time
ncalls tottime percall cumtime percall filename:lineno(function)
7370737 27.547 0.000 30.501 0.000 SymCalc.py:62(ParticleCalculate)
1 3.626 3.626 34.909 34.909 SymCalc.py:8(__init__)
11101110 1.715 0.000 1.715 0.000 {math.pow}
160016 0.502 0.000 0.502 0.000 SymCalc.py:128(ChangeOfPosition)
4440476 0.498 0.000 0.498 0.000 {method 'append' of 'list' objects}
4440444 0.454 0.000 0.454 0.000 {math.fabs}
2250226 0.287 0.000 0.287 0.000 {math.sqrt}
500154 0.280 0.000 0.280 0.000 {range}
1 0.001 0.001 0.002 0.002 random.py:40(<module>)
1 0.001 0.001 0.001 0.001 hashlib.py:55(<module>)
1 0.001 0.001 0.003 0.003 SymCalc.py:2(<module>)
1616 0.000 0.000 0.000 0.000 SymCalc.py:151(Print)
1 0.000 0.000 34.912 34.912 SymMain.py:3(<module>)
32 0.000 0.000 0.000 0.000 SymParticle.py:4(__init__)
---8<--- snip ---8<---
在您的情况下,大部分时间都花在SymCalc.py:62(ParticleCalculate)函数中.看看那个函数,它似乎主要是内存访问和计算.这是查看PyPy的好例子！
PyPy is a fast, compliant alternative implementation of the Python
language (2.7.3 and 3.2.3). It has several advantages and distinct
features:
- Speed: thanks to its Just-in-Time compiler, Python programs often run faster on PyPy. (What is a JIT compiler?)
- Memory usage: large, memory-hungry Python programs might end up taking less space than they do in CPython.
- Compatibility: PyPy is highly compatible with existing python code. It supports cffi and can run popular python libraries like
twisted and django.- Sandboxing: PyPy provides the ability to run untrusted code in a fully secure way.
- Stackless: PyPy comes by default with support for stackless mode, providing micro-threads for massive concurrency.
- As well as other features.
$time pypy SymMain.py
real 0m1.071s
user 0m1.025s
sys 0m0.042s
1.071s – 好多了！将其与我系统上的原始执行时间进行比较：
$time python SymMain.py
real 0m29.766s
user 0m29.646s
sys 0m0.103s
标签：acceleration,python,physics,particles,simulation 来源： https://codeday.me/bug/20190831/1773299.html