第11卷第2期Vol.11No.2
吕梁学院学报
Joornal of Liiliang University
2021年4月
Apn2221
-基础理论研究-
杨艳，范示示，王艳永
(吕梁学院数学系，山西离石233021)
摘要：用二维水波方程，从时间、空间两个角度对连续的水波方程进行离散化处理，推导出二维水波偏微分方程的差分格式，同时从稳定性的角度讨论的参数约束条件，最后用MATLAB模拟水面波动动画。 关键词：水波方程；差分法；参数约束；稳定性
中图分类号：0221.2文献标识码：A文章编号：2095-185X(2221)22-0021-24
水波动画的模拟主要分为三类[1]：第一类是基于波的分析方法，其中最简单的一种是用正弦函数来进行仿真，模拟水波表面•这种方法只是从形态上模拟了水波的特性，并没有考虑水波在运动过程中外部与内部环境以及其他现实因素的影响•第二类是基于粒子的分析，将流体看作是离散的粒子集，每个粒子按照一定的规则不断的运动•第三类是基于物理模型的方法，研究流体的运动，最典型的是通过求解粘性不可压缩流体的动量守恒方程:Navier-stokes方程，得到流体完整的速度场.
本文采用二维水波方程，设定水滴落下的初值条件，得到方程的差分解，并且通过稳定性讨论，得到参数的取值范围，模拟水面荡漾动画1 1模型的建立
本文选用二维水波方程⑵:
T T=«2(昙+昙卜“警，(”，y)e°=[ a,X[c,d],t e[°,+8)⑴
(1)式中Z表示水面高度为空间变量t，为时间变量，其中C、“均是非负的常数,c表示单位时间内水波位移的距离，“表示流体的黏度.
由于在模拟水滴湖面波动的过程中，波面是从滴落位置(％』2)开始呈直径为W的圆环扩散的，在垂直方向的最大位移为h.为了能够产生水滴湖面造成的湖面的波动这一状态，先获取湖面下落的初始圆环的高度人和宽度设定:
z((,2,0)
2
he-gf"22(%_%o)2+(y_2)§(号)
2其他
⑵
为方程的初始条件，表示圆环内部曲面向下凹陷，同时圆环外水面的垂直位移都为零，如图1-同时初值条彳包括:
(3)
边界条件选取
收稿日期：2221-00-01
基金项目：山西省自然科学基金(221821D121212)；吕梁学院优质课程建设项目(YZKC221912)；吕梁学院教学改革创新项目(JXGG222232,JXGG221943)
作者简介：杨艳(1983-)，女，山西柳林人，副教授，研究方向为分数阶微分方程数值解、计算几何。
竺 二o，其中n 为方形区域的边界. （4）
d ” dQ
2差分法求解模型
区间［^b ］和［c,d ］分别兀等分，步长分别为山为了
讨论方便，取方形区域，此时记相同步长为Ad ），彳
结点
（，力）（0wjW e,0w/W 兀）.令 A 为 步长，并用 z （%（, ,J, %）
表示%时刻点（（，力）的水面高度.2.1 方程的离散
将方程（1）中各个导数用如下二阶格式离散,
霭，力,％）
d 2(
d 2Z （，力，％ + 1） — Z （，力，％_i ） *（虫）2
3
力，％）（2+i ，力，％） — 2（2，力，％） +（2-1，力，％）
/人八 2Ad
舟山到普陀山最佳方案
+ （Ad ）
d 2(
d 2力,％）
红河谷歌曲d 2(
d%
2，%）
((，2+i ,%) — 2((，2，%) +((，2_i ,%) / 人八:
A +⑷)(2,2，%+1 ) — 2(2,2，%) +(2,2，%- 1)+ (A% 2
在上述各式中舍去误差项，用（表示（2，力，％）的近似值，代入原模型（1 ）中，得到离散方程,
A %2(
k ck]k k ck]k\
(+1, 一 2(, + (+-1, + (+1 一 2 + (1 )—仏 A d 2 A d 2 丿
+1 k —1z - - — Z --+,
+,2(5)
+ 1,2,…- 1,，二 1,2,-1, k 二 1,2,…，记网比厂二 A d
4 — 8c 2 厂2 _ 2c 2 厂2lAt + 2 ' °2 iAt + 2 '代入(5)式整理，得方程的差分解为：湖南冬季旅游最佳去处
Ql
(— +。2( 1 , + ( + 1 , + _ 1 + 球,+ 1、+a i
a 3_ “A% — 2
fiAt + 2
k+1
2.2初始条件的离散
将
佝(—1 k _ 0,1,2,…
d （ c 、 2（2,2，%）— （2,2，%1d （2,2,0）_ -------------------------------舍去误差项，代入 （3），得：
2A
)+ ( A ) 2
-1
—z --也二o ,
即
（7）式中取k 二0，得：
1
（_
将（）式代入整理，得：
1
2A/
—1
1
(，广
(+。2( Z 0—, + f+i , + (—+ (+1)+
1
1 — tt 3［ +。2（ （1, + （+1，综合（7）式和（9）式，得方程的差分解为：
+ (，1 +(,+1 )]
1———[
+ a 2(Z k 1 +，+ ( • 1 + Z k -+1 )] , k_ 0
(+i 二丿 1 -a 3
， z ，
， +- 1 ++1
bl (—
+。2(, + (i , + (_ i + (+i 、+ 偽彳,，• k _ 1，2,…
⑹
(7)
⑻
(9)
( o )
3解的稳定性分析
阻尼力口为与温度有关的量，可以通过经验得到，数量级为1 0-3.为了讨论稳定性，需确定波速c 和
A t的关系.
首先讨论||Z1||8•假设忖||8=Zj l，由⑼式，得：
=1_1^[%2,+a2（22_1，+，+Z2,,_1+2°，+）]，
故
1121118=1°,，=11a[a2,+a2（°-1，+2+1，+2，1+2+1）]■（11由（6）式，得°2>2,3<2（由于丛为1。-0,A t取充分小，即假设山A W2），故：
141
1-a0=zzZTI>1
此时，（11）式可写为：
II218=1-a+a2（2_!，+Z,1，，+2，-1++1）
假设a»，即：
C W；■（12）
又由于a>。，上式可写为：
1118=1-a[a1+a（2_1，+|Z++1，|+Z2，-+2，+1）]
W丄；[a|22||8+°2（||22||8+II22|8+II22||8+II22||8）]
a1+4a”2”
=^^3112
又a^a2=1，故上式写为||21|8II8，即在假设Cw斗的条件下，第一层稳定，同理可以验证, 3
假设C>；时，不稳定.故以下讨论均假设片w1
（13）同理，%=2,1,2,…时，由（7）式可得：
故：
I11I8w（a1+4a2）||8-a o||21||
II8w（a+4a2）
I2I8=||8■咸阳旅游景点一日游 排名
II8-ao||2-||西安法门寺一日游
记匕+1=（||1+1|8,I1"||8)T，由上式有：
啓w（a+4a7）f
2
记，
4=（a+4a-a，2
要使(13)式稳定，需矩阵A的谱半径p(幻W1+C°A，其中C。为常数.
进一步计算A的特征值为：
入二亍(a[+4。2±槡(a[+4位2)2_4a o),
于a1>2，a1+4a2
4
fjAt+2
Mg V。，故
p（彳）=i（a+4a_槡（位1+4位2）2_4a）2+J槡--A a2
A A t+21+J槡-/At,2-a A
At,+2'
以下证明存在常数c °，使得p() w 1 + c °A2，即存在常数C 。，使
丿8 - -A A+ — a A
a A + 2
w c 2a
+ ,
用反证法，假设对于认识常数c °，都有
丿8 - -A A 22 — a A
a A + 2
> Co △/，
成立，按C 。降轄整理得:
A 22 (“A, + 2)2 C q 2 + 2m A+((A , + 2)C o + A 2 — 8 < 2，
而上式中C 2的系数大于2，即上式对任取的C 。均不成立.故稳定条件为Cr 2 w *，m A w 22在代码的运行中，根据已知的值以及画图的稳定性，取参数A =r= 2.01, = 4
4模拟结果
3、
0 00 0
图3波面荡漾开画面
图2水滴落入波面反弹
历史航班查询0 0
图4整个水域波动
0 0
图5水面趋于平静
以上四幅图是随着时间增加，水面的各个不同时间运动图，从水滴反弹到波面趋于平静的状态.
5结束语
限差分法求解 的水波方程，模拟了随机的水 从而引 面波动的动画2
的参数取值 ，并分析 值和边界值的 •结果表明，该方法可以真实的模拟水 面波动的
现象.本文是作为本科生的数值分析课程的 实验项目.2
[ 2 ] [ 2 ]中均有类似
，创新点是从数学
对 的 进行推导•下一步工作类似于文献[2] [3]，可以结合计算机图形学课程 模型，使得动
画具有真实感，提高学生实验能力.
参考文献：
[1] 李付鹏•基于有限差分法和有限体积法的水流动画模拟[D].合肥：安徽大学,2024.[2] 庞明勇.基于离散模型的 水波实时动态模拟方法[]•水利学报,2027(11).
[3] 程甜甜•面向
水域的动态水面
建模方法研究[D].苏州:
学,2008.
