![bcc24d96f12fb69a5c76c9b9f5316ad1.png](https://i-blog.csdnimg.cn/blog_migrate/cbd2915afe86e80bccdf38cfad8e273a.jpeg)
相信大家都玩过桌球游戏吧,不过大家有没有思考过桌球游戏是如何模拟出真实球与球的撞击的呢?这里我通过我所学的知识给大家讲述一种实现方法,当然方法很多各有优缺点,这里提供一种思路。
![653b93f74f07873fda92b7ff83ed616a.png](https://i-blog.csdnimg.cn/blog_migrate/0017afc671dc6b71433e10bc48bc5d55.png)
如Figure 1所示蓝红两球非对心碰撞,利用高中所学的知识我们可以解决两个小球完全弹性的对心碰撞,对与非对心的弹性碰撞我们这里做一个简化,认为小球没有转动的自由度。首先我们在两球碰撞的界面上建立一个局部坐标系,称之为对心坐标系,如下图所示,对心坐标系的n轴方向为两球球心确定的方向,我们将两个小球的速度投影到对心坐标系,每个小球的速度用va与vn两个方向的速度合成,就可以认为两个小球在对心n的方向上发生对心弹性碰撞,碰撞后两球的vn需要通过动能与动量守恒定律进行求解。
![a793f3465959724d0634c22539c3e3a4.png](https://i-blog.csdnimg.cn/blog_migrate/1ab0dbbd7b9665676c4a474c0b202daa.png)
蓝球速度在对心坐标系下的速度投影为(红球的形式相同):
![23b7cf4bb213d685a92c9f0c5ca6b59f.png](https://i-blog.csdnimg.cn/blog_migrate/a8399935cc4101f53422e5d13a7e794b.png)
其中theta1与theta2如下图,根据两个向量的夹角很容易计算出来
![dd3e82132edc9a6717826a0af850ef71.png](https://i-blog.csdnimg.cn/blog_migrate/9595d5ec586a36db19b633e1ad96d82e.jpeg)
根据对心碰撞计算公式我们可以得到在两球碰撞后在对心坐标系下各自的n轴方向速度如下:
![15c349c03c127affcc6e59080a4df62b.png](https://i-blog.csdnimg.cn/blog_migrate/518c234109bc23e2ba3cf881d0891d50.png)
其中b与r分别代表蓝球与红球的速度。
将碰撞后的速度转化为全局坐标系,同样以蓝球为例,红球的形式相同:
![7f4418b1633348904c70eb4e9b4cbb75.png](https://i-blog.csdnimg.cn/blog_migrate/14a22fcd15270b859386e201c6130bb8.png)
这种方法就可以计算两个小球非对心碰撞的情况,在程序里面每次循环计算两个小球球心的距离就可以判断两个小球是否发生了碰撞。这里贴出了我的代码,欢迎大家交流讨论。这里用的python版本为3.6,小球的信息我用txt进行编辑,每一行为小球的编号,质量,坐标x,坐标y, 速度x,速度y,小球半径。txt文件名称为ball_info2.txt,如下图所示为3个小球的信息txt文档。
![6c64fda40bf07171b7f7b7398b2dd3e8.png](https://i-blog.csdnimg.cn/blog_migrate/917ab0369692748a66e2ed373af689d6.png)
运行结果:
![8c1832cc162f49472a20e512cc2d8814.gif](https://i-blog.csdnimg.cn/blog_migrate/d2afe26ebc4acbdb08e9878492e0f926.gif)
###代码如下######
# -*- coding: utf-8 -*-
"""
Created on Sun Feb 9 12:44:40 2020
@author: Lei Kang
"""
import numpy as np
#from matplotlib import pyplot as plt
import matplotlib.pyplot as plt
from PIL import Image
from matplotlib.patches import Ellipse, Circle
ball_num=2
f=open('ball_info2.txt','rb')
ball_mat=np.loadtxt(f,delimiter=',')
length=2
width=2
ball1_v=np.zeros([2,1])
ball2_v=np.zeros([2,1])
ball_vn=np.zeros([2,1])
dis_mat=np.zeros([ball_num+1,ball_num+1])
v_distb=np.zeros([2,2])
tt=5000
tim=0
dt=0.0001
count=1
tr=0
c=0.00001
def new_vx(length,width,x,y,r,vx,vy,t,c):
bv1=(x>r and x<length-r)
if bv1:
force_x=0
f=vx*np.exp(-c*t)
else:
force_x=0.99
f=-2*vx*force_x+vx
return f
def new_vy(length,width,x,y,r,vx,vy,t,c):
bv2=(y>r and y<width-r)
if bv2:
force_y=0
f=vy*np.exp(-c*t)
else:
force_y=0.99
f=-2*vy*force_y+vy
return f
for ti in range(0,tt):
tim=tim+dt
for i in range(0,ball_num+1):
ball_mat[i][4]=new_vx(length,width,ball_mat[i][2],ball_mat[i][3],ball_mat[i][6],ball_mat[i][4],ball_mat[i][5],tim,c)
ball_mat[i][5]=new_vy(length,width,ball_mat[i][2],ball_mat[i][3],ball_mat[i][6],ball_mat[i][4],ball_mat[i][5],tim,c)
for i in range(0,ball_num+1):
for j in range(i,ball_num+1):
if i!=j:
dis_mat[i][j]=np.sqrt((ball_mat[i][2]-ball_mat[j][2])**2+(ball_mat[i][3]-ball_mat[j][3])**2)
else:
dis_mat[i][j]=ball_mat[i][6]*2+0.01
####利用相对位置矩阵判断是否发生碰撞,并对碰撞的小球的速度进行修改
for i in range(0,ball_num+1):
for j in range(i,ball_num+1):
if dis_mat[i][j]<ball_mat[i][6]+ball_mat[j][6]:
tr=1
#速度分配矩阵
v_distb[0][0]=(ball_mat[i][1]-ball_mat[j][1])/(ball_mat[i][1]+ball_mat[j][1])
v_distb[0][1]=(ball_mat[j][1]+ball_mat[j][1])/(ball_mat[i][1]+ball_mat[j][1])
v_distb[1][0]=(ball_mat[j][1]+ball_mat[j][1])/(ball_mat[i][1]+ball_mat[j][1])
v_distb[1][1]=(-ball_mat[i][1]+ball_mat[j][1])/(ball_mat[i][1]+ball_mat[j][1])
#小球速度在对心坐标上投影
xn=ball_mat[j][2]-ball_mat[i][2]
yn=ball_mat[j][3]-ball_mat[i][3]
cos_theta1=np.abs(xn)/np.sqrt(xn**2+yn**2)
sin_theta1=np.abs(yn)/np.sqrt(xn**2+yn**2)
cos_theta2=sin_theta1
sin_theta2=cos_theta1
ball1_v[0]=ball_mat[i][4]*cos_theta1+ball_mat[i][5]*cos_theta2
ball1_v[1]=-ball_mat[i][4]*sin_theta1+ball_mat[i][5]*sin_theta2
ball2_v[0]=ball_mat[j][4]*cos_theta1+ball_mat[j][5]*cos_theta2
ball2_v[1]=-ball_mat[j][4]*sin_theta1+ball_mat[j][5]*sin_theta2
##根据小球完全弹性对心碰撞计算碰撞后的速度
ball_vn[0]=ball1_v[0]
ball_vn[1]=ball2_v[0]
ball_vn_new=np.dot(v_distb,ball_vn)
ball1_v[0]=ball_vn_new[0]
ball2_v[0]=ball_vn_new[1]
ball_mat[i][4]=ball1_v[0]*sin_theta2-ball1_v[1]*sin_theta1
ball_mat[i][5]=ball1_v[0]*cos_theta2+ball1_v[1]*cos_theta1
ball_mat[j][4]=ball2_v[0]*sin_theta2-ball2_v[1]*sin_theta1
ball_mat[j][5]=ball2_v[0]*cos_theta2+ball2_v[1]*cos_theta1
##时间步推进并更新每个小球的位置
count=count+1
for bi in range(0,ball_num+1):
ball_mat[bi][2]=ball_mat[bi][2]+ball_mat[bi][4]*dt
ball_mat[bi][3]=ball_mat[bi][3]+ball_mat[bi][5]*dt
if count%15==0:
fig = plt.figure(figsize=(5, 5),dpi=80)
axes = fig.add_subplot(111)
for mm in range(0,ball_num+1):
theta=np.arange(0,2*np.pi,0.01)
xcl = ball_mat[mm,2] + ball_mat[mm,6] * np.cos(theta)
ycl = ball_mat[mm,3] + ball_mat[mm,6] * np.sin(theta)
axes.plot(xcl, ycl)
axes.axis([0,2,0,2])
newName = str(count)+".png"
plt.savefig(newName)
imgs = []
for n in range(15,count,15):
img=Image.open(''+str(n)+'.png')
imgs.append(img)
img.save('ball.gif',save_all=True,append_images=imgs,duration=2,)