• 大小: 11KB
    文件类型: .rar
    金币: 1
    下载: 0 次
    发布日期: 2021-06-06
  • 语言: Matlab
  • 标签: Anderson  程序  代码  

资源简介

Anderson计算流体力学入门第7章,四个代码程序,分别对应四个visualize画图程序。参数可以任意修改。matlab编程。

资源截图

代码片段和文件信息

%亚音速-超音速,利用U1,U2,U3守恒形式,激波捕捉_

%Agammadtdx
clearclc
NUM=61;
C=0.5;%Courant number
dx=3/(NUM-1);
step=1400;%时间模拟的步数
Time=0;
gamma=1.4;
p_e=0.6784; %出口对总压的比值,此程序模拟激波
x=0:dx:3;
Cx=0.2;%人工粘性常数 0.01~0.3


A=1+2.2*(x-1.5).^2;
dlnA_dx=4.4*(x-1.5)./A;%对A求导
% dlnA_dx=zeros(1NUM);
% for i=2:1:NUM
% dlnA_dx(i)=(log(A(i))-log(A(i-1)))/dx;
% end
% dlnA_dx(1)=dlnA_dx(2);


rou_ar=zeros(step+1NUM);V_ar=zeros(step+1NUM);T_ar=zeros(step+1NUM);
Ma_ar=zeros(step+1NUM);flux_m=zeros(step+1NUM);P_ar=zeros(step+1NUM);
% drou_t=zeros(1NUM);dV_t=zeros(1NUM);dT_t=zeros(1NUM);
% drou_t2=zeros(1NUM);dV_t2=zeros(1NUM);dT_t2=zeros(1NUM);
% drou_t3=zeros(1NUM);dV_t3=zeros(1NUM);dT_t3=zeros(1NUM);
% rou_e=zeros(1NUM);V_e=zeros(1NUM);T_e=zeros(1NUM);
rou_r=zeros(1NUM);V_r=zeros(1NUM);T_r=zeros(1NUM);
dt_ar=zeros(1step);
rou=zeros(1NUM);T=zeros(1NUM);
dU1_t=zeros(1NUM);dU2_t=zeros(1NUM);dU3_t=zeros(1NUM);
U1_e=zeros(1NUM);U2_e=zeros(1NUM);U3_e=zeros(1NUM);
dU1_t2=zeros(1NUM);dU2_t2=zeros(1NUM);dU3_t2=zeros(1NUM);
dU1_t3=zeros(1NUM);dU2_t3=zeros(1NUM);dU3_t3=zeros(1NUM);
U1_r=zeros(1NUM);U2_r=zeros(1NUM);U3_r=zeros(1NUM);
U1_ar=zeros(step+1NUM);U2_ar=zeros(step+1NUM);U3_ar=zeros(step+1NUM);
S1_e=zeros(1NUM);S2_e=zeros(1NUM);S3_e=zeros(1NUM);
S1_r=zeros(1NUM);S2_r=zeros(1NUM);S3_r=zeros(1NUM);

%初始化
for i=1:1:NUM
    if x(i)>=0 && x(i)<=0.5
       rou(i)=1;
       T(i)=1;
    elseif x(i)>0.5 && x(i)<=1.5
        rou(i)=1-0.366*(x(i)-0.5);
        T(i)=1-0.167*(x(i)-0.5);
    elseif x(i)>1.5 && x(i)<=2.1
        rou(i)=0.634-0.702*(x(i)-1.5);
        T(i)=0.833-0.4908*(x(i)-1.5);
    elseif x(i)>2.1 && x(i)<=3.0
        rou(i)=0.5892+0.10228*(x(i)-2.1);
        T(i)=0.93968+0.0622*(x(i)-2.1);
    end
end
V=0.59./(rou.*A);  %0.59是更配合稳态解的数值。

%
rou_ar(1:)=rou;
V_ar(1:)=V;
T_ar(1:)=T;
P_ar(1:)=rou.*T;
Ma_ar(1:)=V./sqrt(T);
flux_m(1:)=rou.*V.*A;

for k=1:1:step   %模拟的时间步数

rou=rou_ar(k:);
V=V_ar(k:);
T=T_ar(k:);
P=P_ar(k:);

U1=rou.*A;
U2=rou.*V.*A;
U3=rou.*A.*((T/(gamma-1))+(gamma/2)*V.^2);
F1=U2;
F2=(U2.^2)./U1+(gamma-1)/gamma*(U3-(gamma/2)*(U2.^2)./U1);
F3=gamma*U2.*U3./U1-(gamma*(gamma-1)/2)*((U2.^3)./(U1.^2));
J2=(gamma-1)/gamma*(U3-gamma/2*(U2.^2)./U1).*(dlnA_dx);


%确定本次模拟的dt时间间隔
dt_list=C*dx./(T.^0.5+V);
dt=min(dt_list);
%
%预测步
for i=2:1:(NUM-1)

    dU1_t(i)=-(F1(i+1)-F1(i))/dx;
    dU2_t(i)=-(F2(i+1)-F2(i))/dx+J2(i);
    dU3_t(i)=-(F3(i+1)-F3(i))/dx;
    
    S1_e(i)=Cx*((abs(P(i+1)-2*P(i)+P(i-1)))/(P(i+1)+2*P(i)+P(i-1)))*(U1(i+1)-2*U1(i)+U1(i-1));
    S2_e(i)=Cx*((abs(P(i+1)-2*P(i)+P(i-1)))/(P(i+1)+2*P(i)+P(i-1)))*(U2(i+1)-2*U2(i)+U2(i-1));
    S3_e(i)=Cx*((abs(P(i+1)-2*P(i)+P(i-1)))/(P(i+1)+2*P(i)+P(i-1)))*(U3(i+1)-2*U3(i)+U3(i-1));
    
    U1_e(i)=U1(i)+dU1_t(i)*dt+S1_e(i);
    U2_e(i)=U2(i)+dU2_t(i)*dt+S2_e(i);
    U3_e(i)=U3(i)+dU3_t(i)*dt+S3_e(i);
    
end

%BC
U1_e(1)=A(1);
U2_e(1)=2*U2_e(2)-U2_e(3);

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       6109  2012-12-09 18:40  Chapter7\one_D_shockwave_4.m

     文件       3934  2012-12-08 18:12  Chapter7\one_D_subsonic_all_2.m

     文件       3402  2012-12-09 00:46  Chapter7\one_D_supersonic_all_1.m

     文件       5056  2012-12-09 15:50  Chapter7\one_D_supersonic_all_3.m

     文件        388  2012-12-09 15:43  Chapter7\Untitled.m

     文件       5033  2012-12-09 00:23  Chapter7\Untitled3.m

     文件       1208  2012-12-08 11:19  Chapter7\visualize_1.m

     文件       1450  2012-12-08 18:04  Chapter7\visualize_2.m

     文件       1213  2012-12-09 00:28  Chapter7\visualize_3.m

     文件       1218  2012-12-09 18:43  Chapter7\visualize_4.m

     目录          0  2012-12-10 23:15  Chapter7

----------- ---------  ---------- -----  ----

                29011                    11


评论

共有 条评论