[SciPy-user] odeint

LOPEZ GARCIA DE LOMANA, ADRIAN alopez at imim.es
Tue Oct 4 10:31:20 CDT 2005


Hi people, 

I'm using odeint from integrate to solve a system of ODEs. Although the output seems correct a strange error pop up at the terminal:

[alopez at thymus scipy]$ python single_global.integ.py
Traceback (most recent call last):
  File "single_global.integ.py", line 23, in func
    xdot[7] =  + (v_NACT * x[5]**m_NACT) / (x[5]**m_NACT + k_NACT**m_NACT) - k_deg_h * x[7]
ValueError: negative number cannot be raised to a fractional power
odepack.error: Error occured while calling the Python function named func
Traceback (most recent call last):
  File "single_global.integ.py", line 23, in func
    xdot[7] =  + (v_NACT * x[5]**m_NACT) / (x[5]**m_NACT + k_NACT**m_NACT) - k_deg_h * x[7]
ValueError: negative number cannot be raised to a fractional power
odepack.error: Error occured while calling the Python function named func
Traceback (most recent call last):
  File "single_global.integ.py", line 23, in func
    xdot[7] =  + (v_NACT * x[5]**m_NACT) / (x[5]**m_NACT + k_NACT**m_NACT) - k_deg_h * x[7]
ValueError: negative number cannot be raised to a fractional power
odepack.error: Error occured while calling the Python function named func
Traceback (most recent call last):
  File "single_global.integ.py", line 23, in func
    xdot[7] =  + (v_NACT * x[5]**m_NACT) / (x[5]**m_NACT + k_NACT**m_NACT) - k_deg_h * x[7]
ValueError: negative number cannot be raised to a fractional power
odepack.error: Error occured while calling the Python function named func

The input file looks like this:

from Numeric import *

from scipy.integrate import odeint

def func(x, t, *args):

	xdot = [0.46752920416900001, 0.0, 0.00496614094561, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.45751357598800002, 0.0, 0.035314464084199998, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

	xdot[0] =  + R_n - k_deg_n * x[0]

	xdot[1] =  - k_deg_N * x[1] + k_tln_n * x[0] - k_bind_ND * (1/64.0 * x[1] * x[12])

	xdot[2] =  - k_deg_d * x[2] + R_d / (1 + k_inhib_H * x[8]**s_H)

	xdot[3] =  - k_deg_D * x[3] - k_bind_ND * (1/64.0 * x[3] * x[10]) + k_tln_d * x[2]

	xdot[4] =  - k_deg_ND * x[4] + k_bind_ND * (1/64.0 * x[1] * x[12]) - k_diss_ND * (1/8.0 * x[4])

	xdot[5] =  + k_diss_ND * (1/8.0 * x[4]) - k_deg_NACT * x[5]

	xdot[6] =  - k_deg_DACT * x[6] + k_diss_ND * (1/8.0 * x[13])

	xdot[7] =  + (v_NACT * x[5]**m_NACT) / (x[5]**m_NACT + k_NACT**m_NACT) - k_deg_h * x[7]

	xdot[8] =  - k_deg_H * x[8] + k_tln_h * x[7]

	xdot[9] =  + R_n - k_deg_n * x[9]

	xdot[10] =  - k_deg_N * x[10] + k_tln_n * x[9] - k_bind_ND * (1/64.0 * x[10] * x[3])

	xdot[11] =  - k_deg_d * x[11] + R_d / (1 + k_inhib_H * x[17]**s_H)

	xdot[12] =  - k_deg_D * x[12] - k_bind_ND * (1/64.0 * x[12] * x[1]) + k_tln_d * x[11]

	xdot[13] =  - k_deg_ND * x[13] + k_bind_ND * (1/64.0 * x[10] * x[3]) - k_diss_ND * (1/8.0 * x[13])

	xdot[14] =  + k_diss_ND * (1/8.0 * x[13]) - k_deg_NACT * x[14]

	xdot[15] =  - k_deg_DACT * x[15] + k_diss_ND * (1/8.0 * x[4])

	xdot[16] =  + (v_NACT * x[14]**m_NACT) / (x[14]**m_NACT + k_NACT**m_NACT) - k_deg_h * x[16]

	xdot[17] =  - k_deg_H * x[17] + k_tln_h * x[16]

	return xdot

t = arange(0, 780.1, 1.0)

k_deg_d = 0.0437670312482
k_tln_d = 0.198782788861
v_NACT = 0.127802051983
k_tln_n = 0.164704331687
R_d = 0.0156035391386
k_deg_H = 0.180839336779
k_deg_h = 0.136197706451
k_tln_h = 0.177097298407
k_deg_n = 0.0306548917163
k_deg_D = 0.156002070979
k_deg_DACT = 0.136154162928
k_inhib_H = 41625892.9717
k_NACT = 5.53010217661
R_n = 0.0190046778533
k_bind_ND = 38207555.8983
s_H = 2.75203671214
k_deg_N = 0.0882819939502
k_deg_ND = 0.108276907855
m_NACT = 5.58940252216
k_deg_NACT = 0.0335928388267
k_diss_ND = 66522866.1189

parameters = [k_deg_d, k_tln_d, v_NACT, k_tln_n, R_d, k_deg_H, k_deg_h, k_tln_h, k_deg_n, k_deg_D, k_deg_DACT, k_inhib_H, k_NACT, R_n, k_bind_ND, s_H, k_deg_N, k_deg_ND, m_NACT, k_deg_NACT, k_diss_ND]

x_0 = [0.467529204169, 0.0, 0.00496614094561, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.457513575988, 0.0, 0.0353144640842, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

g = open('single_global.veloc.out', 'w')

args = (parameters, g)

x = odeint(func, x_0, t, args)

g.close()

f = open('single_global.out', 'w')

f.write('#\n# generated by ByoDyn\n#\n')

for i in range(len(t)):

	f.write('%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\n'%(t[i], x[i][0], x[i][1], x[i][2], x[i][3], x[i][4], x[i][5], x[i][6], x[i][7], x[i][8], x[i][9], x[i][10], x[i][11], x[i][12], x[i][13], x[i][14], x[i][15], x[i][16], x[i][17]))

f.close()

And the output file, looks OK, there is no single negative number. 

Can anyone reproduce the error? Can anyone tell me what's going on?

Thanks in advance, 

Adrián.



More information about the SciPy-user mailing list