本發(fā)明涉及一種基于晶體滑移機(jī)制的各向異性線彈性本構(gòu)的建立方法，屬于材料力學(xué)性能有限元數(shù)值計(jì)算領(lǐng)域。

背景技術(shù)：

在利用基于真實(shí)微觀組織的建模方法模擬多晶材料宏、微觀力學(xué)性能方面，材料本構(gòu)模型的選用和參數(shù)獲取對(duì)模擬結(jié)果具有決定性的影響，同時(shí)也是一個(gè)亟待解決的薄弱環(huán)節(jié)。國內(nèi)外在研究材料本構(gòu)的問題上，通常是采用各向同性本構(gòu)模型，即用唯象的描述方法考慮介觀層次的相組織模型。獲得相本構(gòu)的方法主要是認(rèn)為材料微觀變形特征與宏觀近似相同，進(jìn)而將具有均一組織結(jié)構(gòu)材料的整體性能賦予該相；而對(duì)于多相材料中的各個(gè)組成相，則是采用具有近似成分的單相材料的宏觀性能進(jìn)行替代。大型商用有限元軟件LS-DYNA是一種通用顯式動(dòng)力分析軟件，在工程應(yīng)用領(lǐng)域具有較高的計(jì)算可靠性。其具有顯/隱式求解功能，求解穩(wěn)定；同時(shí)具有在復(fù)雜邊界加載條件下的非線性動(dòng)力分析和靜力分析功能。后處理功能更是可以實(shí)現(xiàn)內(nèi)部多種求解參數(shù)的可視化，但目前其應(yīng)用主要是基于唯象學(xué)，鮮有考慮材料變形的內(nèi)在機(jī)制。

然而為了考慮材料變形過程中存在的內(nèi)部物理機(jī)制，各向同性的相組織模型已無法滿足要求，因此所建立的模型已逐漸發(fā)展為考慮滑移機(jī)制的晶粒組織模型，即由唯象描述方發(fā)展為晶體塑性描述方法。其中用于研究材料加載變形內(nèi)部機(jī)制的彈塑性自洽模擬方法已得到了廣泛的應(yīng)用，該方法可以用于模擬材料整體的宏觀應(yīng)力-應(yīng)變響應(yīng)；研究材料彈塑性轉(zhuǎn)變過程及過程中表征參量的擬合；擬合晶格演變過程；變形過程中的織構(gòu)演化等。自洽模擬方法的主要理論基礎(chǔ)和假設(shè)有：(1)采用Hill-Hutchinson多晶體彈塑性變形機(jī)制、Vocé硬化定律及硬化系數(shù)矩陣；(2)晶粒的各向異性的實(shí)現(xiàn)方法是通過將晶粒設(shè)置成具有各向異性彈性常數(shù)矩陣的橢圓球體，并且給每個(gè)晶粒賦予對(duì)應(yīng)其晶體結(jié)構(gòu)類型的滑移系及特征參數(shù)；(3)晶粒間的相互作用是通過彈塑性Eshelby張量建立起每一個(gè)晶體與假定的均質(zhì)基體之間的相互作用函數(shù)；(4)采用增量加載模式，每一加載步為小變形加載，即認(rèn)為應(yīng)力-應(yīng)變關(guān)系為瞬時(shí)線彈性。

該方法不僅考慮了多晶材料中晶體的實(shí)際取向，還可以考慮多相材料不同的晶體結(jié)構(gòu)類型和塑性變形滑移機(jī)制、硬化效應(yīng)和取向變化等。但其存在一些缺陷，如：需要配合原位拉伸同步輻射X射線衍射、中子衍射等實(shí)驗(yàn)來獲得與之相映證的數(shù)據(jù)，而這些實(shí)驗(yàn)不僅周期長而且費(fèi)用昂貴；其相互作用是通過“限制張量”在算法內(nèi)部實(shí)現(xiàn)的，無法進(jìn)行可視化，并且無法考慮真實(shí)微觀組織中相鄰晶粒間的相互作用關(guān)系；可設(shè)置的加載邊界條件比較簡(jiǎn)單，通常是由應(yīng)力或應(yīng)變控制的以一定速率的單軸拉伸或壓縮等。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明的目的是提供一種基于晶體滑移機(jī)制的各向異性線彈性本構(gòu)的建立方法，所述方法延續(xù)了自洽模型中采用的描述晶體塑性變形的理論機(jī)制，同時(shí)與有限元方法相結(jié)合，不僅考慮了材料的微觀組織形貌，同時(shí)用有限元網(wǎng)格的單元和節(jié)點(diǎn)化來實(shí)現(xiàn)晶粒間的相互作用，彌補(bǔ)了自洽模型中的不足，并可充分發(fā)揮有限元模擬計(jì)算的優(yōu)勢(shì)，實(shí)現(xiàn)復(fù)雜的邊界加載條件和多樣的可視化后處理操作。

