I have a project where I need ODE solver without dependencies to libraries like Scipy. I decide to implement ODE45. According to tutorials from internet and from what I remember from classes I implement it somehow.The code below contains also example function (pendulum) with exactly same values as in the scipy ODE tutorial https://docs.scipy.org/doc/scipy-0.18.1/reference/generated/scipy.integrate.odeint.html. The resulting plot from my code looks the same.However, I am not really sure if I done it in the good way. Especially I would really appreciate any hints how to improve the functions ode45_step and ode45 and verification whether my implementation is correct. I do not need to check the input variables, it will be done somewhere else.import numpy as npimport matplotlib.pylab as pltdef ode45_step(f, x, t, dt, *args):        One step of 4th Order Runge-Kutta method        k = dt    k1 = k * f(t, x, *args)    k2 = k * f(t + 0.5*k, x + 0.5*k1, *args)    k3 = k * f(t + 0.5*k, x + 0.5*k2, *args)    k4 = k * f(t + dt, x + k3, *args)    return x + 1/6. * (k1 + k2 + k3 + k4)def ode45(f, t, x0, *args):        4th Order Runge-Kutta method        n = len(t)    x = np.zeros((n, len(x0)))    x[0] = x0    for i in range(n-1):        dt = t[i+1] - t[i]         x[i+1] = ode45_step(f, x[i], t[i], dt, *args)    return xdef f(t, y, b, c):        Pendulum example function.        theta = y[0]    omega = y[1]    dydt = [omega, -b*omega - c*np.sin(theta)]    return np.array(dydt)b = 0.25c = 5.0    N = 101x0 = np.array([np.pi - 0.1, 0.0])t = np.linspace(0, 10, N)x = ode45(f, t, x0, b, c)plt.plot(x)plt.show()
ODE45 solver implementation in Python
python;numerical methods