瞿同明 馮云田 ,2) 王孟琦 趙婷婷 狄少丞

* (斯旺西大學(xué)辛克維奇工程計(jì)算中心,英國(guó)斯旺西 SA1 8EP)

? (太原理工大學(xué)機(jī)械與運(yùn)載工程學(xué)院,太原 030024)

** (哈爾濱工程大學(xué)船舶工程學(xué)院,哈爾濱 150001)

引言

顆粒材料廣泛存在于各類(lèi)自然環(huán)境和工程活動(dòng)中,人類(lèi)幾乎所有的基礎(chǔ)設(shè)施都建造于巖土材料這類(lèi)典型的顆粒材料上,理解和預(yù)測(cè)此類(lèi)材料在外載作用下的力學(xué)響應(yīng)具有重要的意義.巖土材料通常由粒徑不一,形狀各異的礦物顆粒以及顆粒之間的空隙組成.顆粒間的互相滑動(dòng)會(huì)引起不可恢復(fù)的塑性變形,因而,顆粒材料在剪切作用下極易發(fā)生體積收縮或膨脹(剪脹性或者剪縮性).在外載荷作用下,顆粒材料主要通過(guò)顆粒之間的接觸力來(lái)傳遞和平衡外力,內(nèi)部常常形成高度復(fù)雜的力鏈網(wǎng)絡(luò)結(jié)構(gòu)(非均勻性)[1].這種力鏈網(wǎng)絡(luò)的穩(wěn)定性來(lái)源于顆粒間的摩擦強(qiáng)度和形狀引起的互鎖效應(yīng),宏觀上體現(xiàn)出顆粒材料的抗剪強(qiáng)度和剛度隨著壓應(yīng)力的增加而增大(壓硬性和摩擦性).除此以外,顆粒材料的力學(xué)響應(yīng)還和應(yīng)力路徑,材料初始各向異性等因素相關(guān)[2-3].這些復(fù)雜的特性使得預(yù)測(cè)顆粒材料的應(yīng)力?應(yīng)變關(guān)系成為一個(gè)極具挑戰(zhàn)的科學(xué)問(wèn)題[4-7].

在過(guò)去幾十年間,基于唯象假設(shè)建立應(yīng)力?應(yīng)變?cè)隽恐g的解析關(guān)系一直是巖土等顆粒材料本構(gòu)問(wèn)題中最主流的研究方式,然而發(fā)展一個(gè)統(tǒng)一的應(yīng)力?應(yīng)變理論模型仍然存在許多困難.并且,許多現(xiàn)有的本構(gòu)模型數(shù)學(xué)形式復(fù)雜,需要借助精細(xì)的室內(nèi)實(shí)驗(yàn)來(lái)標(biāo)定數(shù)十個(gè)自由參數(shù),對(duì)一般工程師而言學(xué)習(xí)和使用成本過(guò)高.由于理論發(fā)展和實(shí)踐應(yīng)用的不相適應(yīng),盡管各類(lèi)高級(jí)模型層出不窮,大量的巖土從業(yè)人員仍然使用基于線彈性理想塑性的摩爾?庫(kù)倫等簡(jiǎn)單本構(gòu)模型開(kāi)展工程問(wèn)題分析.

除唯象模型以外,一些學(xué)者也在尋找新的研究范式來(lái)研究本構(gòu)問(wèn)題.在層級(jí)多尺度模擬方法中,有限元等連續(xù)計(jì)算模型不再使用唯象的本構(gòu)模型,而是通過(guò)低尺度的離散元模擬來(lái)得到材料的應(yīng)力?應(yīng)變關(guān)系[8].在有限元的計(jì)算過(guò)程中,每個(gè)高斯點(diǎn)的位移梯度或者應(yīng)變信息實(shí)時(shí)傳遞給離散元作為邊界條件;離散元計(jì)算后,又將對(duì)應(yīng)的應(yīng)力響應(yīng)或者切向剛度矩陣實(shí)時(shí)傳遞給有限元.然而,此種多尺度計(jì)算方法中每個(gè)高斯點(diǎn)均需要有對(duì)應(yīng)的離散元模型,計(jì)算成本過(guò)高,很難真正應(yīng)用到工程尺度問(wèn)題的研究中.

另外一種代表性的方案是采用神經(jīng)網(wǎng)絡(luò)從大量實(shí)驗(yàn)數(shù)據(jù)中訓(xùn)練代理本構(gòu)模型.這種模型無(wú)需標(biāo)定自由參數(shù),在有限元等宏觀計(jì)算方法中使用方便,并且隨著新數(shù)據(jù)的不斷增加,其預(yù)測(cè)精度可以不斷提高,相比傳統(tǒng)本構(gòu)模型有獨(dú)特的優(yōu)勢(shì).使用神經(jīng)網(wǎng)絡(luò)來(lái)訓(xùn)練巖土材料本構(gòu)模型的思想可以追溯到20 世紀(jì)90 年代[9-12].近年來(lái),隨著深度學(xué)習(xí)技術(shù)的突破性進(jìn)展,基于神經(jīng)網(wǎng)絡(luò)的本構(gòu)模型再次得到廣泛的重視.Wang 等[13-15]結(jié)合博弈論和強(qiáng)化學(xué)習(xí)模型研究了顆粒材料的界面張力和位移的本構(gòu)關(guān)系以及應(yīng)力?應(yīng)變關(guān)系;Zhang 等[16]使用長(zhǎng)短期記憶模型(LSTM)預(yù)測(cè)了顆粒材料的循環(huán)剪切行為;Karapiperis 等[17]在有限元模型中實(shí)施了基于離散元數(shù)據(jù)訓(xùn)練得到的神經(jīng)網(wǎng)絡(luò)本構(gòu)模型.Qu 等[18]基于離散元數(shù)據(jù)探索了3 類(lèi)不同的顆粒材料本構(gòu)訓(xùn)練策略.盡管深度學(xué)習(xí)模型已經(jīng)被證明能有效捕捉特定工況下的顆粒材料力學(xué)響應(yīng),但是要將此類(lèi)純粹的數(shù)據(jù)驅(qū)動(dòng)模型發(fā)展成為可靠可用的,適用于多種加載工況的本構(gòu)模型仍有較遠(yuǎn)的距離.

