形狀記憶合金數(shù)值模型的不確定性分析

李金朋，陳城，侯和濤

(1.山東大學 a.土建與水利學院；b.山東省綠色建筑智能建造工程技術(shù)研究中心;c.山東省工業(yè)技術(shù)研究院 工業(yè)化建筑智能建造協(xié)同創(chuàng)新中心，濟南 250061;2.舊金山州立大學 工程學院，舊金山 94132)

筆者提出對模型參數(shù)進行概率建模的方法，基于SMA棒材的實驗數(shù)據(jù)，采用改進的Graesser和Cozzarelli模型與MCMC算法的組合來分析模型本身固有的不確定性，將模型參數(shù)視為隨機變量，采用Metropolis-Hastings算法來生成樣本參數(shù)集，揭示了參數(shù)的概率特性和參數(shù)之間的潛在相關(guān)性，并從模型參數(shù)的角度研究了SMA模型中固有的不確定性及其對材料能量耗散能力預測值的影響。

1 SMA的數(shù)值模型與材性試驗

1.1 改進的Graesser & Cozzarelli 模型

(1)

(2)

式中：σ為一維應力；β為一維背應力；ε為一維的應變；E為材料的彈性模量；Y為材料的屈服應力；參數(shù)n控制SMA由彈性階段向塑性階段過渡時滯回曲線的“尖銳”程度，取值范圍設計為任意正奇數(shù)；參數(shù)α控制塑性階段滯回曲線斜率，由α=Ey/(E-Ey)得到，其中Ey為屈服后的彈性模量；fT是控制滯回曲線類型和大小的參數(shù)，當fT=0時，模型代表純馬氏體狀態(tài)下的SMA，當fT>0時，模型可模擬SMA的“超彈性”行為；a是控制卸載過程中材料彈性恢復量的參數(shù)；c是控制卸載過程中應力平臺段斜率的參數(shù)；(·)表示對時間的一階導數(shù)；u(x)和erf(x)分別為單位階躍函數(shù)和誤差函數(shù)，數(shù)學表達式為

(3)

但是，該模型無法模擬SMA由奧氏體相向馬氏體相轉(zhuǎn)化完成后出現(xiàn)的“硬化”現(xiàn)象，即材料的彈性模量突然增大，稱為SMA的馬氏體硬化特性。為了描述這一特性，Qian等[10]提出了改進的Graesser & Cozzarelli模型，其數(shù)學表達式為

(4)

(5)

改進后，相較于式(1)，式(4)中參數(shù)n的取值范圍變?yōu)槿我庹龑崝?shù)，更方便該參數(shù)的不確定性分析。另外，新引入的fM和m為控制馬氏體硬化階段曲線的參數(shù)；εMf為SMA轉(zhuǎn)化為全馬氏體時的應變；sgn(x)為符號函數(shù)，數(shù)學表達式為

(6)

相較于式(2)，式(5)添加了描述馬氏體硬化的表達式。在SMA進入馬氏體硬化階段前，這部分值為0，此時式(5)與式(2)完全相同。

1.2 SMA棒材的循環(huán)拉伸試驗

表1 SMA試件原材料的化學成分Table 1 Chemical composition of the testing SMA bar

圖1 SMA循環(huán)拉伸試驗Fig.1 Cyclic tests of SMA bar

將熱處理后的試件進行循環(huán)拉伸試驗，其加載制度如圖2(a)所示。圖2(b)為循環(huán)拉伸試驗得到的滯回曲線，其數(shù)據(jù)用于后續(xù)的不確定性分析。在進行不確定性分析之前，基于峰值應變?yōu)?.08的兩條滯回曲線，采用PSO方法得到使模型具備良好擬合效果的參數(shù)值作為初始參數(shù)組，有助于提高馬爾可夫鏈蒙特卡洛(MCMC)模擬時馬爾可夫鏈的收斂速度。由PSO方法得到的初始參數(shù)具體數(shù)值如表2所示。從圖2(b)可以看到，PSO優(yōu)化方法得到的參數(shù)雖然使得模型具有良好的擬合效果，但模擬結(jié)果與試驗數(shù)據(jù)仍舊存在偏差，因此，有必要研究模型中存在的固有不確定性。

