一、Matlab代码片

%亚声速-超声速等熵喷管流动 非守恒型麦考马克方法数值求解
clear; %清理内存变量
clc; %清理工作窗中的所有显示内容
r=1.4; %比热比
L=3; %喷管长度
i=31; %网格数目
C=0.5; %科朗数
x=linspace(0,L,i); %网格点横坐标
N=1401; %时间步
dx=L/(i-1); %空间步长
dt(1:N)=0; %时间步长
A=1+2.2*(x-1.5).*(x-1.5); %喷管面积
rho(N,i)=0; %流场密度赋值
T(N,i)=0; %流场温度赋值
V(N,i)=0; %流场速度赋值
%预分配内存
rhotav2(1:1400,1:30)=0;  Vtav2(1:1400,1:30)=0;
%初始条件
rho(1,:)=1-0.3146*x; %流场密度初值
T(1,:)=1-0.2314*x; %流场温度初值
V(1,:)=(0.1+1.09*x).*(1-0.2314*x).^0.5; %流场速度初值
%按时间步长推进
for k=1:N-1%计算预估步偏导数rhot(1:i-1)=-V(k,1:i-1).*(rho(k,2:i)-rho(k,1:i-1))/dx-rho(k,1:i-1)....*(V(k,2:i)-V(k,1:i-1))/dx-rho(k,1:i-1).*V(k,1:i-1).*(log(A(2:i))-log(A(1:i-1)))/dx;Vt(1:i-1)=-V(k,1:i-1).*(V(k,2:i)-V(k,1:i-1))/dx-1/r.*((T(k,2:i)...-T(k,1:i-1))/dx+T(k,1:i-1)./rho(k,1:i-1).*(rho(k,2:i)-rho(k,1:i-1))/dx);Tt(1:i-1)=-V(k,1:i-1).*(T(k,2:i)-T(k,1:i-1))/dx-(r-1).*T(k,1:i-1)....*((V(k,2:i)-V(k,1:i-1))/dx+V(k,1:i-1).*(log(A(2:i))-log(A(1:i-1)))/dx);%确定最小时间步长t=C*dx./(V(k,:)+sqrt(T(k,:)));dt(k)=min(t);%计算预估值rho_(1:i-1)=rho(k,1:i-1)+rhot(1:i-1).*dt(k);V_(1:i-1)=V(k,1:i-1)+Vt(1:i-1).*dt(k);T_(1:i-1)=T(k,1:i-1)+Tt(1:i-1).*dt(k);%校正偏导数rhot_(2:i-1)=-V_(2:i-1).*(rho_(2:i-1)-rho_(1:i-2))/dx-rho_(2:i-1)....*(V_(2:i-1)-V_(1:i-2))/dx-rho_(2:i-1).*V_(2:i-1).*(log(A(2:i-1))-log(A(1:i-2)))/dx;Vt_(2:i-1)=-V_(2:i-1).*(V_(2:i-1)-V_(1:i-2))/dx-1/r.*((T_(2:i-1)...-T_(1:i-2))/dx+T_(2:i-1)./rho_(2:i-1).*(rho_(2:i-1)-rho_(1:i-2))/dx);Tt_(2:i-1)=-V_(2:i-1).*(T_(2:i-1)-T_(1:i-2))/dx-(r-1).*T_(2:i-1)....*((V_(2:i-1)-V_(1:i-2))/dx+V_(2:i-1).*(log(A(2:i-1))-log(A(1:i-2)))/dx);%时间导数平均值rhotav(2:i-1)=0.5*(rhot(2:i-1)+rhot_(2:i-1)); rhotav2(k,2:i-1)=abs(rhotav(2:i-1)); Vtav(2:i-1)=0.5*(Vt(2:i-1)+Vt_(2:i-1));Vtav2(k,2:i-1)=abs(Vtav(2:i-1));Ttav(2:i-1)=0.5*(Tt(2:i-1)+Tt_(2:i-1));%流场变量校正值rho(k+1,2:i-1)=rho(k,2:i-1)+rhotav(2:i-1)*dt(k);V(k+1,2:i-1)=V(k,2:i-1)+Vtav(2:i-1)*dt(k);T(k+1,2:i-1)=T(k,2:i-1)+Ttav(2:i-1)*dt(k);%入流边界值V(k+1,1)=2*V(k+1,2)-V(k+1,3);rho(k+1,1)=1;T(k+1,1)=1;%出流边界值rho(k+1,i)=2*rho(k+1,i-1)-rho(k+1,i-2);V(k+1,i)=2*V(k+1,i-1)-V(k+1,i-2);T(k+1,i)=2*T(k+1,i-1)-T(k+1,i-2);%马赫数Ma=V(1:k+1,1:i)./(sqrt(T(1:k+1,1:i)));%流场压强P=rho(1:k+1,1:i).*T(1:k+1,1:i);
end
%绘图 喷管喉道处密度、温度、压力和马赫数的变化
figure;
subplot(2,2,1),plot(1:1401,rho(1:1401,16),'b-');
ylabel('\rho/\rho_0'),xlabel('Timestep'); title('无量纲密度变化');
subplot(2,2,2),plot(1:1401,T(1:1401,16),'r-');
ylabel('T/T_0'),xlabel('Timestep');title('无量纲温度变化');
subplot(2,2,3),plot(1:1401,P(1:1401,16),'k-');
ylabel('P/P_0'),xlabel('Timestep');title('无量纲总压变化');
subplot(2,2,4),plot(1:1401,Ma(1:1401,16),'m-');
ylabel('Ma'),xlabel('Timestep');title('马赫数变化');
%绘图 喷管喉道处无量纲密度和速度时间导数的变化
figure;
plot(1:1400,rhotav2(1:1400,16),'k-');
ylabel('残差'),xlabel('Timestep');
title('喷管喉道处无量纲密度和速度时间导数的变化');
hold on;
plot(1:1400,Vtav2(1:1400,16),'k--');
%绘图 无量纲质量流量在六个不同时刻的瞬时分布
figure;
plot(x,rho(1,1:i).*V(1,1:i).*A(1:i),'k--');
hold on;
plot(x,rho(51,1:i).*V(51,1:i).*A(1:i),'r-');
hold on;
plot(x,rho(101,1:i).*V(101,1:i).*A(1:i),'b-');
hold on;
plot(x,rho(151,1:i).*V(151,1:i).*A(1:i),'g-');
hold on;
plot(x,rho(201,1:i).*V(201,1:i).*A(1:i),'y-');
hold on;
plot(x,rho(701,1:i).*V(701,1:i).*A(1:i),'m-');
title('不同时刻质量流量的变化');
ylabel('无量纲质量流量'),xlabel('无量纲轴向距离');
text(2.4,0.63,'700\Deltat');
text(2.4,0.75,'150\Deltat');
text(2.4,0.53,'200\Deltat');
text(2.4,0.9,'50\Deltat');
text(2.4,1.05,'100\Deltat');
text(2.4,1.3,'0\Deltat');