解析模型和數(shù)據(jù)模型作為兩種不同的科學(xué)研究范式具有不同的優(yōu)缺點(diǎn).解析模型基于一些假設(shè)條件,利用理論方法構(gòu)造整個(gè)預(yù)測(cè)模型.此種范式的優(yōu)點(diǎn)是模型在自由參數(shù)得到標(biāo)定后具有一定的預(yù)測(cè)能力,但是受限于假設(shè)條件的約束,面對(duì)復(fù)雜荷載工況時(shí)往往預(yù)測(cè)與實(shí)際存在不少出入.數(shù)據(jù)模型不需要任何假設(shè),完全利用經(jīng)驗(yàn)數(shù)據(jù)來(lái)構(gòu)造預(yù)測(cè)模型.然而一個(gè)純粹的數(shù)據(jù)驅(qū)動(dòng)模型需要大量且覆蓋全面的訓(xùn)練樣本,對(duì)于本構(gòu)問(wèn)題來(lái)說(shuō),獲取大量的材料應(yīng)力?應(yīng)變數(shù)據(jù)并不容易.另外,當(dāng)數(shù)據(jù)模型遇到超出經(jīng)驗(yàn)樣本以外的情況時(shí),往往預(yù)測(cè)能力有限.本文擬嘗試結(jié)合理論模型和數(shù)據(jù)驅(qū)動(dòng)模型兩種不同的研究范式,通過(guò)將客觀的力學(xué)規(guī)律與應(yīng)力?應(yīng)變數(shù)據(jù)經(jīng)驗(yàn)結(jié)合起來(lái),發(fā)展顆粒材料的數(shù)據(jù)驅(qū)動(dòng)本構(gòu)模型訓(xùn)練方法.

本文首先基于Vogit 假設(shè)推導(dǎo)小應(yīng)變條件下顆粒材料的應(yīng)力?應(yīng)變解析關(guān)系,利用此關(guān)系識(shí)別出一組與顆粒材料本構(gòu)響應(yīng)相關(guān)的重要變量.通過(guò)有向圖的方式,將所有與本構(gòu)響應(yīng)相關(guān)的重要變量包含到從應(yīng)變到應(yīng)力的預(yù)測(cè)中,變量間的映射關(guān)系通過(guò)時(shí)序的深度學(xué)習(xí)模型來(lái)描述.采用離散元常規(guī)三軸以及真三軸數(shù)值試驗(yàn)生成顆粒材料應(yīng)力?應(yīng)變數(shù)據(jù)用于訓(xùn)練,驗(yàn)證和測(cè)試本文提出的深度學(xué)習(xí)模型.

1 基于深度學(xué)習(xí)的數(shù)據(jù)驅(qū)動(dòng)本構(gòu)模型

1.1 基本原理

深度學(xué)習(xí)是機(jī)器學(xué)習(xí)算法中的一類(lèi),主要指基于多層人工神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)的代理模型方法.數(shù)學(xué)形式上,神經(jīng)網(wǎng)絡(luò)是由線性的矩陣運(yùn)算和非線性的激活函數(shù)構(gòu)成的復(fù)合函數(shù),這類(lèi)代理函數(shù)具有強(qiáng)大的表示能力,已經(jīng)被證明能夠以任意的精度來(lái)估計(jì)任何一個(gè)定義在實(shí)數(shù)空間中的有界閉集函數(shù)[19-20].

一個(gè)問(wèn)題是否適合深度學(xué)習(xí)方法,主要有3 個(gè)基本前提:(1)問(wèn)題本身包含某種確定的模式或者規(guī)律;(2)無(wú)法在數(shù)學(xué)上找到精確的解析方法來(lái)描述這種模式;(3)有包含這個(gè)問(wèn)題所有自變量和因變量的完備數(shù)據(jù).對(duì)于顆粒材料而言,本構(gòu)關(guān)系由材料自身物理性質(zhì)所決定,這種關(guān)系必然是受物理材料規(guī)律支配的;并且,已有的大量研究表明:建立一個(gè)精確、統(tǒng)一的解析方程來(lái)描述這種關(guān)系極為困難;再者,通過(guò)各類(lèi)顆粒材料加載試驗(yàn),我們可以得到反映本構(gòu)關(guān)系的大量應(yīng)力?應(yīng)變數(shù)據(jù).這些因素構(gòu)成了發(fā)展基于深度學(xué)習(xí)框架的數(shù)據(jù)驅(qū)動(dòng)材料本構(gòu)關(guān)系的前提.

一般而言,在準(zhǔn)靜態(tài)變形情況下,材料內(nèi)部一點(diǎn)的應(yīng)力張量σij和應(yīng)變張量εmn可以用6 個(gè)獨(dú)立的分量來(lái)表示,即

如圖1 所示,基于深度學(xué)習(xí)的本構(gòu)模型即是利用神經(jīng)網(wǎng)絡(luò)來(lái)建立材料應(yīng)力和應(yīng)變之間的映射關(guān)系.通常輸入?yún)?shù)為應(yīng)變張量,輸出參數(shù)為應(yīng)力張量,中間隱藏層的激活函數(shù)可按如下方式計(jì)算

圖1 基于深度學(xué)習(xí)的本構(gòu)模型示意圖Fig.1 Diagram of deep learning-based constitutive models

假設(shè)第j層為最后一個(gè)隱藏層,第(j+1)層為輸出層,那么應(yīng)力的預(yù)測(cè)值可按照如下公式計(jì)算

