第五次作业

---

题目2.6

Use the Euler method to calculate cannon shell trajectories ignoring both air drag and effect of air density (actually,ignoring the former automatically rules out the latter).Compare your results with those in Figure 2.4,and with the exact solution.

一、摘要

用Euler法解二元二阶常微分方程组。计算一个斜抛运动。
忽略空气阻力和空气密度的影响，计算精确解并与图2.4比较。
图2.4：

二、背景介绍

对于二阶的牛顿方程：

$$\frac{\mathrm{d}^2 x}{\mathrm{d} t^2}=0, \frac{\mathrm{d}^2 y}{\mathrm{d} t^2}=-g$$
设速度$v_{x}=\frac{\mathrm{d} x}{\mathrm{d} t},v_{y}=\frac{\mathrm{d} y}{\mathrm{d} t}$
则可以把方程组降阶为4元一阶方程组
常微分方程组： $$\frac{\mathrm{d} x}{\mathrm{d} t}=v_{x}$$
$$\frac{\mathrm{d} v_{x}}{\mathrm{d} t}=0$$
$$\frac{\mathrm{d} y}{\mathrm{d} t}=v_{y}$$
$$\frac{\mathrm{d} v_{y}}{\mathrm{d} t}=-g$$
然后便可用Euler法解这个方程组

三、正文

import pylab
from math import sqrt
from math import cos
from math import sin
from math import pi
class Euler2Dsimulation:
    def __init__(self,speed,theta,timestep=0.01):
        self._timestep=timestep
        self._x=[0]
        self._y=[0]
        self._vx=[speed*cos(theta)]
        self._vy=[speed*sin(theta)]
        self.g=9.8
        self.calculate()
        self.plot()
    def _accelerate(self,x,y,vx,vy):
        ansax=0
        ansay=-self.g
        return (ansax,ansay)
    def calculate(self):
        while self._y[-1]>=0:
            self._x.append(self._x[-1]+self._vx[-1]*self._timestep)
            self._y.append(self._y[-1]+self._vy[-1]*self._timestep)
            (ax,ay)=self._accelerate(self._x[-1],self._y[-1],self._vx[-1],self._vy[-1])
            self._vx.append(self._vx[-1]+self._vx[-1]*self._timestep*ax)
            self._vy.append(self._vy[-1]+self._timestep*ay)
    def plot(self):
        pylab.plot(self._x, self._y)
        pylab.xlabel('x ($m$)')
        pylab.ylabel('y ($m$)')
Euler2Dsimulation(700,pi/4)
Euler2Dsimulation(700,pi/6)
Euler2Dsimulation(700,pi/3)
pylab.show()        

四、结论

如图

五、致谢

感谢杜威同学在写ptyhon程序方面的帮助。感谢蔡浩老师。