本發(fā)明的目的有以下技術(shù)方案實(shí)現(xiàn)：

一種基于晶體滑移機(jī)制的各向異性線彈性本構(gòu)的建立方法，所述方法包括以下步驟：

步驟一：

基于待求各向異性線彈性本構(gòu)材料的晶粒組織圖像，采用有限元網(wǎng)格剖分方法，建立所述晶粒組織的有限元網(wǎng)格模型，并確定所述有限元網(wǎng)格模型的邊界條件和加載條件；

其中，所述加載條件包括載荷和最大加載步數(shù)，其中，初始加載步數(shù)為0；

步驟二：

采用投影運(yùn)算方法，根據(jù)晶粒組織的本征參數(shù)、初始應(yīng)力應(yīng)變狀態(tài)分量和晶粒取向的初始?xì)W拉角求解晶粒取向歐拉角增量并將所述歐拉角增量與初始?xì)W拉角相加，得到由于應(yīng)變引起晶粒取向變化后的歐拉角根據(jù)變化后的歐拉角獲得有限元網(wǎng)格模型的坐標(biāo)系與晶粒坐標(biāo)系之間的轉(zhuǎn)化矩陣；

其中，所述本征參數(shù)包括晶體結(jié)構(gòu)類型參數(shù)、單晶體彈性剛度矩陣、晶粒內(nèi)的滑移系參數(shù)；

所述滑移系參數(shù)包括滑移系的名稱、等效滑移系的數(shù)量和滑移系的Vocé硬化指數(shù)；所述Vocé硬化指數(shù)包括滑移系開動(dòng)的臨界分切應(yīng)力值τ₀、滑移系從開動(dòng)到累積切應(yīng)變不再增加時(shí)的應(yīng)力增量τ₁、滑移系的初始硬化率θ₀以及滑移系的飽和硬化率θ₁；

所述晶粒取向歐拉角采用Bunge方法定義，即將晶粒坐標(biāo)系繞其z軸，在xy平面內(nèi)轉(zhuǎn)過角；再繞其x軸，在yz平面內(nèi)轉(zhuǎn)過φ角；最后再一次繞z軸，在xy平面內(nèi)轉(zhuǎn)過角，最終與有限元網(wǎng)格模型坐標(biāo)系完全重合，由此得到每個(gè)晶粒的三個(gè)空間取向歐拉角晶粒取向的初始?xì)W拉角利用電子背散射衍射(EBSD)實(shí)驗(yàn)測(cè)量得到；

所述應(yīng)力應(yīng)變狀態(tài)分量分別按照σ_xx，σ_yy，σ_zz，σ_yz，σ_xz，σ_xy和ε_xx，ε_yy，ε_zz，ε_yz，ε_xz，ε_xy的順序給出，初始值均設(shè)為0；

其中，所述σ_xx，σ_yy和σ_zz分別表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)沿X軸、Y軸和Z軸方向的正應(yīng)力；

所述ε_xx，ε_yy和ε_zz分別表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)沿X軸、Y軸和Z軸方向的正應(yīng)變；

所述σ_yz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在Y軸法向平面內(nèi)平行于X軸的切應(yīng)力；

所述σ_xz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸法向平面內(nèi)平行于Z軸的切應(yīng)力；

所述σ_xy表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸法向平面內(nèi)平行于Y軸的切應(yīng)力；

所述ε_yz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在Y軸和Z軸兩方向上微小線段夾角的切應(yīng)變；

所述ε_xz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸和Z軸兩方向上微小線段夾角的切應(yīng)變；

所述ε_xy表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸和Y軸兩方向上微小線段夾角的切應(yīng)變；

步驟三：

基于步驟二所述轉(zhuǎn)化矩陣，通過張量運(yùn)算，得到有限元網(wǎng)格模型坐標(biāo)系下的單晶體彈性剛度矩陣和施密特因子α，其中，所述α＝cosλcosΦ，Φ為滑移面與外力F中心軸的夾角，所述λ為滑移方向與外力F的夾角；

