1. 引言
在众多机器学习算法取得快步发展之前,工程领域通常采用的都是传统的物理模型。这些领域的大部分问题都以偏微分方程(Partial Differential Equation, PDE)的形式刻画或描述。直接求解物理模型可进行精准预测,但复杂的场景下的偏微分方程的解通常难以显式地表示。传统数值求解方法在高维度的情况下,计算代价相当大。相较于基于网格的方法,深度学习方法能弥补传统数值方法求解偏微分方程的先天缺陷 [1] 。因此,在理论研究和工程应用中开发设计出高效精确的基于深度学习的偏微分方程数值解法对降低计算成本和提高计算精度具有重要意义。
Figure 1. Three numerical methods for solving partial differential equations
图1. 求解偏微分方程的三种数值方法
2017年,Raissi教授和他的合作者们提出了一套深度学习算法框架,并将其命名为“Physics-informed Neural Networks”(物理信息神经网络) [2] 。2021年,Karniadakis,Lu等研究人员在著名的Nature子刊Nature Reviews Physics上发表了一篇综述性文章 [3] ,值得注意的是,作者将物理约束和数据驱动结合的方法,归纳统称为物理信息约束的深度学习方法(Physics-informed Deep Learning, PIDL),见图1。
在深度学习方法中,增加网络的深度和复杂性可能会导致损失函数的参数空间中具有梯度爆炸或梯度消失的区域。如果出现梯度爆炸或消失,模型就会变得非常不稳定,并使先前的训练进展变得徒劳无功。传统深度学习方法以链式求导法则和自动微分技术为基础,能够在网络结构中前馈和反向传播梯度信息。但基于自动微分的PINN只能在大量训练配点的情况下实现高精度求解,否则它们很容易偏向于完全的数据驱动方案,进而优化到不满足物理条件的解 [4] 。为了提升PINN在稀疏样本下的训练效率,可以采用数值微分来定义损失函数。基于数值微分的损失函数能够强连接相邻的搭配点,以便在稀疏样本状态下实现有效的训练。
原始的PINN方法需要根据问题的定解条件在一个给定的矩形区域上进行求解。但实际情况中常常因为解函数存在区域内的某些部分无法进行直接观测,从而无法获得该部分的观测信息,进而导致只能在非矩形区域上求解。2019年提出的DeepXDE软件库中提供了7种基本几何结构 [5] ,可以帮助使用者利用基于自动微分的PINN在一些更复杂的问题区域上进行求解和模拟。
受此启发,本文结合区域分解的思想,将基于数值微分的PINN方法应用到非矩形区域上,以扩展原有方法的适用范围。文章第一节简要叙述了PINN的发展历程;第二节给出PINN的一般模型、基于数值微分的PINN的基本原理、区域分解技术和KdV方程的问题模型;第三节展开数值实验,使用基于数值微分的PINN方法对一个非矩形区域上的KdV方程进行求解,以验证方法的有效性;第四节总结了方法的优缺点,并提出了潜在的应用前景。
2. 基本原理
2.1. 物理信息神经网络
在物理信息神经网络的框架下,我们考虑具有如下形式的偏微分方程:
(1)
其中,ut代表了方程的隐藏解,
是一个非线性微分算子,而Ω是Rn的子集,T是终止时刻,下标表示对时间或空间的微分运算。在给定初始、边界条件和物理参数且无法求得解析解时,可以利用有限差分、有限元等传统数值方法求解。PINN则是考虑建立一个神经网络来逼近偏微分方程的解。假设
即为用于逼近解函数
的神经网络,偏微分方程的残差即为:
(2)
并定义物理信息神经网络的损失函数为:
(3)
其中,Lossu表示损失函数中的数据驱动部分。数据驱动部分以在初始、边界条件上,根据实验观测或物理仿真取得的训练数据进行训练。这部分损失函数依赖观测数据作为标签指导模型进行训练,属于监督学习过程。另一部分损失函数Lossr则表示损失函数中的物理约束,也即物理约束部分 [6] 。由于自动微分技术可以高效地获得方程的残差,这部分损失函数在训练时不需要额外的标签数据,属于无监督学习过程。
这种结合监督学习和无监督学习的设计使得损失函数不但可以保证和已知数据的误差尽可能小,而且可以尽可能地满足偏微分方程隐含的物理约束。在输入需要的时空数据后,首先使用一个全连接神经网络来逼近解函数,然后利用自动微分技术求出残差放入损失函数,最后使用梯度下降等优化算法获得全连接神经网络的权重参数,最终训练出可以满足偏微分方程的神经网络参数以逼近解函数。
2.2. 嵌入数值微分技术
作为基于神经网络的无网格方法,PINN相较传统基于网格的数值算法可以节省大量计算开销,这就使其在高维问题上有了用武之地,并且已经证明这种方法在求解椭圆型方程时可以依范数收敛到问题的真实解 [7] 。但自动微分技术需要使用大量的训练配点,增加网络深度又可能导致神经网络在训练过程中出现梯度爆炸或梯度消失的现象 [8] ,对此问题,采用基于数值微分的神经网络技术更为恰当。
在基于数值微分的PINN方法中,使用差分方法而非自动微分来实现数值微分过程。通过选择合适的点并利用泰勒级数展开消除误差项,可以得到数值微分导数项。将函数
,
进行泰勒展开,联立得到
的一阶导数
的二阶中心差分公式为:
(4)
类似的,得到
的二阶导数
的二阶中心差分公式为:
(5)
的三阶导数
的二阶中心差分公式为:
(6)
本文使用上述差分公式进行数值微分,以替代神经网络框架中的自动微分过程。如图2所示,上述差分格式能够以不同权重强耦合相邻的训练配点,以便在稀疏样本状态下的模型仍能获得足够信息,进而实现有效的训练。
Figure 2. Schematic diagram of the central difference format of the second and third derivatives
图2. 二阶导数和三阶导数的中心差分格式示意图
2.3. 区域分解求解思想
区域分解技术基于分治思想,其本身是一类迭代方法。具体来说,区域分解技术基于计算区域的划分:
(7)
然后将计算过程和计算量分配到子区域上进行。对子区域进行划分能够极大地提高算法的可并行性,使得利用并行计算技术加速求解成为可能;同时可以结合各子区域的特点,在各子区域中使用不同的数学模型、计算方法和求解格式,以加快运算速度。受迭代方法中的区域分解思想启发,软件库DeepXDE中的PINN方法首先定义了一些基本几何图形,通过将复杂区域分解为几种这些基本几何图形来实现将原问题在复杂区域进行求解。但这种基于区域分解的求解思想尚未在有数值微分的PINN方法中有过具体应用。
3. 模型建立
KdV方程是一个描述浅水孤立波运动的物理模型,KdV方程与孤子解在流体力学、光学等许多领域得到了广泛的研究和应用。考虑KdV方程:
(8)
给定不同的定解条件,该方程有不同的解。其一种孤立波解的解析形式为:
(9)
该解析解由Stephen Anco博士给出。其中,c代表波速。取c = 2,x0= −1,定义解空间为:
(10)
将该L型区域分解为两个子区域Ω1和Ω2,其中Ω1为
,
的矩形,Ω2为
,
的矩形。在两个子区域上分别定义一个物理信息神经网络进行求解。
根据方程的具体形式,令:
(11)
由式(3),定义物理信息神经网络的损失函数为:
(12)
其中,Nu为观测数据点的个数,Nr为根据残差计算点的个数。设差分步长h = 0.01,使用式(4)和式(6)对方程(8)中的一阶导数和三阶导数进行近似。
4. 数值实验
为了验证上述方法的可行性,本文设计了数值实验进行验证。
实验中的KdV方程的参数c = 2,x0= −1,求解的区域范围如式(12)所示,区域分解方案如图3所示。实验中的物理信息神经网络的拓扑结构为全连接网络,输入层有2个神经元、输出层1个神经元,隐藏层共2层,每层50个神经元,隐藏层间以双曲正切函数为激活函数。神经网络使用Adam优化器以0.01为学习率迭代10,000次,对两个子区域都分别采用这样的神经网络进行求解。实验中的数值微分方法采用式(6)和式(8)的中心差分公式,步长设置为0.01。所有实验均在Intel(R) Core(TM) i7-10750H CPU @ 2.60 GHz上进行,编程语言采用Python。
Figure 4. Solution of PINN based on numerical differentiation
图4. 基于数值微分的PINN的求解结果
图4展示了基于数值微分的物理信息神经网络在L型区域上求解KdV方程的求解结果。总体来看,基于数值微分的PINN在L型区域上比较好地描述了孤立波解的基本性态,可以清晰看到孤立波的运动过程。为了分析在进行区域分解后的PINN的求解效果,下面分别单独讨论PINN在子区域Ω1,Ω2上的求解效果。
Figure 5. Comparison of solution on Ω1
图5. 区域Ω1上的求解结果对比
单独考虑子区域Ω1,图5将PINN的求解结果和真实解的图像进行了对比。图5(a)为PINN求解结果示意图,图5(b)为真实解示意图。直接观察图像,两者的解函数性态并无较大差异,PINN较为清晰地展现了孤立波在该区域上的运动过程。观察PINN在区域Ω1上的绝对误差,见图6(a),在t = 0到t = 0.8秒的区间上,PINN的绝对误差都保持在较低的水平且分布比较平均,在接近t = 1秒时,PINN的误差开始变大。如图6(b),从收敛曲线来看,基于数值的PINN的收敛非常迅速,在迭代初期就迅速下降到了10−2的数量级上,且迭代过程整体保持平稳,并未出现剧烈震荡。
Figure 6. Convergence curve of absolute error and loss function on Ω1
图6. 区域Ω1上的绝对误差和损失函数收敛曲线
单独考虑子区域Ω2,从解出的图像来看,PINN的结果仍然比较清晰。但与PINN在区域Ω1上的绝对误差表现略有不同,基于数值微分的PINN方法在Ω2上的绝对误差分布并不平均,且随着t的增大,绝对误差也有增大的趋势。从收敛曲线来看,基于数值的PINN的损失函数在区域Ω2上仍然迅速下降到了10−2的数量级,但在2000步前出现了一些震荡,2000步后才保持平稳。
Figure 7. Solution and convergence process of PINN on Ω2
图7. PINN在区域Ω2上的求解结果和收敛过程示意图
图4~7展示了PINN求解KdV方程得到的图像,可以直接观察PINN的求解效果。为了展示更为详细的求解结果,分别以
和
为粒度进行采样,将使用上述物理信息神经网络在L型区域上求解KdV方程的部分结果数据汇总如下。
,
时的求解结果在表1中给出;
,
时的求解结果在表2中给出,二者共同组成L型区域上的求解结果。例如,x = 0,t = 0.08时,神经网络预测值
,以此类推。神经网络的求解结果可以通过设置不同的粒度进行输出,极大方便了该方法在各种实际场景下的应用。
Table 1. Display of data sampling of PINN’s solution
表1. 物理信息神经网络求解结果数据采样展示
Table 2. Display of data sampling of PINN’s solution
表2. 物理信息神经网络求解结果数据采样展示
5. 结论
偏微分方程的数值求解方法始终是一类重要的研究热点。传统数值求解方法包括有限差分、有限元、有限体积等。这些方法需要使用网格来近似偏微分方程的定义空间,在高维度的情况下,其计算代价相当大。与之不同的是,经典的机器学习算法都以纯数据驱动为主。相较于基于网格的方法,深度学习方法能使参数数量极大地降低,从而带来更少的计算消耗,使得其对高维问题有更好的扩展性,这就能弥补传统数值方法求解偏微分方程的先天缺陷。
基于数值微分的PINN方法可以看作是传统数值方法与数据驱动算法的结合,这类PINN方法既保留了物理信息神经网络的原有优势,又避免了原始框架中存在的梯度爆炸等隐患,并且能够在稀疏样本下进行训练。既往研究中求解的问题都是在矩形区域上进行的,本文将这一方法应用到一个非矩形区域内。数值实验表明,基于数值微分的物理信息神经网络方法能够有效应对非矩形区域的KdV方程的求解任务,且其求解效果较好、收敛速度较快。针对其他类型的方程的求解还需要进一步研究。
基金项目
四川省大学生创新创业训练计划创新训练项目(S202210619108)。