nihequxian.m

这是我帮一个正在读博士的老师编的关于爆破的m文件

clear all;
U=600;
h=40*10^(-6);
h1=400*10^(-6);
R=2*10^2;
A2=0.332;
b=1.2*10^(-5);
T=40*10^(-6);
for y=0
    x=0;
for n=1:1:100
x=x+12*(n/n);
if x<=1200 
    t=x/U;
    u1=0;
    for i=1:1:3
        u1=u1+(((1/i)*exp(-((i*i*pi*pi*b)*t/(h*h)))*sin((i*pi*y)/h)));               
    end
      u=U*(1-(y/h))-(2*U/pi)*u1;   
    v1=0;
    for w=0:1:3
        c=[1 1 11 375];
        q=[2 120  40320  39916800];        
        v1=v1+(((-0.5)^w)*c(w+1)*(A2^(w+1))/(q(w+1)));    
end 
 v=0.5*u*y/x-0.5*sqrt(b*U/x)*v1;  
 V1=sqrt(u*u+v*v);
 V=V1*T*cos((8*pi/h1)*x*(10^(-6)));
end 
 f(n)=x;
 g(n)=V; 
end
end
for y=10*10^(-6)
    x=0;
for n=1:1:100
x=x+12*(n/n);
if x<=1200 
    t=x/U;
    u1=0;
    for i=1:1:3
        u1=u1+(((1/i)*exp(-((i*i*pi*pi*b)*t/(h*h)))*sin((i*pi*y)/h)));               
    end
      u=U*(1-(y/h))-(2*U/pi)*u1;   
    v1=0;
    for w=0:1:3
        c=[1 1 11 375];
        q=[2 120  40320  39916800];        
        v1=v1+(((-0.5)^w)*c(w+1)*(A2^(w+1))/(q(w+1)));    
end 
 v=0.5*u*y/x-0.5*sqrt(b*U/x)*v1;  
 V1=sqrt(u*u+v*v);
 V=V1*T*cos((8*pi/h1)*x*(10^(-6)));
end 
 f(n)=x;
 g1(n)=V; 
end
end
for y=20*10^(-6)
    x=0;
for n=1:1:100
x=x+12*(n/n);
if x<=1200 
    t=x/U;
    u1=0;
    for i=1:1:3
        u1=u1+(((1/i)*exp(-((i*i*pi*pi*b)*t/(h*h)))*sin((i*pi*y)/h)));               
    end
      u=U*(1-(y/h))-(2*U/pi)*u1;   
    v1=0;
    for w=0:1:3
        c=[1 1 11 375];
        q=[2 120  40320  39916800];        
        v1=v1+(((-0.5)^w)*c(w+1)*(A2^(w+1))/(q(w+1)));    
end 
 v=0.5*u*y/x-0.5*sqrt(b*U/x)*v1;  
 V1=sqrt(u*u+v*v);
 V=V1*T*cos((8*pi/h1)*x*(10^(-6)));
end 
 f(n)=x;
 g2(n)=V; 
end
end
for y=30*10^(-6)
    x=0;
for n=1:1:100
x=x+12*(n/n);
if x<=1200 
    t=x/U;
    u1=0;
    for i=1:1:3
        u1=u1+(((1/i)*exp(-((i*i*pi*pi*b)*t/(h*h)))*sin((i*pi*y)/h)));               
    end
      u=U*(1-(y/h))-(2*U/pi)*u1;   
    v1=0;
    for w=0:1:3
        c=[1 1 11 375];
        q=[2 120  40320  39916800];        
        v1=v1+(((-0.5)^w)*c(w+1)*(A2^(w+1))/(q(w+1)));    
end 
 v=0.5*u*y/x-0.5*sqrt(b*U/x)*v1;  
 V1=sqrt(u*u+v*v);
 V=V1*T*cos((8*pi/h1)*x*(10^(-6)));
end 
 f(n)=x;
 g3(n)=V; 
end
end
for y=40*10^(-6)
    x=0;
for n=1:1:100
x=x+12*(n/n);
if x<=1200 
    t=x/U;
    u1=0;
    for i=1:1:3
        u1=u1+(((1/i)*exp(-((i*i*pi*pi*b)*t/(h*h)))*sin((i*pi*y)/h)));               
    end
      u=U*(1-(y/h))-(2*U/pi)*u1;   
    v1=0;
    for w=0:1:3
        c=[1 1 11 375];
        q=[2 120  40320  39916800];        
        v1=v1+(((-0.5)^w)*c(w+1)*(A2^(w+1))/(q(w+1)));    
end 
 v=0.5*u*y/x-0.5*sqrt(b*U/x)*v1;  
 V1=sqrt(u*u+v*v);
 V=V1*T*cos((8*pi/h1)*x*(10^(-6)));
end 
 f(n)=x;
 g4(n)=V; 
end
end
plot(f,g,'r',f,g1,'k',f,g2,'g',f,g3,'b',f,g4,'y');