熱熔式地膜打孔裝置成孔器-地膜瞬態(tài)熱結構耦合模擬

柯年峰，孫偉，劉小龍，張華，李輝

（730070 甘肅省 蘭州市 甘肅農業(yè)大學 機電工程學院）

0 引言

隨著農業(yè)種植技術的發(fā)展，地膜覆土這一項技術得到了廣泛應用，而膜上覆土技術會造成覆蓋的泥土發(fā)生板結，引起缺苗或者燒苗的情況。使用打孔地膜會有效改善膜上土壤板結，針對地膜打孔問題設計了一種熱熔式成孔器。熱熔式成孔器是基于工具與工件之間傳熱現(xiàn)象，熱量從高溫熱熔式成孔器傳遞到地膜上，達到對地膜形狀和尺寸的加工。有熱量傳遞就會引起溫度的變化。地膜結構由于溫度變化導致其產生熱變形并形成熱應力。同時熱變形也會改變其邊界，又導致地膜發(fā)生破裂[1]。處于熱環(huán)境下的加工，熱應力主要受兩方面影響：一是其內部溫度變化，二是其本身的材料參數(shù)的變化。通過ANSYS Workbench 進行瞬態(tài)熱結構耦合不僅能夠直觀地得到熱熔式成孔器對地膜的力學特性，還能觀測地膜熱應力的變化情況。

學者對薄膜熱應力進行了大量研究。何志童等[2]研究了熱應力對FeCo 薄膜磁晶各向異性的影響，分析得出較大熱應力的引入對于提高磁晶各向異性能有顯著的效果；何濤[3]利用有限元方法研究了不銹鋼316L/Ti/ta-C 膜系從沉積溫度冷卻到室溫時的熱應力分布規(guī)律，分析得出不同沉積溫度、鈦過渡層厚度、不銹鋼基底厚度和ta-C 薄膜厚度情況下，ta-C 薄膜和鈦過渡層所受熱應力均為壓應力；吳靚臻等[4]采用有限元模型對SiO2薄膜熱應力進行模擬薄膜熱應力的大小和分布，分析了不同鍍膜溫度、不同膜厚和不同基底厚度生長環(huán)境下熱應力的大小。

從現(xiàn)有研究看，對不同薄膜的熱應力研究都比較廣泛，但是對地膜的熱應力研究還未開始，因此研究地膜熱應力在不同情況下的變化，包括不同溫度熱熔式成孔器、不同高度熱熔式成孔器，對于解析地膜的物理性質有重要意義。

本文采集熱熔式成孔器以及地膜特性的相關數(shù)據(jù)，建立熱熔式成孔器-地膜的三維模型，設置熱熔式成孔器與地膜熱結構耦合模型，基于ANSYS Workbench 有限元分析軟件對將要打孔的0.01 mm 黑膜進行數(shù)值模擬，獲得不同情況下的熱熔式成孔器對地膜作用時地膜所產生的熱應力。為研究地膜熱應力提供了一定的參考依據(jù)。

1 地膜熱應力問題分析

1.1 瞬態(tài)溫度場有限元求解方程

求解瞬態(tài)溫度場問題的核心是利用相應的數(shù)值方法求解線性常微分方程組，需要通過求解線性常微分方程組，接著通過推導三維溫度場得到一個以時間t 為獨立變量的線性常微分方程組[5]：

式中：C——熱容矩陣；φ——節(jié)點溫度對時間的導數(shù)列陣；K——熱傳導矩陣；——節(jié)點溫度列陣；P——溫度載荷列陣矩陣。C，K，P 的元素由單元相應的矩陣元素集成[6]，為

通過計算，將時間域和空間域的偏微分方程問題在空間域內轉變成N 個節(jié)點溫度常微分方程的初值問題。

1.2 熱-結構耦合分析理論

熱應力是由溫度場和應力場相互作用的結果，溫度場和應力場互相作用為熱固耦合問題。針對這些耦合問題可以從物理方程出發(fā)，熱-結構分析計算的有限元方程為[7]