步驟四：

判斷所述歐拉角增量是否為(0，0，0)，若進(jìn)入步驟五；若進(jìn)入步驟六；

步驟五：

根據(jù)步驟一所述有限元網(wǎng)格模型的邊界條件和加載載荷對(duì)當(dāng)前有限元網(wǎng)格模型進(jìn)行加載，然后針對(duì)當(dāng)前有限元網(wǎng)格模型坐標(biāo)系下的單晶體彈性剛度矩陣，采用有限元計(jì)算方法計(jì)算得到加載后晶粒的應(yīng)力應(yīng)變狀態(tài)分量；將所述加載后晶粒的應(yīng)力應(yīng)變狀態(tài)分量替換步驟二中的初始應(yīng)力應(yīng)變狀態(tài)分量，重復(fù)步驟二至步驟四；其中，每次加載后加載步數(shù)加一；步驟六：

判斷晶粒內(nèi)的潛在要開動(dòng)和已經(jīng)開動(dòng)滑移系的總個(gè)數(shù)n是否大于零；若n＝0，重復(fù)步驟五；若n＞0，進(jìn)入步驟七；

潛在要開動(dòng)和已經(jīng)開動(dòng)滑移系的判斷方法：

若則該滑移系被標(biāo)記為潛在要開動(dòng)，若則該滑移系被標(biāo)記為已經(jīng)開動(dòng)；其中，σ_c表示任意滑移系所在晶粒當(dāng)前的應(yīng)力值，表示第i個(gè)滑移系的臨界分切應(yīng)力值；所述αⁱ表示第i個(gè)滑移系的施密特因子；

步驟七：

根據(jù)Hill-Hutchinson多晶體彈塑性變形理論，先求解n個(gè)滑移系應(yīng)變線性方程組的基向量矩陣；再根據(jù)所述基向量矩陣，求解分配給每個(gè)滑移系的切變量，進(jìn)而得到晶粒在塑性階段的瞬時(shí)剛度矩陣；

步驟八：

判斷當(dāng)前加載步數(shù)是否等于設(shè)定的最大加載步數(shù)，如果等于最大加載步數(shù)，則步驟七得到的瞬時(shí)剛度矩陣即為晶粒的各向異性線彈性本構(gòu)；如果小于最大加載步數(shù)，則將步驟七的瞬時(shí)剛度矩陣替換步驟二中的單晶體彈性剛度矩陣，重復(fù)步驟二至步驟七，直至加載步數(shù)等于步驟一中設(shè)定的最大加載步數(shù)。

有益效果

(1)本發(fā)明所述方法延續(xù)了自洽模型中采用的描述晶體塑性變形的理論機(jī)制，同時(shí)與有限元方法相結(jié)合；采用基于晶粒組織的網(wǎng)格模型，不僅考慮了材料的微觀組織形貌，同時(shí)用有限元網(wǎng)格的單元和節(jié)點(diǎn)化來實(shí)現(xiàn)晶粒間的相互作用，彌補(bǔ)了自洽模型中的不足。

(2)本發(fā)明所述方法利用有限元求解器中的重啟動(dòng)技術(shù)，基于晶粒組織建立的有限元網(wǎng)格模型以一定邊界條件和加載條件進(jìn)行增量式分步重啟動(dòng)加載，每一加載步都繼承了上一加載步結(jié)束時(shí)各單元和節(jié)點(diǎn)的應(yīng)力應(yīng)變狀態(tài)，從而實(shí)現(xiàn)分步增量加載；通過采用有限元求解軟件中提供的各向異性線彈性本構(gòu)模型來體現(xiàn)多晶體材料中各個(gè)晶粒的取向性，同時(shí)結(jié)合晶體塑性變形理論，可以根據(jù)每一步加載后晶粒的應(yīng)力應(yīng)變狀態(tài)，從而達(dá)到在塑性變形機(jī)制層面對(duì)本構(gòu)參數(shù)(即瞬時(shí)模量)實(shí)現(xiàn)實(shí)時(shí)的更新改進(jìn)；充分發(fā)揮有限元模擬計(jì)算的優(yōu)勢(shì)，實(shí)現(xiàn)復(fù)雜的邊界加載條件和多樣的可視化后處理操作。

