Exercise_09: Chapter 3 problem 3.31

一.摘要：

Study the behavior for other types of tables.

二.背景介绍（公式推导）：

1.Between collisions: $$\dot{x}=v_x,(9-1)$$ $$\dot{y}=v_y,(9-2)$$

2.During collisions: $$\boldsymbol{v_{i,perpendicular}=(\boldsymbol{v_i}*\boldsymbol{n})\boldsymbol{n}},(9-3)$$ $$\boldsymbol{v_{i,//}}=\boldsymbol{v_i}-\boldsymbol{v_{i,perpendicular}},(9-4)$$

3.After reflection: $$\boldsymbol{v_{f,perpendicular}}=-\boldsymbol{v_{i.perpendicular}},(9-5)$$ $$\boldsymbol{v_{f,//}}=\boldsymbol{v_{i,//}},(9-6)$$

三.正文：

1.The behavior for the triangle table:

import numpy as np
from pylab import *
x=[0]
y=1
vx,vy=2.3,3
time,i,dt=0,0,0.001
sdt=0.00001
velocity=4
while time<=300:
X=x[i]+vx*dt
Y=y[i]+vy*dt
if Y>0:
if X>=0:
if Y<-X*4.0/3.0+4:
x.append(X)
y.append(Y)
velocity.append(vx)
else:
x1=x[i]+vx*sdt
y1=y[i]+vy*sdt
if y1<-x1*4.0/3.0+4:
x.append(x1)
y.append(y1)
velocity.append(vx)
else:
char1=-vx*0.28-vy*0.96
char2=-vx*0.96+vy*0.28
vx=char1
vy=char2
x.append(x[i]+vx*sdt)
y.append(y[i]+vy*sdt)
velocity.append(vx)
if X<0:
if Y x.append(X)
y.append(Y)
velocity.append(vx)
else:
x2=x[i]+vx*sdt
y2=y[i]+vy*sdt
if y2 x.append(x2)
y.append(y2)
velocity.append(vx)
else:
char3=-vx*0.28+vy*0.96
char4=vx*0.96+vy*0.28
vx=char3
vy=char4
x.append(x[i]+vx*sdt)
y.append(y[i]+vy*sdt)
velocity.append(vx)
if Y<=0:
x3=x[i]+vx*sdt
y3=y[i]+vy*sdt
if y3>0:
x.append(x3)
y.append(y3)
velocity.append(vx)
else:
vy=-vy
x.append(x[i]+vx*sdt)
y.append(y[i]+vy*sdt)
velocity.append(vx)
time=time+dt
i=i+1
plt.figure(figsize=(16,5.5))
subplot(1,2,1)
plt.title("vx0=2.3,vy0=3")
plt.xlabel("x")
plt.ylabel("y")
plt.xticks([-3,-2,-1,0,1,2,3])
plt.yticks([0,1,2,3,4,5])
plt.xlim(-3,3)
plt.ylim(0,5)
plt.plot([0,3],[4,0],color="blue",label="triangle",linewidth=2)
plt.plot([-3,0],[0,4],color="blue",linewidth=2)
plt.plot([-3,3],[0,0],color="blue",linewidth=2)
plt.plot(x,y,label="trajectory",color="red")
plt.scatter(0,1,color="black",alpha=1,linewidth=4,label="initial")
plt.legend()
subplot(1,2,2)
plt.xlabel("x")
plt.ylabel("vx")
for i in range(1000):
if 1000*i<=len(x):
plt.scatter(x[1000*i],velocity[1000*i])
plt.show()

①.when the initial velocity is perpendicular to one side of the triangle:

②.when the initial velocity is not perpendicular to one side of the triangle:

import numpy as np
subplot(1,2,1)
subplot(1,2,2)

2.The behavior for the stadium billiard:

if x[i+1]**2+y[i+1]**2>4:
print x[i+1]**2+y[i+1]**2
while abs(x[i+1]**2+y[i+1]**2-4)>0.001:
#print abs(x[i+1]**2+y[i+1]**2-4)
X=(x[i]+x[i+1])/2
Y=(y[i]+y[i+1])/2
a=X/(X**2+Y**2)**0.5
b=Y/(X**2+Y**2)**0.5
v_x[i+1]=(1-2*a**2)*v_x[i]-2*a*b*v_y[i]
v_y[i+1]=(1-2*b**2)*v_y[i]-2*a*b*v_x[i]
#v_y[i+1]=-((2.0*a**3-2.0*a)/b)*v_x[i+1]-(2*a**2-1)*v_y[i+1]
if X**2+Y**2>4:
x[i+1]=X
y[i+1]=Y
continue
else:
x[i]=X
y[i]=Y
continue
#if abs(abs(X)+abs(Y)-1)<0.1:
# break

3.The behavior for the concentric circle:

四.结论：

1.when the initial velocity is perpendicular to one side of the triangle，the system of triangle is not chaotic；

2.when the initial position is different,the system of stadium billiard is in different states;

3.when changing the initial velocity and position,the plot of $v_x versus x$ is symmetric about x-axis and y-axis in the system of concencircle ,while $v_y versus x$ is not.

五.致谢：

Thanks to Liuxinchen and Dingdongdong senior for reference.And best wishes to you.

未分类

在此输入正文