式中：M——質量矩陣；u——位移；T——溫度載荷；C——結構阻尼矩陣；Ci——比熱矩陣；K——結構剛度矩陣；Kt——熱傳導矩陣；F——總等效結點力列陣；Q——總等效結點熱流率向量。

1.3 耦合熱應力有限元分析

仿真采用順序耦合的方法進行分析。首先建立模型的結構，通過ANSYS Workbench 自帶的模塊建立模型，但是由于其操作復雜，可以嘗試通過SolidWorks 建立模型導入到ANSYS Workbench；其次定義結構的材料屬性，使模型屬性符合實際；接著對建立的模型劃分合適的網格，網格不合適會影響計算精度[8]。在進行瞬態(tài)熱分析時，需要添加一系列條件，對2 個結構的初始溫度、換熱系數(shù)以及換熱方式進行設置。在得到瞬態(tài)溫度場分布以后，將得到的溫度載荷輸入到結構力學分析，利用給定的溫度載荷，添加固定支撐條件以及z 軸方向的位移條件，從而完成耦合應力的計算。最后對得到的結果進行分析和檢驗。

2 熱熔式打孔機構對地膜熱應力數(shù)值計算

2.1 模型建立及網格的劃分

在SolidWorks 中建立模型導入到ANSYS Workbench。熱熔式成孔器長、寬、高分別為10，30，30 mm，其下半部分的圓柱半徑為20 mm。地膜的長、寬、厚分別為20，10，0.01 mm。如圖1 所示。

圖1 熱熔式成孔器-地膜的三維模型Fig.1 Three-dimensional model of hot-melt hole former-mulching film

模型建立以后對模型進行網格劃分。網格劃分的好壞直接影響到計算的精度和速度。熱熔式成孔器與地膜屬于不均勻接觸過程，因此對地膜處的網格需要設置較密。由于熱熔式成孔器是一個不規(guī)則形狀的模型[9]，需要對其網格進行特殊處理，同時為了提高計算精度且節(jié)約計算資源，需要分別劃分網格。地膜處的網格采用Face Meshing 的方法，將選中面生成映射網格。Face Meshing 的優(yōu)點是能夠快速生成高質量的六面體網格。由于熱熔式成孔器下方呈現(xiàn)圓環(huán)狀，系統(tǒng)默認采用掃掠網格劃分法，但是無法生成映射網格，需要在設置此部分網格時點擊mesh—insert—Face Meshing，在Face Meshing中設置Internal Number of Divisions 為6，選中圓環(huán)的環(huán)面生成網格。最后得到有限元模型中包含節(jié)點個數(shù)為30 048，單元個數(shù)為6 999。網格生成圖如圖2 所示。

圖2 網格生成圖Fig.2 Grid generation diagram

2.2 材料熱物性參數(shù)

溫度場是在x，y，z 直角坐標系中連續(xù)介質在各個地方同一時刻的溫度分布。溫度場模擬所需要的材料熱物性參數(shù)有導熱系數(shù)λ（W/m·℃）、比熱容c（J/kg·℃）、密度ρ（kg/m3）、熱膨脹系數(shù)α（1/℃）等。熱熔成孔器主要由結構鋼構成，參考現(xiàn)有的材料參數(shù)[10]，其導熱系數(shù)為44.19 kW/(m·c)，比熱容約為600 J/kg·℃，密度7 850 kg/m3，熱膨脹系數(shù)為1.1×10-5/℃。地膜參考聚乙烯材料參數(shù)[11]，導熱系數(shù)為0.42 kW/(m·c)，比熱容約為0.002 J/kg·℃，密度為950 kg/m3，熱膨脹系數(shù)為0.23×10-5/℃，熔點為105～115℃。將其材料參數(shù)在瞬態(tài)熱分析的Engineering Data 內選定后，進行下一步驟。

2.3 熱-結構耦合模塊設置

