本章会介绍如何用数值模拟的方法解决导弹追踪问题
一、问题提出
二、建立示意图
三、模型建立
1.建立坐标轴
（1）建立B船坐标
以导弹的初始位置为（0，0），那么经过时间t(t>0)时，B船的坐标为
（2）建立导弹坐标
下方是建立了一个速度的微分方程，瞬时速度Vx=dx/dt，Vy=dy/dt。
但是在数值模拟的代码中，并没有体现微分方程，而是在上一个t上加上delta—_d即可
2.设置delta_t，进行模拟
首先将时间T划分为几个间隔，每一个间隔的时间都是delta_t，为了让模拟效果尽量好，我们选取的delta_t越小越好，这里我们不妨取delta_t=0.0000001
之后每经过一次delta_t，船和导弹的坐标将进行更新，当两者的距离足够小时，即可认为他们相撞。当导弹与原点坐标的欧氏距离超过射程时，结束循环，认为不能追上
四、代码求解
1.预备知识
（1）mod(m,n)
mod(m,n)表示求m/n的余数
mod(8,3) ans=2
（2）axis([m n p q])
axis([0 3 0 10]) % 设置横坐标范围为[0, 3] 纵坐标范围为[0, 10]
（3）text(m,n,'xxx')
text(2,4,'清风') % 在坐标为(2,4)的点上标上字符串：清风
2.变量初始化
v=200; % 任意给定B船的速度（后期我们可以再改的）
dt=0.0000001; % 定义时间间隔
x=[0,20]; % 定义导弹和B船的横坐标分别为x(1)和x(2)
y=[0,0]; % 定义导弹和B船的纵坐标分别为y(1)和y(2)
t=0; % 初始化导弹击落B船的时间
d=0; % 初始化导弹飞行的距离
m=sqrt(2)/2; % 将sqrt(2)/2定义为一个常量，使后面看起来很简洁
dd=pdist2(x,y); % 导弹与B船的欧式距离
3.初始化画图参数
(1)其中，固定坐标轴的操作，在我们初次实验中可以不做，因为我们并不知道最后的坐标会落在哪个范围
但是在初次试验之后，就可以做了，便于我们观察走向
(2)定义k的原因时为了控制画图速度，因为我们定义的时间间隔delta_t非常小，我们不肯能每迭代一次就画一个点，这样的画图速度非常慢。
因此，我们选择迭代500次再画一次图（在后续代码中体现）
for i=1:2
plot(x(i),y(i),'.k','MarkerSize',1); % 画出导弹和B船所在的坐标，点的大小为1，颜色为黑色(k)，用小点表示
grid on; % 打开网格线
hold on; % 不关闭图形，继续画图
end
axis([0 30 0 10]) % 固定x轴的范围为0-30 固定y轴的范围为0-10
k = 0; % 引入一个变量 为了控制画图的速度（因为Matlab中画图的速度超级慢）
4.进入循环，开始模拟
我们分x个部分来分别讨论此循环
（1）进入循环，迭代数值
此步骤就是根据我们前面建立的模型，进行数值的迭代，设置循环跳出条件为（d>0.01)
while(dd>=0.001) % 只要两者的距离足够大，就一直循环下去。（两者距离足够小时表示导弹击中，这里的临界值要结合dt来取，否则可能导致错过交界处的情况）
t=t+dt; % 更新导弹击落B船的时间
d=d+3*v*dt; % 更新导弹飞行的距离
x(2)=20+t*v*m; y(2)=t*v*m; % 计算新的B船的位置 （注：m=sqrt(2)/2）
dd=sqrt((x(2)-x(1))^2+(y(2)-y(1))^2); % 更新导弹与B船的距离
tan_alpha=(y(2)-y(1))/(x(2)-x(1)); % 计算斜率，即tan(α)
cos_alpha=sqrt(1/(1+tan_alpha^2)); % 利用公式：sec(α)^2 = (1+tan(α)^2) 计算出cos(α)
sin_alpha=sqrt(1-cos_alpha^2); % 利用公式： sin(α)^2 +cos(α)^2 = 1 计算出sin(α)
x(1)=x(1)+3*v*dt*cos_alpha; y(1)=y(1)+3*v*dt*sin_alpha; % 计算新的导弹的位置
（2）绘制路径图像
我们之前初始化了k，就是在这里使用，利用余数函数mod，设置每迭代500次画一次图
并且，利用 pause(0.001)可以看到整个迭代的过程，如果不加这串命令，则只能看到最终效果，而不能看到整个过程
k = k +1 ;
if mod(k,500) == 0 % 每刷新500次时间就画出下一个导弹和B船所在的坐标 mod(m,n)表示求m/n的余数
for i=1:2
plot(x(i),y(i),'.k','MarkerSize',1);
hold on; % 不关闭图形，继续画图
end
pause(0.001); % 暂停0.001s后再继续下面的操作
end
（3）设置另一个跳出循环的条件
整个过程有两种情况，要么在射程中追上，要么超过射程追不上。此端命令则是判断是否超过射程
if d>50 % 导弹的有效射程为50个单位
disp('导弹没有击中B船');
break; % 退出循环
end
（4）输出最终结果
if d<=50 & dd<0.001 % 导弹飞行的距离小于50个单位且导弹和B船的距离小于0.001（表示击中）
disp(['导弹飞行',num2str(d),'个单位后击中B船'])
disp(['导弹飞行的时间为',num2str(t*60),'分钟'])
end
end
（5）完整循环代码
clear;clc
v=200; % 任意给定B船的速度（后期我们可以再改的）
dt=0.0000001; % 定义时间间隔
x=[0,20]; % 定义导弹和B船的横坐标分别为x(1)和x(2)
y=[0,0]; % 定义导弹和B船的纵坐标分别为y(1)和y(2)
t=0; % 初始化导弹击落B船的时间
d=0; % 初始化导弹飞行的距离
m=sqrt(2)/2; % 将sqrt(2)/2定义为一个常量，使后面看起来很简洁
dd=sqrt((x(2)-x(1))^2+(y(2)-y(1))^2); % 导弹与B船的距离
for i=1:2
plot(x(i),y(i),'.k','MarkerSize',1); % 画出导弹和B船所在的坐标，点的大小为1，颜色为黑色(k)，用小点表示
grid on; % 打开网格线
hold on; % 不关闭图形，继续画图
end
axis([0 30 0 10]) % 固定x轴的范围为0-30 固定y轴的范围为0-10
k = 0; % 引入一个变量 为了控制画图的速度（因为Matlab中画图的速度超级慢）
while(dd>=0.001) % 只要两者的距离足够大，就一直循环下去。（两者距离足够小时表示导弹击中，这里的临界值要结合dt来取，否则可能导致错过交界处的情况）
t=t+dt; % 更新导弹击落B船的时间
d=d+3*v*dt; % 更新导弹飞行的距离
x(2)=20+t*v*m; y(2)=t*v*m; % 计算新的B船的位置 （注：m=sqrt(2)/2）
dd=sqrt((x(2)-x(1))^2+(y(2)-y(1))^2); % 更新导弹与B船的距离
tan_alpha=(y(2)-y(1))/(x(2)-x(1)); % 计算斜率，即tan(α)
cos_alpha=sqrt(1/(1+tan_alpha^2)); % 利用公式：sec(α)^2 = (1+tan(α)^2) 计算出cos(α)
sin_alpha=sqrt(1-cos_alpha^2); % 利用公式： sin(α)^2 +cos(α)^2 = 1 计算出sin(α)
x(1)=x(1)+3*v*dt*cos_alpha; y(1)=y(1)+3*v*dt*sin_alpha; % 计算新的导弹的位置
k = k +1 ;
if mod(k,500) == 0 % 每刷新500次时间就画出下一个导弹和B船所在的坐标 mod(m,n)表示求m/n的余数
for i=1:2
plot(x(i),y(i),'.k','MarkerSize',1);
hold on; % 不关闭图形，继续画图
end
pause(0.001); % 暂停0.001s后再继续下面的操作
end
if d>50 % 导弹的有效射程为50个单位
disp('导弹没有击中B船');
break; % 退出循环
end
if d<=50 & dd<0.001 % 导弹飞行的距离小于50个单位且导弹和B船的距离小于0.001（表示击中）
disp(['导弹飞行',num2str(d),'个单位后击中B船'])
disp(['导弹飞行的时间为',num2str(t*60),'分钟'])
end
end