基于有限体积HWENO格式的一维溃坝流模拟

文章编号,&:8A68%8%!!%&%"%#6%A886%A
基于有限体积Z Y<);格式的一维溃坝流模拟
翁磊&"程国胜&"卢长娜&
摘要
基于浅水波方程组建立一维溃坝流模型&并给出数值模拟结果'其中&空间离散采用E GL M N!E4*V13G41[H34+L26 24,31S T T Y M),6N251T T S3)*Y"格式&时间离散采用四步K U B!K)3S T U S*1S31),B1V1,12H6 1,["@0,[46n033S方法&模拟堤坝溃决时洪水演进过程'模拟结果表明,较采用GL M N格式所得数值解更精确-同时&相比GL M N格式的相应算法&该算法解决一维溃坝流问题能更有效地减弱振荡&对间断具有更高的分辨率'
关键词
溃坝流-有限体积法-E GL M N格式中图分类号N&8#e!
文献标志码>
收稿日期!%&%6%:6%!
资助项目国家青年基金!A%$%:%A7"
作者简介
""翁磊&男&硕士生&主要研究水动力数值模拟'T t(4,[/&:9'5)V'=>引言
,3*)+0531),
""拦河筑坝&既可防洪又可用于灌溉%发电&对社会发展起积极作用&但坝体一旦溃决&形成的洪水波会对下游地区造成严重灾害'洪水波的速度影响着下游的预警和紧急疏散时间&决定受灾范围和程度&因此&溃坝流数值模拟具有非常重要的研究意义'
溃坝流常以浅水波方程组作为模型'溃坝流具有间断性&对浅水波方程组离散时&要求数值方法具有稳定%高分辨率%能较好抑制振荡等特性'近年来&用高分辨率数值格式建立一维溃坝流模型已有一些研究'如,O H)1等*&+基于迎风K U B&迎风K U D和L M N格式讨论了一维溃坝问题并指出K U D格式精度较高&但O<`时间消耗较多-a1,等*!+提出了一种二阶混合K U B有限体积格式&由于K U B格式本身固有的限制&精度至多只能达到二阶&极值点附近降至一阶*9+-鲍远林等*A+基于D)T3Q V S,,方程&采用有限体积n;U C方法&成功模拟了下游为干河道的一维溃坝流-^S V.*4331等*#+采用^5O)*V S5\6I S V42),激波捕捉法成功建立了一维带污染物两相流溃坝模型-文献*:+应用G L6 M N格式模拟了闸门瞬间开启后的水流&水位%流速与解析解拟合较好-文献*8+建立了基于G L M N格式的浅水间断流数值模型&得到了有五阶精度的一维溃坝模拟算法'文献*7+将一种基于E G L M N!E4*6 V13G41[H34+L224,31S T T Y M),6N251T T S3)*Y"格式的限制器引入双曲守恒律方程的求解中&E G L M N格式在重构边界点值时&选用E4*V13插值多项式&其优点是选用较少的模板就能达到较高的精度&同时减少计算量&边界处理也变得简洁'本文将E G L M N格式应用于一维溃坝流进行数值模拟&所得结果与基于G L M N格式的数值结果*8+及精确解*$+比较表明&本文建立的溃坝流模型&格式稳定&分辨率高&能更好地处理类似溃坝这类浅水间断流问题'
>浅水波方程组及时空离散
CHS T T)((S34*(S X44c0S31),2S,+2]S54631V4+125*431Q S31),
K?>有限体积法积分方程
在笛卡尔直角坐标系下&一维浅水波方程组守恒形式如下,
55
其中&2是守恒型向量&7是通量向量&表达式为
2P &
()
5&&"7P 5&5!
&R 1&!
!T 对计算域进行等步长网格划分&记第*个单元
/*[P I *Q &!&I *R ]&!
&单元步长+
I P I *R &
!
Q I *Q &!
&单元中点I *P
(&
!
I
*R &
!
R I *Q )&!
T
对式!&"关于I 求导&得方程,55M 8R 5
5I
&P %T !!"
其中&8P &I
!5&"(
)
I &&P !5&"I
5!
&R 1&!
()
!
I
T 记单元平均值,*P &
+
I -/*
2!I &M "+I &Y *
P &
+
I -
/*8!I &M "+I T 在单元/*
上&对式!&"%!!"关于I 在第*
个单元上积分平均&得+,*+M PQ &+(I 7!2!I *R &!&M ""Q 7!2!I *Q &!&M )""&+Y *
+M PQ &+
(I &!2!I *R &!&M "&8!I *R &!&M ""Q """7!2!I *Q &!&M "&8!I *Q &!&M )""
T 对通量7
!2"&&!2&8"取近似,(72!I *Q &!&M )"6C (y ,T *Q &!&,**Q )&!
-(&2!I *Q &!&M "&8!I *Q &!&M )"63(y ,"T *Q &!&,"**Q &!-Y "T *Q &!&Y "**Q )&!
T
其中&,T *Q &!&,*
*Q &!分别表示2在点!I *Q &!
&M "的左右函数值-Y T
*Q &!&Y *
*Q &!
分别表示8在点!I *Q &
!
&M "的左右函数值T
本文采用a )5S T a S -6;*14+*15H !a a ;"数值通量*7+
7!,*
&,T
"P
(&
!
7
!,*"R 7!,T
"Q !!,*
Q,T
)"&
&!,"*&Y "*-,"T &Y "T
"P
(&
!
&
!,"*&Y "*
"R &!,"T &Y "T
"Q !!Y "*
Q Y "T
)"T 其中&!是C 的I
S 5).1S ,矩阵线性化替换矩阵的谱半径T
K @>Z Y<);重构
为了得到五阶精度&G L M N 格式需要覆盖#个单元&而采用E G L M N 格式&只需要9个单元&具体重构过程如下,
&e !e &"重构,*
*h &!
{分别在模板
/*Q &
&/
*&/*R }&
{&/*Q &
&/}*
{&/*&/*R }&
!I "T 为了使格式达到五阶&{在/*Q &&/*&/*R }&上还需构造一个E 4*V 134四次插值多项式G A !I
"T 重构的多项式分别满足如下条件,
-/*Q &
G L
!I
"+I +I *
P +,*Q &&"L
P &&!&A --/*Q &
G ]
L
!I "+I +I *
P /Y *Q &&"L
P !&A --/*R &
G L
!I
"+I +I *
P +,*R &&"L
P &&9&A --/*R &
G ]
L
!I "+I +I *
P /Y *R &&"L
P 9&A --/*
G L
!I "+I +I *
P,*
&"L P &&!&9&A T ""计算得G L !I *Q &!"&L
P &&!&9&A 的表达式,G &!I *Q &!"P
&9+,*Q &R #:+,*Q &:+
,*R &-G !!I *Q &!"P
#:+,*Q &R &:+,*R +I *
9/Y *
Q &-G 9!I *Q &!"P &9:+,*Q 8:+,*R &R !+I *9/Y *
R &-G A !I *Q &!"P :8&!%/,*Q &R &$9%/,*Q !9&!%/,*R &R 8+I *A %/Y *Q &R 9+I *A %/Y *
R &T 由G A !I *Q &!"P "9
L P &4L G L !I *Q &!"可得线性权,4&P !$7%&4!P !&A%&49P $
7%
T 为减弱振荡&使用非线性权,3L P 73L
"K
73K &其中&73L P 4L !)R "L "&"L P "^-/*
+I !^Q &
5^5I ^
G L !I ()
"!+I &L P &&!&9T )
是为了防止分母为零而引入的微小量&本文取)P &%Q :
-"L 称为光滑因子&取值越小&表示相应
模板越光滑T
重构后,
*
*Q &!
的表达式是,
*
*Q &!P "9
L P &
3L G L !I *Q &!"T
&e !e !"重构Y *
*h &!
{分别在模板
/*Q &
&/
*&/*R }&
{&/*Q &
&/}*
{&/*&/*R }&
上构造E 4*V 134二次重构多项式G &!I "&G !!I
"&G 9!I "T 为使格式达到五阶&{仍需在/*Q &&/*&/*R }&上构造一个E 4*V 134五次插值多项式G A !I
"T 重构的7
8A 翁磊&等'基于有限体积E
GL M N 格式的一维溃坝流模拟'G L M =a 41&43S T 'C1V 0T S 31),)Z ),46+1V 4,21),S T +S V.*4S \.S 24+),Z 1,134X )T 0V 4E G L M N25H4V 4'
-/*Q &
G L
!I
"+I +I *
P +,*Q &&"L
P &&!&A --/*Q &
G ]
L
!I "+I +I *
P /Y *Q &&"L
P &&!&A --/*R &
G L
!I
"+I +I *
P +,*R &&"L
P &&9&A --/*R &
G ]
L
!I "+I +I *
P /Y *R &&"L
P &&9&A --/*
G L
!I "+I +I *
P +,*
&"L P &&!&9&A T 计算得G ]L !I *Q &!"&L P &&!&9&A 的表达式,G ]&!I *Q &!"PQ 9A +I *+,*Q &R &+I *+,*Q &A +I *+,*R &R &!/Y *-G ]!!I *Q &!"PQ !+I *+,*Q &R !+I *
+,*Q &!/Y *Q &Q &!/Y *-G ]9!I *Q &!"P
A +I *+,*Q A +I *
+,*R &R 8!/Y *R 9!/Y *R &-G ]A !I *Q &!"PQ !&&!+I *+,*Q &R !+I *+,*Q &A +I *
+,*R &Q "#&!
/
Y *Q &Q &:/Y *R &&!/Y *R &T 同上&计算可得线性权,4]&P &$&4]!P #:&4]9P &&7
T 此处仍然使用非线性权&构造方式与前文类似T 重构
后的最终表达式为Y
*
*Q &!P "9
L P &
L 3]L G ]L !I *Q &!"T
,T
*h &!
&
Y T
*h &!类似重构可得T
K B>时间离散
本文采用四步K U B @0,[46n 033S 时间离散格
式*8+&记,,
!&"P ,0
R +M !
+!,0"-,!!"P ,0
R +M !
+!,!&""-,!9"P,0R +M +!,!!""-,0R &P &9!,0
R
,!&"R !,
!!"
R ,!9"
"R +M :
+!,!9""T 其中+!,
"表示空间离散化算子T
@>数值实验
M 0V 4*15S T 3423
""用E GL M N 格式模拟溃坝流'在平底光滑
矩形断面的明渠上&堤坝溃决&如图&所示'渠
&V &Q &'I '%-%e &V &%'I '&{
&
将计算域*h
&&&+分成!%%等分&结束时间%e
!2
'图&"溃坝示意
;1['&"B 1S [*S V)Z +S V.*4S \
@K ?>数值试验方法
本文用;)*3*S ,语言对模型进行编程实现&具体计算流程如下,
&"对计算域等步长剖分后&记录单元中点%单
元步长等信息-!"根据初始状态&给定变量初值&设定计算终止时间H
&时间步长+M -9"在每个时间步长内&采用E G L M N 格式由+,0
*
重构单元边界点的左右函数值&计算边界点的数值通量&根据四步@0,[46n 033S
时间离散公式计算下一时间点的变量值+,0F &
*-A "令M PM R +M &若M ZH &则将+,0R &
*的值赋给+,0
*
&返回到9"&继续执行程序&若M XH &则结束循环&输出结果'@K @>数值结果及分析
大坝溃决后&下游水位陡涨&流速骤增&形成溃坝波&向下游传播'%e !2时水位和流量分布如图!和图9所示&图中还给出了基于G L M N 格式的数值解*8+及解析解*$+'图!%图9表明,本文结果较基于G L
M N 格式的数值解更好地拟合精确解-本文所用格式较好地抑制了数值耗散&具有高阶%高分辨率'
B>结论
O ),5T 021),
""本文基于E
G L M N 有限体积格式及四步K U B @0,[46n 033S 法数值离散一维浅水波方程组&建立了
一维溃坝流数学模型'
结果表明&该模型能准确捕捉间断&抑制间断处产生的振荡&具有高阶%高分辨率&
$
8A 学报,自然科学版&!%&%&!!#",A886
A7%I )0*,S T )Z M S ,b 1,[`,1X 4*213Y )Z ?,Z )*V S 31),C514,54S ,+K 45H,)T )[Y ,M S 30*S T C514,54L +131),&!%&%&!!#",A886A7%
图!"水位对比
;1['!"O )V ]S *12),)Z 20*Z S 54T 4X 4T
2
图9"流量对比
;1['9"O )V ]S *12),)Z +125HS *[42
参考文献
@4Z 4*4,542
*&+"O H)1C &<S 1\I '<4*Z )*V S ,543423)Z H1[H *42)T 031),25H4V 42Z )*&B
+S V.*4S \]*).T 4V *I +'n CO LI )0*,S T )Z O 1X 1T L ,[1,44*1,[&!%%&&#!9",!896!7%
*!+"a 1,=;&a S 1I C &=0)G B '<4*Z )*V S ,54)Z H1[H *42)T 031),K U B
25H4V 42Z )*&B+S V.*4S \21V 0T S 31),2
*I +'I )0*,S T )Z 3H4O H1,424,2313034)Z L ,[1,44*2
&!%%#&!7!#",88&687!*9+"刘儒勋&舒其望'计算流体力学的若干新方法*^+'北京,科学
出版社&!%%9a ?`@0-0,&CE `W 1(S ,['C4X 4*S T ,4(V 43H)+2Z )*5)V ]03S 31),S T
Z T 01++Y ,S V 152
*^+'D 41b 1,[,C514,54<*422&!%%9*A +"鲍远林&周晓阳'移动边界的有限体积n ;U C 方法在一维溃坝
波计算中的应用*I +'水利学报&!%%#&9:!&!",&A8%6&A8A D >Nm 0S ,T 1,&l
E N `J 1S )Y S ,['>]]T 15S 31),)Z \1,4315Z T 0-6X 453)*2]T 1331,[V 43H)+1,5S T 50T S 31),)Z &B+S V 6.*4S \(S X 4*I +'I )0*,S T )Z E Y +*S 0T 15L ,[1,44*1,[&!%%#&9:!&!",&A8%6&A8A
*#+"^S V .*4331C &a S *5S ,L &G *S 5H14,BB '&BV )+4T T 1,[)Z +S V 6.*4S \
20*[42(13H Z T )S 31,[+4.*12*I +'D 1)2Y 234V 2L ,[1,44*1,[&!%%7&&%%!!",!$869%7
*:+"刘玉玲&王玲玲&周孝德&等'闸门瞬间开启水流的G L M N 格
式数值计算*I +'武汉大学学报,工学版&!%%$&A!!9",!7&6!79
a ?`m 0T 1,[
&G >M =a 1,[T 1,[&l E N `J 1S )+4&43S T 'M 0V 4*15S T 5S T 650T S 31),Z )*&BZ T )(2.Y *S ]1+)]4,1,[)Z S 2T 0154[S 34(13H G L M N
25H4V 4
*I +'L ,[1,44*1,[I )0*,S T )Z G 0HS ,`,1X 4*213Y &!%%$&A!!9",!7&6!79
*8+"卢长娜'浅水间断流的高精度G L M N 有限体积数值模型研究
*B +'南京,河海大学数理学院&!%%7a `O HS ,[,S '@424S *5H ),H1[H6)*+4*G L M NZ 1,134X )T 0V 425H4V 42
Z )*+125),31,0)022HS T T )((S 34*Z T )(),0,23*0530*4+V 42H2
*B +'M S ,b 1,[
,C5H))T )Z ^S 3H4V S 3152S ,+<HY 2152&E )HS 1`,1X 4*213Y &!%%7*7+"W 10I J &CH0OG'E 4*V 134G L M N25H4V 42S ,+3H41*S ]]T 15S 31),
S 2T 1V 134*2Z )*@0,[46n 033S+125),31,0)02=S T 4*\1,V 43H)+,N
,4+1V 4,21),S T 5S 24*I +'I )0*,S T )Z O )V ]03S 31),S T <HY 2152&!%%A &&$9!&",&&#6&9#
*$+"C3)\4*I I 'G S 34*(S X 42,K H4V S 3H4V S 315S T 3H4)*Y (13H S ]]T 15S 31),2
*^+'M 4(m )*\,G 1T 4Y ?,34*6C514,54&&$$!
5+4#'&1+"%"("%/J 7+4/%0+"%&'7&4H $/&FH &0/7
"%(+%+1/."'#4/Z Y<);069/4/
G L M =a 41&"O E L M ==0)2H4,[&"a `O HS ,[,S
&&C5H))T )Z ^S 3H4V S 3152S ,+<HY 2152&M S ,b 1,[`,1X 4*213Y )Z ?,Z )*V S 31),C514,54_K 45H,)T )[Y &M S ,b 1,["!&%%AA
G H 01$&61"K H12]S ]4*]*424,32S ,4(V )+4T )Z ),46+1V 4,21),S T +S V.*4S \Z T )(2.S 24+),2HS T T )((S 34*4c0S 31),2&1,(H15H E G L M N25H4V 4124V ]T )Y 4+Z )*2]S 51S T +125*431Q S 31),S ,+Z )0*234]2K U B@0,[46n 033S V 43H)+12024+Z )*31V 4+125*431Q S 31),'K H4,0V 4*15S T *420T 32+4V ),23*S 343HS 33H4YS *4V )*4S 550*S 343HS ,3H)24.S 24+),G L M N 25H4V 4-21V 0T 3S ,4)02T Y &3H4]*)[*S V1,3H12]S ]4*5S ,(4S \4,)251T T S 31),Z )*),46+1V 4,21),S T +S V 6.*4S \Z T )(2V )*44Z Z 4531X 4T Y &S ,+])2242242H1[H4**42)T 031),3HS ,3HS 3)Z G L M N25H4V 4'I /2E "$70"+S V 6.*4S \Z T )(2-Z 1,134X )T 0V 4V 43H)+2-E G L M N25H4V 4
%
7A 翁磊&等'基于有限体积E
GL M N 格式的一维溃坝流模拟'G L M =a 41&43S T 'C1V 0T S 31),)Z ),46+1V 4,21),S T +S V.*4S \.S 24+),Z 1,134X )T 0V 4E G L M N25H4V 4'。