0

im trying to solve ODE with python 'solve_ivp' m(x''- u'') + k * (x - u) = 0 where u'' - sine function depended on time, m and k - state scalars Here is my code:

m = 100
k = 0.35
def dSdx(x, S, d2u_dt2, u):
    x, v = S
    return [v,
           d2u_dt2(Amp1, T1, t1) - k*x/m + k*u(Amp1, T1, v0, d0, t1)/m]
v_0 = 0
x_0 = 0
S_0 = (x_0, v_0)
sol = solve_ivp(dSdx, y0=S_0, t_span=t, args=(d2u_dt2, u))
print(sol)

d2u_dt2(Amp1, T1, t1) and u(Amp1, T1, v0, d0, t1) return array

Im getting error:

ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape after 1 dimensions. The detected shape was (2,) + inhomogeneous part.

First part of code:

import numpy  as np
import matplotlib.pyplot as plt
from scipy import integrate
import sympy as sp
from scipy.integrate import odeint
from scipy.integrate import ode
from scipy.integrate import quad
from scipy.integrate import solve_ivp
Amp1 = 80
Amp2 = -70
T1 = 5
T2 = 6
N = 50
v0=0
d0=0
t = np.linspace(0, T1+T2, N)
ax1 = plt.subplot(212)
ax2 = plt.subplot(221)
ax3 = plt.subplot(222)

def d2u_dt2(Amp, T, t):
    return Amp*np.sin((np.pi*t)/(T))

def du_dt(Amp, T, v0, t):
    return -(Amp*T/np.pi)*np.cos(np.pi*t/T) + (v0+(Amp*T/np.pi))

def u(Amp, T, v0, d0, t):
    return -Amp*((T/np.pi)**2)*np.sin(np.pi*t/T) + (v0 + (Amp*T/np.pi))*t + d0

t1 = np.linspace(0, 5, 50)
t2 = np.linspace(0, 6, 50)
t3 = np.linspace(5,11, 50)
ax1.plot(t1, d2u_dt2(Amp1, T1, t1), t3, d2u_dt2(Amp2, T2, t2))
ax2.plot(t1, du_dt(Amp1, T1, v0, t1), t3, du_dt(Amp2, T2, du_dt(Amp1, T1, v0, T1), t2))
ax3.plot(t1, u(Amp1, T1, v0, d0, t1), t3, u(Amp2, T2, du_dt(Amp1, T1, v0, T1), u(Amp1, T1, v0, d0, T1), t2))
plt.show()
bedermau5
  • 11
  • 3
  • Why is there a `t1` in `dSdx` and should it not be `dSdt`? Try replacing `t1` with `x` there. – Lutz Lehmann Mar 13 '22 at 10:34
  • tried, getting just zero array – bedermau5 Mar 13 '22 at 10:46
  • It works when u'' and u just a float num, but with func - no. – bedermau5 Mar 13 '22 at 10:47
  • You should check all your names, especially what you want `x` to be and not to be. Better stay with the convention that the independent variable is named `t`. – Lutz Lehmann Mar 13 '22 at 10:54
  • I understand, but the problem that 'odeint', 'solve_ivp' dont want to eat 'arrays', not about x and t. It works fine, the main question is how to present dynamic var – bedermau5 Mar 13 '22 at 11:01
  • This is also a frequent question, and the answer is as always "interpolation", be it piecewise linear, polynomial or a zero-order hold. https://stackoverflow.com/questions/54266330/how-to-insert-list-in-a-ode, https://stackoverflow.com/questions/32987253/solving-a-system-of-odes-with-changing-constant – Lutz Lehmann Mar 13 '22 at 11:05
  • In the general case you might want to set `y=x-u` and solve the ODE for `y`, reconstructing `x` after the integration. For a non-trivial example with derivatives of the control input and how to avoid them see https://stackoverflow.com/questions/58810439/is-there-a-better-way-to-simulate-pid-control or https://math.stackexchange.com/questions/2829678/confusion-about-heuns-method – Lutz Lehmann Mar 13 '22 at 11:16
  • @LutzLehmann tried interpolation, `f(x) : interp1d(x,80*np.sin((np.pi*x)/(5)), kind='linear')` and got error `TypeError: object of type 'numpy.float64' has no len()` – bedermau5 Mar 13 '22 at 11:50
  • You should make your situation more realistic. I assume that you want in the end to be able to use any input given as function table in the place of `u`. So construct this function table explicitly and remove any on-the-fly evaluations of the example function. This would include that you do not have a reliable source for the second derivative. Or that the first and second derivatives are also part of the given data (which in some sense might be redundant). – Lutz Lehmann Mar 13 '22 at 12:04
  • Please check my [answer](https://stackoverflow.com/questions/71424729/how-to-pass-function-array-into-odeint/71464960#71464960) on your another somewhat similar problem. – Muhammad Mohsin Khan Mar 14 '22 at 10:20

0 Answers0