How to include derivative inside scipy odeint?

Question

I am keep getting a long list of errors when I include a derivative function in a simple differential equation as illustrated in the following code. What am I doing wrong?

import numpy as np
from scipy.misc import derivative
from scipy.integrate import odeint

Imag = 16000.
w = 2*np.pi*60
tau = .05
theta = 1.52
phi = theta - np.radians(90)
t = np.linspace(0,.1,10000)
Ip = np.sqrt(2)*Imag*(np.sin(w*t+phi-theta)-np.exp(-t/tau)*np.sin(phi-theta))

def dI(t):
    return derivative(Ip,t)

def f(y,t):
    Rb = 8.
    N = 240.
    Is = y[0]
    f0 = (dI(t)/N)-Rb*y[0]
    return [f0]

yinit = [0]
sol = odeint(f,yinit,t)
print sol[:,0]

A small fragment of the error output that seems to be repeated several times is as follows:

val += weights[k]*func(x0+(k-ho)*dx,*args)
TypeError: 'numpy.ndarray' object is not callable
odepack.error: Error occurred while calling the Python function named f

Edit: The above code runs without errors in scipy 0.9.0 but gives the above error in 0.12.0

Edit2: Now that the error has been taken care of, I need to add a second function into the integrator as follows:

B = lambda Ip: Ip/(53.05+0.55*abs(Ip))
L = lambda Ip: derivative(B,Ip)*377.2

def f(y,t):
    Rb = 8.
    N = 240.
    Is = y[0]
    f0 = (1/(L+0.002))*((dI(t)*L/N)-Rb*y[0])
    return [f0]

How do I do it?

An example is given here: [http://stackoverflow.com/questions/17533751/method-signature-for-jacobian-func-in-scipy-odeint] — Holger, Nov 14 '13 at 20:09
That does not help me at all. I am not interested in the method as yet. At this point I want the script to run without errors. — amaity, Nov 15 '13 at 01:35

score 1 · Accepted Answer · answered Nov 15 '13 at 07:14

1

The first argument of derivative has to be a function not an ndarray. So you have to replace the ndarray Ip with a function Ip.

Following exapmple works with Python3:

import numpy as np
from scipy.misc import derivative
from scipy.integrate import odeint

def Ip(t):
    return np.sqrt(2)*Imag*(np.sin(w*t+phi-theta)-np.exp(-t/tau)*np.sin(phi-theta))

Imag = 16000.
w = 2*np.pi*60
tau = .05
theta = 1.52
phi = theta - np.radians(90)
t = np.linspace(0,.1,10000)

def dI(t):
    return derivative(Ip,t)

def f(y,t):
    Rb = 8.
    N = 240.
    Is = y[0]
    f0 = (dI(t)/N)-Rb*y[0]
    return [f0]

yinit = [0]
sol = odeint(f,yinit,t)
print(sol[:,0])

answered Nov 15 '13 at 07:14

Holger

2,125
2
19
30

Thanks for the response. This is not working for me (python 2.7, scipy 0.9.0). I am getting the following error: Traceback (most recent call last): File "C:\Windows\system32\config\systemprofile\Desktop\ct_example.py", line 53, in f f0 = (dI(t)/N)-Rb*y[0] TypeError: unsupported operand type(s) for /: 'NoneType' and 'float' odepack.error: Error occured while calling the Python function named f – amaity Nov 15 '13 at 08:50
Strange, on my machine it is also running with python 2.7 and scipy 0.9. What is your numpy version? Maybe, you should update both numpy and scipy to recent versions. – Holger Nov 15 '13 at 09:51
I downloaded the latest version of WinPython and the code runs fine. Thanks. – amaity Nov 15 '13 at 15:53
Sorry, I don't understand, what you want to achieve. I guess you should create a new question for your second edit, because this one is already marked as solved and noone will look at it any longer. – Holger Nov 18 '13 at 10:50
As advised @Holger, I have already posted it as a new question [here](http://stackoverflow.com/questions/20019427/how-to-solve-this-differential-equation-using-scipy-odeint). – amaity Nov 18 '13 at 12:18

How to include derivative inside scipy odeint?

1 Answers1