[發明專利]一種適于中國古建筑木結構用材本構關系的數值模擬方法

在審

• 申請號: 202011256618.6
• 申請日: 2020-11-11
• 公開（公告）號: CN112541285A
• 公開（公告）日: 2021-03-23
• 發明（設計）人: 楊娜;李天昊
• 申請（專利權）人: 北京交通大學
• 主分類號: G06F30/23
• 分類號: G06F30/23;G06F119/14
• 代理公司: 北京衛平智業專利代理事務所(普通合伙) 11392
• 代理人: 張新利;謝建玲

1.一種適于中國古建筑木結構用材本構關系的數值模擬方法，其特征在于，包括以下步驟：

(1)彈性階段：中國古建筑木結構用材彈性階段用以下公式描述：

其中D_ijkl為用材彈性剛度矩陣，i，j，k，l＝1，2，3；σ₁₁，σ₂₂，σ₃₃為用材的主軸正應力，σ₁₂，σ₂₃，σ₁₃為用材的主面剪應力；ε₁₁，ε₂₂，ε₃₃為用材的主軸正應變，ε₁₂，ε₂₃，ε₁₃為用材的主面剪應變；在后續描述中，LR平面內的應力使用σ進行表達，應變使用ε進行表達；RT平面內的應力使用Δ進行表達，應變使用ξ進行表達；TL平面內的應力使用α進行表達，應變使用β進行表達；

用材彈性剛度矩陣D_ijkl中的參數的表達具體如公式(2)所示；

其中E₁，E₂，E₃為用材主軸彈性模量，G₁₂，G₂₃，G₁₃為用材主面剪切模量，△表示中間過渡計算量，ν₁₂，ν₂₃，ν₁₃，ν₂₁，ν₃₂，ν₃₁為用材主面泊松比；

(2)屈服方程：

A.LR平面屈服方程：

將SunChen模型的衍生模型與能夠描述用材LR平面內受力響應的力學假設相結合，得到LR平面屈服方程：

其中a₁₁，a₂₂，b₁₁，b₂₂為LR平面內用材塑性特征參數，μ為LR平面內用材彈性特征參數，k₁是與用材應力有關的參數，a₁₁，a₂₂，b₁₁，b₂₂通過中國古建筑木結構修繕替換件用材材料偏軸受拉和受壓應力-應變曲線獲取，將應力-應變曲線轉換為應力-塑性應變曲線，使用如下公式獲取有效應力-有效塑性應變曲線；

其中σ_x，分別為用材LR平面內受拉或受壓試驗中沿加載方向的試驗應力和試驗塑性應變，和分別為該模型中的有效應力和有效塑性應變，θ為用材偏軸受拉或受壓試驗中的偏軸角度，偏軸角度為L軸與受力軸間的夾角，h(θ)表示中間過渡計算量，通過試值法確定a₁₁，a₂₂，b₁₁，b₂₂的值，將上述有效應力-有效塑性應變曲線擬合成為一條冪函數主曲線求得最優解，擬合過程中的決定系數R²最大，冪函數主曲線形式見公式(5)；

B.RT平面屈服方程：

將SunChen模型與能夠描述用材RT平面內受力響應的力學假設相結合，得到RT平面屈服方程；

其中a₄₄為RT平面內用材塑性特征參數，k₂是與用材應力有關的參數，a₄₄通過中國古建筑木結構修繕替換件用材材料偏軸受壓應力-應變曲線獲取，將應力-應變曲線轉換為應力-塑性應變曲線；使用如下公式獲取有效應力-有效塑性應變曲線；

其中Δ_x，為用材RT平面內受壓試驗中沿加載方向的試驗應力和試驗塑性應變，和分別為該模型中的有效應力和有效塑性應變，φ為用材偏軸受壓試驗中的偏軸角度，偏軸角度為R軸與受力軸間的夾角；i(φ)表示中間過渡計算量，通過試值法確定a₄₄的值，將上述有效應力-有效塑性應變曲線擬合成為一條冪函數主曲線求得最優解，擬合過程中的決定系數R²最大，冪函數主曲線形式見公式(8)；

C.TL平面屈服方程：

將SunChen模型的衍生模型與能夠描述用材TL平面內受力響應的力學假設相結合，得到TL平面屈服方程；

其中c₁₁，c₂₂，d₁₁，d₂₂為TL平面內用材塑性特征參數，ν為TL平面內用材彈性特征參數，k₃是與用材應力有關的參數，c₁₁，c₂₂，d₁₁，d₂₂通過中國古建筑木結構修繕替換件用材材料偏軸受拉和受壓應力-應變曲線獲取，將應力-應變關系轉換為應力-塑性應變曲線，使用如下公式獲取有效應力-有效塑性應變曲線：

