python – 如何更好地使用Cython来更快地求解微分方程?
发布时间:2020-12-20 11:58:37 所属栏目:Python 来源:网络整理
导读:我想降低Scipy的odeint解决差异所需的时间 方程. 为了练习,我使用Python in scientific computations中涵盖的示例作为模板.因为odeint将函数f作为参数,所以我将此函数编写为静态类型的Cython版本并希望如此 odeint的运行时间会明显减少. 函数f包含在名为ode.
我想降低Scipy的odeint解决差异所需的时间
方程. 为了练习,我使用Python in scientific computations中涵盖的示例作为模板.因为odeint将函数f作为参数,所以我将此函数编写为静态类型的Cython版本并希望如此 函数f包含在名为ode.pyx的文件中,如下所示: import numpy as np cimport numpy as np from libc.math cimport sin,cos def f(y,t,params): cdef double theta = y[0],omega = y[1] cdef double Q = params[0],d = params[1],Omega = params[2] cdef double derivs[2] derivs[0] = omega derivs[1] = -omega/Q + np.sin(theta) + d*np.cos(Omega*t) return derivs def fCMath(y,double t,Omega = params[2] cdef double derivs[2] derivs[0] = omega derivs[1] = -omega/Q + sin(theta) + d*cos(Omega*t) return derivs 然后我创建一个文件setup.py来编译函数: from distutils.core import setup from Cython.Build import cythonize setup(ext_modules=cythonize('ode.pyx')) 解决微分方程的脚本(也包含Python import ode import numpy as np from scipy.integrate import odeint import time def f(y,params): theta,omega = y Q,d,Omega = params derivs = [omega,-omega/Q + np.sin(theta) + d*np.cos(Omega*t)] return derivs params = np.array([2.0,1.5,0.65]) y0 = np.array([0.0,0.0]) t = np.arange(0.,200.,0.05) start_time = time.time() odeint(f,y0,args=(params,)) print("The Python Code took: %.6s seconds" % (time.time() - start_time)) start_time = time.time() odeint(ode.f,)) print("The Cython Code took: %.6s seconds ---" % (time.time() - start_time)) start_time = time.time() odeint(ode.fCMath,)) print("The Cython Code incorpoarting two of DavidW_s suggestions took: %.6s seconds ---" % (time.time() - start_time)) 然后我跑: python setup.py build_ext --inplace python solveODE.py 在终端. python版本的时间大约是0.055秒, 是否有人建议改进我的解决方案 编辑 我在两个文件ode.pyx和solveODE.py中加入了DavidW的建议.使用这些建议运行代码大约需要0.015秒. 解决方法
最简单的改变(可能会获得很多)是使用C数学库sin和cos来对单个数字而不是数字进行操作.对numpy的调用和计算它不是一个数组的时间相当昂贵.
from libc.math cimport sin,cos # later -omega/Q + sin(theta) + d*cos(Omega*t) 我很想为输入d分配一个类型(在不改变界面的情况下,没有其他输入可以轻松输入): def f(y,params): 我想我也会像你在Python版本中那样返回一个列表.我不认为你通过使用C数组获得了很多. (编辑:李大同) 【声明】本站内容均来自网络,其相关言论仅代表作者个人观点,不代表本站立场。若无意侵犯到您的权利,请及时与联系站长删除相关内容! |