Solve a ODE for arrays of parameters (Python)

Question

*I am aware that this question is quite simple but I would like to know the best way to set up such a for loop in Python.

I have written a program already to calculate and plot the solution to a 2nd order differential equation (this code is given below).

I would like to know best methods for repeating this calculation for an array of f parameters (hence f_array). I.e. so that the plot shows 20 sets of data referring to solutions as a function of t each with a different value of f.

Cheers for any ideas.

from pylab import *
from scipy.integrate import odeint

#Arrays.
tmax = 100
t = linspace(0, tmax, 4000)
fmax = 100
f_array = linspace(0.0, fmax, 20)

#Parameters
l = 2.5
w0 = 0.75
f = 5.0
gamma = w0 + 0.05
m = 1.0
alpha = 0.15
beta = 2.5

def rhs(c,t):
    c0dot = c[1]
    c1dot = -2*l*c[1] - w0*w0*c[0] + (f/m)*cos((gamma)*t)-alpha*c[0] - beta*c[0]*c[0]*c[0]
    return [c0dot, c1dot]

init_x = 15.0
init_v = 0.0
init_cond = [init_x,init_v]
ces = odeint(rhs, init_cond, t)

s_no = 1
subplot(s_no,1,1)
xlabel("Time, t")
ylabel("Position, x")
grid('on')
plot(t,ces[:,0],'-b')
title("Position x vs. time t for a Duffing oscillator.")
show()

Here is a plot showing the solution to this equation regarding a single value of f for an array of t values. I would like a quick way to repeat this plot for an array of f values.

Warren Weckesser · Accepted Answer · 2014-09-19 19:46:16Z

4

Here's one approach:

Modify rhs to accept a third argument, the parameter f. The definition of rhs should begin

def rhs(c, t, f):
    ...

Iterate over the f_array with a for loop. In the loop, call odeint with the args argument so that odeint gives the value of f as the third argument to rhs. Save the results of each call to odeint in a list. Basically, replace

ces = odeint(rhs, init_cond, t)

with

solutions = []
for f in f_array:
    ces = odeint(rhs, init_cond, t, args=(f,))
    solutions.append(ces)

For each value of f in f_array, you now have a solution in the solutions list.

To plot these, you could put your plot call in another for loop:

for ces in solutions:
    plot(t, ces[:, 0], 'b-')

edited Sep 19, 2014 at 19:46

answered Sep 17, 2014 at 22:12

Warren Weckesser

116k20 gold badges207 silver badges224 bronze badges

Sign up to request clarification or add additional context in comments.

2 Comments

user2992169 Over a year ago

Good suggestion, thank you. I have managed How do I then go onto separate each solution from the list of arrays for plotting?

Warren Weckesser Over a year ago

You can put the call to plot in a loop over solutions. See the updated answer.

Collectives™ on Stack Overflow

Solve a ODE for arrays of parameters (Python)

1 Answer 1

2 Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

2 Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related