1 引言

本部分介紹來自于《斷裂相場法》書籍。

“1998年Francfort和Marigo根據Griffith脆性斷裂理論，提出了一種斷裂力學變分原理，他們以結構內可能的位移場和裂紋面作為自變量，將變形能與斷裂面之和定義為結構總能量，并且認為真實的位移場與裂紋面使得該總能量最小。然而在數值模擬中將離散的裂紋面作為未知量來求解是非常困難的。因此2000年Bourdin等提出了一種相場模型，其中引入了一個連續的標量場，即相場，來近似地描述裂紋。相場值為1和0分別代表材料完全破壞和完好兩種極限狀態，而它們之間的值代表了一種損傷狀態，并且裂紋的彌散程度由相場特征寬度來控制，其值越大彌散寬度越大，反之則越小。然后通過一個與相場相關的裂紋面密度泛函來重構結構內的斷裂能，并將因損傷而退化的變形能與重構的斷裂能代入Francfort-Marigo變分原理就得到了相場模型的基本列式。相場模型中的自變量為兩個連續變化的場，即位移場和相場，因此它可以很方便地由不同數值方法實現。直觀來看，相場模型將一個結構內裂紋萌生與演化問題，轉化為了一個多場耦合情況下求最小能量的優化問題，因此它可以用于直接求解（例如分叉、交叉、融合、扭結等）復雜斷裂問題，而不需要額外的裂紋路徑追蹤方法。”

2 理論

將系統的總勢能表示為如下兩項：

式中第一項能量為：

考慮損傷帶來的退化，彈性能的表達式為：

式中

k為一個小值，用于防止數值不穩定現象。另一項斷裂能為：

因此代入具體表達式可將系統總勢能表達為：

對上述能量進行一階變分可得：

即可得弱形式方程為：

具體外力虛功為：

式中本構方程為：

該弱形式方程是后續推導有限元方程的基礎。同時，通過弱形式方程也可推導得到強形式的控制方程，即位移場和相場的控制方程。對上述弱形式進行分部積分可得：

因次位移場和相場的強形式控制方程為：

以及相應的邊界條件為：

3 有限元離散

為推導有限元離散方程，對位移場和相場控制方程的弱形式進行處理：

對位移場和相場進行插值可得：

m指單元節點的個數。因此相應的梯度場可以插值為：

B矩陣的是由形函數對物理坐標的導數組成的。同理有：

代入到弱形式方程中可得殘值方程；

使用牛頓迭代法求解上述非線性系統。更新格式為:

剛度矩陣為：

為了保證損傷不能愈合，即：

需要做出一些修改，即取歷史上最大的彈性應變能，即：

4 代碼

《斷裂相場法》書中提供了傳統脆性斷裂相場模型的二維UEL代碼。本文將其拓展為三維情況。

UEL需要更新單元的剛度矩陣和右端項，公式在理論部分已詳細給出。

5 測試

5.1 一個單元拉伸破壞

對一個單元施加的邊界條件為，x=0面約束所有位移自由度，即u=0,v=0,w=0;在x=1面進行位移加載，即u=0.1，v=0,w=0。單個單元拉伸破壞時具有解析解。因為相場在單元內是均勻的，因此相場梯度為0，因此可以直接求解相場控制方程得到：

因此可以得到真實應力的表達式為：

數值計算得到的應力應變曲線與解析解的對比結果如下：

取單元上一個積分點，繪制其相場值隨加載時間的變化曲線如下：

5.2 單邊裂紋拉伸

對于單邊裂紋板的拉伸案例，我們選取兩種相場特征裂紋寬度的情況，來展示特征寬度對于結果的影響。

第一個取l=0.015，相場分布圖如下：

第二個取l=0.0075場分布圖如下：

6 參考文獻

[1]. 胡小飛, 張鵬, 姚偉岸. 斷裂相場法. 北京: 科學出版社; 2022.

[2] Kristensen PK, Martínez-Pa?eda E. Phase field fracture modelling using quasi-Newton methods and a new adaptive step scheme. Theoretical and applied fracture mechanics. 2020;107:102446.