基于有限元法的大豆種子碰撞模擬分析

頓國強， 李海山， 于春玲， 楊永振， 高志勇， 毛寧， 張時雨

(1.東北林業(yè)大學機電工程學院，黑龍江 哈爾濱 150040；2.東北林業(yè)大學工程技術(shù)學院，黑龍江 哈爾濱 150040)

大豆作為重要的豆類作物，在其收獲、脫粒、分離以及清選等過程中始終存在大豆與機械部件的接觸作用[1-3]。在脫粒的過程中，大豆植株進入脫?？臻g后，在脫粒部件的作用下從豆莢內(nèi)脫出而實現(xiàn)脫粒，脫出的大豆會與脫粒滾筒上的釘齒、滾筒及凹板和罩殼無規(guī)律的碰撞彈射后進入分離、清選裝置，大豆與脫粒部件的碰撞是大豆機械損傷的主要原因，外部損傷的大豆可由種子分選裝置清除，而具有內(nèi)部損傷的種子易發(fā)霉、生蟲，出芽率低，影響種子的出苗[4-5]。

以大豆種子有限元碰撞過程的應力及位移為主要形式以表征種子內(nèi)部損傷的創(chuàng)新研究方法，因較少應用而未受重視。目前，關(guān)于研究種子的內(nèi)部損傷及顆粒接觸碰撞研究的主要方法有損傷機理分析、試驗研究、力學分析、有限元模擬和離散單元法等[6-8]。MAMLOUK等[9]用有限元法研究了早餐燕麥片的壓縮過程；XU等[10]對水稻籽粒受脫粒齒斜沖擊的過程進行理論分析與有限元模擬；PETRU等[11]建立麻瘋樹種子在壓縮載荷作用下的力學行為的有限元方法模型；CELIK等[12]用有限元法模擬了蘋果在跌落情況下的變形行為；XU等[13]采用離散元法對大豆種子顆粒進行建模和驗證；YOUSEFI等[14]進行了梨果實跌落試驗與有限元模擬；AHMADI等[15]研究了在不同沖擊載荷作用下，蘋果動態(tài)行為的粘彈性有限元分析。上述研究皆為單一內(nèi)容的研究，結(jié)合上述研究成果，采用綜合的方法模擬分析大豆種子與脫粒部件的碰撞過程，揭示大豆種子碰撞的損傷機理。

本研究在合農(nóng)60品種參數(shù)統(tǒng)計結(jié)果的基礎(chǔ)上，利用有限元法模擬大豆種子的碰撞過程，分析其內(nèi)部應力分布與位移形變隨時間的變化，并通過全因素試驗分析碰撞速度及接觸半徑對最大應力和最大位移試驗指標的影響。此外，以大豆種子模型在長度、寬度和厚度3個方向的模擬撞擊過程，分析了碰撞過程中大豆種子體積變化對最大應力和最大位移的影響,由于大豆種子在脫粒機中工作的實際情況非常復雜，分析大豆種子在脫粒碰撞過程中的接觸面視為剛性平面，本文的研究是大豆脫粒過程的一基礎(chǔ)性研究，以期為大豆種子脫粒機的仿真設(shè)計和產(chǎn)品研發(fā)提供理論參考。

1 大豆種子的碰撞模擬分析

1.1 大豆種子物理參數(shù)測定

選用合農(nóng)60大豆種子為試驗材料，隨機選取完整的大豆種子100粒，含水率為13.1%，利用電子分析天平(量程100 g，精度0.001 g)測量單個大豆種子質(zhì)量；利用精度0.01 mm的數(shù)顯游標卡尺測量大豆的三軸尺寸(長度L、寬度W、厚度T)，測量結(jié)果見表1。

由測量結(jié)果可知，大豆種子的三軸尺寸均近似服從正態(tài)分布，其三軸尺寸的正態(tài)分布如圖1所示。

表1 大豆種子測量結(jié)果

圖1 大豆種子的三軸尺寸的分布Fig.1 Distribution of triaxial dimensions of the soybean seed

由圖1可知，大豆種子在長度方向分布相對集中，均值為6.601 mm；在寬度方向和厚度方向分布相對更加分散，均值分別為6.645mm和6.021 mm。

1.2 大豆種子有限元模型的建立

由于大豆種子的幾何形狀近似于橢球形，故大豆種子三維建模時，可將大豆種子模型簡化為橢球體，以實際測量的大豆種子的三軸尺寸的平均長為6.635 mm、平均寬為6.510 mm、平均厚為5.985 mm作為橢球的三軸尺寸。大豆種子的幾何尺寸如圖2所示，有限元三維網(wǎng)格模型如圖3所示。