其中α_x，為用材TL平面內受拉或受壓試驗中沿加載方向的試驗應力和試驗塑性應變，和分別為該模型中的有效應力和有效塑性應變，γ為用材偏軸受拉或受壓試驗中的偏軸角度，偏軸角度為T軸與受力軸間的夾角，j(γ)表示中間過渡計算量，通過試值法確定c₁₁，c₂₂，d₁₁，d₂₂的值，將上述有效應力-有效塑性應變曲線擬合成為一條冪函數主曲線求得最優解，擬合過程中的決定系數R²最大，冪函數主曲線形式見公式(11)；

(3)用材受拉：中國古建筑木結構用材LR平面、RT平面和TL平面受拉時均出現脆性斷裂，整體應力-應變曲線呈現線性，LR平面屈服方程、RT平面屈服方程和TL平面屈服方程描述的屈服面上的點為應力-應變曲線非線性階段起始點或線性階段終點，因此用材受拉狀態下的應力水平達到屈服面時，用材出現斷裂，斷裂時的應變增量對應的應力變為零，采用如下方程進行描述：

假設第n步迭代得到的應力張量為n+1步中給出一個新的應變增量dε_kl，計算得到試探應力：

將試探應力代入屈服方程，屈服方程根據變形出現的平面擬定，分別為LR平面屈服方程、RT平面屈服方程和TL平面屈服方程；若則說明用材仍處于彈性階段，其中n'＝1，2，3，n'＝1時對應LR平面屈服方程中的k₁，n'＝2時對應RT平面屈服方程中的k₂，n'＝3時對應TL平面屈服方程中的k₃，令作為新的應力張量，按照公式(12)進行循環；若則說明用材已經進入屈服階段：用材處于受拉狀態時發生斷裂，保留該應變增量，應變增量步后對應的應力為零；

(4)用材受剪：中國古建筑木結構用材LR平面、RT平面和TL平面受剪時均出現脆性斷裂，整體應力-應變曲線呈現線性，LR平面屈服方程、RT平面屈服方程和TL平面屈服方程描述的屈服面上的點為應力-應變曲線非線性階段起始點或線性階段終點，因此用材受剪狀態下的應力水平達到屈服面時，用材出現斷裂，斷裂時的應變增量對應的應力變為零，同樣采用公式(12)進行描述，同樣將試探應力代入屈服方程，屈服方程根據變形出現的平面擬定，分別為LR平面屈服方程、RT平面屈服方程和TL平面屈服方程，若則說明用材仍處于彈性階段，令作為新的應力張量，按照公式(12)進行循環，若則說明用材已經進入屈服階段：用材處于受剪狀態時發生斷裂，保留該應變增量，應變增量步后對應的應力為零；

(5)用材LR平面受壓：用材LR平面受壓，簡化為3個階段，彈性階段由公式(1)和公式(2)確定；用材LR平面受壓時的應力狀態達到公式(3)、公式(4)和公式(5)確定的屈服面后，用材進入理想塑性流動階段，過程如下：

假設第n步迭代得到的應力張量為n+1步中給出一個新的應變增量dε_ij，根據公式(12)計算得到試探應力若在LR平面中，k_n'＝k₁，則說明用材已經進入屈服階段，該試探應力不是真實用材應力，需要將超出屈服面的部分減去：將總應變增量中的塑性應變增量引起的虛假應力從試探應力中減去，從而使應力值返回屈服面已實現理想塑性流動狀態，具體方程如下：

總應變增量dε_ij由彈性應變增量和塑性應變增量之和進行表達，

若n+1步的真實應力為：

其中為塑性應變增量，公式(14)能夠使用材應力始終處于屈服面上，實現理想塑性流動狀態，其中塑性應變增量根據相關聯的流動法則求解：

其中△λ為一致性系數，根據公式(15)，LR平面上的塑性應變表示為：

其中為用材沿L軸方向的塑性正應變增量，為用材位于LR平面內的塑性剪應變增量，為用材沿R軸方向的塑性正應變增量，△λ_LR求解是上述過程的關鍵，假設用材在第n步屈服，那么用材在第n+1步同樣屈服，應力狀態均在屈服面上，即：

其中為第n步試探應力代入公式(3)的結果，而為第n+1步試探應力代入公式(3)的結果；

故此

其中為與的差值；

將上式按照LR平面屈服方程進行一階泰勒展開，得到：

其中△σ₁₁為沿著用材L軸方向的正應力增量，△σ₂₂為沿著用材R軸方向的正應力增量，△σ₁₂為位于用材LR平面內的剪應力增量；

其中M、X、Y、Z為運算過程中的過渡量；根據徑向返回算法，公式(19)中的應力增量表述為：

