python - How do I improve the speed of this code? (Solving ODEs with scipy.integrate.odeint) -

python - How do I improve the speed of this code? (Solving ODEs with scipy.integrate.odeint) -

i trying solve relatively big ode scheme scipy.integrate.odeint module. implemented code , can solve equation correctly. process slow. profile code , realized of computing time spent on computing or generating ode scheme itself, sigmoid function expensive have take guess. here piece of code i'm using:

def __sigmoid(self, u): # homecoming .5 * ( u / np.sqrt(u**2 + 1) + 1) homecoming 0.5 + np.arctan(u) / np.pi def __connectionistmodel(self, g, t): """ returning ode scheme """ g_ia_s = np.zeros(self.ngenes * self.ncells) in xrange(0, self.ncells): in xrange(0, self.ngenes): g_ia = self.params["r"][a] *\ self.__sigmoid( sum([ self.params["w"][b + a*self.ngenes]*g[self.ngenes*i + b] b in xrange(0, self.ngenes) ]) +\ self.params["wm"][a]*self.mdata[i] +\ self.params["h"][a] ) -\ self.params["l"][a] * g[self.ngenes*i + a] # uncomment part including diffusion if == 0: g_ia += self.params["d"][a] * ( - 2*g[self.ngenes*i + a] + g[self.ngenes*(i+1) + a] ) elif == self.ncells-1: g_ia += self.params["d"][a] * ( g[self.ngenes*(i-1) + a] - 2*g[self.ngenes*i + a] ) else: g_ia += self.params["d"][a] * ( g[self.ngenes*(i-1) + a] - 2*g[self.ngenes*i + a] + g[self.ngenes*(i+1) + a] ) g_ia_s[self.ngenes*i + a] = g_ia homecoming g_ia_s def solve(self, inp): g0 = np.zeros(self.ngenes * self.ncells) t = np.arange(self.t0, self.t1, self.dt) self.integratedexpression = odeint(self.__connectionistmodel, g0, t, full_output=0) homecoming self.integratedexpression

as see in each iteration, should generate ncells*ngenes (100*3=300) equations , pass odeint. although i'm not sure guess generating equations expensive in comparing solving them. in experiment, solving whole scheme takes 7sec consists of 1sec of odeint , 6secs of __connectionistmodel.

i wondering if there way improve or not? tried utilize sympy define symbolic ode scheme , pass symbolic equations odeint didn't works since couldn't define array of symbols later access array.

in worst case, have deal or utilize cython speed solving process, wanted create sure i'm doing right , there no way improve it.

thanks help in advance.

[update]: profiling result,

ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 7.915 7.915 grntest.py:1(<module>) 1 0.000 0.000 7.554 7.554 grn.py:83(solve) 1 0.000 0.000 7.554 7.554 odepack.py:18(odeint) 1 0.027 0.027 7.554 7.554 {scipy.integrate._odepack.odeint} 1597 5.506 0.003 7.527 0.005 grn.py:55(__connectionistmodel) 479100 1.434 0.000 1.434 0.000 grn.py:48(__sigmoid) 479102 0.585 0.000 0.585 0.000 {sum} 1 0.001 0.001 0.358 0.358 grn.py:4(<module>) 2 0.001 0.001 0.207 0.104 __init__.py:10(<module>) 27 0.014 0.001 0.185 0.007 __init__.py:1(<module>) 7 0.006 0.001 0.106 0.015 __init__.py:2(<module>)

[update 2]: made code publicly available: pystgrn

vectorize, vectorize, vectorize more. , utilize info structures facilitate vectorization.

the function __connectionistmodel uses lot of access pattern a[i*m+j], equivalent access row i , column j in 2d array total of m columns. suggests 2d array right way store data. can eliminate loop on i function using numpy slicing notation , vectorizing follows:

def __connectionistmodel_vec(self, g, t): """ returning ode scheme """ g_ia_s = np.zeros(self.ngenes * self.ncells) g_2d = g.reshape((self.ncells, self.ngenes)) w = np.array(self.params["w"]) mdata = np.array(self.mdata) g_ia_s = np.zeros((self.ncells, self.ngenes)) in xrange(0, self.ngenes): g_ia = self.params["r"][a] *\ self.__sigmoid( (w[a*self.ngenes:(a+1)*self.ngenes]*g_2d).sum(axis=1) +\ self.params["wm"][a]*mdata +\ self.params["h"][a] ) -\ self.params["l"][a] * g_2d[:,a] g_ia[0] += self.params["d"][a] * ( - 2*g_2d[0,a] + g_2d[1,a] ) g_ia[-1] += self.params["d"][a] * ( g_2d[-2,a] - 2*g_2d[-1,a] ) g_ia[1:-1] += self.params["d"][a] * ( g_2d[:-2,a] - 2*g_2d[1:-1,a] + g_2d[2:,a] ) g_ia_s[:,a] = g_ia homecoming g_ia_s.ravel()

as far can see, returns same values original __connectionistmodel. bonus, function more compact. optimized function has same inputs , outputs original, performance, might want organize info in numpy arrays instead of lists avoid conversion lists arrays on every call. , i'm sure there other minor performance tweaks well.

anyway, original code gave me these profiling results (insert mandatory "my computer faster yours" brag here):

ncalls tottime percall cumtime percall filename:lineno(function) 1597 3.648 0.002 5.250 0.003 grn.py:52(__connectionistmodel) 479100 0.965 0.000 0.965 0.000 grn.py:48(__sigmoid) 479100 0.635 0.000 0.635 0.000 {sum} 1 0.017 0.017 5.267 5.267 {scipy.integrate._odepack.odeint} 1598 0.002 0.000 0.002 0.000 {numpy.core.multiarray.zeros}

with __connectionistmodel_vec, get:

ncalls tottime percall cumtime percall filename:lineno(function) 1597 0.175 0.000 0.247 0.000 grn.py:79(__connectionistmodel_vec) 4791 0.031 0.000 0.031 0.000 grn.py:48(__sigmoid) 4800 0.021 0.000 0.021 0.000 {method 'reduce' of 'numpy.ufunc' objects} 1 0.018 0.018 0.265 0.265 {scipy.integrate._odepack.odeint} 3197 0.013 0.000 0.013 0.000 {numpy.core.multiarray.array}

python scipy ode odeint