优化下列MATLAB代码：Mx=1000; Mt=1000; T=0.2; %计算 0.2 时刻的值 L=2; gama=1.4; deltax=L/Mx; deltat=T/Mt;%在 x，t 上划分网格并设置参数 n=1; U1=[0;0;0]; U2=[0;0;0]; F1=[0;0;0]; F2=[0;0;0]; p=zeros(ceil(T/deltat)+1,Mx+1); u=zeros(ceil(T/deltat)+1,Mx+1); P=zeros(ceil(T/deltat)+1,Mx+1); E=zeros(ceil(T/deltat)+1,Mx+1);%设定初始矩阵 p(1,:)=[ones(1,Mx/2) 1 0.125*ones(1,Mx/2)]; u(1,:)=[0.75*ones(1,Mx/2) 0.75 zeros(1,Mx/2)]; P(1,:)=[1*ones(1,Mx/2) 1 0.1*ones(1,Mx/2)]; E(1,:)=P(1,:)./(p(1,:)*(gama-1))+0.5*u(1,:).^2;%赋初始值 while n*deltat<=T for j=1:1:Mx+1 if j==1||j==2||j==Mx+1||j==Mx u(n+1,j)=u(n,j); P(n+1,j)=P(n,j); p(n+1,j)=p(n,j); else U1=U(n,j,u,p,P,gama); F2=L_W(deltax,deltat,u,P,p,n,j,gama); F1=L_W(deltax,deltat,u,P,p,n,j-1,gama); U2=U1-deltat/deltax*(F2-F1); p(n+1,j)=U2(1); u(n+1,j)=U2(2)/U2(1); P(n+1,j)=(gama-1)*(U2(3)-0.5*p(n+1,j)*u(n+1,j)^2); end end n=n+1; end x=-1:deltax:1; y1=p(n,:);y2=u(n,:);y3=P(n,:); plot(x,y1,'r',x,y2,'g',x,y3,'k'); xlabel('X'); ylabel('Y(0.02)');%设置横，纵坐标并画图 title('row/u/p--t'); legend('row','u','p');

可以尝试如下优化：

1. 将重复计算的值提取出来，避免重复计算，如：

ceil_T_deltat = ceil(T/deltat)+1;
Mx_plus_1 = Mx+1;
half_Mx = Mx/2;

2. 将循环中的条件判断提取出来，避免每次循环都进行判断，如：

for j=1:Mx_plus_1
    if j~=1 && j~=2 && j~=Mx_plus_1 && j~=Mx
        ...
    else
        u(n+1,j)=u(n,j);
        P(n+1,j)=P(n,j);
        p(n+1,j)=p(n,j);
    end
end

3. 将计算过程中用到的计算公式提取出来，避免代码重复，提高代码可读性，如：

U1 = U(n,j,u,p,P,gama);
F2 = L_W(deltax,deltat,u,P,p,n,j,gama);
F1 = L_W(deltax,deltat,u,P,p,n,j-1,gama);
U2 = U1-deltat/deltax*(F2-F1);
p(n+1,j) = U2(1);
u(n+1,j) = U2(2)/U2(1);
P(n+1,j) = (gama-1)*(U2(3)-0.5*p(n+1,j)*u(n+1,j)^2);

4. 可以使用向量化计算，避免循环，提高代码效率，如：

x = -1:deltax:1;
n = 1:ceil_T_deltat;
j = 3:Mx;
while n(end)*deltat<=T
    U1 = U(n,j,u,p,P,gama);
    F2 = L_W(deltax,deltat,u,P,p,n,j,gama);
    F1 = L_W(deltax,deltat,u,P,p,n,j-1,gama);
    U2 = U1-deltat/deltax*(F2-F1);
    p(n+1,j) = U2(1);
    u(n+1,j) = U2(2)./U2(1);
    P(n+1,j) = (gama-1)*(U2(3)-0.5*p(n+1,j).*u(n+1,j).^2);
    p(n+1,[1 2 Mx_plus_1]) = p(n,[1 2 Mx_plus_1]);
    u(n+1,[1 2 Mx_plus_1]) = u(n,[1 2 Mx_plus_1]);
    P(n+1,[1 2 Mx_plus_1]) = P(n,[1 2 Mx_plus_1]);
    n = n+1;
end
y1 = p(n,:);
y2 = u(n,:);
y3 = P(n,:);
plot(x,y1,'r',x,y2,'g',x,y3,'k');
xlabel('X');
ylabel('Y(0.02)');
title('row/u/p--t');
legend('row','u','p');