2'nd order ode python - cfd - sphere motion in fluid flow

Question

phisical problem is given by 2'nd order ODE: m*x''= s*(v - x')^2 mathematical solution is rewriting to 2x 1'nd order ODE: u = x' , u' = s(v - u)^2

boundary conditions: u0 = 0, x0 = 0

t0=0, tmax=10., dt=1.

but i think theres some bug, and cant find it

code:

import numpy as np
from scipy.integrate import odeint
import pylab as pl

rho=10e3
Cd=1
r=0.05
m=5
s=rho*Cd*(r**2)*np.pi/(2*m)
v=1

t_max, dt = 10., 1.
time_vec = np.linspace(0, t_max, num=np.round(t_max/dt)+1)
x_vec=[0.,0.]

def d(x_vec,t):
    return np.array([s*((v-x_vec[0])**2), x_vec[0]]) # u'[0]=s*((v-u_vec[0])**2), u[1]'=u_vec[0]

x,u = odeint(func=d, y0=x_vec, t=time_vec).T

pl.plot(time_vec, u[:], label='u[m/s]=x/t')
pl.plot(time_vec, x[:], label="x[m]")
pl.legend()
pl.show()

Answer 1

It looks to me like you may want to use x and u instead of u and u'. See examples here (page 2) and here and also odeint docs. With the system x'=u and u' = s(v - u)^2, d can take [x,u] as input and return [x',u'] as odeint expects.