黃鐘民 謝 臻 張易申 彭林欣 ,?,2)

* (廣西大學土木建筑工程學院,南寧 530004)

? (廣西大學廣西防災減災與工程安全重點實驗室,工程防災與結構安全教育部重點實驗室,南寧 530004)

引言

面內功能梯度材料薄板結構在土木工程、海洋工程等領域應用廣泛,該結構由功能梯度材料所組成(functionally graded material,FGM),其材料特性隨著空間位置的變化而表現(xiàn)為梯度性的變化[1].在現(xiàn)有研究中,于天崇等[2]研究了面內變剛度薄板在特定邊界下彎曲問題的Levy 解.朱竑禎等[3]研究了周邊固支圓形面內變剛度薄板軸對稱彎曲問題的級數(shù)解答.何建璋等[4]研究了面內變剛度矩形薄板自由振動問題的辛彈性解.以上的理論解答僅針對特定的功能梯度函數(shù)及特定邊界才成立,一般的情況下難以得出理論解答,而在數(shù)值解法上仍以有限元為主.Santare 和Lambros[5]發(fā)展了一種針對材料屬性為指數(shù)分布的梯度有限元求解格式.Kim 和Paulino[6]研究了梯度單元以及分層單元在不同荷載下的計算性能.黃立新等[7]基于分層法思想分析了功能梯度材料的平面應力問題.田云德和秦世倫[8]采用分層法研究了功能梯度厚板的熱應力問題.對于面內變剛度功能梯度薄板,采用分層法,薄板求解域采用有限元網(wǎng)格劃分,每個單元的材料參數(shù)為常數(shù),而其材料參數(shù)則根據(jù)功能梯度函數(shù)由單元內特定點進行計算.有限元網(wǎng)格劃分越密其計算結果越精確,而在實際計算中,越精細的網(wǎng)格會導致總體剛度矩陣規(guī)模巨大,需要耗費大量的計算機內存.無論采用何種數(shù)值方法,其最終目的均是求得面內變剛度薄板彎曲控制偏微分方程的近似解答,為進一步豐富該類研究,本文擬結合深度學習技術并發(fā)展求解該類問題的新解法.

在早期就有研究[9-10]將人工神經(jīng)網(wǎng)絡作為一類偏微分方程的求解器用于求解偏微分方程,但由于其對計算機計算能力的要求過高以及優(yōu)化算法中存在的問題,這一解法在當時并未得到很好的發(fā)展.而如今,自深度學習在計算機視覺、語音文字識別取得成功的應用后,深度學習技術也在各個學科領域加速發(fā)展.在力學領域,Weinan 和Yu[11]提出深度Ritz 法,該方法采用變分求解形式對偏微分方程進行求解.Raissi 等[12]提出了用于求解高階非線性偏微分方程的物理驅動的神經(jīng)網(wǎng)絡(physics-informed neural networks,PINNs).Sirignano 和Spiliopoulos[13]則提出求解高階微分方程的深度伽遼金法(deep galerkin method,DGM).Samaniego 等[14]建立了深度能量法并將其應用于求解彈性、超彈性等力學問題.瞿同明等[15]基于深度學習技術,研究了細觀力學中的顆粒本構關系.謝晨月等[16]發(fā)展了一種模擬湍流大渦的神經(jīng)網(wǎng)絡方法.劉宇翔等[17]基于卷積神經(jīng)網(wǎng)絡研究了無網(wǎng)格方法中影響域的優(yōu)化問題.郭宏偉和莊曉瑩[18]采用深度配點法以及深度能量法求解了薄板彎曲問題.陳豪龍和柳占立[19]基于數(shù)據(jù)驅動的神經(jīng)網(wǎng)絡模型求解了熱傳導反問題.

在上述研究中,神經(jīng)網(wǎng)絡解法[11-14]并不像有限元解法一樣可以輕松施加邊界條件,早期的研究采取根據(jù)邊界條件構造滿足偏微分方程特解試函數(shù)的形式來處理邊界條件,但采用該方法會使得簡支邊、自由邊試函數(shù)的表達式變得復雜,導致程序的實現(xiàn)較為復雜.近期的研究則采用罰函數(shù)的方法將邊界處的誤差納入神經(jīng)網(wǎng)絡的訓練誤差中,從而將原問題轉換為無約束優(yōu)化問題,在實際計算中,也會存在著由于邊界誤差項難以收斂而影響求解精度的情況[20].