在進行耦合之前，需要對兩物體接觸進行分析。由于熱熔式成孔器與地膜表面相切，并且會沿z 軸方向相對移動，因此兩物體屬于非線性接觸。將模型導入Mechanical 后，系統(tǒng)自動將二者設置為Frictionless。Frictionless 接觸類型能夠迭代多次且允許法向分離。在Frictionless 設置中還需要選擇接觸面與目標面，接觸面選擇熱熔式成孔器的圓弧面，目標面選擇地膜的上表面。其中接觸面顯示紅色，目標面顯示藍色。為了防止在接觸過程中熱熔式成孔器穿透地膜表面，需要選擇非對稱接觸選項。生成接觸對以后，還需要對接觸行為、接觸算法等方式進行設置。

完成接觸設置后，進行熱結構耦合相關設置。仿真分為2 部分：第1 部分為分析熱熔成孔器與地膜換熱時的溫度場；第2 部分為將所求的溫度場作為載荷輸入到結構力學進行熱應力分析。

第1 部分分析結構隨時間的變化而變化的瞬態(tài)溫度場。在瞬態(tài)熱結構內設置熱熔式成孔器模型為恒定熱源，表面溫度一直維持在設定值，同時設置地膜溫度為22 ℃。計算時，對于復雜的分析需要設置一些控制器選項。在Transient Thermal 中的Anslysis Settings 進行熱分析設置。對瞬態(tài)熱的載荷步控制采用時間作為跟蹤參數(shù)，表示真實時間歷程的變量在變化，計算得到的結果是時間的相關函數(shù)。設置結束的時間為10 s，初始時步為0.1 s，最大時步為1 s，最小時步為0.001 s，同時需要開啟時間積分設置。在求解器的設置上，系統(tǒng)有直接求解器Direct 和迭代求解器Iterative 供選擇。此處選擇程序控制，即在計算過程中由系統(tǒng)自動選擇兩者中的一種求解方式，有利于得到較為精確的結果。

第2 部分是將所得的溫度載荷輸入到瞬態(tài)結構力學模塊計算得到應力結果。為了約束實體在x，y，z 方向的不規(guī)則移動，需要添加固定支撐條件。在對模型進行固定支撐約束時，選取地膜的4 個面作為固定約束。由于熱熔成孔器是自上而下移動與地膜接觸傳熱，因此需要添加位移約束，即約束實體在x，y，z 方向的移動，在z 軸方向輸入一個具體的值，表示在該方向強制移動此距離，其他2 個方向的位移約束均設置為0，便能模擬出熱熔成孔器的運動過程。在Transient 中的Anslysis Settings 進行結構分析設置，載荷步控制同樣采用時間作為跟蹤參數(shù)，設置結束的時間為10 s，初始時步為0.1 s，最大時步為1 s，最小時步為0.001 s，與Transient Thermal 的設置保持一致。在求解器的設置上，與Transient Thermal 設置有所不同，需要將求解器控制中的Weak springs 打開。因為在有限元軟件中，經過網格離散化，幾何體可能處于微小的不平衡狀態(tài)，打開后便能使得模型在一個平衡狀態(tài)下計算。由于地膜的材料屬性受熱易變形，還需要將Large Deflection 打開。

Workbench 結果后處理功能很強大，能得到用戶所需的多種結果。如變形、應力應變、接觸工具、疲勞工具等。由于本仿真探求的是關于地膜的應力問題，需要在Solution 工具條Strain 中的Stress 查看應力命令。等效應力/應變又稱馮米斯(von-Mises)等效應力/應變，是塑形材料計算時最常用的結果，計算結果需滿足第四強度理論σVon≤σs[12]，所以只需查看等效應力的命令便能得到地膜應力的結果。

由于熱熔成孔器與地膜接觸時具有速度快、接觸時間短、加熱面積小等特點，為了更直觀地觀測地膜溫度變化，可從地膜表面溫度場進行分析。

在熱熔成孔器的表面溫度分別設置為200，210，220℃，地膜表面溫度22℃，位移條件為0.008 m。取2.5，5，7.5，10 s 四個時間研究溫度場的分布。通過ANSYS Workbench 后處理分別得到4 個時間段的溫度場，如圖3 所示。

