1、1五类黎曼间断问题数值模拟题目:1对于一维流 Euler 方程: 0, 0txtxUFUA( 当 地 线 性 化 后 为 )试用通量的矩阵分裂格式的思想,构造此方程的数值解法。左右边界不考虑均使用 Neumann 边界条件: 10x2编写 MatLab 程序计算本章第六节中所述的 5 类黎曼间断问题中的任意三类问题。3模拟初场如下的非定常流,并尝试分析常数 k 取何数值可使流场最大马赫数大于或小于 0.3:(0,)max(0,)1,1,.,pkuConst可 设计算过程如下:1.通量的矩阵分裂格式:将方程组: 0UFtx离散后得: 11/2/()niiii2方程空间导数项的差分计算转化为网格单
2、元界面处的通量计算。而双曲型问题的主要特点是存在特征线,并且扰动以特征速度沿特征线传播。不同的特征分量,传播方向不同。为体现迎风的思想,需要考虑这种方向性。以一维情况为例,由于影响来自上游,对于向右传播的分量,应该采用左边的值计算,对于向左传播的分量,采用右边的值计算。按照这种思想,需要进行特征分裂。对于方程组0UFtx,当地线性化后为:0UAtx。可将 A写成 1T的形式,其中 12diag(,)n 。方程组左乘1,有: 11 0Ttx 定义 chU,使得:1chU,则有:ChChUt将 分解成正部和负部,即 ,其中:1()2,1()2对这两种方向不同的特征波,考虑其影响域,可以导出界面处通
3、量 F的表达式。由于: F U()RLRLA也可以写成:F11Ch ()FhTTT对于特征形式的通量 ,同样有:ChhFURChhCh=()()()LRLUChhF3也即有: ChRh()LChFF界面处的通量为: hhC()=()+LR界 面也可以直接由左侧或右侧得到: ChhRhhCLC()=()()()-LLR FFF界 面界 面上式左乘 ,得:TRL=hCTF界 面界 面对从左右两侧得到的通量进行平均,则得: h1()2RLU“”表示代入 Roe 平均值进行计算。Roe平均值:, ,1LRuD1LRHRLD,LR2()au前次作业已经导出矩阵 A的形式,并求出其特征值:123,uau,
4、且相应的矩阵 T及向量 chU的形式为:211TMPauH4212()chpauUTa因此,通量的矩阵分裂格式: 11/2/Ch()F()2niiiiRLtUFxTU2.黎曼间断问题:一维流动,初始条件: ,(,),01(,),0LRupxupx即所有物理量在 t=0 时刻的值在 x=0 处有间断,左右两侧常数分布。这种初始条件的间断一般是不满足间断关系的,因而是不稳定的。在 后,t这种初始间断立即分解为若干满足间断关系的间断及中心稀疏波或激波。由于初始情况不同,一般该问题可以分解为下图五种类型中的一种。(虚线表示接触间断,单条线表示激波,一束线表示稀疏波)Riemann 问题解的类型的判别规
5、则:5因此对于该问题,用 roe 格式离散,初始条件按上面的方式给定,边界条件采用 Neumann 边界条件: 10xU即可编写程序进行计算,现计算出三类结果,罗列如下:(1)左、右都是激波 , ,0,.1RRup,3,.2LLup6(2)左稀疏波,右激波 ,1,3LLup1,0,.1RRup(3)中间间断面,两边稀疏波,2,0,1.5LLup,0,1RRup7上面三幅图中,速度、压力、密度同时出现突变,对应于激波间断;速度、压力、密度同时出现连续变化,对应于稀疏波;而速度、压力不变,密度突变,对应接触间断。3.带压力脉冲的初场计算:按条件: (0,)max(0,)1,1,.,pkuConst
6、可 设在程序中设定初场,即: ax(0)1pkA, ,0u,在每个时间步计算各节点的马赫数分布,实际计算中发现,const在时间推进的过程中,马赫数分布逐渐在两侧出现尖峰,之后尖峰增大,直到达到某个最大值,然后再减小,并且尖峰的位置向两侧传播。因此在每个时间步求出各节点马赫数的极大值,看其是否超过 0.3。计算结果如下:CFL=0.5,N=100 时, ,最大马赫数大于 0.3。 时,某一17.3k1.3k8时刻的马赫数分布如下,最大马赫数 0.30039经过计算发现,不同的 CFL 数和节点个数 N 的设置对计算结果有影响。节点个数增多,局部马赫数大于 0.3 需要的 k 值减小;CFL 数减小,需要的 k 值则增大。具体如下:(1)CFL=0.5,N=200 时, 时即出现了马赫数大于 0.314.8k9(2)CFL=0.3,N=100 时, 时出现了马赫数大于 0.318.3k