同時由于彎曲剛度函數(shù)是面內坐標的連續(xù)函數(shù),面內變剛度薄板彎曲問題的控制方程為一包含了彎曲剛度導數(shù)項的復雜4 階偏微分方程,在實際計算中采用DGM 和PINN 等方法對其求解時,會存在由于彎曲剛度偏導數(shù)在域內不收斂而導致網(wǎng)絡擬合不佳的問題.

基于上述原因,本文針對薄板彎曲問題求解的特點,結合前面所述的兩種邊界處理方案,建立了一種針對面內變剛度薄板彎曲問題的非全連接前饋神經(jīng)網(wǎng)絡模型,該模型包含撓度網(wǎng)絡與彎矩網(wǎng)絡:撓度網(wǎng)絡用于預測薄板的撓度,彎矩網(wǎng)絡用于預測薄板的彎矩,進而將問題轉換為求解4 個二階偏微分方程組.在邊界條件的處理上,本文仍采用罰函數(shù)方法,不同之處在于本文模型的輸出為撓度、彎矩,因而可根據(jù)位移邊界條件對撓度網(wǎng)絡構造試函數(shù),根據(jù)廣義應力邊界條件對彎矩網(wǎng)絡構造試函數(shù),這使得本文模型對于常見的邊界條件的施加更為簡便,進而減小邊界誤差項帶來的影響,同時計算效率也得到提高.本文采用Pytorch 深度學習框架編寫求解程序,選取不同邊界條件的面內變剛度薄板算例,在Ubuntu Kylin 操作系統(tǒng)上進行計算,計算機的CPU配置為Intel(R) Core(TM) i7-8700 CPU @ 3.20GHz,8GB 內存,并將計算所得結果與理論解、有限元解進行對比分析,以驗證本文方法的有效性.

1 面內變剛度薄板彎曲理論

本文研究變厚度薄板或彈性模量參數(shù)在面內變化的薄板的彈性小變形彎曲,設薄板的厚度函數(shù)為h(x,y),材料的泊松比ν為常數(shù),彈性模量函數(shù)為E(x,y).

根據(jù)Kirchhoff 板理論基本假定,幾何方程為

物理方程為

其中n表示邊界的外法線方向,s表示邊界的切線方向.

將式(4)～ 式(6)代入平衡方程(3)即可得面內變剛度薄板彎曲偏微分控制方程

2 神經(jīng)網(wǎng)絡模型

2.1 神經(jīng)網(wǎng)絡模型的構建

本文方法并非直接設計網(wǎng)絡來求解方程(12),而是采用兩個神經(jīng)網(wǎng)絡模型來進行求解,如圖1 所示,將待求解的4 階偏微分控制方程轉換為求解4 個二階偏微分方程組,該解法本質上仍屬于強形式的求解方案.如果僅以撓度作為預測解,在試函數(shù)的構造上對于不同形狀的求解域以及簡支、自由邊界條件的構造會出現(xiàn)困難.本文采用撓度網(wǎng)絡預測薄板撓度{彎矩網(wǎng)絡預測薄板彎}矩,這樣的做法可以使得位移邊界條件由撓度網(wǎng)絡施加,廣義應力邊界條件由彎矩網(wǎng)絡施加.

圖1 本文神經(jīng)網(wǎng)絡模型示意圖Fig.1 The schematic diagram of neural network model in this paper

2.2 誤差函數(shù)

根據(jù)廣義應力?應變關系式(6)～ 式(8)可求得撓度二階偏導的彎矩表達式

誤差函數(shù)的構造是神經(jīng)網(wǎng)絡訓練的核心,由于本文方法引入兩個網(wǎng)絡進行計算,故在訓練中需要考慮兩者之間的耦合誤差.若采用無約束優(yōu)化方案,本文誤差函數(shù)主要根據(jù)撓度與彎矩網(wǎng)絡在邊界處的誤差、彎矩網(wǎng)絡在力平衡方程(3)中的誤差、預測的撓度與彎矩通過式(13)～ 式(15)建立的耦合誤差來構造.

采用均方誤差(mean square error,MSE)來衡量神經(jīng)網(wǎng)絡的擬合誤差,記撓度網(wǎng)絡與彎矩網(wǎng)絡的內部參數(shù)分別為 θw,θm,本文模型的誤差函數(shù)可構造為

