晶體塑性模型的案例
熱-彈-黏塑性晶體塑性模型文章推薦
傳統室溫本構模型通常需要依賴大量不同溫度、不同加載路徑下的實驗數據進行擬合，很難真正解釋“溫度如何影響晶體滑移和多晶塑性響應”。
Cyr 等人針對這一問題提出了一個三維熱-彈-黏塑性晶體塑性模型，即 TEV 模型，用于描述 FCC 多晶材料，特別是 AA5754 鋁合金在升溫條件下的力學行為。該模型的核心思想是：材料變形不僅包含彈性變形和晶體塑性滑移，還需要顯式考慮熱膨脹變形。因此，總變形梯度被分解為彈性/剛體轉動部分、熱變形部分和塑性變形部分。
在本構層面，作者保留了 FCC 晶體的 12 個 {111}<110> 滑移系，并采用冪律型滑移率方程描述率相關塑性流動。與常規晶體塑性模型不同的是，該模型把溫度效應系統地引入到多個關鍵物理量中：首先，單晶彈性常數 C11、C12、C44 隨溫度變化；其次，滑移阻力引入熱軟化函數，用來描述溫度升高后滑移更容易發生的現象；再次，單滑移硬化參數也被寫成溫度函數，包括參考臨界分切應力、初始硬化率和硬化指數。
這個模型的優勢在于，它不是簡單地給宏觀應力-應變曲線加一個溫度修正系數，而是從晶體滑移層面描述溫度對材料響應的影響。換句話說，它可以同時分析宏觀應力變化、微觀滑移活動、織構演化、局部應變集中和熱軟化機制。因此，它比普通經驗型熱塑性模型更適合用于多晶材料溫成形模擬。
作者首先利用 AA5754 鋁合金在 25 ℃、148 ℃、204 ℃ 和 232 ℃ 下的單軸拉伸實驗數據標定溫度相關硬化參數。隨后，又預測了 177 ℃ 和 260 ℃ 下的拉伸響應。
展開 黃永剛晶體塑性模型耦合相場方法模擬多晶斷裂
模擬包含200個晶粒的多晶模型，使用平面應變簡化，拉伸變形30%，模擬拉伸過程中裂紋的產生和發展，其中斷裂總能量包含彈性變形能和塑性耗散功兩部分，模擬結果如下
初始多晶模型：
網格劃分（CPE4網格）：
相場分布（0：材料完好，1：材料完全失效）：
退化程度分布：
歐拉角（phi）分布：
可以看到耦合相場的晶體塑性模型具有潛在的預測裂紋萌生和發展的能力，其準確程度取決于斷裂能參數的選擇，與更精細的實驗對比，如原位的ebsd拉伸將成為良好的校核手段。這可能成為介觀尺度下斷裂力學的應用提供良好的參考
晶體塑性每日文章推薦（三）
文章doi：10.1016/j.ijplas.2009.11.004
推薦理由：作者基于Orowan硬化方程提出了考慮基于位錯的晶體塑性模型（SCCE-D）這通過修改傳統的模型中硬化實現，并于廣泛使用的唯象單晶本構模型（SCCE-T）模型進行了對比，通過實驗拉伸曲線的數據值反演得到適合兩種模型的最佳參數，通過兩種模型與實驗中單晶變形進行比較，結果表明，基于位錯的本構模型在單晶變形預測中精度更高，作者分析認為，基于位錯密度的模型更高反應力單晶變形的物理過程，而傳統唯象模型反應的更接近于多晶的平均特征
作者的研究思路
一，從經典的唯象模型出發，根據Orowan硬化理論，推導基于位錯的晶體塑性模型，兩類模型采用相同的流動方程，修改體現在滑移系的硬化方程中，昨日推薦的文獻則在硬化和流動都進行了對應的修改，兩類模型的硬化方程分別為：
(1)唯象模型（SCCE-T）：
(2)位錯密度模型（SCCE-D）：
二，研究FCC-Cu和BCC-Fe兩種材料，基于特定方向擬合兩組本構方程得到最適合的材料參數，并應用于不同取向單晶變形行為預測,發現基于位錯密度的模型預測結果更接近于實驗結果
擬合圖（Cu）：
預測圖：
擬合圖（Fe）
預測圖
三，作者使用RVE比較了基于單晶擬合獲得的參數在多晶預測結果的情況，法向唯象的模型預測結果顯著高于位錯密度的模型，說明宏觀變形預測對硬化模型十分敏感，并進一步比較兩種方法在織構演化預測方面的差異，發現結果無明顯差異，因此可以認為織構演化對硬化模型不敏感
對作者分析感興趣的可對huang程序進行簡單修改，并根據對應參數復現作者的研究結果
展開 一個好用的Abaqus晶體塑性模型生成插件-Voronoi模型V8.0
V7.0版本介紹：
一個好用的Abaqus晶體塑性模型生成插件-Voronoi模型
https://zhuanlan.zhihu.com/p/611427546
2. V8.0版本新增功能：
2.1 二維核殼晶體模塊
圖2.1 二維圓形核殼晶體模塊
圖2.2 二維多邊形核殼晶體模塊
2.2 三維核殼晶體模塊
圖2.3 三維球形核殼晶體模塊
圖2.4 三維多面體核殼晶體模塊
2.3 桁架模型模塊
圖2.5 桁架結構模型生成模塊
2.4 圓形和圓柱邊界加權晶體模塊
圖2.6 二維圓形邊界加權晶體模塊
圖2.7 三維圓柱邊界加權晶體模塊
2.5 二維梯度晶體模塊
圖2.8 二維梯度晶體模塊
2.6 三維圓柱邊界梯度模塊
圖2.9 三維圓柱邊界梯度晶體模塊
展開 
基于HUANG晶體塑性模型的三點彎曲案例------案例十
?基于HUANG晶體塑性模型的三點彎曲案例
案例教學如下
1，建立幾何模型：0.6*0.3*0.05（mm）的平板，半徑0.05，長度0.3的下壓棍和支撐并完成裝配
三點彎曲模型圖
2，分配材料屬性：分別賦予板材晶體屬性參數，棍子純彈性參數，并假設下壓輥和支撐的剛度遠遠大于平板，即變形過程中約束為剛體。
模型的材料屬性分配
3，建立接觸條件：板材與棍子支撐之間法向硬接觸，切向摩擦系數設置為0.2
建立接觸屬性
4，網格劃分對板材和棍子，支撐進行網格劃分，板材網格類型C3D4（保留晶界形狀），棍子支撐C3D8
5，載荷施加：固定下邊支撐，同時對上面的下壓輥施加向下的0.04的位移，其他自由度全部限制為0。
6，提交作業與后處理
三點彎曲位移分布情況
累計塑性應變分布情況
模型的應力分布情況
模型的變形流動方向分布情況
展開 考慮晶界影響的晶體塑性模型
這個思想非常適合和晶體塑性模型結合，因為晶體塑性本來就是逐滑移系計算的。
在作者的模型中，晶界通透性可以進一步轉化為晶界障礙應力。通透性越高，障礙應力越小；通透性越低，障礙應力越大。這樣一來，晶界對塑性滑移的影響就可以直接進入滑移率方程：只有當有效分切應力足夠克服晶界障礙時，晶界附近的滑移才能繼續發展。
這篇文章的另一個重要部分是位錯重分配。作者把有限元單元中的位錯內容等效成“超位錯”，再根據塑性滑移活動對可動位錯進行重新分布，并計算由這些位錯分布產生的背應力。這個處理很有啟發性：它不是直接追蹤每一根真實位錯，那樣計算量太大；但它也不是完全經驗化地加一個強化項，而是在連續體模擬和位錯物理之間做了一個折中。作者的模型概念圖：
積分流程圖：
從結果來看，作者的模型能夠再現單個位錯塞積問題中的位錯密度分布；在二維和三維規則晶粒陣列中，也能預測出與實驗同量級的 Hall–Petch 斜率。對于粗晶鐵多晶拉伸響應，這個兩尺度模型比傳統 CP-FEM 或 Taylor 類模型給出了更好的預測。此外，作者還比較了 Fe-3%Si 柱狀晶樣品中的晶格曲率，模型預測的晶界附近曲率峰值與實驗結果基本一致。
作者發現模型可以非常準確的預測晶粒尺寸效應：
我認為這篇文章的價值不只是“提出了一個更復雜的模型”，而是提供了一種很清楚的建模思路：晶界強化不一定只能通過經驗晶粒尺寸項來描述，也可以從滑移傳遞、位錯通量和局部障礙應力出發，逐步把晶界的物理作用放進晶體塑性框架中。
我們可以作者提出的模型完整的構建一個考慮晶界多尺度模型，演示如何計算每個滑移系對應的晶界通透系數，并將其轉化為晶界障礙應力引入晶體塑性本構中。通過對比是否考慮晶界障礙項，可以觀察晶界附近滑移活動、位錯密度分布以及應變局部化特征的變化。
展開 晶體塑性：構建Dream3D pipeline用于將EBSD模型制作成Abaqus可執行文件 ￥180
晶體塑性：構建Dream3D pipeline用于將EBSD模型制作成Abaqus可執行文件
案例實操
用于生成模型的Dream3D pipeline文件，只需要你設置EBSD數據的路徑和導出路徑即可，可以直接生成abaqus的晶體塑性模型，提供原始文件！
包含老版本Dream3D 6.5的管道文件，并且根據官方的使用說明文件已經成功移植到最新版Dream3D 7.4版本了。
晶體塑性每日文章推薦（二十三）
發現考慮傳統模型可以捕捉FCC結構的應力應變行為，但無法很好的預測低層錯能金屬的織構特征，尤其是平面應變壓縮或冷軋變形（低層錯能金屬通常獲得的織構為brass型織構，而中高層錯金屬通?？棙嫗閏oppor型織構），作者分析認為對于低層錯能的FCC金屬（如黃銅等），考慮微剪切帶的建模對于織構預測的準確性是至關重要的，即黃銅型織構和銅型織構之間的主要差異是由于剪切帶而非變形孿晶引起的。然而，變形孿晶在影響應變硬化響應中起著重要作用，也可能在微尺度剪切帶的萌生中起到重要作用。因此提出了一個考慮滑移+孿生＋微剪切帶的FCC結構的晶體塑性模型，并成功預測了低層錯能金屬的織構特征，該方案被大量引用是目前FCC結構考慮孿生+剪切帶的主流晶體塑性模型，并已被集成到damask軟件平臺。作者使用的模型介紹如下：
晶體滑移區域與孿晶區域的應力分別為
其中L表示四階彈性剛度張量，E是柯西格林應變，T是PK2應力
兩個區域對應的四階彈性剛度張量的關系為
晶體架構下的柯西應力可以簡單的表示為：
其中f是每個孿生系統對應的體積分數
速度梯度的貢獻包含三部分：
對于FCC結構，滑移最有可能局限于共面或近共面滑移系統。因此，預計Ns-tw將遠小于Ns，故右側的第三項通常不考慮
考慮孿生的極性，孿生體積分數的演化表示為：
此外，滑移和孿生的應變硬化表示為：
滑移：
孿晶：
使用作者提出的滑移+孿晶模型預測的平面應變壓縮織構為：
考慮剪切帶效應的修正模型：
考慮滑移+孿晶+剪切帶效應模型預測的織構結果：
作者最終分析得到的結論為：，銅型和黃銅型織構之間差異的主要成分是由于微觀尺度剪切帶，而不是由于變形孿晶。
展開 晶體塑性每日文章推薦（五）
文章名稱：《Comparison of the implicit and explicit finite element methods using crystal plasticity》
doi：10.1016/j.commatsci.2006.08.002
推薦理由：作為顯式晶體塑性與隱式晶體塑性模型比較，以及適用性討論最經典的文章之一，詳細介紹了顯式與隱式求解器的區別，其研究結果表明：對于更簡單的加載條件，隱式方法的求解時間更短。在涉及接觸的載荷條件下，顯式方法被證明是優選的選擇。與使用多個處理器解決分析的隱式方法相比，顯式方法顯示出持續高水平的并行化效率。
眾所周知：在隱式方法中通常使用牛頓迭代方法，解涉及迭代，直到每個增量都滿足收斂標準。因此計算收斂于精確解（隱式：基于t+Δt時刻確定t+Δt時刻的狀態變量），而顯式方法中的有限元方程使用的通常是歐拉向前方法，并且在這種形式下，它們可以直接求解，以在增量結束時確定解，當增量步小于穩定時間增量時，總是可以保證計算的進行，但結果不一定收斂于精確解，（顯式：基于t時刻確定t+Δt時刻的狀態變量）因此準靜態問題使用隱式往往更加高效，而涉及到接觸則使用顯式可以保證計算穩定
由于晶體塑性模型考慮了介觀尺度的塑性變形的真實物理過程，因此被廣泛用于已在模擬金屬和金屬基材料中的介觀尺度下大變形和應變局部化，然而其高度非線性的積分過程，相較于傳統模型，其數值成本往往很高，因此基于復雜的晶體塑性模型更能體現顯式于隱式積分的差異，這里作者的討論顯式程序和隱式程序使用huang的亞彈性框架進行
由于隱式程序在大量的博士論文可以找到詳細內容，這里不做贅述，這里主要提到文獻關于顯式的一些內容。
顯式計算存在臨界穩定時間的概念，當最大增量步大于該時刻時會導致計算結果發散，從而使得計算結果失去意義。
展開 晶體塑性每日文章推薦（十八）
文章doi：10.1016/j.ijplas.2009.10.004
上一篇推文介紹了基于L2范數最小化計算GND的推文，基于類似的思想可以實現率無關數值模型的構建，作者的創作思想就是利用塑性變形過程中最大能量耗散原理，將單晶屈服問題視作約束優化問題，其中約束就是每一個滑移系統的屈服函數，并將傳統率無關晶體塑性模型中的數值奇異問題，通過創建約束條件的組合進行優化分析得到率無關的晶體塑性數值解。該率無關的本構計算方法計算效率相對較高，且數值穩定性很好，與以往研究和實驗結果具有良好的一致性。此外也有很多率相關模型使用奇異值分解進行數值求解。
作者的思想
彈塑性問題通常被定義為約束優化問題，旨在尋找給定應變增量的最佳應力張量和內部變量。在這樣的問題中，目標函數是基于最大耗散原理定義的，約束是屈服函數。示意圖和公式為：
以塑性變形率方程為切入點：
λ為一致性參數
在單晶中，整體塑性變形是多個滑移系統上滑移的結果。在晶體塑性問題中，變形由多個屈服面定義，這些屈服面不平滑相交，屈服函數的數量取決于晶體中滑移系的數量。
展開 梯度晶體塑性模型對應的umat子程序 ￥1200
文獻一：《Gradient plasticity in gradient nano-grained metals》
文獻二：《Grain rotations during uniaxial deformation of gradient nano-grained metals using crystal plasticity finite element simulations》
推薦理由：兩篇文章使用了類似的研究方法，通過構建具有梯度分布的晶粒模型，基于原始的唯象晶體塑性模型進行修改，將初始屈服，硬化模量，飽和強度，以及率相關系數構造為晶粒尺寸的函數，實現建立具有尺寸效應的多晶本構模型，這對目前金屬梯度結構介觀尺度下力學性能的表征具有一定的啟發性
文獻一的研究使用Voronoi鑲嵌方法構建梯度納米晶結構，使用的本構模型如下：
流動方程：
硬化方程為：
通過假設：單晶水平上的所有抗滑移參數與局部晶粒尺寸D的平方根成反比
修正對應的參數為：
其中彈性參數對應Cu的參數
有限元模型為：
研究了平面應變條件下簡單拉伸不同區域的應力應變分布特征
CPFE結果揭示了GNG-Cu橫截面中的梯度應力和梯度塑性應變。這些空間梯度是由于在具有梯度尺寸的晶粒中逐漸達到屈服點以及相應的梯度滑動阻力而產生的。
CPFE結果還揭示了梯度應力和梯度塑性應變的非均勻空間分布，這是隨機晶粒取向和梯度晶粒尺寸共同作用的結果。
展開 
多尺度晶體塑性模擬文章推薦
這篇文章對我們的啟發在于：晶體塑性并不只能用于單晶拉伸、RVE 或微觀變形分析，也可以嵌入顯式動力學框架，用于研究真實工程結構中的局部變形、吸能和織構演化。對于高溫合金、鋁合金薄壁件、微尺度構件等問題，如果材料存在明顯織構或晶粒尺度效應，將晶體塑性與結構有限元耦合，能夠提供比傳統本構更豐富的物理信息。
我們可以將我之前推文提到的umat-taylor模型轉化為vumat子程序，進一步使用晶體塑性模型模擬大變形結構尺度材料變形行為。案例展示如下：
初始模型參考文章的設置（上下兩層鋼板，中間為薄殼結構）：
使用通用接觸，摩擦系數設置為0.5，共4000個單元，每個單元包含50個具有不同初始取向晶粒。共20萬晶粒。
邊界條件設置為下端鋼板固定，上端下壓。
模擬結果如下：
應力分布結果：
晶粒1的剪切滑移：
晶粒2的剪切滑移：
晶粒50的剪切滑移：
單元標號5變形結束后的50個歐拉角分布：
展開 塑性工程學報：Custom450鋼拉伸的晶體塑性有限元分析
而研究微細觀尺度的變形不均勻性是新材料開發及優選的重要準則，晶體塑性有限元方法將晶體塑性理論和有限元軟件進行了恰當的融合，成為研究細觀層次塑性變形行為的強有力工具。來自華東理工大學機械與動力工程學院的艾鑫團隊，基于Voronoi方法建立了Custom450 鋼拉伸的二維晶體塑性模型，分析了初始硬化模量、參考剪切應變率、應變率敏感系數、初始屈服應力以及飽和流動應力對材料應力——應變曲線的影響，并對晶體塑性參數進行了標定。
在文獻中，作者所建立的單晶本構模型參考了HUANGY的單晶體模型的子程序UMAT，此率相關硬化晶體塑性模型需要確定的參數包括初始硬化模量h0、初始屈服應力τ0、參考剪切應變率γ，應變率敏感系數n和飽和流動應力τs，其他參數通過計算和查找文獻獲得。基于Voronoi方法，作者在有限元軟件Abaqus中建立了Custom450材料的多晶體二維幾何模型并將本構關系嵌入軟件中，進行拉伸過程的模擬。
圖1所示是微結構模型及其網格劃分，幾何模型尺寸長度為0. 2 mm，寬度為0. 5mm，共包含100個晶粒，大小和形狀隨機，且晶粒取向隨機分布。
圖1包含100個晶粒的微結構模型及其網格劃分
圖2是邊界條件的約束情況，模型的上端面和下端面的所有節點在y方向上具有均勻的位移，左側所有節點在x方向上設置約束，使其不能橫向移動，y方向自由，在右邊界施加載荷，右側的所有節點x方向上經受同等應變載荷，而在y方向上是自由的。
圖2邊界條件示意圖
對于體心立方晶體來說，3個滑移系包括1個主滑移系和2個次滑移系。分別對包含1、2、3組滑移系開動的情形進行模擬，結果如圖3所示，只有主滑移系 { 110} < 111 >啟動時，應力——應變曲線在彈塑性區間過渡的位置存在明顯拐點，并與試驗曲線吻合良好。
展開 基于密西西比州立大學晶體塑性模型預測不同變形下織構演化案例教學------案例八 ￥99
基于密西西比州立大學晶體塑性模型預測不同變形下織構演化
官方使用原始案例
案例一，單向壓縮75%（FCC）
加載條件
織構演化結果
1， 案例二，單向拉伸75%（FCC）
加載條件
織構演化結果
密西西比州立大學晶體塑性有限元代碼和黃永剛院士的程序一樣，均是開源代碼，可免費獲得，并且同時可以考慮FCC，BCC，HCP的滑移和孿晶變形，有著廣泛的應用，目前該代碼已經集成到FEPX計算軟件中，支持并行運算，計算效率很高。本案例采用該代碼，研究FCC，BCC兩種結構在單向拉伸，壓縮，平面應變壓縮等75%的變形量下織構的演變（需要注意的是，這個代碼的輸出使用的Kocks輸出表示取向，為了使用方便，已經在程序中修改bunge標號，可用MTEX直接繪制極圖）
使用包含500個隨機取向的單元預測取向演化
初始隨機取向
一：FCC織構演化
單向壓縮75%的取向分布
單向拉伸75%的取向分布
平面應變壓縮75%取向分布
一：BCC織構演化
單向壓縮75%取向分布
單向拉伸75%取向分布
平面應變壓縮75%取向分布
其中FCC和壓縮和拉伸與官網所提供的案例保持一致，FCC，BCC的平面應變壓縮與已有文獻的典型織構一致，完全正確。
展開 低率相關性的晶體塑性模型，不同迭代方案計算時間的差異性
常規的唯象晶體塑性模型的流動方程通常使用冪律形式：
其中m為率相關系數，對于較小的m值，如≤0.01，整體的響應結果被認為接近率無關響應，然而該參數顯著影響積分效率，對于不同的迭代方案，其對穩定性的影響也不僅相同，這里嘗試進行簡單對比，對比指標和總計算增量步數和計算時間（所有程序均使用單核計算）：
所有迭代方案使用相同的硬化模型和相同材料參數，并對包含200個晶粒的多晶模型進行20%的拉伸變形模擬。如下圖所示：
（1）對于以彈性變形梯度和塑性速度梯度為迭代變量的寫法：
計算時間：
增量步數：
（2）以PK2應力和滑移系統當前強度為迭代變量的雙重迭代全隱式迭代方案：
計算時間：
增量步數：
（3）以PK2應力和滑移系統當前強度為迭代變量的雙重迭代半隱式迭代方案：
計算時間：
增量步數：
（4）以PK2應力和滑移系統當前強度為迭代變量的單次迭代求解方程組全隱式迭代方案：
計算時間：
增量步數：
（5）以滑移系剪切應變為迭代變量的迭代方案：
無法完成模擬的收斂！??！
增量步數：
（6）以柯西應力為迭代變量的迭代方案：
計算時間：
增量步數：
（7）以偏應力為迭代變量的迭代方案：
計算時間：
增量步數：
模擬得到的效果圖：
展開 