其中,e為最后一個(gè)隱藏層的單元數(shù)量,W和b分別為神經(jīng)網(wǎng)絡(luò)的權(quán)重和偏置系數(shù).

在一個(gè)確定的神經(jīng)網(wǎng)絡(luò)架構(gòu)下尋找描述對(duì)應(yīng)映射關(guān)系的網(wǎng)絡(luò)權(quán)重和偏置的過(guò)程稱為訓(xùn)練或者學(xué)習(xí).如圖2 所示,一個(gè)神經(jīng)網(wǎng)絡(luò)的訓(xùn)練過(guò)程主要包括前向傳播、反向傳播和迭代更新網(wǎng)絡(luò)參數(shù)(權(quán)重和偏置).在訓(xùn)練前,網(wǎng)絡(luò)參數(shù)通常是隨機(jī)初始化的,此時(shí)通過(guò)前向傳播計(jì)算得到的預(yù)測(cè)值必然是偏離實(shí)際值的.這個(gè)預(yù)測(cè)值和實(shí)際值的差別通常用損失函數(shù)L來(lái)衡量.訓(xùn)練一個(gè)可靠的預(yù)測(cè)模型的過(guò)程可以轉(zhuǎn)化為一個(gè)以最小化損失函數(shù)L為目標(biāo)的優(yōu)化問(wèn)題,即

圖2 人工神經(jīng)網(wǎng)絡(luò)的訓(xùn)練過(guò)程Fig.2 Basic procedures of training artificial neural networks

其中,argmin 是指使得目標(biāo)函數(shù)取最小值時(shí)的自變量值.

常規(guī)神經(jīng)網(wǎng)絡(luò)一般通過(guò)基于一階梯度的優(yōu)化算法,如自適應(yīng)矩估計(jì)算法(Adam)等來(lái)不斷迭代,從而找到使得L達(dá)到極小值時(shí)的網(wǎng)絡(luò)參數(shù).現(xiàn)代神經(jīng)網(wǎng)絡(luò)中,這些梯度信息都使用反向傳播(BP)算法通過(guò)自動(dòng)微分和鏈?zhǔn)椒▌t計(jì)算得到.對(duì)于一個(gè)已經(jīng)訓(xùn)練好的模型,也即網(wǎng)絡(luò)中的權(quán)重和偏置參數(shù)能很好地描述映射關(guān)系時(shí),我們可以通過(guò)向前傳播,計(jì)算得到較為可靠的應(yīng)力預(yù)測(cè)值.

1.2 材料本構(gòu)模型的時(shí)序本質(zhì)與深度學(xué)習(xí)模型的選用

顆粒材料應(yīng)力響應(yīng)具有明顯的應(yīng)變歷史依賴和應(yīng)變路徑依賴特征,一個(gè)能反映復(fù)雜應(yīng)力?應(yīng)變路徑條件下的材料本構(gòu)響應(yīng)模型,必須具備“記憶”應(yīng)力或者應(yīng)變歷史的能力[21],因此本構(gòu)模型本質(zhì)上是一種時(shí)序模型或者過(guò)程模型.在當(dāng)前深度學(xué)習(xí)模型中,LSTM 和GRU 等循環(huán)神經(jīng)網(wǎng)絡(luò)(RNN)模型主要適用于包含時(shí)間序列的預(yù)測(cè)問(wèn)題.由于幾種RNN 模型預(yù)測(cè)能力近似,同等網(wǎng)絡(luò)規(guī)模條件下,GRU 相比LSTM 需要更少的訓(xùn)練參數(shù),因此本文主要采用GRU 作為基本的深度學(xué)習(xí)模型來(lái)訓(xùn)練材料的本構(gòu)行為.

2 基于Vogit 均質(zhì)化假設(shè)的顆粒材料小應(yīng)變應(yīng)力?應(yīng)變關(guān)系

本文嘗試引入顆粒物質(zhì)力學(xué)理論來(lái)指導(dǎo)深度學(xué)習(xí),以期降低深度學(xué)習(xí)在訓(xùn)練階段對(duì)數(shù)據(jù)量的需求,提高模型在預(yù)測(cè)階段的泛化能力.對(duì)于顆粒材料而言,目前主要有兩種唯象本構(gòu)模型:連續(xù)模型和細(xì)觀力學(xué)模型.基于連續(xù)假設(shè)的彈塑性模型采用屈服面來(lái)描述材料的塑性,借助流動(dòng)法則來(lái)描述塑性的發(fā)展.細(xì)觀力學(xué)模型雖然不再需要假想的屈服面,但仍需要通過(guò)預(yù)測(cè)細(xì)觀組構(gòu)的演化來(lái)描述應(yīng)力?應(yīng)變關(guān)系.不管是基于連續(xù)假設(shè)還是細(xì)觀力學(xué)的唯象模型,在預(yù)測(cè)材料隨外部應(yīng)變引起的應(yīng)力響應(yīng)時(shí),常常需要使用一些從試驗(yàn)數(shù)據(jù)中擬合得來(lái)的參數(shù)(如硬化參數(shù)等)來(lái)描述主應(yīng)力?主應(yīng)變關(guān)系.為避免使用此類(lèi)擬合參數(shù),在準(zhǔn)靜態(tài)加載條件下,本研究將時(shí)序的應(yīng)力?應(yīng)變過(guò)程簡(jiǎn)化成單個(gè)時(shí)刻的應(yīng)力?應(yīng)變響應(yīng)來(lái)考慮,即小應(yīng)變條件下的顆粒材料應(yīng)力?應(yīng)變關(guān)系.

散粒體在小變形情況下,顆粒體系內(nèi)的彈性勢(shì)能U可按照如下方式計(jì)算

當(dāng)顆粒接觸為線性彈性接觸時(shí),儲(chǔ)存在接觸k內(nèi)的彈性勢(shì)能Uk為