圖2 大豆種子模型幾何尺寸Fig.2 Geometry of the soybean seed model

圖3 大豆種子的FEM網(wǎng)格劃分Fig.3 FEM mesh of the soybean seed

利用SolidWorks有限元分析模塊創(chuàng)建大豆種子的碰撞模型算例。查閱相關(guān)文獻[15-19]確定大豆種子的材料屬性參數(shù)見表2。

表2 大豆種子材料屬性Table 2 Material parameters of soybean

1.3 大豆種子碰撞模擬及結(jié)果分析

在本研究中，利用跌落試驗研究模擬撞擊過程，使用Solidworks Simulation FEM code 對大豆種子撞擊情況進行了仿真。為了模擬物體掉落并撞擊剛性表面的過程有限元程序必須迭代求解以下一般方程，應為隨著物體掉落和撞擊剛性平面，力和剛度不斷變化。

在接觸過程中，如果接觸壓力小于屈服壓力，該接觸被視為彈性且由赫茲公式?jīng)Q定[16-18]。內(nèi)應力F為

(1)

(2)

(3)

式中:E*為有效楊氏模量，MPa；R*為曲率半徑，mm；φ為接觸面積，mm2；E1為大豆種子模型的彈性模量，MPa；E2為接觸剛性平面的彈性模量，MPa；σ1為大豆種子模型的泊松比；σ2為接觸剛性平面的泊松比；py為大豆種子模型的屈服力，MPa。

在本研究中，大豆種子模型與剛性平面發(fā)生初始碰撞之前不考慮轉(zhuǎn)動，假定大豆種子是線性各向同性材料，根據(jù)表2材料參數(shù)，建立大豆種子的碰撞模型，輔助創(chuàng)建撞擊平面及標定速度方向，如圖4所示。

圖4 撞擊過程方案示意圖Fig.4 Schematic diagram of collision process

根據(jù)建立好的大豆種子模型，對其撞擊過程進行仿真，設(shè)定顆粒與接觸面的碰撞速度為2 m·s-1，摩擦系數(shù)為0.4。沖擊后的求解時間設(shè)定為120 μs，圖解數(shù)50。運行撞擊仿真算例，提取顆粒碰撞過程中應力最大時應力云圖及位移云圖的xoy面剪裁如圖5所示。圖5以真實比例顯示模型的應力和位移效果。

從大豆種子模型的應力分布xoy剪裁可以看出，大豆的應力最大主要發(fā)生在模型的中下部，向四周呈輻射狀擴散，最大應力值為19.5 MPa;從模型的位移云圖xoy面剪裁可見，因軟件在計算模型的變形時，以接觸面為參考基準，這時，模型的變形情況會與實際變形相反，最大位移及最小位移發(fā)生在模型的上部及下部，呈輻射狀擴散，位移分布關(guān)于y軸對稱，最大位移值為0.063 73 mm。為更加直觀形象的表達變形大豆種子的應力分布，如圖6a所示，在變形大豆種子模型的截面隨機探測3點，探測其隨時間變化的曲線如圖6b所示。

圖5 模型的應力和位移Fig.5 Stress and displacement of the model

由圖6b可以看出，在模型內(nèi)部y軸向應力分布情況均為先增大后減小，在90 μs到120 μs減小到0，模型與剛性平面接觸點處應力較大，探測點的應力呈現(xiàn)先增大后減小直至消失。

圖6 模型探測點和曲線Fig.6 Detection point and curve of the model

為分析大豆種子的應力分布及變形隨時間的變化，提取大豆種子單個碰撞周期不同時間內(nèi)部應力分布如圖7所示,單個碰撞周期不同時間內(nèi)部應力分布如圖8所示。

由圖7可以看出，在2.98 μs時，大豆碰撞開始，在模型的底部接觸區(qū)開始產(chǎn)生應力，應力最大位置位于接觸點，應力值為2.5 MPa，隨著時間的增加，應力區(qū)域上移并呈輻射狀逐漸增大，應力值也同步增加，在t=44.97 μs時，應力值達到最大，此時應力最大值為19.5 MPa，此后因大豆回彈現(xiàn)象的發(fā)生，大豆反向運動逐漸脫離接觸面，應力區(qū)域逐漸縮小，應力值也相應減小，直至89.99 μs，大豆脫離碰撞接觸面，但應力區(qū)域沒有隨大豆的脫離而馬上消失，而在大豆的內(nèi)部存在一定的殘余應力，以應力波的形式傳輸于大豆內(nèi)部，且應力波逐漸衰減。