圖2 加載制度和試驗及優(yōu)化模型模擬結(jié)果示意圖Fig.2 Loading protocol and comparison between experiment results and PSO methods

表2 基于循環(huán)拉伸試驗數(shù)據(jù)的PSO優(yōu)化參數(shù)值Table 2 Model parameter values from cyclic tests of SMA bars using PSO method

2 MCMC不確定性分析

2.1 MCMC方法

2.1.1 DRAM方法 采用MCMC方法中兩個重要方法的結(jié)合，即延遲拒絕法(DR法，Delaying Rejection)和自適應采樣(AM法，Adaptive Metropolis Samplers)，簡稱為DRAM方法[13-14]。DR法通過適當調(diào)整標準Metropolis-Hastings(MH)算法中馬爾可夫鏈每一步中的提議分布來提高MCMC算子的效率。AM法則是基于馬爾可夫鏈的歷史來調(diào)整提議分布。當提議分布的方差非常小時，標準MH算法傾向于以小步長“遍歷”目標分布，而無法有效地探索狀態(tài)空間，且樣本點分布偏移樣本空間中心，產(chǎn)生偏差。Haario等[14]的研究證明了AM方法能夠解決“探索范圍未覆蓋整個樣本空間”的問題，DR方法能夠解決“樣本集中位置偏移樣本空間中心”的問題，而DR和AM的組合，即DRAM方法可以同時解決這兩個問題。

2.1.2 似然函數(shù) 似然函數(shù)定義為給定的一組參數(shù)值下，模型模擬結(jié)果與實驗數(shù)據(jù)一致的概率，也可以將其視為模型預測和實驗測量之間的誤差概率。ss函數(shù)是似然函數(shù)的一部分，用于描述歸一化后的誤差[15]。似然函數(shù)及ss函數(shù)公式為

(7)

(8)

式中：Lj為第j步的似然函數(shù)值；ssj為第j步的歸一化誤差值；n為實驗數(shù)據(jù)的數(shù)據(jù)點個數(shù)；σ2為設定的容許誤差，用以限制ss的大??；F*(ti)為在時間點ti的實驗測量結(jié)果；fj(ti)為第j步時在時間點ti的模型輸出結(jié)果。

DRAM模擬方法的具體流程如圖3所示，采用MATLAB編程并模擬。模擬次數(shù)(Nsimu)設置為2×106次，為了研究合理范圍內(nèi)的不確定性，容許誤差σ2設為0.022。

圖3 DRAM模擬方法的流程示意圖Fig.3 Flow chart of simulation process

2.1.3 參數(shù)的先驗概率 在模擬開始之前，首先要建立參數(shù)的先驗概率，包括參數(shù)的范圍及其在該范圍內(nèi)的分布。研究表明，先驗分布并不是MCMC模擬得到參數(shù)所收斂的后驗分布的決定性因素，而是影響收斂速度的關(guān)鍵因素[16]。因此，模擬中僅根據(jù)數(shù)學或物理要求對參數(shù)的范圍進行合理規(guī)定，不指定其分布方式。參數(shù)的先驗概率設置如表3所示。

表3 參數(shù)的先驗概率設置Table 3 Prior of model parameters

2.2 MCMC分析結(jié)果

2.2.1 模型參數(shù) 基于峰值應變?yōu)?%的加載循環(huán)試驗數(shù)據(jù)得到的參數(shù)不確定性分析結(jié)果如表4、圖4所示。表4列出了10個參數(shù)的概率特性，包括均值、方差和偏度。圖4為所有參數(shù)的頻率分布直方圖，圖4中大多數(shù)參數(shù)的頻率分布呈現(xiàn)出單峰結(jié)構(gòu)，表明對參數(shù)采樣的馬爾可夫鏈收斂服從該參數(shù)

表4 參數(shù)分布的均值、標準差和偏度Table 4 Moments of model parameters

圖4 參數(shù)頻率分布直方圖Fig.4 Histogram of model parameters

