一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法
2023-10-19 15:03:57 1
一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法
【專利摘要】本發明涉及一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，步驟為：在摩擦熱-應力-磨損耦合過程的時間軸上將該過程離散為N個增量步，瞬態熱傳導分析、應力分析和磨損分析採用相同時間軸和增量步；建立幹滑動摩擦的有限元模型，進行瞬態熱傳導分析，獲得輸出溫度場；設置有限元模型應力分析的單元類型，進行應力分析，獲得輸出接觸壓力場；根據溫度場、接觸壓力場及兩接觸面接觸節點的相對滑移速率，計算兩接觸面接觸節點的磨損量並確定其空間方向；根據磨損量和方向，修正接觸節點的位移，更新有限元模型；判斷應力-磨損耦合過程的模擬是否完成。本發明計算精度損失小且計算效率較高，能廣泛應用於各類幹滑動摩擦的摩擦熱-應力-磨損弱耦合問題的模擬。
【專利說明】一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法
【技術領域】
[0001]本發明涉及一種幹滑動摩擦耦合模擬方法，特別是關於一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法。
【背景技術】
[0002]幹滑動摩擦界面存在摩擦熱-應力-磨損耦合現象，其中，摩擦熱主要產生於兩個接觸面之間的摩擦接觸表層，摩擦熱的熱流密度與切向摩擦應力和相對滑移速度有關。摩擦熱會直接導致接觸面的溫度迅速升高，迅速升高的溫度致使接觸面發生熱變形，從而最終使得接觸壓力在接觸面上重新不均勻的分布，在接觸壓力和相對滑移速度的綜合作用下，接觸面逐漸磨損。由於接觸面上接觸壓力分布不均以及接觸面不同位置處相對滑移速度存在差異，因此，接觸面不同位置處的磨損量不同。這種接觸面上磨損量的不均勻性導致接觸壓力的再次重新分布，進而改變切向摩擦應力的分布和熱流密度的大小。另一方面，接觸面溫度的改變還會顯著改變兩接觸材料的摩擦特性，比如摩擦係數和磨損係數,除此之夕卜，還會改變材料的本構特性，比如屈服極限和硬度等。因此，接觸面的摩擦熱(溫度場)、應力(接觸壓力場)和磨損(磨損量)之間是相互耦合的，這種耦合行為幾乎不可能求得解析解，而採用數值方法進行求解數值解時，如果考慮摩擦熱-應力-磨損三者之間的雙向耦合效應，那麼計算的效率極低，無法求解較長時間的摩擦熱-應力-磨損耦合問題。
【發明內容】
[0003]針對上述問題，本發明的目的是提供一種計算精度損失小且計算效率高的幹滑動摩擦熱-應力-磨損順序耦合的模擬方法。
[0004]為實現上述目的，本發明採取以下技術方案:一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，包括以下步驟:1)時域離散:假設摩擦熱-應力-磨損耦合過程的時間軸為(V tN)，在該時間軸上將摩擦熱-應力-磨損耦合過程離散為N個增量步，瞬態熱傳導分析、應力分析和磨損分析採用相同的時間軸和增量步；2)瞬態熱傳導分析:建立幹滑動摩擦的有限元模型，在步驟I)的時間軸Uci^n)上採用瞬態熱傳導分析程序進行瞬態熱傳導分析，其步驟為:①初始化有限元模型的溫度，選定瞬態熱傳導分析程序，設置有限元模型中兩構件的單元類型、生熱熱流邊界和散熱邊界對每個增量步進行瞬態熱傳導分析；
③獲得並輸出溫度場；3)應力分析:設置有限元模型應力分析的單元類型，定義邊界條件、載荷和接觸條件，將步驟2)輸出的溫度場作為熱載荷施加在有限元模型中，並對有限元模型進行應力分析，獲得並輸出接觸壓力場；4)確定磨損量:根據熱傳導分析得到的溫度場、應力分析得到的接觸壓力場以及兩接觸面接觸節點的相對滑移速率，計算兩接觸面接觸節點的磨損量，並且確定其空間方向；5)更新網格:根據兩接觸面接觸節點的磨損量和方向，修正該接觸節點的位移，更新有限元模型；6)判斷應力-磨損耦合過程的模擬是否完成，若增量步達到N，則模擬過程結束；否則，進入下一個增量步，回到步驟2)。
[0005]所述步驟I)中的瞬態熱傳導分析中的數值積分方法採用向後歐拉方法，應力分析中的數值積分方法採用Newmark- β方法，磨損分析中的數值積分方法採用顯式歐拉算法。
[0006]所述步驟2)中，瞬態熱傳導分析程序中的生熱熱流邊界由摩擦功率引起，其熱流密度表達式為.Λ = μΤ,式中，f為兩接觸面接觸節點的摩擦熱流；μ為兩接觸面接觸節點的摩擦係數為接觸面的平均接觸壓力；f為兩接觸表面在接觸點的相對滑動速率；摩擦係數μ為兩接觸面接觸節點的溫度Τ、平均接觸壓力P等因素的函數，其表達式為://二/I,...):平均接觸壓力$和兩接觸表面在接觸點的相對滑動速率^根據有限元模型的具體結構和相對運動求解；摩擦熱流f通過熱流自定義子程序DFLUX施加在構件的底面上，兩接觸面接觸節點溫度T、接觸節點到構件邊沿的距離均從計算結果文件中讀取；用戶子程序DFLUX包括以下步驟:①從瞬態熱傳導分析程序中跳轉至用戶子程序DFLUX ;②初始化常量和係數；③根據兩接觸面接觸節點坐標計算迴轉半徑？，根據兩接觸面接觸節點溫度T和平均接觸壓力]H十算摩擦係數μ ;④根據迴轉半徑？和摩擦係數μ計算兩接觸面接觸節點處的摩擦熱流^ ；⑤返回瞬態熱傳導分析程序。
[0007]所述步驟4)中的磨損量大小的計算採用廣義非線性Archard磨損公式:Ah =K.ρ.As,式中，Ah為兩接觸面接觸節點的磨損增量；K為廣義Archard磨損公式的磨損係數；Ρ為兩接觸面接觸節點的接觸壓力；△ s為兩接觸表面在接觸點的相對滑動位移增量；上式的左右兩邊同時除以增量步時間，則得到如下公式:A = ，式中，A為兩接觸面接觸節點的磨損速率；對於參數K，為溫度和接觸壓力等因素的函數，其表達式為:
=，由實驗測量得到；1)和As從有限元分析結果的文件中讀取彳從下式計算得到:
'—tsS
[0_] r=-_
[0009]式中，At為增量步時間大小。
[0010]所述步驟4)中的磨損量的方向的確定方法如下:兩接觸面接觸節點分為兩類，一類是位於接觸面邊緣的節點，簡稱邊界節點；另一類是位於接觸面內部的節點，簡稱非邊界節點；對於非邊界節點的磨損量方向是由節點法方向確定，節點法方向是直接從有限元分析結果文件中讀取。
[0011]所述步驟4)中的磨損量大小和方向的確定通過網格移動用戶子程序UMESHMOT1N實現，包括以下步驟:①初始化常量和係數；②讀取兩接觸面接觸節點的接觸壓力P和溫度T ;③計算兩接觸面接觸節點的磨損量，確定兩接觸面接觸節點的磨損方向；
④輸出兩接觸面接觸節點磨損信息，用於繪製磨損量雲圖。
[0012]所述步驟5)中，兩接觸面接觸節點位移的修正和有限元模型的更新採用任意拉格朗日歐拉法實現；需要修正位移的兩接觸面接觸節點為接觸壓力與相對切向滑移均不為零的接觸節點。
[0013]本發明由於採取以上技術方案，其具有以下優點:1、本發明由於在實際工程應用，如制動過程中，摩擦熱-應力屬於弱耦合問題，可以忽略應力和磨損對摩擦熱傳導單方向的耦合項，本發明採用的順序耦合方法模擬摩擦熱-應力-磨損耦合行為，即將摩擦熱-應力-磨損弱耦合問題分解為瞬態熱傳導分析和應力-磨損耦合分析，因此，在計算精度損失較小的情況下極大地提高計算效率。2、本發明由於基於現有商用的有限元軟體平臺，因此設置簡單、容易上手且便於推廣。3、本發明由於採用的瞬態熱傳導分析可以直接使用商用軟體的算法，應力-磨損耦合分析在商用有限元軟體的應力分析算法的基礎上二次開發完成，磨損的計算完全嵌入有限元軟體之中，免去了數據傳輸和交互，因此，使得計算效率大大提高。本發明適用範圍廣，可以廣泛應用於各類幹滑動摩擦的摩擦熱-應力-磨損弱耦合問題的模擬。
【專利附圖】
【附圖說明】
[0014]圖1是本發明的流程示意圖；
[0015]圖2是本發明的用戶子程序UMASFL流程示意圖；
[0016]圖3是本發明的用戶子程序DFLUX流程示意圖；
[0017]圖4是本發明的用戶子程序UMESHMOT1N流程示意圖；
[0018]圖5是本發明的銷-盤的有限元模型示意圖；
[0019]圖6是本發明的銷底部的局部放大示意圖；
[0020]圖7是圖2中I處局部放大示意圖；
[0021]圖8是採用順序耦合方法和完全耦合方法得到的銷-盤磨損試驗的溫度分布對比示意圖，其中圖(a)為採用完全耦合方法得到的銷-盤磨損試驗的溫度分布，圖(b)為採用順序耦合方法得到的銷-盤磨損試驗的溫度分布；
[0022]圖9是本發明採用順序耦合方法和完全耦合方法得到的銷-盤磨損試驗的接觸節點溫度分布對比示意圖，其中，表示節點N101101採用完全耦合方法得到的銷-盤磨損試驗的接觸節點溫度，O表示節點N101105採用完全耦合方法得到的銷-盤磨損試驗的接觸節點溫度，表示節點N101109採用完全耦合方法得到的銷-盤磨損試驗的接觸節點溫度，表示節點N101113採用完全耦合方法得到的銷-盤磨損試驗的接觸節點溫度，
表示節點N101501採用完全耦合方法得到的銷-盤磨損試驗的接觸節點溫度，* 表示節點N101101採用順序稱合方法得到的銷-盤磨損試驗的接觸節點溫度， 表示節點NlOl 105採用順序耦合方法得到的銷-盤磨損試驗的接觸節點溫度，表示節點NlOl 109採用順序耦合方法得到的銷-盤磨損試驗的接觸節點溫度， 表示節點N101113採用順序耦合方法得到的銷-盤磨損試驗的接觸節點溫度，表示節點N101501採用順序耦合方法得到的銷-盤磨損試驗的接觸節點溫度；
[0023]圖10是本發明採用順序耦合方法和完全耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力對比示意圖，其中， 表示節點N101101採用完全耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，-表示節點NlOl 105採用完全耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，表示節點NlOl 109採用完全耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力， ' 表示節點N101113採用完全耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，表示節點N101501採用完全耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，表示節點Ν101101採用順序耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，O*表示節點Ν101105採用順序耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，表示節點Ν101109採用順序稱合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，表示節點Ν101113採用順序耦合方法得到的銷-盤磨損試驗的接觸節點接觸壓力，' ^表示節點Ν101501採用順序稱合方法得到的銷-盤磨損試驗的接觸節點接觸壓力；
[0024]圖11是本發明採用順序耦合方法和完全耦合方法得到的銷-盤磨損試驗的接觸節點磨損量對比示意圖，其中，?~表示節點Ν101101採用完全稱合方法得到的銷-盤磨損試驗的接觸節點磨損量， 表示節點Ν101105採用完全耦合方法得到的銷-盤磨損試驗的接觸節點磨損量，表示節點Ν101109採用完全耦合方法得到的銷-盤磨損試驗的接觸節點磨損量，表示節點Ν101113採用完全耦合方法得到的銷-盤磨損試驗的接觸節點磨損量， 表示節點Ν101501採用完全耦合方法得到的銷-盤磨損試驗的接觸節點磨損量，.0表示節點Ν101101採用順序耦合方法得到的銷-盤磨損試驗的接觸節點磨損量，
表示節點Ν101105採用順序耦合方法得到的銷-盤磨損試驗的接觸節點磨損量，
表示節點Ν101109採用順序耦合方法得到的銷-盤磨損試驗的接觸節點磨損量，表示節點Ν101113採用順序稱合方法得到的銷-盤磨損試驗的接觸節點磨損量，表示節點Ν101501採用順序耦合方法得到的銷-盤磨損試驗的接觸節點磨損量。
【具體實施方式】
[0025]下面結合附圖和實施例對本發明進行詳細的描述。
[0026]如圖1所示，本發明提供一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，包括以下步驟:
[0027]I)時域離散:假設摩擦熱-應力-磨損耦合過程的時間軸為(V tN)，在該時間軸上將摩擦熱-應力-磨損耦合過程離散為N個增量步，對瞬態熱傳導分析、應力分析和磨損分析採用相同的時間軸和增量步，增量步的設定需要滿足:
【權利要求】
1.一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，包括以下步驟: 1)時域離散:假設摩擦熱-應力-磨損耦合過程的時間軸為tN)，在該時間軸上將摩擦熱-應力-磨損耦合過程離散為N個增量步，瞬態熱傳導分析、應力分析和磨損分析採用相同的時間軸和增量步； 2)瞬態熱傳導分析:建立幹滑動摩擦的有限元模型，在步驟I)的時間軸(^t,)上採用瞬態熱傳導分析程序進行瞬態熱傳導分析，其步驟為: ①初始化有限元模型的溫度，選定瞬態熱傳導分析程序，設置有限元模型中兩構件的單元類型、生熱熱流邊界和散熱邊界； ②對每個增量步進行瞬態熱傳導分析； ③獲得並輸出溫度場； 3)應力分析:設置有限元模型應力分析的單元類型，定義邊界條件、載荷和接觸條件，將步驟2)輸出的溫度場作為熱載荷施加在有限元模型中，並對有限元模型進行應力分析，獲得並輸出接觸壓力場； 4)確定磨損量:根據熱傳導分析得到的溫度場、應力分析得到的接觸壓力場以及兩接觸面接觸節點的相對滑移速率，計算兩接觸面接觸節點的磨損量，並且確定其空間方向； 5)更新網格:根據兩接觸面接觸節點的磨損量和方向，修正該接觸節點的位移，更新有限元模型； 6)判斷應力-磨損耦合過程的模擬是否完成，若增量步達到N，則模擬過程結束；否則，進入下一個增量步，回到步驟2)。
2.如權利要求1所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟I)中的瞬態熱傳導分析中的數值積分方法採用向後歐拉方法，應力分析中的數值積分方法採用Newmark-β方法，磨損分析中的數值積分方法採用顯式歐拉算法。
3.如權利要求1所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟2)中，瞬態熱傳導分析程序中的生熱熱流邊界由摩擦功率引起，其熱流密度表達式為:
Q = MPY， 式中，q為兩接觸面接觸節點的摩擦熱流；μ為兩接觸面接觸節點的摩擦係數；尹為接觸面的平均接觸壓力^為兩接觸表面在接觸點的相對滑動速率； 摩擦係數μ為兩接觸面接觸節點的溫度Τ、平均接觸壓力P等因素的函數，其表達式為:
μ = μ(Τ,ρ,...)., 平均接觸壓力P和兩接觸表面在接觸點的相對滑動速率^根據有限元模型的具體結構和相對運動求解； 摩擦熱流通過熱流自定義子程序DFLUX施加在構件的底面上，兩接觸面接觸節點溫度Τ、接觸節點到構件邊沿的距離均從計算結果文件中讀取；用戶子程序DFLUX包括以下步驟: ①從瞬態熱傳導分析程序中跳轉至用戶子程序DFLUX ； ②初始化常量和係數； ③根據兩接觸面接觸節點坐標計算迴轉半徑？，根據兩接觸面接觸節點溫度T和平均接觸壓力聲計算摩擦係數μ ； ④根據迴轉半徑？和摩擦係數μ計算兩接觸面接觸節點處的摩擦熱^ ⑤返回瞬態熱傳導分析程序。
4.如權利要求2所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟2)中，瞬態熱傳導分析程序中的生熱熱流邊界由摩擦功率引起，其熱流密度表達式為:
q = μΤπ， 式中，H為兩接觸面接觸節點的摩擦熱流；μ為兩接觸面接觸節點的摩擦係數；戶為接觸面的平均接觸壓力；^為兩接觸表面在接觸點的相對滑動速率； 摩擦係數μ為兩接觸面接觸節點的溫度Τ、平均接觸壓力#等因素的函數，其表達式為:
μ = μ{Τ,ρ,...)., 平均接觸壓力P和兩接觸表面在接觸點的相對滑動速率#根據有限元模型的具體結構和相對運動求解； 摩擦熱流^通過熱流自定義子程序DFLUX施加在構件的底面上，兩接觸面接觸節點溫度Τ、接觸節點到構件邊沿的距離均從計算結果文件中讀取；用戶子程序DFLUX包括以下步驟: ①從瞬態熱傳導分析程序中跳轉至用戶子程序DFLUX； ②初始化常量和係數； ③根據兩接觸面接觸節點坐標計算迴轉半徑？，根據兩接觸面接觸節點溫度T和平均接觸壓力多計算摩擦係數μ ； ④根據迴轉半徑？和摩擦係數μ計算兩接觸面接觸節點處的摩擦熱流歹； ⑤返回瞬態熱傳導分析程序。
5.如權利要求1或2或3或4所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟4)中的磨損量大小的計算採用廣義非線性Archard磨損公式: Ah= K.P.As， 式中，Ah為兩接觸面接觸節點的磨損增量；K為廣義Archard磨損公式的磨損係數；P為兩接觸面接觸節點的接觸壓力；△ s為兩接觸表面在接觸點的相對滑動位移增量； 上式的左右兩邊同時除以增量步時間，則得到如下公式:
h-K' p-f j 式中，《為兩接觸面接觸節點的磨損速率； 對於參數K，為溫度和接觸壓力等因素的函數，其表達式為= 六…)，由實驗測量得到；p和△ S從有限元分析結果的文件中讀取 >從下式計算得到:.Av Y=—
At.1 式中，At為增量步時間大小。
6.如權利要求1或2或3或4所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟4)中的磨損量的方向的確定方法如下:兩接觸面接觸節點分為兩類，一類是位於接觸面邊緣的節點，簡稱邊界節點；另一類是位於接觸面內部的節點，簡稱非邊界節點；對於非邊界節點的磨損量方向是由節點法方向確定，節點法方向是直接從有限元分析結果文件中讀取。
7.如權利要求5所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟4)中的磨損量的方向的確定方法如下:兩接觸面接觸節點分為兩類，一類是位於接觸面邊緣的節點，簡稱邊界節點；另一類是位於接觸面內部的節點，簡稱非邊界節點；對於非邊界節點的磨損量方向是由節點法方向確定，節點法方向是直接從有限元分析結果文件中讀取。
8.如權利要求1到7任一項所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟4)中的磨損量大小和方向的確定通過網格移動用戶子程序UMESHMOT1N實現，包括以下步驟: ①初始化常量和係數； ②讀取兩接觸面接觸節點的接觸壓力P和溫度T； ③計算兩接觸面接觸節點的磨損量，確定兩接觸面接觸節點的磨損方向； ④輸出兩接觸面接觸節點磨損信息，用於繪製磨損量雲圖。
9.如權利要求1到8任一項所述的一種幹滑動摩擦熱-應力-磨損順序耦合的模擬方法，其特徵在於:所述步驟5)中，兩接觸面接觸節點位移的修正和有限元模型的更新採用任意拉格朗日歐拉法實現；需要修正位移的兩接觸面接觸節點為接觸壓力與相對切向滑移均不為零的接觸節點。
【文檔編號】G06F17/50GK104200034SQ201410457959
【公開日】2014年12月10日 申請日期:2014年9月10日 優先權日:2014年9月10日
【發明者】桂良進, 張方宇, 範子傑 申請人:清華大學