This article

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.

I want to solve the equation of motion in Python

What is life

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.

The equation of motion you want to solve this time

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.

code

Overall flow

The code flow is as follows

Load the required library
Describe the differential equation you want to solve
Decide the initial conditions and ask odeint to solve the equation.
Draw the result

Preparation

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

Write differential equations

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

Decide the initial conditions and have the differential equation solved

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))

Draw the result in a diagram

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.

Next is

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.

Solving the equation of motion in Python (odeint)