二、计算结果展示

1.喷管喉道处密度、温度、压力和马赫数的变化

2.喷管喉道处无量纲密度和速度时间导数的变化

3.无量纲质量流量在六个不同时刻的瞬时分布

拟一维喷管流动的数值解——亚声速-超声速等熵喷管流动的非守恒型CFD解法（MacCormack方法）相关推荐

拟一维喷管流动的数值解——亚声速-超声速等熵喷管流动的守恒型CFD解法（MacCormack方法）
一.Matlab代码片 %亚声速-超声速等熵喷管流动守恒形CFD解法 MacCormack方法 clear; %清理内存变量 clc; %清理工作窗中的所有显示内容 r=1.4; %比热比 L=3; ...
亚声速-超声速等熵喷管拟一维流动的CFD解法（附完整代码）
入门CFD,主要参考书目<计算流体力学基础及其应用>(John D.Anderson 著,吴颂平等译) 实现了第 7.3 节亚声速-超声速等熵喷管流动的CFD解法的代码,采用的是 ...
亚声速 – 超声速等熵喷管流动数值模拟(文字）
亚声速 – 超声速等熵喷管流动问题描述流过拉伐尔喷管的定常等熵流动在喷管,流动经过等熵地加速已达到超声速喷管出口处流体的压力.温度.速度和马赫数分别记为在喷管收缩段,流动是亚声速的:在喷管的 ...
喷管流动的守恒型CFD解法及激波捕捉（附完整代码）
入门CFD,主要参考书目<计算流体力学基础及其应用>(John D.Anderson 著,吴颂平等译) 实现了第 7.6 节激波捕捉的代码,采用的是 MacCormack 方法, ...
拟一维喷管流动的数值解——全亚声速等熵喷管流动的非守恒型CFD解法（MacCormack方法）
一.Matlab代码片 %全亚声速等熵喷管流动非守恒型麦考马克方法数值求解 clear; %清理内存变量 clc; %清理工作窗中的所有显示内容 r=1.4; %比热比 L=3; %喷管长度 i=3 ...
Matlab 编程《计算流体力学基础及应用（约翰D安德森）》全亚声速等熵喷管流动CFD解法拉瓦尔喷管非守恒形式方程解法
Matlab 编程 <计算流体力学基础及应用(约翰D安德森)> 全亚声速等熵喷管流动CFD解法拉瓦尔喷管非守恒形式方程解法问题之全亚声速等熵喷管流动的CFD解法初始化参数迭代过 ...
持续不定期更新：CFDC++之拟一维喷管流动的数值解（1）
前言: 学习openfoam,或者准确来说,CFD,断断续续也有5个多月.3月底在本科学校写下的第一篇openfoam学习的博文,到现在在研一的学校里,尝试自己写代码计算cfd问题,慢慢一点一点地构建 ...
持续不定期更新：CFDC++之拟一维喷管流动的数值解（2）
这篇博文将剩下的问题解决完.第一篇在: https://mp.csdn.net/postedit/101038218 在初始化步骤之后,就到了计算下一时间步的步骤了.计算之前先讲一讲这里用到的计算方法 ...
用Python实现拟一维喷管流动的数值解
import math import decimal import numpy as np from decimal import * wnx=31#最后一个网格记数31 dx=3/(wnx-1)#网 ...

拟一维喷管流动的数值解——亚声速-超声速等熵喷管流动的非守恒型CFD解法（MacCormack方法）

一、Matlab代码片

二、计算结果展示

1.喷管喉道处密度、温度、压力和马赫数的变化

2.喷管喉道处无量纲密度和速度时间导数的变化

3.无量纲质量流量在六个不同时刻的瞬时分布

拟一维喷管流动的数值解——亚声速-超声速等熵喷管流动的非守恒型CFD解法（MacCormack方法）相关推荐

最新文章

热门文章