errors while solving ode's python

1._

odeint(deriv_x, xinit, t)

uses xinit as its initial guess for x. This value for x is used when evaluating deriv_x.

deriv_x(xinit, t)

It looks like you are trying to solve the second-order ODE

r '' = -C rhat
   -- -- -- -- -
   |
   r | ** 2

Suggestion : 2

Posted February 27, 2013 at 02:31 PM | categories: ode | tags:

import numpy as np
import matplotlib.pyplot as plt

Tr = 0.9
Vr = np.linspace(0.34, 4, 1000)

#analytical equation
for Pr
Prfh = lambda Vr: 8.0 / 3.0 * Tr / (Vr - 1.0 / 3.0) - 3.0 / (Vr ** 2)
Pr = Prfh(Vr) # evaluated on our reduced volume vector.

# Plot the EOS
plt.plot(Vr, Pr)
plt.ylim([0, 2])
plt.xlabel('$V_R$')
plt.ylabel('$P_R$')
plt.savefig('images/ode-vw-1.png')
plt.show()

1._

import numpy as np
import matplotlib.pyplot as plt

Tr = 0.9
Vr = np.linspace(0.34, 4, 1000)

#analytical equation
for Pr
Prfh = lambda Vr: 8.0 / 3.0 * Tr / (Vr - 1.0 / 3.0) - 3.0 / (Vr ** 2)
Pr = Prfh(Vr) # evaluated on our reduced volume vector.

# Plot the EOS
plt.plot(Vr, Pr)
plt.ylim([0, 2])
plt.xlabel('$V_R$')
plt.ylabel('$P_R$')
plt.savefig('images/ode-vw-1.png')
plt.show()

2._

>>> >>> >>> >>> >>> ... >>> >>> >>> ... [<matplotlib.lines.Line2D object at 0x1c5a3550>]
   (0, 2)
   <matplotlib.text.Text object at 0x1c22f750>
      <matplotlib.text.Text object at 0x1d4e0750>

3._

from sympy
import diff, Symbol
Vrs = Symbol('Vrs')

Prs = 8.0 / 3.0 * Tr / (Vrs - 1.0 / 3.0) - 3.0 / (Vrs ** 2)
print diff(Prs, Vrs)

5._

from scipy.integrate
import odeint

def myode(Pr, Vr):
   dPrdVr = -2.4 / (Vr - 0.333333333333333) ** 2 + 6.0 / Vr ** 3
return dPrdVr

Vspan = np.linspace(0.334, 4)
Po = Prfh(Vspan[0])
P = odeint(myode, Po, Vspan, rtol = 1e-4)

# Plot the EOS
plt.plot(Vr, Pr) # analytical solution
plt.plot(Vspan, P[: , 0], 'r.')
plt.ylim([0, 2])
plt.xlabel('$V_R$')
plt.ylabel('$P_R$')
plt.savefig('images/ode-vw-2.png')
plt.show()

6._

... >>> >>> >>> >>> ... [<matplotlib.lines.Line2D object at 0x1d4f3b90>]
   [<matplotlib.lines.Line2D object at 0x2ac47518e710>]
      (0, 2)
      <matplotlib.text.Text object at 0x1c238fd0>
         <matplotlib.text.Text object at 0x1c22af10>

from sympy
import diff, Symbol
Vrs = Symbol('Vrs')

Prs = 8.0 / 3.0 * Tr / (Vrs - 1.0 / 3.0) - 3.0 / (Vrs ** 2)
print diff(Prs, Vrs)

from scipy.integrate
import odeint

def myode(Pr, Vr):
   dPrdVr = -2.4 / (Vr - 0.333333333333333) ** 2 + 6.0 / Vr ** 3
return dPrdVr

Vspan = np.linspace(0.334, 4)
Po = Prfh(Vspan[0])
P = odeint(myode, Po, Vspan, rtol = 1e-4)

# Plot the EOS
plt.plot(Vr, Pr) # analytical solution
plt.plot(Vspan, P[: , 0], 'r.')
plt.ylim([0, 2])
plt.xlabel('$V_R$')
plt.ylabel('$P_R$')
plt.savefig('images/ode-vw-2.png')
plt.show()

