首页
学习
活动
专区
圈层
工具
发布
社区首页 >问答首页 >加速odeint计算

加速odeint计算
EN

Stack Overflow用户
提问于 2017-03-30 12:12:12
回答 1查看 1.6K关注 0票数 0

我正在使用scipy的odeint运行一个相当大的ODEs系统。为了保持系统的连贯性、易读性和简洁性,我将其大部分划分为类,以便相互调用。从可读性的角度来看,它似乎工作得很好,但是由于对python对象的重复调用,它大大减慢了实际的ODE解算器的速度。我研究了一些优化代码的方法,以保持一定程度的可读性,同时使其在计算时更有效(在将其转换为numpy感知的lambda函数之前,使用Numba或Sympy ),但是到目前为止,我还不是很成功。我想知道在这种情况下有什么策略会有帮助。我已经提供了下面代码的一个简单版本。提前谢谢你。

代码语言:javascript
复制
class NaL():
    g = 0
    def I(self, v, E): return self.g * (v - E) #Leak current

class Neuron():
    C_m = 1.0 #membrane capacitance
    V = 0  # Voltage
    m, h = 0, 0  #activating variables
    def I_Na(self): #Sodium Currents
        return self.NaT.I(self.V, self.m, self.h)+ self.NaL.I(self.V)  
    def __init__(self):
        self.NaT = NaT()
        self.NaL = NaL()
        self.NaT.g = 3000
        self.NaL.g = 20

I1 = Neuron()

def Run(X,t):
    I1.V, I1.m, I1.h = X 
    dydt = [(0 - I1.I_Na()) / I1.C_m, #dV/dt for neuron
            I1.NaT.dmdt(I1.V, I1.m), #dm/dt for sodium channel
            I1.NaT.dhdt(I1.V, I1.h) #dh/dt for sodium channel
           ]
    return dydt
X = odeint(Run, [-70, 0.0050, 0.9961], t)
EN

回答 1

Stack Overflow用户

发布于 2017-03-30 16:16:43

集成的瓶颈几乎肯定是计算派生的Run,因为集成的每个步骤都必须多次调用它。现在,由于所有的类成员和函数调用,Python的开销很大。

为此,一种解决方案是真正使用模块JiTCODE,其中不是定义积分器要使用的函数,而是符号地定义微分方程(使用SymPy)。然后,它们被转换成C代码,C代码被编译成Python函数,而Python函数又用于集成。这样,在实际集成之前,您的类和类函数只会被调用几次,并且集成的效率与定义没有类或类似类的派生时一样高。(此外,您可以通过编译获得相当大的速度提升,特别是如果您有许多神经元。)

按照现在的情况,将示例转换为JiTCODE很简单:

代码语言:javascript
复制
from jitcode import jitcode, provide_basic_symbols

t, y = provide_basic_symbols()

class NaL():
    […]

class Neuron():
    […]

I1 = Neuron()

I1.V, I1.m, I1.h = y(0), y(1), y(2)

Run = [
    (0 - I1.I_Na()) / I1.C_m,
    I1.NaT.dmdt(I1.V, I1.m),
    I1.NaT.dhdt(I1.V, I1.h)
]

t = range(100)

ODE = jitcode(Run)
ODE.set_initial_value([-70,0.0050,0.9961], t[0])
X = np.vstack(ODE.integrate for time in t[1:])

这段代码现在不能工作,原因与示例代码不能工作的原因相同。请注意,类中涉及动态变量的所有算术都是符号化的。对于加法和乘法,这是开箱即用的,但是如果你有像math.sinsympy.sin这样的东西,你必须用sympy.sin代替它。

票数 0
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/43108429

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档