（计算数学专业论文）多介质流体计算的laxfriedrichs格式(1).pdf

摘要 多介质流体( m u l t i f l u i df l o w s ) 的计算是近年来计算流体力学中的个热点相 较于单介质的流体，它的主要难点在于介质界面的处理上很多对于单介质流体十 分成功的计算格式应用到多介质流体上时，在介质界面处都产生了伪振荡如何消 除这些不必要的振荡，也就成了算法研究的重点 在这篇论文中，我们设计了一个基于非平衡状态的l a x - f r i e d r i c h s 格式，用以 对多介质流体进行数值模拟我们从对于两种介质流体标准的l a x - l m e d r i c h s 格式 入手进行研究，发现不同介质间的数值热传导是造成压强上伪振荡的主要原因基 于这个发现。我们通过使两种介质绝热从而阻止它们之间不应该有的热传导。从而 清f 象压强上的伪振荡这么傲，流体在介质交界面处的网格内将不再处于热平衡状 态，伴随而来我们作了几个必要的前提假设由此我们所设计的格式对于总质量 各介质的分质量、动量、以及总能量是完全守恒的在计算中，它保持总密度及各 种介质的分密度始终非负更重要的是在流体通过介质边界处时速度和压强能保持 一致我们用具体数值例子来验证算法的有效性 关键字：多介质流，平衡态，非平衡态，数值热传导 a b s t r a c t n u m e r i c a ls i m u l a t i o no fm u l t i f l u i df l o w si sah o ts p o ti nc o m p u t a t i o n a lf l u i d d y n a m i e s c o m p a r e dw i t ht h a to fs i n g l ef l u i dc o m p u t a t i o n ，t h em a i nd i 伍c u l t yo f m u l t i f u i l ds i m u l a t i o ni 8t ot r e a tt h em a t e r i a li n t e f f a e e s s c h e m e st h a ta r es u c c e s s f u l i ns o l v i n g 咖西ef i l i l df l o ww i l lp r o d u c es p u r i o u so s c f f i a t m u sn e a rm a t e r i a li n t e r - f a c e s h o wt oe l i m i n a t et h e s es p l l r i o a so s c i l l a t i o n st h e nb e c o m e st h em a i np o i n ti n d e s i g n i n gg o o dn u m e r i c a lm e t h o d s i nt h i sp a p e rw ep r e s e n tal a x - f r i e d r i c h ss c h e m eb a s e do nn o n e q u i l i b r i u m t h e r m a ls t a r e sf o rn u m e r i c a ls i m u l a t i o no fm u l t i f l u i df l o w s w js t a r tw i t ha ni n - v e s t i g a t i o no ft h es t a n d a r dl a x - f r i e d r i c h ss c h e m ef o rt w o - f l u i df l o w sa n df i n dt h a t t h en u m e r i c a lh e a tc o n d u c t i o nb e t w e e nt h ed i f f e r e n ts p e c i e si st h em a i nc a u s ef o r p r o d u c i n gs p l l r i o a so s c i l l a t i o n si np r e s s u r ef i e l d u p o nt h i sf i n d i n gw es t o pt h en u o m e r i c a lh e a tc o n d u c t i o nb e t w e e nd i f f e r e n ts p e c i e sb yk e e p i n gt h e ma d i a b a t i cw i t h e a c ho t h e ri nt h es c h e m e i nt h i sw a y , t h ef l o wn o wi ne a c hg r i dc e l li si nan o n e q u i - l i b r i u ms t a t e o nw h i c hs e v e r a la s s u m p t i o n sa r em a d e t h es c h e m ed e s i g n e di nt h i s a p p r o a c hi sf u l l yc o n s e r v a t i v ei nt h es e n s et h a ti tc o n s e r v e st h et o t a la 8w e l la 8i n d i - v i d u a lm a s s e 8 m o m e n t u ma n dt o t a le n e r g y , i tm a i n t a i a st h ep o s i t i v i t yo ft h et o t a l a n di n d i v i d u a ld e n s i t i e s ，a n dm o r ei m p o r t a n ti tm a i n t a i n st h eu n i f o r m n e s so ft h e v e l o c i t ya n dp r e s s u r ea c r o s sm a t e r i a li n t e r f a c e s n u m e r i c a le x a m p l e sa r ed i s p l a y e d t os h o wt h ee f f i e n c yo ft h es c h e m ei ns i m u l a t i n gt w o - c o m p o n e n tf l o w s k e yw o r d s ：m u l t i f l u i df l o w ，e q u i l i b r i u m ，n o n e q u i l i b r i u m ，n u m e r i c a lh e a t c o n d u c t i o n 原创性声明 本人声明：所呈交的论文是本人在导师指导下进行的研 究工作除了文中特别加以标注和致谢的地方外，论文中不包 含其他人已发表或撰写过的研究成果参与同一工作的其他同 志对本研究所做的任何贡献均已在论文中作了明确的说明并表 示了谢意。 签名：日期： 本论文使用授权说明 本人完全了解上海大学有关保留、使用学位论文的规 定，e p = 学校有权保留论文及送交论文复印件，允许论文被查 阅和借阅；学校可以公布论文的全部或部分内容 ( 保密的论文在解密后应遵守此规定) 签名：导师签名。日期： 第一章绪论 数值模拟多介质流体流动的复杂现象是近些年来计算流体力学中的大研究热 点在这一章中，我们将主要介绍多介质流的同题描述，困难所在以及这领域的 发展现状 1 1 多介质气体动力学问题的介绍以及它的发展现状 多介质流体计算的的主要难点在于数值模拟可能造成速度和压强在介质界面附 近的伪振荡，因为即使是最好的具有一阶精度的格式，例如l a x - f r i e d r i c h s 格式，以 及g o d u n o v 格式，在计算多介质流体介质分界面时已经产生了伪振荡，而从这些 格式扩展而来的更高阶的格式则更难消除这些振荡关于这一点在文献 1 】、1 1 3 中有较详细的论述在接触面以外的位置，流体实际上是单介质流，对于单介质流 已有很多有效的数值计算方法而如何处理多介质流介质边界的间断问题，则成了 主要的研究对象处理问断的数值方法大致可以分为两类，跟踪法和捕捉法，它们 各有利弊 1 2 问题的数学描述 首先，我们给出本文所研究的问题的数学描述一般的无粘可压缩气体运动都 可以用e u l e r 方程来描述，它由流体的三个守恒律所组成，即质量守恒，动量守恒， 以及能量守恒具体方程如下， 风+ ( 班) 。= 0 ( 刖) f + ( _ p t 2 + p k = 0 ( 1 2 1 ) 最+ ( u ( e + 力k = 0 ， 其中p 是密度，“是速度，e = i 1 艘2 + e 是总能量，e 是内能如果再加上一个 状态方程，那么整个方程组就是封闭的，并且在初始值和边界条件给定的情形下可 以求得唯一解 对于单介质的理想气体，我们有它的状态方程 p = ( ，y 一1 ) e ，( 1 2 2 ) 其中7 是气体的热容比 】 2 0 0 6 上海大学硕士学位论文 2 如果流体是有多种介质构成的，并且在开始状态各种介质是互相分开，没有发 生混合。那么由于它们没有粘性，以及互相绝热的特性，在以后的时间里他们将仍 然互相分开，不发生混合在不同的介质处我们将有不同的状态方程，例如所有的 理想气体都有( 1 2 2 ) 形式的状态方程，唯一的区别是不同的理想气体会有不同的 热容比 当然所有的流体介质均在介质界面处与别的介质相连接，而事实上这个 连接处也正是大多数应用同题主要关注的 然而在对多介质流的大多数有限体积法的数值模拟中，迄今为止所有的数值格 式都将引入一定量的数值耗散，不同介质的流体在介质边界处会发生一定量的数值 混合因此人们提出了一些改进过的模型替代( 1 2 1 ) 方程组来描述多介质流体， 以应付必然要发生的数值混合至今所有的改进模型都是基于混合区域达到热平衡 这一假设之上的，由于在单介质流体的时候，所有区域都处于热平衡状态，而且所 有的热力学量如压强、温度、内能熵和焓等等，只有在热平衡的时候才有定义， 医此这么做似乎无可厚非于是基于这个假设，混合区域的状态方程也能够被推导 出来如果所有介质均为理想气体，那么混合气体也是理想气体，因而有它自己的 热容比仉它是所有组成它的单种介质气体的热容比的某种意义上的加权平均在 本文中我们将考虑仅含两种理想气体的流体，而模型的多种流体介质可以简单地予 以推广由两种气体混合而成的气体热容比满足如下公式， 7 = 碟譬篙笋， z p l g ，l t + p 2 g 2 t ” 其中g l ，c a 分别是两种介质的等容热容比，第二章中我们将对它作具体推导 个常用的描述两种介质流体的改进模型叫做质量分配模型【1 】或热力学模型 c 3 3 l ，它的具体形式如下， n t + ( p l u ) 。= 0 彬+ ( 助“) z2 0 ( 1 2 4 ) ( 膨h + ( p u 2 + p ) 。= 0 皿+ m ( f + p ) b = 0 ， 其中m ，砌是两种介质各自的密度，他们满足p = p l + m ，t 是速度，e = ! p u 2 + e 是总能量，e 是内能而热容比，y 现在作为( 茁，t ) 的函数由( 1 2 3 ) 给出当然还有 其他很多改进模型，在【1 1 中有简单的介绍当气体没有粘性，不发生热传导的时 候。这些模型都相互等价，但是如果有粘性或者热传导的时候，即使很小，它们将 不再相互等价 2 0 0 6 上海大学硕士学位论文 3 个好的数值模拟应该是守恒形的。对于总质量以及各种介质的分质量、动量、 总能量均守恒但是所有那些对于单介质流体的守恒形格式推广到多介质流得时候 都产生了伪振荡这些伪振荡首先出现在混合区域流体的压强上，然后再影响到该 区域内所有其他物理量在第三章中我们将具体叙述这些振荡的确是由于数值耗散 所引起，特别是其中两种不同气体介质之闻伪热传导所造成的但是由( 1 2 1 ) e u l e r 方程组所描述的流体是无牯且绝热的，因此不会有任何混合发生，在介质交界面处 也就不应该有任何振荡即使是由n a v i e r - s t o k e s 方程组所描述的实际气体，虽然 不同气体介质问存在粘性和热传导，但是混合区域的宽度十分狭小，当粘性和热传 导很弱的时候压强在这个区域的振荡是不可见的问题是对于数值格式，粘性和热 传导都被放大到了网格宽度的程度，这远远大于实际，因而流体压强上的振荡就变 得可见了 为了解决这个问题，很多物理现象的实质在数值模拟的时候不得不被牺牲掉 近年来，个普遍的做法是牺牲解的守恒性，很多非守恒以及半守恒的方法被发展 出来，我们可以在【1 、 2 1 、f 3 】、【1 1 j 、【1 2 】，【13 】、【1 4 】【1 5 1 、【1 6 1 ( 2 6 】、 【2 8 】，【2 9 】、印l 、【3 2 】和 3 3 】等参考文献中找到这些方法但在很多问题中，解 的守恒性其实是不可替代的，他直接关系到间断位置是否正确，因面对于那些对间 断位置十分敏感的实际应用，现有的方法显然难以满足 1 3 本文给出的解决方法 我们从最基本的热力学和流体力学理论出发，发现了造成一般对于单介质流体 十分有效而推广到多介质流体却会在介质界面处产生严重振荡的原因，在于介质边 界所在的网格内发生了数值热传导，从而破坏了介质界面两边流体压强和速度的一 致性，这最终导致了解在介质界面处的振荡于是我们通过假设两种气体介质间互 相绝热，在介质边界处的网格不需达到整体热平衡，而只对每种介质达到局部的热 平衡，杜绝两种介质问的数值热传导这样做使得我们能够在保持格式的守恒性前 提下，消除数值结果在介质界面处伪振荡最终对于两种介质的流体我们提出了改 进过的基于非热平衡状态的l a x - f r i e d r i c h s 格式，它对总质量、各种介质的质量、 动量总能量均守恒同时在数值例子的测试中也取得了良好的计算结果 第二章混合气体的热平衡和热容比 在这一章中我们将讨论一些相关的热力学背景，了解气体在网格内是否热平衡 对计算的影响 2 1 热力学知识背景 为了更深入地理解混合流体问题，有必要先回顾一些相关热力学知识在热力 学中，我们的研究对象为闭合的热力学系统，它可以具备任意量的物质，并且性质 能唯一且完全地被一定的宏观参量所描述这些物质被物理墙壁与周围分开，系统 可以与周围交换能量，但不能和周围交换质量热力学中最重要也是最基本的概念 是热力学平衡，它是个系统经过足够长的时间后自动到达的一种状态，当系统达 到热力学平衡后，它的热力学特征，比如温度，压强、熵等描述系统状态的参数。将 不再随着时间而改变，同时整个系统的热力学状态量也只在平衡状态下才有定义 祖对地，个封闭系统在达到热力学平衡前所有的状态都是不平衡态，它的热力学 状态量无法定义，因为系统将随着时间不断改变从热力学平衡这个概念伴随而来 的个重要结论是，如果两个系统相互之间处于热平衡状态那么他们将有相同的温 度，这也是为什么温度计可以测量温度的原因一个系统若处于热力学平衡状态，那 么他的任意子系统都将处于热力学平衡状态而对于一个处于非平衡态的系统，它 有可能可以被分成若干个处于热平衡状态的子系统例如个温度计和被测量物体 组成的系统，在温度计停止变动前系统处于非平衡态，但是每一时刻温度计和被测 物体本身构成的子系统都处于热平衡态。他们有自已的温度从上面这个椤l 子中我 们也可发现，一个系统达到热平衡的主要手段是自身以及和周围环境进行热交换 更多的热力学概念可以在热力学教课书【1 0 1 中找到 2 2 混合气体得热容比 现在，处于平衡态的气体必须有相同的温度假设我们有一由两种介质构成的 混合气体，分别用下标1 和2 表示，并且有不同的等温比热和等压比热g ，g ，以 及不同的热容比，y 两种介质的密度分别为n ，脚假设混合气体处于热平衡并且 也是一种理想气体，然后我们计算混合气体新的热容比，y 首先，我们知道，当混合气体处于热平衡的时候，混合气体中的每一种介质作 4 2 d 0 6 上海大学硕士学位论文 5 为一个系统处于热平衡状态并且它们热力学性质由他们自身的状态方程所控制 因此。各种介质的分压强可以由( 1 2 2 ) 得到， 鲁铂，丧一 ( z 2 t ) 并且两种介质相同的温度可以分别用下面两个公式计算， 击= z 厩e 2 = z ( 2 2 2 ) 丽j 刮丽i 训 卫j 由d a l t o n 的分压定律，我们得到 p = p 1 + 砌，( 2 2 3 ) 同时由能量守恒定律，我们得到 e = e l + e 2 ( 2 2 4 ) 注意到混合气体处于热平衡状态并且被看作是一种理想气体，因此( 1 2 2 ) 式对它 同样成立将( 2 2 3 ) 和( 2 2 4 ) 代入其中，我们得到 _ = p l _ + 一p 2 = e 1 + 勖 (225)1一 -、。， 然后将( 2 2 ，1 ) 和( 2 2 2 ) 代入( 2 2 3 ) ，可得到 p l + 纯= 陬一1 ) m g ，1 t + 如一1 ) 纯岛，2 t ( 2 2 6 ) 由( 2 2 4 ) 可得到 e l + e , 2 一p l g ，1 t + 戊q ，2 t ( 2 2 7 ) 最后把( 2 2 4 ) ，( 2 2 5 ) 带入( 2 2 7 ) 得到 = 止等警等掣， 偿z 固 经过整理最终得到 ，y = - 7 1 p l 瓦c v d 万t + 鬲7 2 瓦p 2 c v r , 2 t ( 2 2 9 ) 7 p l g 1 7 + 见g ，2 r ” 这一混合气体的热容比的形式被很多计算流体力学的相关研究论文所引用，我们在 这里介绍它的推导只是为了保证我们叙述的完整性接下来我们将要看到。正是由 于使用这一热容比1 的形式，最终导致数值方法中在介质界面处压强的不一致。 第三章基于平衡态的l a x - f r i e d r i c h s 格式 l a x - f r i e d f i c h s 格式是个经典的守恒形格式，如【2 0 】和【2 1 】中所描述的那样， l a x - f r i e d r i c h s 格式可以被看成是g o d u n o v 格式的一种变型两种方法都通过将数 值解作为网格平均的逼近，然后解各个网格交界处的r i e m a n n 问题，并在下一时 刻对这些解取网格平均作为这一层上的数值解唯一的区别在于l a x - f r i e d r i c h s 格 式定义在偏移的网格上，而g o d u n o v 格式定义在普通的正交网格上精确地说， l a x - f r i e d r i c h s 格式在个时间层上定义在以为中心的网格上，然后下个时间 层上定义在以o ，+ l ，2 为中心的网格上其中巧= j h ，巧+ l 2 = 0 + 1 2 ) h 下面我 们将具体描述( 1 2 4 ) 方程组的l a x - f r i e d r i c h s 格式的物理解释这里气体在混合区 域被认为是热平衡的，并且使用( 1 2 2 ) 所描述的状态方程，热容比，y 由( 1 2 3 ) 所 决定我们还将分析通过介质界面的时候压强的一致性是怎样被数值耗散，特别是 两种流体介质间的热传导所破坏的 3 1 格式的推导 假设第n 层上的数值解定义在以巧为中一5 - 的的网格内，作为对真解的网格平 均逼近， w 7 = ( ( p ，垮，( 您骖，( p 垮，( 曰) ) 7 ( 3 1 1 ) 然后我们的算法按下述重构- 演化一取平均的过程进行t 1 重构一个关于x 的分段常数函数 r ( 叫w “) = w ? ， 一l 肛 2 0 0 6 上海大学硕士学位论文 9 0 如图( 3 2 1 ) 所示，开始它位于网格边界+ l ，2 处，在接触面两侧，流体是 两个常状态， 叼= 豫p i , 舢o , p l 椰f t , e 矧| ) t , 3 勺s3 0j o ， ( 3 2 9 ) 其中 蜀= 釉铲+ 南，历= 1 ( p , ) a 2 + 占， ( 3 2 1 0 ) 声为接触面两边共同的压强，而，y 1 和 2 是两种气体各自的热容比在物理上就相 当于取出+ l ，2 处的小隔膜，使得两边的气体互相接触长达 t 的一段时间，接 着在，+ l 处放上小隔膜阻止流体流动，使得气体在偏斜的网格里力学意义上 充分混合，并热力学意义上达到热平衡 因为接触面以速度面运动，可得到 以= 言( 1 + a 霞) p z ，见= 言( 1 一h a ) p , ( 3 2 1 1 ) 以及 p u = 言( 1 + a 豇) 向面+ 言( 1 一施) p r 面 ( 3 2 1 2 ) 注意到户= , o l + 忽，从( 3 2 1 1 ) ，( 3 2 1 2 ) 可得到，当气体在偏移网格内充分混合后有 u = t ( 3 2 1 3 ) 这意味着速度在计算中保持一致 然而，压强却无法保持这性质下面由能量入手分析流体的压强由能( 3 1 8 ) 中最后一式得到 e = 言( 1 + a 面) 局+ 丢( 1 一她) 耳 ( 3 2 1 4 ) 注意到( 3 2 1 0 ) ，( 3 2 1 2 ) ，( 3 2 1 3 ) ，( 3 2 1 4 ) 我们得到h l 层上的内能 e = ；( 糟+ 再1 - a f t j p _( 3 。舫) 8 2 j ( 百= t + ；i 渺 【3 2 1 5 j 现在在这个偏斜的网格内混合气体在t 。+ 1 时刻达到热平衡它拥有( 1 2 2 ) 式的状 态方程，并且7 可以像( 2 。1 9 ) 中那样计算最终我们得到 p鬻筹杀磐(17。+一aa2 1 1 a a ) c , , 2 p ,1 + 糟1 涫 ( 3 z 加) ( + a 面) q l n + ( 一，y 1 一 复一 7 ， 、 在上式中，我们可以看到，p ，万几乎不可能相等，除非是在单介质流的时候饥= 蚀， 这时压强才能保持一致 2 0 0 6 上海大学硕士学位论文 1 0 我们从上面的讨论中可以看到在通过介质界面时压强的一致性被混合气体在混 合区域达到热平衡这一假设所破坏这也最终导致了数值解在介质界面附近的剧烈 振荡因为混合气体的热容比1 与两种气体介质的热容比1 l ，y 2 均不相同，所以由 ( 3 2 1 6 ) 所计算的压强在这个混合区域也和两边的压强均不相同由于( 3 2 1 6 ) 中 压强的振荡在单介质流的情况时并不存在，我们可以断定这个振荡的问题正是由于 两种介质问的热传导所造成的 在以e u l e r 方程组所刻划的流体内，嵇体被认为是无粘性且绝热的，因此在介 质边界附件并不发生混合即使在由n a v i e r - s t o k e s 方程所刻划的实际气体中，虽 然由于粘性和热传导，存在一个气体混合并局部达到热平衡的区域，但在r e y n o l d s 数很大，热传导很弱或气体几乎绝热的情形，这个混合区域将非常之小因此在这 些流体中压强的振荡应该是肉眼所不可见的然而在l a x - f r i e d r i c h 格式的重构步 骤中，在每个小网格内气体充分混合且达到热力学平衡，这使得气体的粘性和热传 导都被扩大到网格长度o 的程度，从而使得压强的振荡变得可见 第四章基于非平衡态的l a x - f r i e d r i c h s 格式 在这一章中将介绍我们所设计的用于计算双介质流基于非平衡态的l a x - l i i e d r i c h a 格式就像我们在绪论中所说的那样，我们将放弃混合气体处于热平衡这一假设， 而不是牺牲解的守恒性，并且这样设计的守恒形格式同样能消除压强上伪振荡 4 1 方法的前提假设 在前一章中，我们看到压强的振荡是由于数值耗散所引起的，最主要的就是其 中两种介质问伪热传导为了保证格式的稳定性我们不能完全消除格式的数值粘 性同样的原因，我们也不能取消单种介质内的数值热传导但是我们能够通过使 两种气体介质之间绝热来阻止两种介质问的数值热传导。现在，两种气体介质在热 力学意义上被分隔开来，从而它们之问没有任何热交换这样在有混合发生的网格 内的气体不再处于平衡态，为此这我们将引入如下的假设， 1 ) 两种气体介质在力学意义上充分混合总密度等于两种介质分密度之和两 种气体的速度保持一致 2 ) 每种单一的介质自身达到热平衡从而仅仅相关于该种气体介质的压强，和 内能具有( 1 2 2 ) 形式的状态方程 3 ) 两种气体介质之间保持绝热，因而它们之间可能没有相同的温度，更没有一 个类似平均的整体热容比1 存在 4 ) 气体压强仍旧满足d a l t o n 分压定律，即混合气体总压强等于所有介质分压 强之和 5 ) 由于压强而对气体所做的功按照气体中每种介质的质量比例分别施加到每 一气体介质上 下面一节中我们将根据这些前提假设推导出我们基于非平衡态的l a x - f r i e d i c h s 格式 4 2 格式的推导 首先，由假设3 ) 我们知道两种介质的内能不像在第三章中那样由一个相同的 1 1 2 0 0 6 上海大学硕士学位论文 7 = ( ( p - ) ；，( p 2 ) ? ，伽) ；，) ；，池) ；) r ， ( 4 2 1 ) ( b ) ；= ；( m ) ? ( 嵋) 2 + ( e 浮，( 易) ；= ；( 砌) ；( 田) 2 + 池) ；， ( 4 2 2 ) h 护箬 ( e 。肛等 ( 4 2 3 ) 蜞= ( p t ) n + l 2 + 她) 渤： 粥2 2 = 五1 吩，n + 略- 嗡- ) 一( 晴- ( 嗡，) 2 + 臻一露( 婶) 2 一赡) 黧。嘶鼢 h 。q ( 蜥班= 丽意转r 一叫 0 1 ) 荔2 = ( 饥一1 ) ( ( 蜀) 并2 一孙1 ，1 ，j n + + l l 2 ( t t 茹马2 ) 2 ) ( 见) ，n + + l l 2 = ( 1 1 ) ( ( 晚) j n + + l l 2 一 ( 以) j n + + l l 2 ( 譬筹办2 ) 2 ) ( 4 2 7 ) 。 p 筋2 = ) ，n + + l l 2 + 慨) 撑 2 2 0 0 6 上海大学硕士学位论文 须要注意的几点是，在我们的计算格式中两种气体介质互相绝热，因而我们也 将不用涉及到取两种气体介质热容比1 l ，他某中间值的混合气体热容比7 在下一 节的理论分析中我们将看到正是由于这特性使得我们的格式在介质界面附近不会 产生任何伪振荡另外( 4 2 5 ) 中最后两个方程并不是以守恒形式出现的，但是这 并不说明我们的格式是非守恒的，将它们相加我们会得到总能量仍旧是守恒的，而 事实上在上一章所描述的最初始的l a x - n i e d r i c h s 格式中两种气体介质的各自能量 也是不守恒的，它们通过各种途径与另一种介质交换能量，而我们所做的就是把这 种交换限制在做功上 4 3 对新l a x - f r i e d r i c h s 格式的理论分析 我们最主要的理论成果归结为下面这一定理 定理4 1 所设计的l a x - f r i e d r i c h s 格式在满足c f l 条件的前提下，有如下3 个 性质， 格式在计算过程中，总质量和各介质的质量，动量，以及总能量均守恒 纠格式在计算过程中，每种介质的密度始终非负 别格式在计算过程中，速度和压强在介质界面两边僻合区域，保持一致 证明，1 ) 由于对于总质量、分质量、动力以及总能量的计算格式就是对守恒形 方程的离散化，因此数值解在c f l 条件得到满足的情况下，将仍旧满足守恒律 2 ) 我们注意到最初的两个方程可以写成 ( n ) j n + + 1 l 2 = ；( 1 + a 哆) ( 风圬4 - ；( 1 一a t 鼻1 ) ( 胁) 务l ， ( 4 3 8 ) 其中 = 1 ，2 ，a = 基是网格的步长比由c f l - 条件- a m 衅 瑶+ l “2 ) 0 3 ) 为了证明第三个结论，我们同样考虑一个纯粹的介质界面问题，并且流体在 介质界面处的传输速度为西 0 但这一次我们需要考虑两种介质混合的一般情 况我们观察在t 。时刻两个以奶。和x j o + 1 为中心的相邻网格点我们假设在这两 个网格里速度和压强相同那么在这两个网格点里相应的数值解为 叼= 戮p z , , p t , 2 , p t u 呱, e t , 1 。, 嚣，篙+ 1 ， 性s 加， i( l ，n 口，p 面，晶i ，蜀j 2 ) r ，j ；如+ ， 。 2 0 0 6 上海大学硕士学位论文1 4 i 算 i p 。r ，, l 嵋。i ；e f ：1 ，e i ；e r 。，i + 膳和 l 鼍芝；墓笔l + 艘和 t i m e = t r h l 置4 3 1 个介质界面以建度霞 0 从左向右传播。相应地在它两边流体的状态是抽，l ，九，2 ，m h 局1 ，固，2 ) r 以及( 九l ，p r ，2 ，l r b ；l ，日，2 尸，其中帆声h ，l + 甩，和嘶= ( m l + p 2 ) 面分别是左右两边流体地动量 交畀面的位置在+ i z 越k 并且在个时坷同隔内向右移动一= 丢埔 其中 n = n ，l4 - n ，2 ，店= 1 + 2 蜀，1 = 内，1 铲4 - e t ，l ，马，2 = 面2 + e f ，2 ( 4 3 1 1 ) 蜀j 1 = 互1 ，k l u - 2 + e l ，b ，2 = 互1 - 2 + e r ，2 从上面两个式子中我们不难看出，在左右两个网格里可能都含有两种气体介质由 假设2 ) 和4 ) ，并且在这两个网格中压强是相同的，我们得到t ，= 器，劬= 嚣， ， p 。a 1 ：? 。珥2 ( 4 3 1 2 ) 2 枯，铆= 昔， 睁。纠 乒= p r , l + 鼽2 ， 其中声是两个网格内流体的共同压强，我们需要证明l a x - f r i e d r i c h s 格式( 4 2 5 ) 能 够维持流体速度和压强不变为了这个目的，方便起见，我们记第如+ - 层上的数值 2 0 0 6 上海大学硕士学位论文 解为， 瞩n + + l l ，2 = h ，以， ，局，岛) 7 ( 4 3 1 3 ) 我们仍旧从流体速度开始研究从l a x - f r i e d r i e h e s 格式( 4 2 5 ) 我们知道两种 介质的密度由下面两式计算 n 5 黔翌 - + 黔攀1 ，。 ( 4 。从) 砌= ；( 1 一a 豆) 肛j 2 + ( 1 + 面) 胁，2 、 同样的我们有流体动量计算式子。 o u = 1 ( 1 一a 缸) 一面+ ；( 1 + a 面) 见霞 ( 4 3 1 5 ) 综合上面两个方程，我们得到，第t 。+ l 层上的流体速度为 缸= 聚端等等器嵩黔掣矾( 4 3 1 6 ， 可见在如+ 1 层上流体速度保持一致 接下去我们再来研究压强从l a x - p d e d f i c h 格式中我们得到两种介质的各自 的总能量能够如下计算， 蜀2 ；( 1 一舾) 毋，l + ( 1 + 1 面) 目1 ( 4 3 1 7 ) 马= ；( 1 一她) 蜀- 2 + ( 1 + a 面) 蜀t 2 、。 因而第一种介质的内能是s e - = 蜀一互1 p - 铲= ；( 1 一a 面) ( 矗z 一 l 霞2 ) + ；( 1 + a 霞) ( 蜀，。一庙，。霞2 ) ， ( 4 3 1 8 ) 或者由( 4 2 9 ) ，( 4 2 1 2 ) ，( 4 2 1 5 ) ，亦可写为 e ，= ；( 1 一a 面) e 。，+ ；( 1 + a 露) e ( 4 3 1 9 ) 同理可得第二种流体介质的内能是 、e 。= 殴一；见舻= ；( 1 一a 霞) ( 民：一“，。豆2 ) + ；( 1 + a 面) ( 局，z 一印，。面2 ) ，( 4 3 2 0 ) 或者由( 4 2 9 ) ，( 4 2 1 2 ) ，( 4 2 i 5 ) ，亦可写为 e 2 = 互1 ( 1 一旭) e 。+ ；( 1 + 旭) 日，。 ( 4 3 2 1 ) 2 0 0 6 上海大学硕士学位论文 1 6 再由两种介质的状态方程，我们得到 p l = h 1 ) e - = ；( 1 - a 面h l + ；( 1 + a 咖，1 ( 4r 3 2 2 ) 和 现兰( 饥一1 ) e 2 = ；( 1 - a 面) 以2 + ；( 1 + a 百) p l f 2 ( 4 3 2 3 ) 最终由( 4 2 1 0 ) 及d a l t o n 定律我们得到， p = p l + p 2 = 多， ( 4 3 2 4 ) 即流体压强在厶+ 1 层上也保持一致，整个定理证明完毕 注4 1 虽然这里推导的l a x - f f i e d f i c h s 格式仅具有一阶精度，但高阶精度的算 法可以用f 2 3 】和 2 7 1 中的方法获得 注4 2 虽然在这里我们仅仅考虑了两种理想气体的情况，但这种使得各种流 体介质之间绝热的思想肯定可以推广至更为一般的情况中去并且即使当实际应用 中的气体状态方程不具有理想气体的模式，只要对单一流l a x - f r i e d r i c h s 格式能够 维持速度和压强通过切向间断后保持一致，那么对于多介质流体将仍然保持一致 4 4 新l a x - f r i e d r i c h s 格式所描述的数学模型 在这一章的最后我们给出一个描绘多介质流的数学模型，它的离散形式恰恰就 是我们新导出的l a x - f r i e d r i c t m 格式 反+ 0 q b = 0 ( p 1 ) f + ( p l u ) ：= 0 熬：2 p u 幺嚣：0 t 筋， ( 肚x + ( 2 + p k = 、 。 ( 蜀h + ( 易b + m ( 叼k = 0 ( 易k + ( u e z k + 蚝( t 纠。= 0 ， 其中 1 毋= ；( p 1 ) “2 + e 1 ，骂2 = ；( p 2 ) “2 + e 2 ( 4 4 2 6 ) 和 m = 忐，y 2 = a p l + p 2 ( 4 4 2 7 ) p l 十m 2 0 0 6 上海大学硕士学位论文 1 7 是两种介质之间的质量比，它们满足归一化条件m + k 兰1 为了使得整个方程 组有唯解，我们需要一个类似状态方程的式子，它可以这样定义。 p = 击+ 者 ( a 嬲) p 2 弓可+ 弓可。 ( 4 4 。2 8 ) 显然，当流体内只有一种介质的时候，这个数学模型就自动退化为e u l e r 方程组 ( 1 2 1 ) 第五章数值例子 这一章中我们将用几个具体的数值例子来说明新的计算格式在模拟两种介质流 体运动时的有效性在下面这些例子中点代表我们改进的l a x - f r i e d i r i c h s 格式的计 算结果，而连续曲线代表某种意义上的精确解由于l a x - y 、 i e d r i c h s 格式比g o d n o v 格式有更多的数值粘性，相对而畜可能需要更多网格点来使得计算足够精确 5 1 基本的切向间断 这是个最基本的切向问断问题，但能很好地对比出我们方法的好坏它的初 值如下 0 $ 0 2 0 2 z 1 0 以= 1 0 ，冉= o 1 2 5 铆= 1 0 ，t l r = 1 0( 5 1 1 ) a ；1 0 ，肼= 1 0 m = 1 4 ，佛= 1 2 边界条件设置为无限远边界条件，即中间流体的变化不会影响到边界，边界点附近 流体的状态始终保持不变 图( 5 1 1 ) 是原基于平衡态的l a x - f r i d r i c h s 格式在t = 0 1 2 时的计算的结果， 取2 0 0 个网格点，c f l 数取0 5 左上圈为密度，右上图为速度，下图为压强图 中虚线为数值解，直线是精确解( 这个精确解是显然的，可以直接写出) 可以看到 在速度和压强上，计算结果出现了非常严重的振荡，密度上也有一些小的振荡 1 8 2 0 0 6 上海大学硕士学位论文 j ；j 。 | 圉5 1 1l = 0 1 2 时切向间断同题基于平衡杏l a x - f r i e d r l d l 格式的教值解取2 0 0 十同格点。c f o 致取0 & 左上圉代表密度右上圈代表建度，中下图代表压强 而图( 5 1 2 ) 是新设计的基于非平衡态的l a x - f r i d r i c h s 格式在t = 0 1 2 时的计 算结果，同样取2 0 0 个网格点，c f l 数取0 5 左上图为密度，右上图为速度，下 图为压强图中虚线为数值解，直线是精确解可以看速度和压强完全没有任何多 余振荡，与精确解完全吻合，在密度上也没有任何振荡，只含有一些由于低阶格式 而造成的磨损 2 0 0 6 上海大学硕士学位论文 圈5 , 1 ，2t ；0 1 2 时切向间断同盟基于非平衡杏k n i e d r 妯目格式的数值解，取2 0 0 个同格点c f l - 数取o 5 左上图代表密度。右上图代表速度中下圈代表压强 5 2s o d 激波管问题 我们考虑初始间断在茹= 0 5 处的r i e m a n n 问题，并且有如下的初值 0 z 0 5 0 5 z 1 0 庙= 1 0 ，办= o 1 2 5 铆= 0 0 ，t r = 0 0( 5 2 2 ) 鼽= 1 0 ，辫= 0 1 m = 1 4 ，* = 1 , 2 边界条件设置为无限远边界条件，中间位置产生的激波，稀疏波以及切向间断，在 计算过程中不会影响到边界 这就是著名的s o d 激波管问题，很多相关的流体计算论文中都用它来测试自己 的计算格式我们可以在第一章中所介绍的大多数参考文献中发现它被用作测试例 子 2 0 0 6 上海大学硕士学位论文 2 1 图( 5 2 1 ) 是原先基于平衡态的l a x - f r i e d r i c h s 格式在t = o 2 5 处的计算结果， 左上图代表流体密度，右上图代表流体速度，下图代表流体压强其中虚线是基于非 平衡态l a x - f f i e d r i c h s 格式的计算结果，使用了2 0 0 个网格节点，c f l 数取0 5 而直线是用迭代法求得的r i e m a n n 问题的真解可以看到速度和密度在介质边界 处都有明显的振荡 囝5 , 1 2t 暑o 2 5 时s o d # 獭波管问魍基于平衡态l a x - f r i e d r l c h 格式的散值解，取2 0 0 个丹格点，c f l - 数取 0 5 左上图代表密度。右上阿代表逮度。中下图代表压强 2 0 0 6 上海大学硕士学位论文 图( 5 2 2 ) 是基于非平衡态的l a x - f r i e d r i c h s 格式在t = 0 2 5 时的计算结果，可 以看到除了在间断，特别是两种流体接触面处，由于l a x - h i e d r i c h s 格式的高耗散 所造成的磨损比较严重以外，在流体接触面两边，流体的速度，压强很好得保持一 致，并且格式具备良好得守恒性，无论是总密度，或者每种流体介质得分密度，以 及压强都始终保持非负 田5 , 2 ，2t = o 2 5 时s o d 5 激被管问题基于非平衡态l a x - f r i e d r i d 暗格式的散值解。取o 个同格点，c f l - 效 取0 5 左上圈代表密度右上图代表速度，中下圉代表压强 为了更好的展示我们的格式在处理强间断以及他们的相互作用的能力。我们用 2 0 0 6 上海大学硕士学位论文 下面的两种介质流体的b l a s t - w a v e 问题来测试我们的格式它有如下的初始值t 0 z s 0 1 ， n = 1 0 ， 蛳= 0 0 ， a = 1 0 0 0 0 ， m = 1 4 ， 0 1 霉0 9 ，0 9 1 0 加= 1 0 ，办= 1 0 u r n = 0 0 ，撕= 0 0( 5 3 3 ) 舫= 0 0 1 ，西= 1 0 0 0 =