式中,kn和ks為顆粒之間的法向和切向剛度.

如圖3 所示,通過(guò)將顆粒位移等效為相應(yīng)連續(xù)體內(nèi)對(duì)應(yīng)點(diǎn)處的位移,則顆粒相對(duì)位移與對(duì)應(yīng)連續(xù)體內(nèi)的局部應(yīng)變(顆粒A 與B 之間的等效應(yīng)變)可表示為如下關(guān)系

圖3 顆粒接觸與接觸位移Fig.3 A contact between particles and contact displacements

根據(jù)彈性力學(xué)理論,連續(xù)體的應(yīng)力張量可通過(guò)應(yīng)變能密度u對(duì)應(yīng)變張量求取偏導(dǎo)得到,即

假設(shè)顆粒試樣內(nèi)的變形是統(tǒng)計(jì)均勻的(Vogit’s 假設(shè)),即在對(duì)應(yīng)連續(xù)體內(nèi)任意位置處的整體應(yīng)變 εij等于材料中沿著顆粒接觸方向的局部應(yīng)變

基于式(13),假設(shè)散粒體和對(duì)應(yīng)連續(xù)體在小應(yīng)變情況下的彈性勢(shì)能相等,則式(12)可表示為

其中,V為顆粒試樣的體積.

將應(yīng)力張量對(duì)應(yīng)變張量求偏導(dǎo),可得到對(duì)應(yīng)的彈性剛度矩陣為

式中,δin為克羅內(nèi)克函數(shù),當(dāng)i=n時(shí),δin=1;當(dāng)i≠n時(shí),δin=0.若令

式(17) 表明Cijmn由顆粒間的法向接觸剛度kn、切向剛度ks以及兩個(gè)單位體積內(nèi)的組構(gòu)張量aijmn和bijmn構(gòu)成.在顆粒試樣承受外界荷載作用時(shí),假設(shè)顆粒未發(fā)生破碎或其他材料劣化,則kn和ks保持不變;隨著試樣在宏觀上發(fā)生不可恢復(fù)的塑性變形(微觀上表現(xiàn)為顆粒間發(fā)生相互滑移)時(shí),組構(gòu)張量aijmn和bijmn隨之不斷演化.

小應(yīng)變情況下的顆粒材料的宏觀應(yīng)力?應(yīng)變關(guān)系可表示為當(dāng)切應(yīng)力和切應(yīng)變?yōu)? 時(shí),3 個(gè)主應(yīng)力主應(yīng)變的增量關(guān)系展開(kāi)后為

其中,C1122=C2211,C1133=C3311,C2233=C3322.因而彈性剛度矩陣中共有6 個(gè)獨(dú)立分量描述主應(yīng)力與主應(yīng)變的關(guān)系.對(duì)于準(zhǔn)靜態(tài)加載而言,式(19)的關(guān)系在試樣加載過(guò)程中總是成立,但隨著顆粒材料組構(gòu)的不斷演化,彈性剛度矩陣Cijmn中各個(gè)分量的值也會(huì)隨之變化.基于細(xì)觀力學(xué)的唯象模型通常需要借助流動(dòng)法則和硬化函數(shù)等來(lái)擬合得到試樣內(nèi)部組構(gòu)的演化過(guò)程,本文中,這些微觀結(jié)構(gòu)的演化過(guò)程則直接由深度學(xué)習(xí)來(lái)預(yù)測(cè)得到.

3 力學(xué)知識(shí)指導(dǎo)下的深度學(xué)習(xí)策略

由于顆粒材料內(nèi)部描述細(xì)觀結(jié)構(gòu)的組構(gòu)張量屬于內(nèi)變量,在有限元等連續(xù)介質(zhì)方法中計(jì)算應(yīng)力?應(yīng)變關(guān)系時(shí),此類(lèi)內(nèi)變量不能直接作為輸入.本文嘗試引入有向圖理論,將微觀結(jié)構(gòu)的演化過(guò)程包含在應(yīng)力?應(yīng)變關(guān)系的描述中.

有向圖是圖論中用于描述物與物之間聯(lián)系的一種方法,主要由結(jié)點(diǎn)和含有方向的邊(連接兩個(gè)結(jié)點(diǎn))構(gòu)成.結(jié)點(diǎn)表示某一個(gè)物理量或者對(duì)象,邊線表示兩個(gè)物理量之間存在某種聯(lián)系并且這種聯(lián)系是從源結(jié)點(diǎn)(尾部)指向目標(biāo)結(jié)點(diǎn)(頭部)的.根據(jù)有向圖理論和小應(yīng)變條件下得到的顆粒材料應(yīng)力?應(yīng)變關(guān)系,本文設(shè)計(jì)的本構(gòu)關(guān)系訓(xùn)練策略如圖4 所示.

圖4 基于有向圖包含組構(gòu)演化的本構(gòu)訓(xùn)練方式Fig.4 A directed graph-based constitutive training approach incorporating fabric evolution

圖4 中,初始結(jié)點(diǎn)為主應(yīng)變,作為輸入?yún)?shù);目標(biāo)結(jié)點(diǎn)為主應(yīng)力,作為輸出參數(shù);中間節(jié)點(diǎn)為6 個(gè)相互獨(dú)立的彈性剛度張量分量.每條邊之間的映射關(guān)系通過(guò)深度神經(jīng)網(wǎng)絡(luò)來(lái)描述.將該有向圖展開(kāi)后,可用兩條子神經(jīng)網(wǎng)絡(luò)來(lái)描述整個(gè)有向圖結(jié)構(gòu):