(3)本發(fā)明所述方法可以更進(jìn)一步地考慮材料的真實(shí)微觀組織形貌，從相組織模型發(fā)展為晶粒組織模型，從而可以更深層次的考慮材料變形過程中存在的內(nèi)部物理機(jī)制；采用本發(fā)明提供的方法可以計(jì)算得到晶粒組織有限元網(wǎng)格模型中每個(gè)晶粒的隨著加載過程不斷變化的各向異性線彈性本構(gòu)，對(duì)再現(xiàn)材料組織演變過程、分析材料局部受力狀態(tài)，以及預(yù)測(cè)材料失效形式等提供了更為準(zhǔn)確的研究方法。

附圖說明

圖1為本發(fā)明所述方法的流程示意圖；

圖2為所述晶粒組織的有限元網(wǎng)格模型的坐標(biāo)系(左)和晶粒的坐標(biāo)系(右)示意圖；

圖3為實(shí)施例1中所述Vocé硬化指數(shù)含義示意圖；

圖4為實(shí)施例1中所述施密特因子定義示意圖。

具體實(shí)施方式

下面結(jié)合附圖和具體實(shí)施例來詳述本發(fā)明，但不限于此。

實(shí)施例1

如圖1所示，一種基于晶體滑移機(jī)制的各向異性線彈性本構(gòu)的建立方法，所述方法包括以下步驟：

步驟一：

基于待求各向異性線彈性本構(gòu)材料的晶粒組織圖像，基于有限元軟件LS-DYNA，采用有限元網(wǎng)格剖分方法，建立所述晶粒組織的有限元網(wǎng)格模型，并確定所述有限元網(wǎng)格模型的邊界條件和加載條件；

所述加載條件包括載荷和加載步數(shù)，其中，初始加載步數(shù)為0；

所述晶粒組織的有限元網(wǎng)格模型的坐標(biāo)系和晶粒的坐標(biāo)系示意圖如圖2所示；

步驟二：

采用Fortran程序讀取晶粒組織的本征參數(shù)、初始應(yīng)力應(yīng)變狀態(tài)分量和晶粒取向的初始?xì)W拉角并根據(jù)應(yīng)變?cè)隽康耐队瓣P(guān)系，所述歐拉角增量由晶粒的應(yīng)變?cè)隽堪匆韵鹿接?jì)算得出：

并將所述歐拉角增量與初始?xì)W拉角相加，得到由于應(yīng)變引起晶粒取向變化后的歐拉角根據(jù)變化后的歐拉角獲得有限元網(wǎng)格模型的坐標(biāo)系OXYZ與晶粒坐標(biāo)系oxyz之間的轉(zhuǎn)化矩陣R；

其中，所述本征參數(shù)包括晶體結(jié)構(gòu)類型參數(shù)、單晶體彈性剛度矩陣、晶粒內(nèi)的滑移系參數(shù)；

所述晶體結(jié)構(gòu)類型參數(shù)包括晶體結(jié)構(gòu)類型和晶胞的軸比c/a；所述晶體結(jié)構(gòu)類型用ihcp表示，ihcp＝1表示為密排六方晶體結(jié)構(gòu)(hexagonal close-packed,HCP)，其晶向指數(shù)和晶面指數(shù)均用四坐標(biāo)表示；ihcp＝0表示為立方晶體結(jié)構(gòu)(CUBIC)，其晶向指數(shù)和晶面指數(shù)用三坐標(biāo)；所述晶胞的軸比c/a用rca表示，以TC6鈦合金中的α和β兩相為例，對(duì)于其中的α相，rca＝1.597；對(duì)于其中的β相，rca＝1.0。

所述單晶體彈性剛度矩陣C_ij，為6*6的矩陣，單位為GPa，對(duì)于密排六方晶體結(jié)構(gòu)具有對(duì)稱性C₁₁＝C₂₂，C₁₃＝C₂₃，C₄₄＝C₅₅＝C₆₆；而對(duì)于立方晶體結(jié)構(gòu)則為C₁₁＝C₂₂＝C₃₃，C₁₂＝C₁₃＝C₂₃，C₄₄＝C₅₅＝C₆₆；其中，TC6鈦合金中的α相，其晶體結(jié)構(gòu)類型參數(shù)為：hexagonal close-packed(HCP)，晶體結(jié)構(gòu)類型ihcp＝1，晶胞的軸比rca＝1.597，單晶體彈性剛度矩陣為：

