这是一个分析子弹入水过程速度、位移和加速度随时间变化关系的代码。目前是可以运行了，得到的结果是矩形，现在想求解出图像（速度v, 位移z, 加速度a与时间t)的关系，图像代码该怎么编呢？）

clear,clc,close all
%frenshxyn
M = 50/1000;%Kg
viscosity = 1.01*10^-3;
den_w = 1000;%Kg/m^3
den_p = 8.9*1000;%Kg/m^3
V = M/den_p;
d = 0.0221;%单位米
h = (V-2/3*pi*(d/2)^3)/(pi*(d/2)^2);%单位米
v(1) = 300;%米/秒
z(1)=0;
t(1)=0;
Fd(1)=0;
a(1)=0;
m(1)=0;
Cd0=0.2601;
%位移步长
dP=0.005;%单位米
%弹丸外形
z=0:0.00005:0.01835;%半球头子弹的高，单位米
r=zeros(size(z));
for i=1:length(z)
if (z(i)>=0)&(z(i)<0.01105)
r(i)=((d/2)^2-(d/2-z(i))^2)^0.5;
elseif (z(i)>=0.01105)&(z(i)<0.01835)
r(i)=d/2;
end
n=2211;%代表交界面
k=1;
while (v(k)>0)&(z(i)<0.0184)
if k<=n-1
j=k;
dz(k)=z(k+1)-z(k);
else j=n-1;
dz(k)=dP;
end 
m(k+1)=2/3.*den_w*pi*[r(j+1)]^3;
kd(k)=M/(M+m(k+1));
A(k+1)=pi*[r(j+1)]^2;
dm(k)=m(k+1)-m(k);
Cdm(k)=2/(den_w*A(k+1))*dm(k)/dz(k);
Cd(k)=kd(k)*[Cd0+Cdm(k)];
Fd(k+1)=1/2*den_w*[v(k)]^2*A(k+1)*Cd(k);
a(k+1)=-Fd(k+1)/M;
v(k+1)=([v(k)]^2+2*a(k+1)*dz(k))^0.5;
dt(k)=2*dz(k)/(v(k)+v(k+1))
z(k+1)=z(k)+dz(k);
t(k+1)=t(k)+dt(k);
k=k+1;
end
%绘图
end

参考图

参考下这个链接
https://b23.tv/iATkV8J

没有语法错误不代表没有逻辑错误，一般我遇到这种情况，会通过多输出一些检查点、或注释一些代码、用debug的方法去定位原因
另外，你写的代码最好是用插入代码的格式，不然不好分析

贴张程序代码的图片看看，提问中的代码被改格式有语法错误。