由圖8可以看出，大豆單個碰撞周期的形變位移的形變過程與大豆應力分布隨時間的變化具有相同的趨勢，位移最大時間為44.97 μs，位移值0.063 73 mm，之后變形逐漸減小，直至89.99 μs脫離接觸面，但因大豆種子內(nèi)部殘余應力的作用，大豆并未恢復原狀，而是沿y軸方向以“伸長-壓縮-伸長”的方式非周期性的衰減變形。

為提取大豆種子模型的應力和位移的時間歷史圖解，在模型表面分別隨機選取4點，以探測其應力時變歷程和位移歷程如圖9所示。由圖9b可見，在單個碰撞周期中，在模型表面分別隨機選取四點的模型表面的應力隨時間變化的趨勢基本相同，呈現(xiàn)先增大后減小，在80 μs時，應力達到穩(wěn)定狀態(tài)趨于0；模型表面的位移曲線與表面應力歷程曲線具有相類似的趨勢，但位移歷程曲線較為光滑，而應力歷程具有一定的波動，且在當大豆種子脫離接觸面后，應力曲線并未消失，而是以小振幅的波動曲線存在，這與碰撞過程的應力理論分析相一致；而模型的位移曲線在大豆種子脫離接觸面后以一定的速度線性增大，這主要是由于在碰撞模型中，位移的計算基準是接觸面，且在模型中不存在外部加速度，因此，位移曲線產(chǎn)生了如圖9d所示的趨勢。

圖7 單個碰撞周期不同時間的應力分布Fig.7 Stress distribution at different times in a single collision cycle

圖8 單個碰撞周期不同時間的位移分布Fig.8 Displacement distribution at different times in a single collision cycle

圖9 大豆模型中應力和位移的時變歷程Fig.9 The time-varying history of stress and displacement in soybean model

2 大豆碰撞仿真試驗

2.1 撞擊形態(tài)的分析

大豆脫粒機械技術(shù)是極為復雜的過程，脫粒轉(zhuǎn)筒內(nèi)大豆種子運動的隨機性，以及轉(zhuǎn)筒轉(zhuǎn)速和裝載量對大豆的隨機運動都有影響，所研究大豆種子在不同的撞擊形態(tài)下與脫粒裝置撞擊瞬間內(nèi)部應力以及位移變化有重要的實際意義。大豆種子模型被簡化為橢球體，當與模擬的剛性平面撞擊時,由于模型的重心與剛性平面的接觸點的接觸半徑的不同，其撞擊形態(tài)也不同，分別為對稱撞擊、非對稱撞擊和極限撞擊。其撞擊狀態(tài)如圖10所示。

圖10 大豆種子模型不同的撞擊形態(tài)Fig.10 The different impact forms of the soybean seed model

為分析大豆種子接觸半徑變化及碰撞速度改變對大豆有限元模型應力分布及形變位移的影響，基于上述大豆種子有限元碰撞分析結(jié)果，以大豆碰撞接觸半徑及碰撞速度為試驗因素，應力最大值及形變位移最大值為試驗指標，進行全因素試驗，分析碰撞參數(shù)變化對大豆應力及變形的影響，為排除其他應因素的干擾，模型材料屬性、模型網(wǎng)格劃分及接觸面靜摩擦設(shè)為定值。全因素試驗因素水平如表3所示。

利用軟件Design-expert 8.0.6對試驗數(shù)據(jù)進行全因素試驗并進行方差分析，試驗結(jié)果的方差分析見表4。

表3 試驗因素水平Table 3 Test factors and levels

表4 方差分析Table 4 Variance analysis

由表4試驗結(jié)果方差分析可知，在信度水平為0.01時，對于最大應力這一試驗指標，接觸半徑和碰撞速度兩個因素對其影響極顯著；對于最大位移這一試驗指標，接觸半徑和碰撞速度兩個因素對其影響極顯著。

在大豆種子撞擊過程中，接觸半徑和碰撞速度對其影響都比較顯著，進行了不同接觸半徑在不同碰撞速度的撞擊仿真過程，獲得了大豆種子撞擊過程最大應力和最大位移與接觸半徑和碰撞速度之間關(guān)系曲線如圖11所示。

圖11 接觸半徑對最大應力和最大位移的影響Fig.11 Effect of contact radius on maximum stress and maximum displacement

根據(jù)圖11a，在上述幾種撞擊狀態(tài)下，隨著接觸半徑增大，撞擊過程中最大應力呈先減小再增大后減小再增大的趨勢。隨著碰撞速度的增大，在各撞擊形態(tài)下，最大應力呈現(xiàn)增大的趨勢。根據(jù)圖11b，隨著接觸半徑增大，撞擊過程中最大位移呈先減小再增大再趨于平緩的趨勢。隨著碰撞速度的增大，在各撞擊形態(tài)下，最大位移呈現(xiàn)增大的趨勢。