對(duì)于TC6鈦合金中的β相，其晶體結(jié)構(gòu)類型參數(shù)為：CUBIC，晶體結(jié)構(gòu)類型ihcp＝1，晶胞的軸比rca＝1.0，單晶體彈性剛度矩陣為：

所述滑移系參數(shù)包括滑移系的名稱、等效滑移系的數(shù)量和滑移系的Vocé硬化指數(shù)；晶粒內(nèi)所有需要考慮的滑移系，包括文件中滑移系的總數(shù)量、需要讀取的滑移系的總數(shù)量(Num_slip)；

滑移系的名稱：對(duì)于密排六方晶體結(jié)構(gòu)來說常見的滑移系有四個(gè)分別是：{0002}<11-20>基面滑移系、{10-10}<11-20>柱面滑移系、{10-11}<11-20>一階錐面滑移系和{10-11}<11-23>二階錐面滑移系，對(duì)于體心立方晶體結(jié)構(gòu)來說常見的滑移系有三個(gè)：{1-10}<111>，{11-2}<111>和{12-3}<111>；

等效滑移系的數(shù)量用(m_slip)表示；

滑移系的Vocé硬化指數(shù)：其包括滑移系開動(dòng)的臨界分切應(yīng)力值τ₀、滑移系從開動(dòng)到累積切應(yīng)變不再增加時(shí)的應(yīng)力增量τ₁、滑移系的初始硬化率θ₀、滑移系的飽和硬化率θ₁、需要讀取的滑移系間的相互作用系數(shù)h^ij(i＝j(luò)表示自硬化系數(shù)，i≠j表示互硬化)；具體含義由應(yīng)力應(yīng)變曲線表示，如圖3所示；在Fortran程序中將所有需要讀取的滑移系的所有等效滑移系順次標(biāo)號(hào)；以滑移系{0001}<11-20>為例，有6個(gè)等效滑移系(包含負(fù)向滑移)，如表1所示；

表1晶體內(nèi)滑移系參數(shù)設(shè)置

所述晶粒取向歐拉角采用Bunge方法定義，即將晶粒坐標(biāo)系繞z軸，在xy平面內(nèi)轉(zhuǎn)過角；再繞x軸，在yz平面內(nèi)轉(zhuǎn)過φ角；最后再一次繞z軸，在xy平面內(nèi)轉(zhuǎn)過角，最終與有限元網(wǎng)格模型坐標(biāo)系完全重合，由此得到每個(gè)晶粒的三個(gè)空間取向歐拉角晶粒取向的初始?xì)W拉角利用EBSD實(shí)驗(yàn)測(cè)量得到；

所述應(yīng)力應(yīng)變狀態(tài)分量分別按照σ_xx，σ_yy，σ_zz，σ_yz，σ_xz，σ_xy和ε_xx，ε_yy，ε_zz，ε_yz，ε_xz，ε_xy的順序給出，初始值均設(shè)為0；該文件中包含晶粒個(gè)數(shù)用變量n_grain存儲(chǔ)；

其中，所述σ_xx，σ_yy和σ_zz分別表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)沿X軸、Y軸和Z軸方向的正應(yīng)力；

所述ε_xx，ε_yy和ε_zz分別表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)沿X軸、Y軸和Z軸方向的正應(yīng)變；

所述σ_yz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在Y軸法向平面內(nèi)平行于X軸的切應(yīng)力；

所述σ_xz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸法向平面內(nèi)平行于Z軸的切應(yīng)力；

所述σ_xy表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸法向平面內(nèi)平行于Y軸的切應(yīng)力；

所述ε_yz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在Y軸和Z軸兩方向上微小線段夾角的切應(yīng)變；

所述ε_xz表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸和Z軸兩方向上微小線段夾角的切應(yīng)變；

所述ε_xy表示在步驟一所述有限元網(wǎng)格模型的坐標(biāo)系中，晶粒所受到的空間內(nèi)在X軸和Y軸兩方向上微小線段夾角的切應(yīng)變；

步驟三：

基于步驟二所述轉(zhuǎn)化矩陣，通過張量運(yùn)算，得到有限元網(wǎng)格模型坐標(biāo)系下的單晶體彈性剛度矩陣Cij′和施密特因子α，其中，所述α＝cosλcosΦ，Φ為滑移面與外力F中心軸的夾角，所述λ為滑移方向與外力F的夾角，如圖4所示；