如果采用構造試函數(shù)的形式使得邊界誤差強制滿足,則誤差函數(shù)無需計算邊界誤差

其中,向量x=(x,y)表示神經(jīng)網(wǎng)絡的輸入;n和s分別為邊界的法線、切線方向;Ω 表示求解域;Γ1,Γ2,Γ3分別為固支、簡支、自由邊界,? Ω=Γ1∪Γ2∪Γ3;為施加于邊界處的彎矩,Mn(x;θm),Qn(x;θm)分別為根據(jù)彎矩網(wǎng)絡輸出求得的彎矩、剪力;kp為網(wǎng)絡耦合系數(shù),取值范圍為1～ 1000,該系數(shù)的選取會影響撓度網(wǎng)絡與彎矩網(wǎng)絡之間的耦合效果;k1,k2,k3為邊界處的罰系數(shù),取值范圍為1～ 10 000.

式(22)中

本文的誤差函數(shù)表達式中包含撓度、彎矩對自變量的二階偏導數(shù)項,對于這些偏導項,一方面可以利用神經(jīng)網(wǎng)絡的輸出構造差分求解格式來近似求解,但采用該方案需要較大的計算量才能得到精確的計算結果;另一方面,基于計算圖的自動微分技術(automatic differentiation,AD)可以高效地處理神經(jīng)網(wǎng)絡對輸入變量求導過程,當前的深度學習框架如Tensoflow,Pytorch,MindsSpore 等均支持自動微分.本文基于Pytorch 提供的自動微分接口實現(xiàn)對上述偏導項及誤差函數(shù)梯度的計算.

在實際計算中,也可靈活采用混合邊界誤差的形式進行求解,如對于部分簡單的邊界條件構造特解,而對于復雜的邊界采用相應的無約束優(yōu)化方案.建立本文的誤差函數(shù)后,在每個訓練批次(epoch)中均需計算其梯度并結合誤差反向傳播算法更新網(wǎng)絡的內部參數(shù),關于該過程,Tang 和Yang[21]對其進行了詳細的討論.

2.3 學習率選取方案

學習率的選取可直接影響神經(jīng)網(wǎng)絡的訓練,目前神經(jīng)網(wǎng)絡學習率的選取仍帶有一定的經(jīng)驗性,但總體而言,在訓練初期選取較大的學習率可以加快誤差收斂速度,在訓練后期,此時神經(jīng)網(wǎng)絡模型已經(jīng)學習到相應的特征,此時往往需要降低學習率,以便對神經(jīng)網(wǎng)絡內部參數(shù)進行微調,使得誤差波動幅度不至于過大.經(jīng)過本文的實踐,本文學習率選取方案如下

其 中,t為訓練次數(shù).

2.4 算法流程

算法1.本文算法

3 算例

3.1 雙向面內變剛度圓形薄板軸對稱彎曲

圖2 所示受橫向均布荷載q(x,y)=?q0作用的周邊固支圓形薄板,其半徑R,ν=0.3,彎曲剛度函數(shù)沿半徑變化D(ρ)=D0e?mρ/a,其中m 為梯度系數(shù),q0,D0為常數(shù).該問題存在理論解[14].

圖2 圓形面內變剛度薄板Fig.2 Circular thin plate with in-plane stiffness gradient

本文選取梯度參數(shù)m分別為0,0.5,1,2 的情況進行計算.本算例在計算過程中僅需要施加位移邊界條件,設w(x,y;θw)為撓度網(wǎng)絡的輸出,考慮到本算例邊界條件較為簡單,構造撓度試函數(shù)為w?=w(x,y;θw)

本算例的撓度模型以及彎矩模型均采用具有6 層隱藏層,每層隱藏層具備30 個神經(jīng)元的網(wǎng)絡結構(記為6 × 30),采用x2作為激活函數(shù),kp=100,NΩ=640 ;采用Adam 優(yōu)化算法,對各個工況所采用的學習率方案均一致.本文算例在不同的梯度下上述基本參數(shù)不變,僅更改梯度系數(shù).為更加詳細地顯示誤差函數(shù)的變化情況,如無特殊說明,本文均取誤差的十進制對數(shù)作為等效誤差并繪制相應的訓練誤差曲線.