其中為沿著用材L軸方向的試探正應力增量，為沿著用材R軸方向的試探正應力增量，為位于用材LR平面內的試探剪應力增量；

將公式(20)和公式(21)代入公式(19)，△λ_LR被求解：

在總應變達到應力軟化閾值時，應力軟化閾值為彈性總應變閾值的1.2倍，出現應力軟化現象；

(6)用材RT平面受壓：用材RT平面受壓，簡化為3個階段，彈性階段由公式(1)和公式(2)確定，用材RT平面受壓時的應力狀態達到公式(6)、公式(7)和公式(8)確定的屈服面后，用材進入理想塑性流動階段，過程如下：

假設第n步迭代得到的應力張量為n+1步中給出一個新的應變增量dξ_ij，根據公式(12)計算得到試探應力若則說明用材已經進入屈服階段，在RT平面中，k_n'＝k₂，該試探應力不是真實用材應力，需要將超出屈服面的部分減去：將總應變增量中的塑性應變引起的虛假應力增量從試探應力中減去，從而使應力值返回屈服面已實現理想塑性流動狀態，總應變增量由彈性應變增量和塑性應變增量之和進行表達，見公式(13)，真實應力計算見公式(14)，相關聯塑性流動法則見公式(15)，RT平面內的塑性應變表示為：

其中用材沿R軸方向的塑性正應變增量，為用材位于RT平面內的塑性剪應變增量，為用材沿T軸方向的塑性正應變增量，△λ_RT的求解是上述過程的關鍵，假設用材在第n步屈服，那么用材在第n+1步同樣屈服，應力狀態均在屈服面上，即：

其中為第n步試探應力代入公式(6)的結果，而為第n+1步試探應力代入公式(6)的結果；故此：

其中為與的差值，將上式按照LR平面屈服方程進行一階泰勒展開，得到：

其中△△₂₂為沿著用材R軸方向的正應力增量，△△₃₃為沿著用材T軸方向的正應力增量，△△₂₃為位于用材RT平面內的剪應力增量；

其中U、V、W為運算過程中的過渡量，公式(29)中的應力增量表述為：

其中為沿著用材R軸方向的試探正應力增量，為沿著用材T軸方向的試探正應力增量，為位于用材RT平面內的試探剪應力增量；

將公式(30)和公式(31)代入公式(29)，△λ_RT被求解：

在總應變達到應力硬化閾值時，應力硬化閾值為彈性總應變閾值的2倍，用材出現應力硬化現象，木材領域統稱為二次強化現象；

(7)用材TL平面受壓：用材TL平面受壓，簡化為3個階段，彈性階段由公式(1)和公式(2)確定；用材RT平面受壓時的應力狀態達到公式(9)、公式(10)和公式(11)確定的屈服面后，用材進入理想塑性流動階段，過程如下：

假設第n步迭代得到的應力張量為n+1步中給出一個新的應變增量dβ_ij，根據公式(12)計算得到試探應力若在TL平面中，k_n'＝k₃，則說明用材已經進入屈服階段，該試探應力不是真實用材應力，需要將超出屈服面的部分減去：將總應變增量中的塑性應變引起的虛假應力增量從試探應力中減去，從而使應力值返回屈服面已實現理想塑性流動狀態，總應變增量由彈性應變增量和塑性應變增量之和進行表達，見公式(13)，真實應力計算見公式(14)，相關聯塑性流動法則見公式(15)，TL平面內的塑性應變表示為：

其中為用材沿T軸方向的塑性正應變增量，為用材位于TL平面內的塑性剪應變增量，為用材沿L軸方向的塑性正應變增量，△λ_TL求解是上述過程的關鍵，假設用材在第n步屈服，那么用材在第n+1步同樣屈服，應力狀態均在屈服面上，即：

其中為第n步試探應力代入公式(9)的結果，而為第n+1步試探應力代入公式(9)的結果；

故此：

其中為與的差值；

將上式按照TL平面屈服方程進行一階泰勒展開，得到：

其中△α₁₁為沿著用材L軸方向的正應力增量，△α₃₃為沿著用材T軸方向的正應力增量，△α₁₃為位于用材TL平面內的剪應力增量；

其中：

其中L、G、H、I為運算過程中的過渡量，公式(40)中的應力增量表述為：

其中為沿著用材T軸方向的試探正應力增量，為沿著用材L軸方向的試探正應力增量，為位于用材TL平面內的試探剪應力增量；

將公式(41)和公式(42)代入公式(40)，△λ_TL被求解

在總應變達到應力硬化閾值時，應力硬化閾值應變為彈性總應變的2倍，應力出現硬化，木材領域統稱為二次強化現象。