所述張量運(yùn)算公式為：Cij′＝R^-1·Cij；

步驟四：

判斷歐拉角增量是否為(0，0，0)，若進(jìn)入步驟五；若進(jìn)入步驟六；

步驟五：

根據(jù)步驟一所述有限元網(wǎng)格模型的邊界條件和加載載荷對(duì)當(dāng)前有限元網(wǎng)格模型進(jìn)行加載，然后針對(duì)當(dāng)前有限元網(wǎng)格模型坐標(biāo)系下的單晶體彈性剛度矩陣，采用有限元計(jì)算方法計(jì)算得到加載后晶粒的應(yīng)力應(yīng)變狀態(tài)分量；將所述加載后晶粒的應(yīng)力應(yīng)變狀態(tài)分量替換步驟二中的初始應(yīng)力應(yīng)變狀態(tài)分量，重復(fù)步驟二至步驟四；其中，每次加載后加載步數(shù)加一；步驟六：

判斷晶粒內(nèi)的潛在要開動(dòng)和已經(jīng)開動(dòng)滑移系的總個(gè)數(shù)n是否大于零；若n＝0，表示該晶粒仍處于彈性變形階段，則重復(fù)步驟五；若n＞0，表示該晶粒進(jìn)入塑性變形階段，則進(jìn)入步驟七；

潛在要開動(dòng)和已經(jīng)開動(dòng)滑移系的判斷方法：

若則該滑移系被標(biāo)記為潛在要開動(dòng)，若則該滑移系被標(biāo)記為已經(jīng)開動(dòng)；其中，σ_c表示滑移系所在晶粒當(dāng)前的應(yīng)力值，表示第i個(gè)滑移系的臨界分切應(yīng)力值；所述αⁱ表示第i個(gè)滑移系的施密特因子；

步驟七：

根據(jù)Hill-Hutchinson多晶體彈塑性變形理論，認(rèn)為每個(gè)晶粒的塑性應(yīng)變是所有滑移系應(yīng)變量(γⁱ)的累加值：

根據(jù)胡克定律，每一個(gè)加載狀態(tài)都對(duì)應(yīng)著一個(gè)聯(lián)系晶粒應(yīng)力與應(yīng)變的瞬時(shí)模量(L_c)，即晶粒的瞬時(shí)應(yīng)力(σ_c)與應(yīng)變(ε_c)關(guān)系滿足：

其中，晶粒的總應(yīng)變(ε_c)為彈性應(yīng)變和塑性應(yīng)變的加和，代入(1)和(2)得：

其中M_c為L_c的逆矩陣。

再根據(jù)Hill(1996)的理論滑移系的臨界開動(dòng)應(yīng)力值(τⁱ)與滑移系應(yīng)變量(γⁱ)存在如下關(guān)系：

將(1)～(4)代入分切應(yīng)力的關(guān)系式：

τⁱ＝σ_cαⁱ (5)

中得：

其中，最終得到的n個(gè)滑移系應(yīng)變線性方程組的基向量矩陣X^ij為：

X^ij＝h^ij+αⁱL_cα^j (9)

將(9)帶入(8)，求解分配給每個(gè)潛在滑移系的切變量γⁱ：

γⁱ＝fⁱε_c (10)

其中，

Y^ij為X^ij的逆。

最終在原有的晶粒應(yīng)力應(yīng)變瞬時(shí)模量(L_c)的基礎(chǔ)上求解出考慮了塑性變形階段的瞬時(shí)剛度矩陣

步驟八：

判斷當(dāng)前加載步數(shù)是否等于設(shè)定的最大加載步數(shù)，如果等于最大加載步數(shù)，則步驟七得到的瞬時(shí)剛度矩陣即為晶粒的各向異性線彈性本構(gòu)；如果小于最大加載步數(shù)，則將步驟七的瞬時(shí)剛度矩陣替換步驟二中的單晶體彈性剛度矩陣，重復(fù)步驟二至步驟七，直至加載步數(shù)等于步驟一中設(shè)定的最大加載步數(shù)。

