Reputation: 258
I'm working with an ODE of the form:
a*dv/dt + (b+k1)*v + c*integral_0->t_(vdt) = k1*v1 + k2*integral_0->t_(v1dt)
I'm trying to implement odeint to get a solution for this system, but I'm not sure how to do that with an integral in the ODE. v1 is a known input, so that integral on the right side isn't a concern.
Upvotes: 0
Views: 525
Reputation: 25992
Set x
as the integral of v
, so that x'=v
, x''=v'
, similarly x1
for v1
, so that your equation reads as second order differential equation
a*x''+(b+k1)*x'+c*x=k1*v1+k2*x1
Given v1
as input, the state vector needs to contain the three integrated variables x, v, x1
which gives the ODE function
def odesys(y,t):
x, v, x1 = y
v1 = eval_v1(t)
return [ v, (k1*v1+k2*x1 - (b+k1)*v-c*x )/a, v1 ]
To use it with odeint you can for example do
t = np.linspace(0,T,2001); # define the end time T before
y0 = [ 0, 0, 0 ] # standard convention is that everything is zero for negative times
y = odeint(odesys, y0, t) # add arguments for higher accuracy if needed
x, v, x1 = y.T # transpose of a list of tuples is a tuple of lists
plt.plot(t,x); plt.show() # as example that should work
Upvotes: 1