在網(wǎng)絡(luò)學(xué)習(xí)階段,兩條子網(wǎng)絡(luò)單獨(dú)訓(xùn)練;在使用階段,兩條子網(wǎng)絡(luò)合并成一條完整的信息流:子網(wǎng)絡(luò)1 的輸出連同3 個(gè)主應(yīng)變序列一起作為子網(wǎng)絡(luò)2 的輸入,共同預(yù)測(cè)材料的應(yīng)力響應(yīng).因此,基于有向圖的神經(jīng)網(wǎng)絡(luò)模型雖然充分利用了顆粒材料內(nèi)部的微觀結(jié)構(gòu)演化信息,但在使用過(guò)程中并不違背本構(gòu)模型僅僅將主應(yīng)變和主應(yīng)力作為最終輸入和輸出信息的基本原則.

4 數(shù)據(jù)準(zhǔn)備與深度學(xué)習(xí)模型的實(shí)施

4.1 基于離散元的數(shù)值三軸試驗(yàn)與模擬工況

細(xì)觀尺度的離散元仿真已經(jīng)被證明能夠有效模擬顆粒材料的本構(gòu)行為[22],并且在層級(jí)多尺度計(jì)算模型中,離散元模擬直接代替了傳統(tǒng)的唯象本構(gòu)關(guān)系進(jìn)行有限元計(jì)算[23].本文采用離散元模擬三軸實(shí)驗(yàn)來(lái)產(chǎn)生顆粒材料的應(yīng)力?應(yīng)變數(shù)據(jù),用于深度學(xué)習(xí)訓(xùn)練.如圖5 所示是本文采用的數(shù)值三軸試樣,每個(gè)試樣共包含4037 個(gè)顆粒,半徑為2～ 4 mm 均勻分布;顆粒之間采用線性接觸模型,法向和切向接觸剛度分別為105N/m 和5 × 104N/m;顆粒密度為2600 kg/m3,局部阻尼系數(shù)為0.5.為生成一個(gè)較密實(shí)的均勻試樣,在指定空間內(nèi)生成顆粒后,采用一個(gè)較小的制樣顆粒摩擦系數(shù)(如 0.05),使用伺服方法通過(guò)移動(dòng)試樣的6 個(gè)邊界來(lái)施加三向等壓壓力,直到試樣穩(wěn)定地達(dá)到200 kPa 的目標(biāo)圍壓值,再將制樣階段的接觸摩擦系數(shù)改成實(shí)際樣本摩擦系數(shù)0.5.

圖5 離散元三軸試驗(yàn)?zāi)Ｐ虵ig.5 Triaxial compression models via discrete element modelling

本文主要考慮了3 類(lèi)典型的三軸試驗(yàn)工況,具體包括:

(1)包含應(yīng)力?應(yīng)變加卸載循環(huán)的常規(guī)三軸壓縮試驗(yàn)291 組;

(2)平均有效應(yīng)力p恒定的真三軸加卸載試驗(yàn)122組(p=σ11+σ22+σ33,σ11=σ22);

(3)主應(yīng)力系數(shù)b恒定的真三軸壓縮試驗(yàn)20 組

所有工況的最大軸向應(yīng)變?yōu)?2%,在每組加卸載工況中,卸載應(yīng)變和再加載應(yīng)變的取值各不相同,在等b真三軸加載中,b以0.05 的間隔在(0,1]之間取值.在440 組所有三軸試驗(yàn)中,288 組用于訓(xùn)練,72 組用于驗(yàn)證(validation),剩下的80 組用于測(cè)試(test).注意驗(yàn)證不同于測(cè)試,驗(yàn)證組用于神經(jīng)網(wǎng)絡(luò)模型在訓(xùn)練過(guò)程中調(diào)整超參數(shù)組合,而測(cè)試集則僅用于測(cè)試目的,測(cè)試數(shù)據(jù)在模型訓(xùn)練過(guò)程中從未接觸到神經(jīng)網(wǎng)絡(luò),這樣的設(shè)計(jì)可避免在測(cè)試神經(jīng)網(wǎng)絡(luò)模型時(shí)發(fā)生數(shù)據(jù)泄露的可能性,從而保證神經(jīng)網(wǎng)絡(luò)測(cè)試結(jié)果的可靠性.

4.2 深度學(xué)習(xí)模型實(shí)施過(guò)程

在構(gòu)建深度學(xué)習(xí)模型之前,離散元模型中計(jì)算得出的數(shù)據(jù)需要預(yù)處理.一方面,現(xiàn)實(shí)問(wèn)題中的輸入變量往往具有不同的數(shù)據(jù)大小、單位以及分布特征,而差異較大的原始輸入數(shù)據(jù)會(huì)降低訓(xùn)練效率,阻礙神經(jīng)網(wǎng)絡(luò)在訓(xùn)練過(guò)程中的收斂.因此,原始數(shù)據(jù)需要通過(guò)標(biāo)準(zhǔn)化將其縮放到平均值為0、標(biāo)準(zhǔn)差為1 的分布空間,這將有利于降低訓(xùn)練過(guò)程中網(wǎng)絡(luò)陷入局部最優(yōu)的可能性.另一方面,循環(huán)神經(jīng)網(wǎng)絡(luò)模型的輸入數(shù)據(jù)需要滿足某些特定特征:所有的輸入數(shù)據(jù)必須是三維矩陣的形式,其中第一維表示樣本,一段應(yīng)力?應(yīng)變響應(yīng)的序列數(shù)據(jù)即是一個(gè)樣本;第二維是時(shí)間步(time steps),也稱移動(dòng)窗口的長(zhǎng)度,這是指循環(huán)神經(jīng)網(wǎng)絡(luò)在給出下一個(gè)預(yù)測(cè)輸出時(shí)需要依據(jù)在此之前的輸入變量的數(shù)量,若時(shí)間步取20,則在預(yù)測(cè)當(dāng)前應(yīng)力時(shí),需要輸入當(dāng)前應(yīng)變值以及在此之前的19 個(gè)應(yīng)變值;第三維則是數(shù)據(jù)特征,比如3 個(gè)主應(yīng)變分量可當(dāng)作3 個(gè)數(shù)據(jù)特征.