為分析大豆種子在對稱狀態(tài)、非對稱狀態(tài)和極限狀態(tài)3種撞擊狀態(tài)，即不同接觸半徑時有限元模型的應變分布及形變位移的影響，選定碰撞速度為2 m·s-1，提取3種撞擊狀態(tài)時相對應的模型應力分布及形變位移云圖，并提取應力最大值及位移最大值隨接觸半徑變化而變化的云圖，如圖12、圖13所示。

圖12 不同撞擊形態(tài)的大豆種子最大應力Fig.12 Maximum stress of the soybean seed model in different impact forms

圖13 不同撞擊形態(tài)的大豆種子最大位移Fig.13 Maximum displacement of the soybean model in different impact forms

由圖12、圖13可見,在對稱撞擊狀態(tài)下撞擊，大豆種子模型應力集中區(qū)域基本呈對稱式分布，大豆種子模型位移呈對稱式分布；在非對稱撞擊狀態(tài)下撞擊，即由圖12d、13d大豆種子模型應力和位移集中區(qū)域基本呈對稱式分布，由圖12e、12f、13e、13f可見,大豆種子模型應力和位移集中區(qū)域基本呈非對稱式分布，最大位移量發(fā)生在模型上部；在極限狀態(tài)下撞擊，由圖12g、13g可知，大豆種子模型應力和位移集中區(qū)域基本呈非對稱式分布，最大位移量發(fā)生在模型上部。

2.2 不同大豆種子體積對撞擊過程的影響

為了研究大豆種子體積變化對撞擊過程的影響，根據(jù)表1測得百粒大豆種子的厚度尺寸分布在5.120～6.940 mm之間,平均長厚比為1.109 mm，平均寬厚比為1.089 mm，根據(jù)三軸尺寸基于一定的關(guān)系，以大豆種子的厚度尺寸為基礎(chǔ)數(shù)據(jù)，即在大豆種子厚度尺寸最小值和最大值區(qū)間段內(nèi)定量選取數(shù)據(jù)，根據(jù)長厚比、寬厚比即可算出大豆種子的長度、寬度尺寸,進行如表5的大豆種子撞擊仿真方案。

表5 大豆種子撞擊仿真方案和結(jié)果Table 5 Simulation scenario and results of the soybean seed impingement

根據(jù)上述方案，分別以大豆種子模型的長度方向、寬度方向和厚度方向進行撞擊仿真過程，得到的數(shù)據(jù)繪制曲線如圖14所示。

圖14 大豆種子體積對最大應力和最大位移的影響Fig.14 Effect of the volume of the soybean seed on maximum stress and maximum displacement

根據(jù)圖14a，大豆種子模型在長度方向、寬度方向和厚度方向的撞擊過程中，大豆種子體積較小時，寬度方向和厚度方向撞擊產(chǎn)生的最大應力大于長度方向撞擊產(chǎn)生的最大應力，隨著大豆種子體積的增大，長度方向撞擊產(chǎn)生的最大應力大于寬度方向和厚度方向撞擊產(chǎn)生的最大應力；無論是長度方向、寬度方向還是厚度方向撞擊，隨著大豆種子體積的增大，最大應力都呈現(xiàn)先增大再減小的趨勢。根據(jù)圖14b，寬度方向和長度方向撞擊產(chǎn)生的最大位移大于厚度方向撞擊產(chǎn)生的最大位移，隨著大豆種子體積的增大，撞擊過程中最大位移都呈現(xiàn)線性增大的趨勢。

3 結(jié)論

1)基于合農(nóng)60品種參數(shù)統(tǒng)計的結(jié)果，利用有限元法模擬大豆種子的碰撞過程，分析其內(nèi)部應力分布與位移形變隨時間的變化，并以碰撞速度及接觸半徑為試驗因素，最大應力及最大位移為試驗指標，進行全因素試驗。結(jié)果表明，在信度水平為0.01時，對于最大應力和最大位移，接觸半徑和碰撞速度2個因素對其影響極顯著。

2)通過大豆種子模型在長度、寬度和厚度3個方向的模擬撞擊過程，分析碰撞過程中大豆種子體積變化對最大應力和最大位移的影響。研究發(fā)現(xiàn)，隨著大豆種子體積的增大，撞擊過程中最大應力都呈現(xiàn)先增大再減小的趨勢；隨著大豆種子體積的增大，3個方向撞擊過程中最大位移都呈現(xiàn)線性增大的趨勢。