楊云磊,侯木舟,羅建書
(1. 中南大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,湖南 長(zhǎng)沙 410083; 2. 湖南交通工程學(xué)院 高科技研究院,湖南 長(zhǎng)沙 421001)
科學(xué)和工程領(lǐng)域中的很多問(wèn)題涉及到偏微分方程(PDE)的數(shù)值解,大量的物理現(xiàn)象,從流體流動(dòng)到電磁場(chǎng)領(lǐng)域,熱擴(kuò)散和海洋聲波的傳播等等廣泛的科學(xué)問(wèn)題,建立的很多數(shù)學(xué)模型都是對(duì)流擴(kuò)散方程,關(guān)于模型的準(zhǔn)確解或近似解的求解,對(duì)研究數(shù)學(xué)模型的意義非常重要,然而很多情況下,求得模型的準(zhǔn)確解或者良好的數(shù)值解是非常困難的或不可行的,因此研究對(duì)流擴(kuò)散方程的數(shù)值解很有必要.
PDE數(shù)值解法的傳統(tǒng)研究包括差分格式的建立,網(wǎng)格剖分,線性方程組的求解及差分格式的相容性、收斂性和穩(wěn)定性的證明[1],圍繞這些研究產(chǎn)生了很多PDE數(shù)值方法,如比較熟悉的差分方法、有限元方法、有限體積法等等[2]. 一維常系數(shù)對(duì)流擴(kuò)散方程是一類特殊的PDE,這些數(shù)值方法也適用于計(jì)算對(duì)流擴(kuò)散方程.
神經(jīng)網(wǎng)絡(luò)的應(yīng)用涉及很多領(lǐng)域,如模式識(shí)別[3]、圖像處理[4]、風(fēng)險(xiǎn)評(píng)估[5]、系統(tǒng)控制[6]、預(yù)測(cè)[7-8]、分類[9]等方面. 隨著神經(jīng)網(wǎng)絡(luò)研究的發(fā)展,關(guān)于神經(jīng)網(wǎng)絡(luò)逼近能力的研究促進(jìn)了神經(jīng)網(wǎng)絡(luò)在微分方程數(shù)值解中的應(yīng)用. 如: 侯木舟等提出構(gòu)造小波神經(jīng)網(wǎng)絡(luò)逼近實(shí)函數(shù)[10],構(gòu)造小波RBF神經(jīng)網(wǎng)絡(luò)[11],L2(R)上的RBF神經(jīng)網(wǎng)絡(luò)[12]及衰減的RBF神經(jīng)網(wǎng)絡(luò)[13]逼近任意多元函數(shù),Guangbin Huang提出極限學(xué)習(xí)機(jī)的概念[14],并在其很多研究[15-18]中給出極限學(xué)習(xí)機(jī)逼近能力的證明.
根據(jù)神經(jīng)網(wǎng)絡(luò)的逼近能力,可以構(gòu)造微分方程問(wèn)題的近似解,從而建立微分方程數(shù)值解的神經(jīng)網(wǎng)絡(luò)模型. 如微分方程兩點(diǎn)邊值問(wèn)題的神經(jīng)網(wǎng)絡(luò)解法[19],傳輸線方程的神經(jīng)網(wǎng)絡(luò)模型算法[20]等,均采用余弦基函數(shù)作為神經(jīng)網(wǎng)絡(luò)隱層激活函數(shù). 本文根據(jù)勒讓德多項(xiàng)式的微分性質(zhì)及矩陣張量積的性質(zhì),構(gòu)造一維對(duì)流擴(kuò)散方程的勒讓德神經(jīng)網(wǎng)絡(luò)模型,采用改進(jìn)的極限學(xué)習(xí)機(jī)(IELM)算法求解網(wǎng)絡(luò)權(quán)值,并進(jìn)行數(shù)值實(shí)驗(yàn)來(lái)分析網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)對(duì)計(jì)算結(jié)果的影響.
勒讓德多項(xiàng)式P0(x),P1(x),…,Pn(x)滿足三項(xiàng)遞推公式
(1)
(2)
(3)
(4)
(5)
由方程 (4)~(5) 可得
(6)
其次, 求(n+1)Pn+1(x)=(2n+1)xPn(x)-nPn-1(x)關(guān)于x的導(dǎo)數(shù)
(7)
由方程 (6) 和 (7) 可得
(8)
由方程 (8), 經(jīng)過(guò)遞推易得P′(x)=P(x)M和矩陣M.
定義1[21]設(shè)A=(aij)∈Cm×n,B=(bij)∈Cp×q, 則如下分塊矩陣
為A與B的張量積, 簡(jiǎn)記為A?B=(aijB).
矩陣張量積性質(zhì):
如果A∈Ck×m,B∈Cp×s,C∈Cm×n,D∈Cs×q, 則(A?B)(C?D)=(AC)?(BD).
一維對(duì)流擴(kuò)散方程第一邊值問(wèn)題一般形式為
(9)
式中: 系數(shù)u,v是常數(shù);x∈[0,a];m是已知常數(shù). 問(wèn)題(9)就稱為常系數(shù)一維對(duì)流擴(kuò)散方程的數(shù)值解問(wèn)題.
定理2對(duì)任意兩元連續(xù)二階可導(dǎo)函數(shù)ρ(x,t)∶D→R,D=[a,b]×[c,d]存在自然數(shù)n, 常數(shù)βij(i=0,1,…,n;j=0,1,…,n), 勒讓德多項(xiàng)式P0(x),P1(x),…,Pn(x),P0(t),P1(t),…,Pn(t), 使得當(dāng)n→∞時(shí),有
(10)
‖ρ(x,t)-ρ*(x,t)‖2=
?D[ρ(x,t)-ρ*(x,t)]2dxdt=
(11)
式中:ρ*(x,t)為ρ(x,t)的最佳逼近, 即‖ρ-ρ*‖ 能夠取得最小值.
證明假設(shè)
β00P0(x)P0(t)+β01P0(x)P1(t)+…+
β0nP0(x)Pn(t)+β10P1(x)P0(t)+….
(12)
令
E=‖ρ(x,t)-ρ*(x,t)‖2=
?D[ρ(x,t)-ρ*(x,t)]2dxdt.
(13)
為了找到能夠使得E取得最小值的βij的值, 令E關(guān)于所有βij的偏導(dǎo)數(shù)等于零,可得如下方程組
(14)
令
P(x)=[P0(x),P1(x),…,Pn(x)],
(15)
P(t)=[P0(t),P1(t),…,Pn(t)],
(16)
P(x,t)=P(x)?P(t)=
[P0(x)P0(t),P0(x)P1(t),…,P0(x)Pn(t),
P1(x)P0(t),…,P1(x)Pn(t),…,Pn(x)Pn(t)],
(17)
β=[β00,β01,…,β0n,β10,…,β1n,…,βn0,…,βnn],
(18)
(19)
方程 (14) 可以改寫成
(20)
由矩陣分析理論可知,方程 (20) 有唯一解, 所以能夠求出βij的值及近似解ρ*(x,t), 使得E取得最小值.
定理3假定ρ(x,t)是問(wèn)題(9)的準(zhǔn)確解,由定理2可知,存在ρ*(x,t)是ρ(x,t)的近似解.
證明由定理1,矩陣張量積性質(zhì)及定理2,很容易構(gòu)造ρ(x,t)的近似解ρ*(x,t).
假定ρ*(x,t)是ρ(x,t)的近似解,且存在自然數(shù)n, 常數(shù)βij(i=0,1,…,n;j=0,1,…,n), 勒讓德多項(xiàng)式P0(x),P1(x),…,Pn(x),P0(t),P1(t),…,Pn(t)使得當(dāng)n→∞時(shí),有
(21)
且
P(x,t)(I?M)β,
(22)
P(x,t)(M?I)β,
(23)
P(x,t)(M2?I)β,
(24)
ρ*(x,t)滿足初邊值條件
(25)
(26)
(27)
假定ρ*(x,t)是ρ(x,t)的近似解,則
(28)
將式(21)~(27)代入式(28),可得
(29)
此時(shí),方程(29)是權(quán)值β的方程組,通過(guò)優(yōu)化算法求出方程組的解,即β的值,可得ρ(x,t)的近似解ρ*(x,t),即是問(wèn)題(9)的近似解.
(30)
可得問(wèn)題(9)的離散形式
Hβ=T.
(31)
此時(shí),問(wèn)題(9)就轉(zhuǎn)化為求解方程組(31).
定理4方程(31)是可解的, 且滿足如下三種情況:
1) 若矩陣H為方陣且是可逆的,則β=H-1T;
證明定理的證明可參考矩陣論[21]中矩陣廣義逆理論及Guang-Bin Huang的文章[14-18]. 算法步驟:
1) 選擇兩個(gè)勒讓德多項(xiàng)式的乘積P(x)P(t)作為激活函數(shù),建立問(wèn)題(9)的神經(jīng)網(wǎng)絡(luò)模型,構(gòu)造近似解ρ*(x,t);
2) 將ρ*(x,t)代入問(wèn)題(9)得到方程(29);
3) 將方程(29)表示成矩陣形式(31);
這里以數(shù)值實(shí)驗(yàn)來(lái)驗(yàn)證勒讓德神經(jīng)網(wǎng)絡(luò)方法的有效性.
例1當(dāng)u=0.01,v=1,a=1時(shí),求對(duì)流擴(kuò)散方程
例2當(dāng)u=0.1,v=0.1,a=2時(shí),求對(duì)流擴(kuò)散方程
表 1 實(shí)驗(yàn)結(jié)果
由表 1 中數(shù)據(jù)可知,隨著隱層神經(jīng)元的增加,計(jì)算精度提高,但絕不是越多越好. 當(dāng)n=20時(shí),例 1 的計(jì)算精度為6.445 4×10-9,此時(shí)計(jì)算時(shí)間僅需要0.210 7 s. 當(dāng)n=10時(shí),例 2 的計(jì)算精度為6.376 1×10-5,此時(shí)計(jì)算時(shí)間僅需要0.037 9 s. 采用所提出的新方法,能夠由簡(jiǎn)單的網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu),用較快的速度達(dá)到一定的計(jì)算精度.
圖1 和圖2 表示隨著n值不同,例 1和例 2 的計(jì)算時(shí)間. 圖3 表示誤差與神經(jīng)元之間的關(guān)系. 圖4 表示計(jì)算時(shí)間與神經(jīng)元之間的關(guān)系.
圖1 例 1 所用時(shí)間Fig.1 Calculation time of example 1
圖2 例 2 所用時(shí)間Fig.2 Calculation time of example 2
圖3 n取不同值的計(jì)算誤差Fig.3 Error with different values of n
圖4 n取不同值時(shí)的平均計(jì)算時(shí)間Fig.4 Average calculation time with different values of n
本文主要介紹求解一維對(duì)流擴(kuò)散方程問(wèn)題的勒讓德神經(jīng)網(wǎng)絡(luò)方法,建立該問(wèn)題的神經(jīng)網(wǎng)絡(luò)模型,當(dāng)近似解滿足偏微分方程及初邊值條件時(shí),根據(jù)矩陣分析理論,可以采用基于廣義逆矩陣的IELM算法求解近似解中的權(quán)值,為一類常系數(shù)一維對(duì)流擴(kuò)散方程問(wèn)題的數(shù)值解提供一種新的方法. 數(shù)值實(shí)驗(yàn)結(jié)果表明,采用IELM算法求解神經(jīng)網(wǎng)絡(luò)模型中的權(quán)值參數(shù),當(dāng)樣本量給定時(shí),計(jì)算精度及速度只與網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)中隱層神經(jīng)元數(shù)量有關(guān),且只需少量的隱層神經(jīng)元,通過(guò)改變神經(jīng)元數(shù)量即可得到相應(yīng)的計(jì)算精度. 本文的實(shí)驗(yàn)結(jié)果可以為該方法在應(yīng)用時(shí)隱層神經(jīng)元的選擇提供參考,未來(lái)可以研究該算法的推廣能力,如隱層神經(jīng)元的最優(yōu)取值區(qū)間等.