(1)假設(shè)取向的8晶粒模型，其中4個(gè)為α相晶粒，4個(gè)為β相晶粒。有限元以10^-3應(yīng)變速率每步加載0.1％的應(yīng)變，經(jīng)過10步加載后可以判斷出α相晶粒中的4個(gè)晶粒全部進(jìn)入塑性變形階段，β相中有3個(gè)晶粒進(jìn)入塑性變形階段，從輸出結(jié)果可首先對(duì)每個(gè)晶粒在當(dāng)前加載步下是否進(jìn)入塑性階段做出判斷：“進(jìn)入塑性階段”或者“沒有進(jìn)入塑性階段”；其次可以得知晶粒中哪些滑移系已開動(dòng)以及每個(gè)滑移系對(duì)應(yīng)的標(biāo)號(hào)及開動(dòng)量等信息，所有晶粒內(nèi)滑移系的開動(dòng)信息如表2和表3所示：

表2假設(shè)取向的8晶粒模型中鈦合金α相晶粒內(nèi)滑移系的開動(dòng)結(jié)果信息

表3假設(shè)取向的8晶粒模型中鈦合金β相晶粒內(nèi)滑移系的開動(dòng)結(jié)果信息

其中以α相的晶粒1為例，其中共有6個(gè)滑移系開動(dòng)，分別是第3個(gè)、第5個(gè)、第7個(gè)、第29個(gè)、第39個(gè)和第40個(gè)，每個(gè)滑移系對(duì)應(yīng)的開動(dòng)量分別為：0.1155E-01、0.4345E-02、0.1872E-01、0.2221E-01、0.3231E-02和0.1375E-01；

同時(shí)還給出了每個(gè)晶粒的各向異性線彈性本構(gòu)，即6*6的Cij矩陣(單位為GPa，見表4和表5)，經(jīng)對(duì)比后可以發(fā)現(xiàn)，滑移系開動(dòng)數(shù)量越多、開動(dòng)量越大的晶粒，各向異性線彈性本構(gòu)變化越大。

表4假設(shè)取向的8晶粒模型相中鈦合金α相每個(gè)晶粒的C_ij矩陣

表5假設(shè)取向的8晶粒模型相中鈦合金α相每個(gè)晶粒的C_ij矩陣

(2)假設(shè)取向的108晶粒模型，其中70個(gè)為α相晶粒，38個(gè)為β相晶粒。有限元以10^-3應(yīng)變速率每步加載0.1％的應(yīng)變，經(jīng)過10步加載后可以判斷出α相晶粒中的61個(gè)晶粒進(jìn)入塑性變形階段，β相中有28個(gè)晶粒進(jìn)入塑性變形階段。與(1)類似，可首先對(duì)每個(gè)晶粒在當(dāng)前加載步下是否進(jìn)入塑性階段做出判斷，模型中的晶粒數(shù)量越多，配分于每個(gè)晶粒上的應(yīng)力越接近于真實(shí)受力情況，此時(shí)晶粒除了有“進(jìn)入塑性階段”和“沒有進(jìn)入塑性階段”兩種情況外，還有可能出現(xiàn)某個(gè)晶粒在上一步加載時(shí)進(jìn)入了塑性階段，但該晶粒中的滑移系硬化效應(yīng)致使當(dāng)前加載步增加的應(yīng)力不足以使其中的任何一個(gè)滑移系再次開動(dòng)，在計(jì)算過程中體現(xiàn)為基向量方程組無解，即在當(dāng)前加載步下該晶粒的狀態(tài)為“由塑性回到彈性”；其次，可以得知晶粒中哪些滑移系已開動(dòng)及每個(gè)滑移系對(duì)應(yīng)的標(biāo)號(hào)及開動(dòng)量等所有晶粒內(nèi)滑移系的開動(dòng)信息；最終還可以給出每個(gè)晶粒的各向異性線彈性本構(gòu)，即6*6的C_ij矩陣(單位為GPa)，列舉α相第27個(gè)晶粒的滑移系開動(dòng)結(jié)果信息和各向異性線彈性本構(gòu)C_ij矩陣信息如表6所示。

表6假設(shè)取向的108晶粒模型中α相第27個(gè)晶粒的滑移系開動(dòng)結(jié)果信息和各向異性線彈性本構(gòu)C_ij矩陣

本發(fā)明包括但不限于以上實(shí)施例，凡是在本發(fā)明精神的原則之下進(jìn)行的任何等同替換或局部改進(jìn)，都將視為在本發(fā)明的保護(hù)范圍之內(nèi)。