每個梯度參數(shù)下神經(jīng)網(wǎng)絡模型的計算誤差如圖3 所示,由此可見,m=0 時為剛度恒定的薄板,此時誤差函數(shù)收斂較快,相比其他工況其誤差最終的收斂值最小;隨著梯度系數(shù)的增大,訓練誤差最終的收斂值出現(xiàn)增大的趨勢.

表1 本文方法計算與理論解對比(無量綱)Table 1 Comparison of dimensionless calculated by neural network method and the theoretical solution

表1 本文方法計算與理論解對比(無量綱)Table 1 Comparison of dimensionless calculated by neural network method and the theoretical solution

注:本文誤差計算公式為 u 為本文方法計算結果,u?為理論解 (Note:The relative error in this paper is calculated by 100%,where u is the calculation results of this paper,u ? is the theoretical solution)

表2 本文方法計算與理論解對比(無量綱)Table 2 Comparison of dimensionless calculated by neural network method and the theoretical solution

表2 本文方法計算與理論解對比(無量綱)Table 2 Comparison of dimensionless calculated by neural network method and the theoretical solution

圖3 神經(jīng)網(wǎng)絡訓練誤差曲線Fig.3 The convergence curve of neural network training error

圖4 PINN 求解圓形面內變剛度薄板彎曲問題的訓練誤差收斂曲線(m=2)Fig.4 Training error convergence curve of PINN (m=2)

為了說明本文方法在求解面內變剛度薄板彎曲問題上的優(yōu)點,本算例也利用PINN 來求解其四階偏微分控制方程(12),采用隱藏層層數(shù)為6,每層隱藏層具有30 個神經(jīng)元的神經(jīng)網(wǎng)絡模型對工況m=2進行求解,激活函數(shù)為Tanh 函數(shù),訓練的數(shù)據(jù)點由求解域中隨機生成,數(shù)據(jù)點的產(chǎn)生有兩種方案:

方案(1)為在整個求解域中隨機生成訓練數(shù)據(jù)點,此時的數(shù)據(jù)點可在原點附近生成;

方案(2)為在求解域中隨機生成的數(shù)據(jù)點但離原點較遠.此時兩方案在訓練過程中的誤差收斂情況如圖4 所示,可見采用相同的模型,而生成的數(shù)據(jù)點不同則會導致模型的訓練出現(xiàn)不同的結果,雖然采用數(shù)據(jù)點生成方案(1)的模型訓練也收斂,但由于其誤差此時收斂于一個較大的值,得不到正確解.經(jīng)過本文分析,這主要是由于本算例的彎曲剛度函數(shù)D的二階導數(shù)在靠近原點區(qū)域出現(xiàn)“爆炸”式變化的原因,即剛度函數(shù)的二階偏導在原點處不收斂,在靠近原點處等的解答急劇增大,這會導致PINN 采用方案(1)訓練時,遇到靠近原點處的點,計算所得域內誤差突然增大,進而導致誤差訓練難以收斂.PINN 最初提出時并未考慮求解域內存在奇異點的情況,對于該情況,一般情況下可在生成的數(shù)據(jù)點中排除掉奇異點,但對于本算例中奇異點處被施予荷載的情況,如果不能很好地處理則會影響求解的精度.對此,本文認為可以弱化相應的偏微分控制方程再利用PINN 求解,也可參考本文思路,結合神經(jīng)網(wǎng)絡解法的特點,根據(jù)具體問題對原偏微分控制方程等效化處理.本文方法在求解時并非直接從方程(12)入手,而是通過求解一系列偏微分方程組來逼近真實解,避開了對彎曲剛度函數(shù)求偏導數(shù),故其求解僅與域內的剛度值有關,其適應性更強,對薄板彎曲問題的求解更具“魯棒性”.

3.2 單向面內變剛度方形薄板受非線性荷載作用

如圖5 所示邊長為a的方形薄板,厚度h,ν=0.3,1,2,3 邊固支,4 邊簡支,其彎曲剛度函數(shù)為為常數(shù),m為梯度系數(shù),受橫向非線性荷載作用,利用本文方法求解m分別為0,1,2 情況下的撓度、內力.

圖5 方形面內變剛度薄板Fig.5 Thin square plates with in-plane stiffness gradient