本文的深度學(xué)習(xí)模型通過(guò)Tensorflow 和keras來(lái)構(gòu)建和訓(xùn)練.神經(jīng)網(wǎng)絡(luò)的大小(層數(shù)和神經(jīng)元數(shù)量)決定了神經(jīng)網(wǎng)絡(luò)所能達(dá)到的上限表示能力、網(wǎng)絡(luò)訓(xùn)練過(guò)程中學(xué)習(xí)速率等超參數(shù)則決定了模型能在當(dāng)前數(shù)據(jù)集下將其預(yù)測(cè)能力發(fā)揮到何種水平.本文采用的網(wǎng)絡(luò)結(jié)構(gòu)和超參數(shù)組合通過(guò)試錯(cuò)法得到:兩個(gè)子網(wǎng)絡(luò)除了輸入層和輸出層以外,還包括兩個(gè)含有120 個(gè)神經(jīng)元的GRU 隱藏層和20 個(gè)神經(jīng)元的全連接隱藏層.全連接層采用Sigmoid 激活函數(shù),輸出層采用linear 激活函數(shù).移動(dòng)窗口的長(zhǎng)度為40,采用Adam 優(yōu)化器訓(xùn)練神經(jīng)網(wǎng)絡(luò);學(xué)習(xí)速率為0.01,批量數(shù)據(jù)的尺寸(batch size)為128.損失函數(shù)為預(yù)測(cè)值與實(shí)際值的平均絕對(duì)誤差函數(shù)(MAE).模型在訓(xùn)練了500 個(gè)epochs 后基本達(dá)到收斂狀態(tài),兩個(gè)神經(jīng)網(wǎng)絡(luò)在訓(xùn)練集和驗(yàn)證集上的預(yù)測(cè)表現(xiàn)分別如圖6 中的學(xué)習(xí)曲線所示.

圖6 學(xué)習(xí)曲線Fig.6 Learning curves

需要說(shuō)明的是,選擇40 的移動(dòng)窗口長(zhǎng)度并不意味著不關(guān)注每條應(yīng)力-應(yīng)變曲線上前39 步應(yīng)變引起的應(yīng)力值,整個(gè)序列中的應(yīng)力響應(yīng)都應(yīng)考慮進(jìn)去.實(shí)踐中為了確保總能得到長(zhǎng)度為40 的應(yīng)變序列,對(duì)于前40 組應(yīng)力值,我們通過(guò)在第一個(gè)應(yīng)變值前面不斷充填零值,使得應(yīng)變序列的長(zhǎng)度總為選定的移動(dòng)窗口長(zhǎng)度.這樣的措施確保了整個(gè)應(yīng)力序列都能夠得到預(yù)測(cè),實(shí)踐表明具有較好的預(yù)測(cè)效果.此外,充零的方式還具有較合理的物理意義,即在荷載作用之前,試樣并不產(chǎn)生應(yīng)變響應(yīng).

4.3 深度學(xué)習(xí)模型的評(píng)價(jià)指標(biāo)

預(yù)測(cè)模型的好壞通過(guò)量化預(yù)測(cè)值與實(shí)際值的平均差異得到.將平均絕對(duì)誤差函數(shù)作為損失函數(shù)是評(píng)價(jià)模型預(yù)測(cè)能力的一種有效方式,但此種誤差評(píng)價(jià)指標(biāo)往往與訓(xùn)練對(duì)象的值有關(guān),且缺乏對(duì)預(yù)測(cè)能力的直觀體現(xiàn).本文除了采用平均絕對(duì)誤差函數(shù)外,還采用另外一種評(píng)分的方式來(lái)對(duì)模型進(jìn)行評(píng)價(jià).該方法首先需要計(jì)算得到第j條應(yīng)力?應(yīng)變曲線上第i個(gè)點(diǎn)的平方誤差損失SSEij

式中,Nj為第j條應(yīng)力?應(yīng)變曲線上的樣本點(diǎn)數(shù)目,將曲線上所有點(diǎn)的SSEij按照升序排列,可以為每條曲線計(jì)算一個(gè)對(duì)應(yīng)的得分

式中,εP%是累計(jì)分布函數(shù)上對(duì)應(yīng)P%時(shí)的SSE值,這個(gè)值被用來(lái)作為一個(gè)代表值來(lái)評(píng)價(jià)整條應(yīng)力?應(yīng)變曲線上的預(yù)測(cè)得分;當(dāng)預(yù)測(cè)應(yīng)力?應(yīng)變序列中εP%對(duì)應(yīng)的SSE小于εcrit時(shí),可認(rèn)為預(yù)測(cè)足夠準(zhǔn)確,得分為1.在本文中,P%取90%,εcrit為0.001.

以上介紹了單條應(yīng)力?應(yīng)變預(yù)測(cè)曲線的打分標(biāo)準(zhǔn),在評(píng)價(jià)每個(gè)神經(jīng)網(wǎng)絡(luò)模型的預(yù)測(cè)能力時(shí),一個(gè)更合理的指標(biāo)是判斷該模型在測(cè)試集上所有樣本的平均預(yù)測(cè)表現(xiàn),即使用已經(jīng)訓(xùn)練好的模型,得到對(duì)所有測(cè)試樣本上的平均MAE 值或者平均得分.

5 預(yù)測(cè)結(jié)果

5.1 總體預(yù)測(cè)表現(xiàn)

在80 組測(cè)試樣本中,該深度學(xué)習(xí)模型在測(cè)試集上的平均MAE 值為0.013;平均得分為0.986;87.5%的測(cè)試樣本(70 組)獲得了滿分,其中最差的預(yù)測(cè)得分為0.639,對(duì)應(yīng)MAE 值為0.051;次差的一組預(yù)測(cè)得分為0.823,對(duì)應(yīng)MAE 值為0.031;在70 組滿分預(yù)測(cè)中,最低的MAE 值為0.004 2 和0.004 5.為直觀體現(xiàn)模型的預(yù)測(cè)能力,兩組最好和最差的預(yù)測(cè)情況如圖7 所示.

