It solves the equation of motion using a library called odeint in python (scipy to be exact). I am using python3.6.3 which I put in anaconda.
It solves the equation of motion (second-order differential equation), but it can also be applied to other differential equations such as the first-order differential equation.
There are times in life when you want to solve a simple equation of motion. In particular, there are times when you want to solve and plot the state of exercise. What to do at that time
This time, the purpose is to plot easily using the python library "odeint".
odeint
Looks like oden it. I think it's ODE integral. ODE = Ordinary Differential Equation
Official reference https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.odeint.html
There is also a library called ode https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.ode.html
Upon examination, it seems that there is a newer package "solve_ivp", which is recommended. I want to write an article about this one day. https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.solve_ivp.html#scipy.integrate.solve_ivp
By the way, odeint and the other two (ode, solve_ivp) have different argument order. Be careful not to get confused.
Any equation of motion is fine, but this time I will solve a simple throw.
m\ddot y=-rv+mg
In odeint, it is necessary to change the second-order differential equation, especially to two first-order differentials. For example, in this case
\dot y=v \\
m\dot v=-rv+mg
Is actually solved by odeint.
The code flow is as follows
First, put in various libraries.
import
import numpy as np #numpy
from scipy.integrate import odeint #odeint
import matplotlib.pyplot as plt #to draw graphs
Then define the differential equation. The air resistance r and the mass m of the mass point are defined in the next step. Since the gravitational constant g is not changed, it is defined here this time.
function
def func(s, t, r, m):
y, v = s #s is a pair of variables y and v
g=9.80665 #m/s^2
dsdt = [v, (-r*v-m*g)/m]
return dsdt
Determine the air resistance r and mass m, determine the initial conditions, and solve the differential equation. Note one point about the last argument of sol = odeint (func, y0, t, args = (r, m)). You pass a list like args = (r, m), so if you have one argument you have to do something like args = (r,). This time it's not a problem, but be careful when solving other equations of motion.
solve
r=12
m=100
y0 = [0,10]#Throw up from position 0 at initial velocity 10
t = np.linspace(0, 7, 201)#Calculate in 201step increments from time 0 to 7
sol = odeint(func, y0, t, args=(r,m))
Let's visualize the result.
visualize
plt.plot(t, sol[:, 0], 'b', label='y')#about y plot
plt.plot(t, sol[:, 1], 'g', label='v')#about v plot
plt.legend(loc='best')#Add a legend
plt.xlabel('t')
plt.grid()#Add a grid
plt.show()
Then you should get this figure.
At first, changing t from 0 to 100 to check the terminal velocity, or solving the problem when throwing it diagonally is a good practice to understand the code.
Recommended Posts