本算例選取6 × 30 的撓度網(wǎng)絡模型,5 × 50 的彎矩網(wǎng)絡模型,激活函數(shù)選擇Tanh 函數(shù),kp=100,NΩ=450,采用Adam 優(yōu)化算法.訓練誤差曲線如圖6 所示.設w(x,y;θw)為撓度網(wǎng)絡的輸出,彎矩網(wǎng)絡的輸出為M(y},考慮本算例?的位移及廣義應力邊界條件,撓度及彎矩試函數(shù)構造為

圖6 神經(jīng)網(wǎng)絡訓練誤差曲線Fig.6 Neural network training error-curve

如果將本算例的應力邊界條件也通過撓度網(wǎng)絡進行構造,則撓度表達式將復雜許多.

將本算例的計算結果與有限元解答對比,有限元計算中每個單元的彎曲剛度根據(jù)單元的形心坐標計算,采用50× 50 的矩形薄板非協(xié)調單元來對求解域進行離散,離散方案通過小片測試,對本文的3 種工況該網(wǎng)格離散方案均收斂.將本文與有限元計算結果的無量綱撓度、彎矩進行對比分析,無量綱計算公式為w?和M?為實際計算所得的撓度、彎矩值.

由撓度計算圖7 可知,在撓度的求解上,本文解法與有限元解法一致.由圖8 的彎矩對比圖可發(fā)現(xiàn),本文彎矩解與有限元解答基本吻合,而當梯度系數(shù)m=2 時,雖然本文解與有限元解在部分點上彎矩的相對誤差增大,但整體上解答吻合.

圖7 不同梯度參數(shù)下本文撓度計算結果與有限元對比(y=0.5 m)Fig.7 Comparison of dimensionless deflection calculation results in this paper with FEM when different gradient parameters (y=0.5 m)

圖8 不同梯度參數(shù)下本文彎矩 計算結果與有限元對比(y=0.5 m)Fig.8 Comparison of dimensionless bending moment calculation results of this paper with FEM when different gradient parameters(y=0.5 m)

3.3 三角形面內變剛度薄板受線性分布荷載作用

圖9 三角形面內變剛度薄板Fig.9 Thin triangular plate with in-plane stiffness gradient

圖10 三角形面內變剛度薄板沿軸線x=0 上的撓度 分布Fig.10 Dimensionless deflection variation of thin triangular plate with in-plane stiffness gradient along axis x=0

圖11 m=0.2 時三角形面內變剛度薄板彎矩的有限元計算結果Fig.11 Finite element calculation of bending moment of thin plate with triangular in-plane variable stiffness (m=0.2)

圖12 m=0.2 時三角形面內變剛度薄板彎矩的本文計算結果Fig.12 Neural network method calculation of bending moment of thin plate with triangular in-plane variable stiffness (m=0.2)

4 本文方法計算效率及其有限元對比分析

以T表示算例各個工況下的平均用時,T′表示各個工況下誤差函數(shù)的構建及其梯度計算的平均耗時,本文各算例求解時迭代所需的數(shù)據(jù)點數(shù)及所需的平均時間、內存如表3 所示,可看出本文誤差函數(shù)的構建及其梯度的求解占據(jù)了神經(jīng)網(wǎng)絡訓練總時長的70%左右.

表3 本文各算例求解所需的數(shù)據(jù)點數(shù)、內存、時間Table 3 The number of training data points,computational memory and computing time of numerical examples in this paper

本文的有限元求解程序采用python 語言進行編寫,選用每個節(jié)點有3 個自由度的薄板彎曲單元計算,剛度矩陣以 compressed sparse column (CSC)格式的稀疏矩陣存儲,利用科學計算庫scipy 中的線性求解器求解剛度方程.根據(jù)有限元解答的最小撓度判斷有限元解答是否收斂,算例2 和算例3 中各個工況下有限元計算收斂時的節(jié)點數(shù)目如表4 所示,本文有限元計算所需的節(jié)點數(shù)與計算所需內存的關系如圖13 所示.

圖13 有限元計算所需內存與節(jié)點數(shù)的關系(薄板單元)Fig.13 The relationship between the number of nodes and the memory needed in finite element calculation using thin plate bending element

表4 算例2、算例3 的有限元求解收斂所需節(jié)點數(shù)目Table 4 The number of nodes needed for the convergence of the finite element solution of numerical example 2 and 3