的后驗概率分布。圖4顯示，除了參數(shù)E、參數(shù)Y和參數(shù)εMf外，其余參數(shù)后驗分布的均值均大于PSO優(yōu)化方法得到的初始參數(shù)值。這意味著PSO的優(yōu)化結(jié)果是否最優(yōu)有待商榷，原因可能在于其算法中粒子游走出現(xiàn)某一方向的偏差而錯失最優(yōu)路徑。此外，絕大部分參數(shù)的分布峰值位置與初始值基本吻合，表明優(yōu)化方法得到的參數(shù)值通常具有普遍性；然而，參數(shù)c的分布中，峰值出現(xiàn)的位置距離初始值有較大偏差，表明PSO優(yōu)化方案應用在SMA模型中可能具有片面性以及偏差性，即優(yōu)化方法提供的參數(shù)值可能并不是該參數(shù)最普遍采用的值。這表明SMA具有本身固有的不確定性，利用MCMC方法研究SMA不確定性是必要的。

圖5顯示了參數(shù)α、εMf、m和fM之間的成對相關(guān)性。兩個參數(shù)的相關(guān)性越強，則對應橢圓的離心率越大；橢圓長軸方向為左下至右上時表示參數(shù)之間存在正相關(guān)關(guān)系，反之呈負相關(guān)關(guān)系。由圖5可見，參數(shù)m和α、a和c、m和fM兩兩之間存在較強正相關(guān)性；fT和α之間存在較強負相關(guān)性；fT和a、c和fM等之間不存在線性相關(guān)性或相關(guān)性較弱。參數(shù)之間的相關(guān)性可以為數(shù)學模型的研究提供參考。

圖5 參數(shù)相關(guān)性示意圖Fig.5 Correlation between the simulated parameters

2.2.2 能量耗散 為了更好地說明不確定性研究的必要性，研究通過模型參數(shù)的概率分布建立材料耗能能力的概率特征，其概率密度示意圖如圖6(a)所示。與實驗數(shù)據(jù)得到的結(jié)果相比，在累積概率密度為15%時，材料的能量耗散能力相對誤差高達20%；累積概率密度為85%時，相對誤差為10%。圖6(b)顯示了對應能量耗散的累積概率密度值分別為0.15、0.5、0.64和0.85時模型的擬合效果。結(jié)果顯示，此時模型具有非常好的擬合效果，但沒有PSO優(yōu)化結(jié)果好。這表明更新模型參數(shù)不會消除或補償模型模擬的偏差，但會將偏差控制在可接受的范圍內(nèi)。

圖6 基于能量耗散能力預測的不確定性分析結(jié)果示意圖Fig.6 The results based on the energy dissipation analysis

能量耗散通常通過等效粘性阻尼(EVD)ξeq來衡量，ξeq是一個與尺寸無關(guān)的指數(shù)，表示為

(7)

式中：WD為每個加載周期SMA的能量耗散值；WE為相應線性系統(tǒng)的應變能。EVD在很大程度上取決于滯回環(huán)的形狀，較大的滯回環(huán)會得到較大的EVD值，滯回曲線越飽滿，EVD值越大。由圖7可知，隨著峰值應變的增加，EVD略微增加，當峰值應變達到8%時，EVD略有下降。這表明當加載峰值應變達到近6%時，材料達到其最佳能量耗散能力。這一結(jié)論對基于SMA的地震控制裝置的設計具有借鑒意義。

圖7 基于不同應變峰值的EVD不確定性分析結(jié)果示意圖Fig.7 Uncertainty analysis of EVD for different peak-strain

3 結(jié)論

1)基于形狀記憶合金棒材循環(huán)拉伸試驗數(shù)據(jù)的DRAM算法采樣得到的馬爾可夫鏈體現(xiàn)出模型參數(shù)的概率特性。樣本的分布特征(均值、方差等)體現(xiàn)出優(yōu)化方法可能存在偏差，部分參數(shù)之間存在線性相關(guān)性，在進行數(shù)值模型研究時應予以重視。

2)數(shù)值模型的不確定性也體現(xiàn)在模型的耗能預測上，在累積概率密度為15%時，材料的能量耗散能力相對誤差高達20%；累積概率密度為85%時，相對誤差為10%。加載應變峰值對材料的耗能性能有明顯影響，等效粘滯阻尼分布顯示，加載峰值應變?yōu)?%時，材料耗能性能較其他對比組更好。