圖7 兩組最佳與最差預(yù)測(cè)Fig.7 Examples of the two best and worst predictions

在得分最低的兩組預(yù)測(cè)中,盡管在加卸載處的應(yīng)力?應(yīng)變預(yù)測(cè)值與實(shí)際值存在差別,但總體的應(yīng)力響應(yīng)趨勢(shì)已經(jīng)被神經(jīng)網(wǎng)絡(luò)模型捕捉到.在損失函數(shù)MAE 最低的兩組預(yù)測(cè)中,深度學(xué)習(xí)模型則幾乎重現(xiàn)了實(shí)際的應(yīng)力?應(yīng)變響應(yīng).

5.2 內(nèi)插與外推預(yù)測(cè)表現(xiàn)

根據(jù)已有數(shù)據(jù)以及模型預(yù)測(cè)未知區(qū)域的目標(biāo)對(duì)象時(shí),若預(yù)測(cè)點(diǎn)在已有數(shù)據(jù)范圍內(nèi)通常稱為內(nèi)插預(yù)測(cè);若數(shù)據(jù)點(diǎn)完全在已有數(shù)據(jù)范圍之外,則為外插預(yù)測(cè).本文所用的訓(xùn)練數(shù)據(jù)集僅在常規(guī)三軸和等p真三軸中包含了一次應(yīng)力?應(yīng)變加卸載循環(huán),在等b真三軸壓縮工況中僅包含單調(diào)加載工況.圖8 分別給出了在測(cè)試集內(nèi)的常規(guī)三軸和等p真三軸一次加卸載,以及等b真三軸加載工況中的幾個(gè)代表性內(nèi)插預(yù)測(cè).

圖8 幾組代表性內(nèi)插預(yù)測(cè)Fig.8 Some representative interpolation predictions

注意到在常規(guī)三軸加載工況下,由于水平圍壓在加載過(guò)程中保持穩(wěn)定,σ22和σ33均保持不變;在平均有效應(yīng)力p恒定的真三軸加載工況中,σ22和σ33保持了相同的變化趨勢(shì);在中主應(yīng)力系數(shù)b恒定的真三軸加載工況中,σ33維持大小不變.經(jīng)過(guò)訓(xùn)練的神經(jīng)網(wǎng)絡(luò)模型均準(zhǔn)確預(yù)測(cè)到了這些趨勢(shì).

本研究在測(cè)試數(shù)據(jù)集中準(zhǔn)備了常規(guī)三軸和等p真三軸單調(diào)加載工況,也準(zhǔn)備了6 組包含2 次及以上加卸載循環(huán)的常規(guī)三軸加載工況.這些工況完全在模型訓(xùn)練過(guò)程中所用的訓(xùn)練集和測(cè)試集范圍之外,因此被用來(lái)測(cè)試神經(jīng)網(wǎng)絡(luò)的外推預(yù)測(cè)能力.8 組測(cè)試樣本中除了一組預(yù)測(cè)得分為0.98 以外,其余預(yù)測(cè)都得到了滿分.幾組代表性預(yù)測(cè)工況如圖9 所示,結(jié)果表明深度學(xué)習(xí)模型能準(zhǔn)確地適應(yīng)單調(diào)加載和多次加卸載的外推工況.

圖9 幾種代表性外推預(yù)測(cè)Fig.9 Some representative extrapolation predictions

6 討論

顆粒材料的多次加卸載等復(fù)雜工況本質(zhì)上都體現(xiàn)為細(xì)觀結(jié)構(gòu)的往復(fù)演化.在材料參數(shù)不發(fā)生較大變化的情況下,預(yù)測(cè)試樣應(yīng)力?應(yīng)變響應(yīng)的關(guān)鍵在于預(yù)測(cè)試樣組構(gòu)的演化.由于采用了基于時(shí)序的循環(huán)神經(jīng)網(wǎng)絡(luò)模型,在預(yù)測(cè)當(dāng)前應(yīng)力值時(shí),包含了足夠長(zhǎng)的應(yīng)變歷史信息來(lái)推斷當(dāng)前的組構(gòu)演化和應(yīng)力響應(yīng),因此深度學(xué)習(xí)模型能有效區(qū)分加載和卸載,做出準(zhǔn)確的應(yīng)力預(yù)測(cè).另外,盡管包含多次加卸載的工況未曾直接出現(xiàn)在模型訓(xùn)練過(guò)程中,但是類(lèi)似的細(xì)觀組構(gòu)演化過(guò)程可能已經(jīng)被訓(xùn)練數(shù)據(jù)所覆蓋,這也是深度學(xué)習(xí)模型對(duì)多次加卸載工況具備較好外推能力的原因.外推預(yù)測(cè)結(jié)果表明當(dāng)訓(xùn)練集中包含豐富的一次加卸載樣本時(shí),模型將具有一定的兩次或三次等加卸載工況的外推預(yù)測(cè)能力,這在實(shí)踐中能有效降低準(zhǔn)備訓(xùn)練樣本的數(shù)據(jù)規(guī)模.