可以發(fā)現(xiàn)有限元解答收斂時所需的節(jié)點(網(wǎng)格)數(shù)目隨著梯度系數(shù)的增大而增大,節(jié)點數(shù)與所需內存并非呈現(xiàn)線性關系,節(jié)點數(shù)的增多會導致所需計算內存的急劇增大,而神經(jīng)網(wǎng)絡方法在梯度系數(shù)增大時,在單次迭代時仍可以較少的數(shù)據(jù)點進行迭代,這使得本文方法在求解時所需內存較小.在時間上,以算例3 中梯度系數(shù)m=0.5 為例,此時有限元求解收斂時,所需內存為291.5 MiB,剛度方程從組裝到求解用時1.798 s,可以看出本文解法的求解速度明顯慢于有限元,一方面是由于本文方法求解的是強形式的偏微分控制方程,與求解弱形式的方程相比,往往需要更多的迭代次數(shù),另一方面,神經(jīng)網(wǎng)絡方法是一類以數(shù)據(jù)為驅動的解法,在求解過程中需要往復迭代,而對于本文研究的線性問題,有限元僅需求解一次剛度方程.

5 網(wǎng)絡耦合系數(shù) kp 、隱藏層層數(shù)與神經(jīng)元個數(shù)、激活函數(shù)對網(wǎng)絡收斂性的影響分析

5.1 網(wǎng)絡耦合系數(shù) kp

本文模型由于引入兩個神經(jīng)網(wǎng)絡模型對面內變剛度薄板彎曲問題進行計算,與采用單網(wǎng)絡的神經(jīng)網(wǎng)絡方法相比,在訓練過程中需考慮網(wǎng)絡之間的耦合誤差,本文引入網(wǎng)絡耦合系數(shù)kp以加強網(wǎng)絡之間的耦合,kp的取值過小會導致網(wǎng)絡耦合不佳,過大則可能會導致訓練不收斂,因此需要對kp的取值進行討論,選取kp=1,10,100,1000 四種情況對算例2 進行計算分析(梯度系數(shù)m=2),以網(wǎng)絡訓練時薄板的變形能、外力功變化情況衡量不同kp下本文模型的訓練效果.

薄板變形能計算公式為

k為根據(jù)彎矩網(wǎng)絡的輸出由式(8)和式(9)計算所得.

外力所作實功為

其則根據(jù)撓度網(wǎng)絡的輸出值計算.

采用高斯積分,由上述公式計算神經(jīng)網(wǎng)絡訓練時薄板的變形能、應變能的變化情況,根據(jù)撓度網(wǎng)絡輸出計算外力實功Wext,根據(jù)彎矩網(wǎng)絡輸出計算薄板變形能Wint.當神經(jīng)網(wǎng)絡解答收斂時,根據(jù)能量原理,其外力實功應與變形能相等.此時各工況下的變形能與外力實功在神經(jīng)網(wǎng)絡訓練過程中的變化情況如圖14 所示.可以看出,隨著訓練的進行,外力實功與薄板的變形能逐漸收斂,在各個工況下,變形能的收斂速度相比于外力功要慢,進而說明撓度網(wǎng)絡收斂速度較快.當kp=1 時,此時薄板變形能與外力功收斂時,兩者相差較大;而隨著kp的增大,神經(jīng)網(wǎng)絡計算的變形能收斂加快,此時當變形能與外力功收斂時,兩者之間的差別減小,而當kp=1000 時,其對訓練的影響與kp=100 時差別不大,在實際計算中,該系數(shù)不能過大,過大會導致神經(jīng)網(wǎng)絡學習不到正確的特征,進而導致訓練不收斂.

圖14 神經(jīng)網(wǎng)絡訓練時不同 kp 下薄板變形能與外力功的變化情況Fig.14 Changes in the deformation energy and the external force work during the training process with different kp

5.2 隱藏層層數(shù)與神經(jīng)元個數(shù)

隱藏層層數(shù)與每層的神經(jīng)元個數(shù)的選取均會影響到神經(jīng)網(wǎng)絡的訓練效率,不失一般性,在其他計算參數(shù)不變的情況下,本節(jié)選取算例2 中m=2 的情況分別討論隱藏層層數(shù)、神經(jīng)元個數(shù)的改變對計算過程中誤差函數(shù)收斂情況的影響,計算結果如圖15所示.

圖15 隱藏層數(shù)、每層神經(jīng)元個數(shù)對計算誤差的收斂影響對比Fig.15 Comparison of the effects of different number of hidden layers and neurons on the convergence of loss function