圖3 不同時刻地膜表面溫度分布Fig.3 Surface temperature distribution of mulch film at different times

剛開始加熱時，熱熔成孔器溫度還不能及時傳遞到地膜上，2.5 s 時中心處的溫度分布不均勻，直到10 s 時中心處的溫度才逐漸呈均勻分布狀態(tài)，隨著時間變化，熱熔成孔器與地膜接觸時間長，與地膜交互的溫度增大。由圖3 可知，地膜表面溫度隨熱熔成孔器表面溫度同步增大，當熱熔成孔器表面溫度達到220 ℃時，地膜表面溫度達到其熔點。

由于熱熔式成孔器與地膜接觸傳熱，導致地膜的溫度發(fā)生變化，這種變化使得地膜內部產生熱應力。將上一步計算得到的溫度載荷導入到瞬態(tài)結構模塊，并結合所給定的固定支撐條件以及位移條件進行地膜的熱應力分析。

在熱熔成孔器的表面溫度設置成為220 ℃，地膜表面溫度22℃，位移條件分別為0.004，0.006，0.008 時的熱應力分布情況。取2.5，5.0，7.5，10 s 四個時間研究熱應力的分布。通過ANSYS Workbench 后處理分別得到4 個時間段的熱應力圖如圖4 所示。

從圖4 可以看出，考慮位移約束以后地膜的應力變化隨高度增大而變大。采用約束條件以后，可以分析出不同情況下應力的變化情況。在位移初期，由于交互位移變化較小，應力的輪廓面比較不規(guī)則。隨著時間的增大，位移向下移動與地膜接觸面積增大，應力輪廓面呈出一個橢圓形的形狀，且邊緣較為平整。在移動過程中，由于熱熔式成孔器下移高度的變化，其應力也是隨之變化的。在下移8 mm時，呈現(xiàn)的應力集中在熱熔式成孔器與地膜接觸處，這可能造成地膜的破裂。輪廓面上的應力圖尺寸符合打孔尺寸。

在熱應力的作用下，地膜會發(fā)生變形或者破裂。孫偉[13]研究發(fā)現(xiàn)地膜在應力8.49e6 Pa 時會發(fā)生破裂。在對輪廓面的處理上，設置Capped Isosurface 的值為8.49e6 Pa，并且點擊Top 選項便能從輪廓面上看到大于選定值的區(qū)域，這樣能夠直觀地看到應力作用于地膜造成破裂的形狀及尺寸大小，如圖5 所示。

3 結論

通過ANSYS Workbench 將瞬態(tài)熱分析與瞬態(tài)結構力學分析相結合的方法，開展針對地膜在200～220℃內的瞬態(tài)溫度場仿真，以及針對熱熔打孔器在-4，-6，-8 mm 的瞬態(tài)力學仿真。結論如下：

（1）熱熔式打孔器在工作過程中與地膜接觸，由于溫度梯度而引起內能交換，會產生受熱不均的現(xiàn)象，如果沒達到合適溫度會造成成孔不完整或未穿透現(xiàn)象。

（2）通過瞬態(tài)熱分析計算，得到了熱熔式成孔器對地膜傳熱的溫度場變化，在其他條件相同的情況下熱熔是成孔器達到220℃時，地膜與熱熔式成孔器接觸處表面最高溫度達到其臨界溫度。溫度的變化主要集中在最初的10 s 內，后期的變化較平緩。且在升溫過程中，通過對地膜熱結構耦合計算，得到了地膜的應力場分布，整個過程中地膜應力場分布均勻，只在兩工件接觸區(qū)域附近應力值較大，其余大部分區(qū)域無應力或者應力變化較小。當下降高度為8 mm 時所產生的應力作用得到的孔洞符合尺寸要求。

（3）對地膜受熱過程中發(fā)生的變形進行分析，引起地膜變形的包括軸向變形和徑向變形，由于設置固定支撐及位移約束，軸向變形隨高度的變化而改變。通過觀察仿真輪廓圖地膜中心處的變形區(qū)域，對實際生產起到指導作用。