1. 引言
作为用来刻画未知变量、时间导数、空间变量和导数之间的关系的偏微分方程的一种——二阶双曲型偏微分方程是最为常见和十分重要的 [1] - [6]。这是因为双曲型偏微分方程的基本的模型是波动方程,而波动现象普遍存在于我们的生活和生产之中。例如空气中传播的声波、绳运动产生的波、水波等,这些波叫做都是波,具体讲是机械波 [7] - [13]。另外,还有工业生产和工程中应用的光波、无线电波等,这些波具体讲作叫做电磁波 [14] [15] [16] [17],故求解二阶双曲型的偏微分方程具有重要意义。
求解二阶双曲型的偏微分方程有多种方法。但是,当二阶双曲型偏微分方程的解析解难以求出时,人们常诉诸于数值模拟方法求其数值解。这其中,人们最常用的方法有显式差分格式和隐式差分格式方法 [18] [19] [20] [21] [22]。然而,这两种方法各有什么特点,他们分别适合解决什么问题呢?尚没有明确的结论。为此,本文对这个问题做一研究。
2. 利用显式差分格式法进行数值求解
2.1. 数值解和误差分析
现给出一个双曲型波动方程问题,如下所示:
,
,
,
,
该问题精确解为
。
首先,利用显式差分格式的方法求解其数值解,并求出此数值解的误差的绝对值。具体做法是:取空间步长h为1/100,时间步长
为1/200,此时的步长比s为1/2,计算固定
时,时间从0.1到1秒的数值解以及数值解和精确解之差的绝对值(误差),结果如表1所示。再次取空间步长和时间步长分别为1/100和1/100,结果如表2所示。第三次取他们的值分别为1/100和1/50,此时计算结果如表3所示。
![](Images/Table_Tmp.jpg)
Table 1. Exact solution, numerical solution and absolute error at some nodes ( h = 1 / 100 , τ = 1 / 200 )
表1. 部分节点处精确解、数值解和误差绝对值
![](Images/Table_Tmp.jpg)
Table 2. Exact solution, numerical solution and absolute error at some nodes ( h = 1 / 100 , τ = 1 / 100 )
表2. 部分节点处精确解、数值解和误差绝对值
![](Images/Table_Tmp.jpg)
Table 3. Exact solution, numerical solution and absolute error at some nodes ( h = 1 / 100 , τ = 1 / 50 )
表3. 部分节点处精确解、数值解和误差绝对值
由表1、表2、表3中的数据可以看出,当取空间步长和时间步长分别为1/100和1/200时,有效数字位数可以达到6位;当取空间步长和时间步长分别为1/100和1/100时,有效数字位数可以达到5位,他们都很好的逼近了精确解。然而,当空间步长和时间步长分别为1/100和1/50时误差非常大,此时显式差分格式方法也不收敛。
进一步,我们绘制出了空间步长和时间步长分别为1/10和1/20,1/20和1/40,1/40和1/80时,在某一时刻的误差曲线和整体的误差曲面图,如图1和图2所示。由这两个图可以看出,当空间步长和时间步长取1/10和1/20时数值解的误差较大,随着步长减小,每个坐标的误差均减小。
![](//html.hanspub.org/file/15-1251203x25_hanspub.png)
Figure 1. The error curves corresponding to the length of asynchrony when t = 1
图1. t = 1时不同步长时对应的误差曲线
![](//html.hanspub.org/file/15-1251203x26_hanspub.png)
Figure 2. Error surfaces corresponding to asynchronous length
图2. 不同步长时对应的误差曲面
为分析不同步长下显式差分格式的准确程度,我们固定步长比为1/2,依次取空间步长为1/10、1/20一直到1/1280,计算在
处的最大误差,以及步长细分时最大误差的比值。结果如表4所示。
![](Images/Table_Tmp.jpg)
Table 4. Maximum error of asynchronous long time numerical solution
表4. 取不同步长时数值解的最大误差
根据表4可知,当时间步长从1/10变化到1/1280时,最大误差从
减小到
,误差大幅度减小,说明步长越小,数值解的误差越小。同时,当时间步长与空间步长同时缩小为原来的1/2时,最大误差缩小为原来1/4。
2.2. 显式差分格式的计算时间
程序的运行时间反应了算法的复杂程度和计算效率,因此下面计算当步长比固定为1/2时不同步长情况下的显式差分格式的程序运行时间。考虑计算机运行的不稳定性,5次重复运行,取时间的平均值,结果如表5,平均时间变化图如图3所示。
从表5可以看出,当空间步长和时间步长为1/10和1/20时,程序运行的时间约为0.001秒,而当空间步长和时间步长为1/5120和1/10,240,运行时间约为0.8秒,时间的变化非常大。从图3可以看出,当步长减小时,程序的运行时间逐渐增加,并且增加速度越来越快。
![](Images/Table_Tmp.jpg)
Table 5. Calculation time of explicit difference scheme with different step size
表5. 不同步长下的显式差分格式计算时间
![](//html.hanspub.org/file/15-1251203x33_hanspub.png)
Figure 3. Time mean of explicit difference scheme under asynchronous long time
图3. 显式差分格式不同步长下的时间均值
3. 利用隐式差分格式进行数值求解
3.1. 数值解与误差分析
下面利用隐式差分格式的方法再次求解上述问题的数值解,并求出其数值解的误差的绝对值。取空间步长为1/100,时间步长为1/100,此时的步长比为1,计算固定
时,时间从0.1到1秒的数值解以及数值解和精确解之差的绝对值(误差)。结果如表6所示。其次再取空间步长和时间步长分别为1/200和1/200,再次求其数值解和误差,结果如表7所示。
进一步,我们也绘制出了当空间步长和时间步长取1/10和1/20,1/20和1/40,1/40和1/80时在某一时间的误差曲线和整体的误差曲面图,如图4和图5所示。由这两个图可以看出,当空间步长和时间步长取1/10和1/20时数值解的误差较大,随着步长减小,每个坐标的误差均减小。
![](Images/Table_Tmp.jpg)
Table 6. Exact solution, numerical solution and absolute error at some nodes ( h = 1 / 100 , τ = 1 / 100 )
表6. 部分节点处精确解、数值解和误差绝对值
![](Images/Table_Tmp.jpg)
Table 7. Exact solution, numerical solution and absolute error at some nodes ( h = 1 / 200 , τ = 1 / 200 )
表7. 部分节点处精确解、数值解和误差绝对值
![](//html.hanspub.org/file/15-1251203x39_hanspub.png)
Figure 4. The error curves corresponding to the length of asynchrony when t = 1
图4. t = 1时不同步长时对应的误差曲线
![](//html.hanspub.org/file/15-1251203x40_hanspub.png)
Figure 5. Error surfaces corresponding to asynchronous length
图5. 不同步长时对应的误差曲面
为分析不同步长下显式差分格式的准确程度,我们固定步长比为1,依次取空间步长为1/10、1/20一直到1/640,计算在
处的最大误差,以及步长细分时最大误差的比值,结果如表9所示。
根据表8中数据可知,当时间步长从1/10变化到1/640时,最大误差从0.003508减小到
,误差大幅度减小,说明步长越小,数值解的误差越小。同时,当时间步长与空间步长同时缩小为原来的1/2时,最大误差缩小为原来1/4,最大误差的变化趋势与显式差分格式的变化趋势相似。
![](Images/Table_Tmp.jpg)
Table 8. Maximum error of asynchronous long time numerical solution
表8. 取不同步长时数值解的最大误差
3.2. 隐式差分格式的计算时间
下面计算当步长比固定为1时不同步长情况下的隐式差分格式的程序运行时间,考虑计算机运行的不稳定性,也重复运行5次,取时间的平均值,结果如表9所示,平均时间变化图如图6所示。
![](Images/Table_Tmp.jpg)
Table 9. Calculation time of implicit difference scheme with different step size
表9. 不同步长下的隐式差分格式计算时间
从表9可以看出,当空间步长和时间步长为1/10和1/10时,程序运行的时间约为0.0045秒,而当空间步长和时间步长为1/5120和1/10240,运行时间约为1.85秒,时间的变化也非常大。从图6可以看出,当步长减小时,程序的运行时间逐渐增加,并且增加速度越来越快。与显式差分格式方法比较,隐式差分格式方法的运行时间较长,这说明隐式差分格式的复杂度较高。
![](//html.hanspub.org/file/15-1251203x46_hanspub.png)
Figure 6. Time mean of implicit difference scheme under asynchronous long time
图6. 隐式差分格式不同步长下的时间均值
4. 两种差分格式的比较和分析
本文对显式差分格式和隐式差分格式方法的特点进行了比较和分析,结论如下:
1) 运用显式差分格式得到数值解的收敛性取决于选取的步长比,当步长比小于等于1时算法是收敛的,因此需要注意设置的时间和空间步长;而隐式差分格式在任何情况下都是收敛的;
2) 在算法收敛的情况下显式差分格式和隐式差分格式得到的数值解的误差都比较小,并且设置的步长越小,划分的格数越多,运用两种方法得到的数值解的误差都越小,说明这两种方法都能很好地逼近精确解;
3) 对于显式差分格式和隐式差分格式方法,当步长比固定时,步长以2的倍数减小时,数值解的最大误差几乎以4的常速率减小;
4) 运用显式差分法和隐式差分法求解运行的时间均随划分的格数的增加而增加,说明格数越多,计算的次数越多,复杂度越高;
5) 当划分的格数相同时,运用显式差分法求解运行的时间要比隐式差分法运行的时间短,说明显式差分格式法比隐式差分格式法计算简单,求解效率更高。
综上所述,隐式差分格式法和显式差分格式法都是求解双曲型偏微分方程数值解的有效方法,但是各有特点。显式差分格式法对设置的步长比有限制,太小的步长比可能导致数值解不收敛,但是其程序运行时间短,计算效率高,而隐式差分格式法对任意的步长都可以求解,但是程序运行的时间较长,计算较复杂。