可以發(fā)現(xiàn),在一定程度內誤差函數(shù)的收斂速度隨著網(wǎng)絡層數(shù)、隱藏層每層神經(jīng)元數(shù)的增大而加快,當兩者都增大到一定程度時,此時誤差的收斂速度會趨向于一個“飽和”狀態(tài).適當?shù)卦黾訉訑?shù)或神經(jīng)元的個數(shù)有利于誤差函數(shù)的收斂,目前在利用神經(jīng)網(wǎng)絡方法求解偏微分方程的研究中,在兩者的選取上仍帶有一定的經(jīng)驗性.結合上述的計算結果,考慮到計算機硬件能力的限制,本文算例的隱藏層層數(shù)在4～ 6 層間選取,每層神經(jīng)元數(shù)目在30～ 50 之間選取.

5.3 激活函數(shù)

在神經(jīng)網(wǎng)絡的訓練過程中,有多種非線性激活函數(shù)可以選擇,非線性的激活函數(shù)是使神經(jīng)網(wǎng)絡具備擬合非線性函數(shù)能力的重要原因,常用的激活函數(shù)有Tanh,ReLU,Sigmoid,Swish 函數(shù)等.為討論激活函數(shù)對神經(jīng)網(wǎng)絡訓練的影響,其余計算參數(shù)不變,本節(jié)選取x2作為激活函數(shù)并與Tanh,Swish 函數(shù)進行對比分析,對算例1、算例3 進行計算,訓練過程中的誤差走向如圖16 所示.在誤差函數(shù)的收斂速度上,x2優(yōu)于Tanh,Swish 函數(shù),同時由于其函數(shù)形式較為簡單,在自動微分計算中其所需計算量較小.本文經(jīng)驗表明,在薄板彎曲問題的求解上,采用多項式函 數(shù)x2作為激活函數(shù)可加快神經(jīng)網(wǎng)絡的收斂.

圖16 不同的激活函數(shù)對計算誤差的收斂影響對比Fig.16 Comparison of the effects of different activation functions on the convergence of loss function

6 結論與展望

本文基于深度學習技術與強形式的求解方案建立了一種直角坐標下求解面內變剛度薄板彎曲問題的神經(jīng)網(wǎng)絡方法,通過幾個算例分析,得出以下結論:

(1)本文解答與理論解、有限元解吻合,證明了本文方法在求解面內變剛度薄板彎曲問題上的正確性,本文的神經(jīng)網(wǎng)絡模型不需要對彎曲剛度函數(shù)求偏導,其適應性更強,同時在薄板的位移邊界條件、應力邊界條件的施加上較為方便.

(2)本文方法屬于強形式的數(shù)值解法,其計算所得結果具備連續(xù)性與可導性.理論上,本文方法可以求解彈性模量以及厚度在面內連續(xù)變化的薄板彎曲問題.彎矩網(wǎng)絡的求解受到梯度系數(shù)的影響,在梯度變化較大處彎矩網(wǎng)絡的求解精度受到一定的影響,但對撓度網(wǎng)絡的求解精度影響不大.

(3)由于神經(jīng)網(wǎng)絡方法為迭代類解法,本文方法在薄板線性彎曲問題求解上的收斂速度較有限元慢,但其計算所需內存較小.通過本文的模型結構可看出,神經(jīng)網(wǎng)絡方法具備相當大的靈活性,根據(jù)這一特點,可進一步發(fā)展求解面內變剛度功能梯度薄板非線性彎曲問題的神經(jīng)網(wǎng)絡方法,神經(jīng)網(wǎng)絡方法在非線性問題的求解中具備潛在優(yōu)勢.

(4)本文模型仍存在優(yōu)化空間:一方面在本文模型的訓練過程中,誤差函數(shù)及其梯度的計算在整個訓練過程中占據(jù)大部分的時間,可以考慮優(yōu)化誤差函數(shù)的構建過程,如引入有限元中形函數(shù)的思想對算法進行優(yōu)化;另一方面為了使得撓度、彎矩網(wǎng)絡具備較強的表達能力,本文模型采用了兩個具有獨立參數(shù)的網(wǎng)絡進行計算,這導致了本文模型的訓練參數(shù)較多,為此后續(xù)優(yōu)化中可將本文的兩個網(wǎng)絡合并為一個網(wǎng)絡(2 個輸入,4 個輸出),對網(wǎng)絡結構進行改良,以減少訓練參數(shù).