print myode(Po, 0.34)

Vspan = np.linspace(0.334, 4)
Po = Prfh(Vspan[0])
P = odeint(myode, Po, Vspan)

# Plot the EOS
plt.plot(Vr, Pr) # analytical solution
plt.plot(Vspan, P[: , 0], 'r.')
plt.ylim([0, 2])
plt.xlabel('$V_R$')
plt.ylabel('$P_R$')
plt.savefig('images/ode-vw-3.png')
plt.show()

Suggestion : 3

In scipy, there are several built-in functions for solving initial value problems. The most common one used is the scipy.integrate.solve_ivp function. The function construction are shown below:,The above figure shows the corresponding numerical results. As in the previous example, the difference between the result of solve_ivp and the evaluation of the analytical solution by Python is very small in comparison to the value of the function.,The way we use the solver to solve the differential equation is: solve_ivp(fun, t_span, s0, method = 'RK45', t_eval=None),for an initial value \(S_0 = 0\). The exact solution to this problem is \(S(t) = \sin(t)\). Use solve_ivp to approximate the solution to this initial value problem over the interval \([0, \pi]\). Plot the approximate solution versus the exact solution and the relative error over time.

import matplotlib.pyplot as plt
import numpy as np
from scipy.integrate
import solve_ivp

plt.style.use('seaborn-poster')

   %
   matplotlib inline

F = lambda t, s: np.cos(t)

t_eval = np.arange(0, np.pi, 0.1)
sol = solve_ivp(F, [0, np.pi], [0], t_eval = t_eval)

plt.figure(figsize = (12, 4))
plt.subplot(121)
plt.plot(sol.t, sol.y[0])
plt.xlabel('t')
plt.ylabel('S(t)')
plt.subplot(122)
plt.plot(sol.t, sol.y[0] - np.sin(sol.t))
plt.xlabel('t')
plt.ylabel('S(t) - sin(t)')
plt.tight_layout()
plt.show()

sol = solve_ivp(F, [0, np.pi], [0], t_eval = t_eval, \
   rtol = 1e-8, atol = 1e-8)

plt.figure(figsize = (12, 4))
plt.subplot(121)
plt.plot(sol.t, sol.y[0])
plt.xlabel('t')
plt.ylabel('S(t)')
plt.subplot(122)
plt.plot(sol.t, sol.y[0] - np.sin(sol.t))
plt.xlabel('t')
plt.ylabel('S(t) - sin(t)')
plt.tight_layout()
plt.show()

F = lambda t, s: -s

t_eval = np.arange(0, 1.01, 0.01)
sol = solve_ivp(F, [0, 1], [1], t_eval = t_eval)

plt.figure(figsize = (12, 4))
plt.subplot(121)
plt.plot(sol.t, sol.y[0])
plt.xlabel('t')
plt.ylabel('S(t)')
plt.subplot(122)
plt.plot(sol.t, sol.y[0] - np.exp(-sol.t))
plt.xlabel('t')
plt.ylabel('S(t) - exp(-t)')
plt.tight_layout()
plt.show()

F = lambda t, s: np.dot(np.array([
   [0, t ** 2],
   [-t, 0]
]), s)

t_eval = np.arange(0, 10.01, 0.01)
sol = solve_ivp(F, [0, 10], [1, 1], t_eval = t_eval)

plt.figure(figsize = (12, 8))
plt.plot(sol.y.T[: , 0], sol.y.T[: , 1])
plt.xlabel('x')
plt.ylabel('y')
plt.show()

Suggestion : 4

An example of using ODEINT is with the following differential equation with parameter k=0.3, the initial condition y0=5 and the following differential equation. ,Differential equations are solved in Python with the Scipy.integrate package using function odeint or solve_ivp. ,Another Python package that solves differential equations is GEKKO. See this link for the same tutorial in GEKKO versus ODEINT. ,y0: Initial conditions of the differential states

 y = odeint(model, y0, t)