相比于顆粒材料室內(nèi)三軸試驗(yàn),基于離散元的虛擬三軸試驗(yàn)所需成本更低,效率更高,且方便模擬各類(lèi)加載工況.采用基于離散元的虛擬三軸試驗(yàn)產(chǎn)生顆粒材料應(yīng)力?應(yīng)變數(shù)據(jù)主要包含兩個(gè)方面的意義:(1) 離散元方法已經(jīng)表明具備定性再現(xiàn)顆粒材料基本響應(yīng)特征的能力[24],通過(guò)低成本的離散元模擬,可為探索數(shù)據(jù)驅(qū)動(dòng)本構(gòu)模型提供經(jīng)驗(yàn).(2) 高質(zhì)量的離散元模擬本身有能力對(duì)實(shí)際顆粒材料的宏觀力學(xué)響應(yīng)做出定量預(yù)測(cè).在離散元的參數(shù)得到合理標(biāo)定[25-26],使用更加真實(shí)合理的接觸模型,并考慮復(fù)雜形狀因素的離散元模擬中[27-28],實(shí)際顆粒材料的本構(gòu)行為能很好的被離散元所捕捉到,此時(shí)基于離散元數(shù)據(jù)基本能代表真實(shí)材料數(shù)據(jù)[22,29-30].

本文的理論推導(dǎo)主要基于線彈性接觸模型,由于從線彈性接觸到非線性彈性接觸的映射關(guān)系并不復(fù)雜,深度神經(jīng)網(wǎng)絡(luò)完全有能力捕捉到這個(gè)映射關(guān)系的轉(zhuǎn)變.因此本文提出的深度學(xué)習(xí)框架也適用于預(yù)測(cè)非線性彈性接觸模型的應(yīng)力?應(yīng)變響應(yīng).此外,若需嚴(yán)格地發(fā)展基于非線性彈性接觸模型的深度學(xué)習(xí)策略,可在本文提出的訓(xùn)練框架下,參考文獻(xiàn)[25]中推導(dǎo)的基于非線性彈性接觸模型的顆粒材料應(yīng)力?應(yīng)變關(guān)系進(jìn)行深度學(xué)習(xí)訓(xùn)練.

數(shù)據(jù)驅(qū)動(dòng)本構(gòu)模型相比于傳統(tǒng)基于唯象假設(shè)的解析模型有幾個(gè)突出優(yōu)點(diǎn):(1)可以根據(jù)新數(shù)據(jù)的補(bǔ)充而不斷改進(jìn)和優(yōu)化,而唯象模型在參數(shù)標(biāo)定好以后,成為一個(gè)封閉的體系,無(wú)法充分利用工程活動(dòng)中產(chǎn)生的海量數(shù)據(jù)資源;(2)復(fù)雜的唯象模型往往需要標(biāo)定多達(dá)十多個(gè),甚至數(shù)十個(gè)自由參數(shù),并且有些模型在有限元中實(shí)施較為復(fù)雜;然而數(shù)據(jù)驅(qū)動(dòng)的本構(gòu)模型不論考慮何種復(fù)雜的荷載條件,在有限元等連續(xù)方法中實(shí)施時(shí)均能夠做到標(biāo)準(zhǔn)化和模塊化.

數(shù)據(jù)驅(qū)動(dòng)模型和理論驅(qū)動(dòng)模型具有不同的優(yōu)點(diǎn)和缺點(diǎn),數(shù)據(jù)驅(qū)動(dòng)本構(gòu)模型并不能替代理論本構(gòu)模型,兩種研究范式是相互補(bǔ)充和促進(jìn)的.一個(gè)典型的例子是,盡管離散元模擬可以得到代表性試樣中的剪應(yīng)力分量,在室內(nèi)試驗(yàn)條件下我們一般僅能測(cè)到試樣中的3 個(gè)主應(yīng)力和主應(yīng)變分量;深度學(xué)習(xí)模型在主應(yīng)力和主應(yīng)變數(shù)據(jù)集下也只能訓(xùn)練得到主應(yīng)力與主應(yīng)變空間的映射關(guān)系.一個(gè)完整的本構(gòu)模型需要進(jìn)一步將主應(yīng)力與主應(yīng)變空間的應(yīng)力關(guān)系變換到一般應(yīng)力?應(yīng)變空間當(dāng)中.由于顆粒材料在加載過(guò)程中主應(yīng)力和主應(yīng)變?cè)隽糠较虼嬖诜枪草S特征[31],因此純粹的坐標(biāo)變換方法將不可避免存在誤差,考慮顆粒材料主應(yīng)力和主應(yīng)變關(guān)系的非共軸特性是唯象模型和數(shù)據(jù)模型共同的挑戰(zhàn)[32].

7 結(jié)論

本文結(jié)合小應(yīng)變應(yīng)力?應(yīng)變理論關(guān)系和深度神經(jīng)網(wǎng)絡(luò)模型,通過(guò)有向圖發(fā)展了一種基于深度循環(huán)神經(jīng)網(wǎng)絡(luò)的顆粒材料應(yīng)力?應(yīng)變數(shù)據(jù)驅(qū)動(dòng)模型訓(xùn)練方法.通過(guò)常規(guī)三軸加卸載、等p真三軸加卸載和等b真三軸加載等多軸加載工況對(duì)模型進(jìn)行訓(xùn)練和測(cè)試.結(jié)果表明使用的深度神經(jīng)網(wǎng)絡(luò)模型能有效預(yù)測(cè)訓(xùn)練和驗(yàn)證過(guò)程從未見(jiàn)過(guò)的測(cè)試數(shù)據(jù),在僅含有常規(guī)三軸和等p真三軸一次加卸載循環(huán)的訓(xùn)練集中,深度學(xué)習(xí)模型有足夠的泛化能力預(yù)測(cè)到包含單調(diào)加載和多個(gè)加卸載循環(huán)的應(yīng)力響應(yīng),具有良好的內(nèi)插和外推預(yù)測(cè)能力.結(jié)果表明通過(guò)結(jié)合先驗(yàn)的應(yīng)力?應(yīng)變解析范式和經(jīng)驗(yàn)的數(shù)據(jù)驅(qū)動(dòng)范式,可以發(fā)展出適應(yīng)多種加載工況的顆粒材料主應(yīng)力?主應(yīng)變預(yù)測(cè)模型.