王 睿

(安徽農(nóng)業(yè)大學(xué) 工學(xué)院，合肥 230036)

SPH- FEM耦合方法在切削加工數(shù)值模擬中的研究

王 睿

(安徽農(nóng)業(yè)大學(xué) 工學(xué)院，合肥 230036)

針對(duì)SPH與FEM的各自特點(diǎn)，采用SPH與FEM耦合的方法進(jìn)行切削加工數(shù)值模擬。在變形比較小處采用FEM，以提高計(jì)算精度和效率；在變形較大處采用SPH，以避免網(wǎng)格畸變?cè)斐捎?jì)算困難。在ABAQUS軟件環(huán)境中，建立了彎頸刨刀對(duì)45號(hào)鋼切槽的SPH- FEM耦合模型，模擬了刨削過程材料的應(yīng)力、切削應(yīng)力，驗(yàn)證了耦合算法的有效性。

切削；數(shù)值模擬；SPH- FEM耦合算法

切削加工在機(jī)械制造工藝中占據(jù)極其重要的地位。國內(nèi)外學(xué)者對(duì)切削加工做了大量的研究，主要研究方法有解析理論法、實(shí)驗(yàn)法和有限元法(finite element method, FEM)等。由于切削加工過程是一個(gè)復(fù)雜的彈塑性變形和斷裂的過程，解析理論法進(jìn)行理論分析和推導(dǎo)比較困難。實(shí)驗(yàn)法則需要精密的實(shí)驗(yàn)儀器和檢測(cè)設(shè)備，并且數(shù)據(jù)結(jié)果受到實(shí)驗(yàn)人員的操作水平的影響。有限元法是建立在彈塑性力學(xué)模型的基礎(chǔ)上，通過變分方法，使得誤差函數(shù)達(dá)到最小值并產(chǎn)生穩(wěn)定解的數(shù)值模擬方法。它可以模擬切削中材料的變形過程，精確的顯示出切削過程的瞬時(shí)應(yīng)力、應(yīng)變、溫度等信息。通過有限元方法可以揭示切削機(jī)理，優(yōu)化切削參數(shù)，對(duì)于降低生產(chǎn)成本，保證產(chǎn)品質(zhì)量，提高生產(chǎn)效率具有重要意義。

但是依賴于網(wǎng)格的傳統(tǒng)有限元法對(duì)于切削模擬并不是十分理想，原因是切削引起的大變形會(huì)導(dǎo)致網(wǎng)格嚴(yán)重畸變，使計(jì)算無法進(jìn)行下去。目前解決的辦法是采取網(wǎng)格重劃分和刪除畸變網(wǎng)格單元的技術(shù)。當(dāng)網(wǎng)格畸變達(dá)到設(shè)定準(zhǔn)則時(shí)，進(jìn)行網(wǎng)格重劃分，將原有的網(wǎng)格單元信息映射到新的網(wǎng)格單元，從而形成新的網(wǎng)格，對(duì)于畸變依然嚴(yán)重的網(wǎng)格單元?jiǎng)t及時(shí)進(jìn)行刪除。這種方法雖然能夠保證切削模擬的順利進(jìn)行，但是頻繁地劃分網(wǎng)格不僅占用了大量的計(jì)算時(shí)間，而且舊網(wǎng)格單元向新網(wǎng)格單元的映射以及刪除網(wǎng)格單元不可避免地會(huì)引入誤差，影響模擬精度，限制了有限元切削模擬的推廣。

光滑粒子動(dòng)力學(xué)(smoothed particle hydrodynamics, SPH)法是一種拉格朗日形式的無網(wǎng)格法。SPH法采用粒子劃分，不依賴于網(wǎng)格，具有很好的自適應(yīng)性，可以避免網(wǎng)格畸變，適合切削引起的大變形問題。然而SPH法對(duì)每個(gè)粒子進(jìn)行計(jì)算時(shí)，需要搜索影響區(qū)域內(nèi)近鄰的粒子信息。粒子物理量計(jì)算和搜索信息都比較費(fèi)時(shí)，

因此計(jì)算效率比普通的有限元法低，對(duì)于三維模型占用計(jì)算機(jī)資源較大。

針對(duì)SPH與FEM的各自特點(diǎn)，為提高計(jì)算效率并消除網(wǎng)格畸變，采用SPH與FEM耦合的方法解決切削數(shù)值模擬問題。在變形大的區(qū)域采用SPH，避免FEM的網(wǎng)格畸變過大造成計(jì)算困難。在變形小的區(qū)域采用FEM，以提高計(jì)算效率。SPH與FEM耦合算法分為固定耦合算法和自適應(yīng)耦合算法。固定耦合算法在計(jì)算之前就已確定SPH區(qū)域和FEM區(qū)域。自適應(yīng)耦合算法則在計(jì)算之前都是FEM網(wǎng)格，在計(jì)算過程中自動(dòng)地將大變形的有限元網(wǎng)格單元轉(zhuǎn)換為光滑粒子，并按SPH法計(jì)算物理量。本文在大型通用有限元軟件ABAQUS 6.14的環(huán)境下，建立彎頸刨刀對(duì)45號(hào)鋼切槽的SPH和FEM自適應(yīng)耦合模型，對(duì)刨削過程的材料的應(yīng)力、切削應(yīng)力和切削機(jī)理進(jìn)行了分析。

1 SPH法的基本原理

SPH法的核心是插值理論，它將連續(xù)介質(zhì)離散成一系列具有質(zhì)量的拉格朗日粒子，用積分近似表示場(chǎng)函數(shù)，再對(duì)域內(nèi)的相鄰粒子的值積分。SPH方程的構(gòu)造，先進(jìn)行核估計(jì)，再進(jìn)行粒子近似。在SPH法中，對(duì)任意一個(gè)函數(shù)f(x)，其核估計(jì)的表達(dá)式為：

=∫Ωf(x')W(x-x',h)dx'

(1)

式中W(x-x′,h)為核函數(shù)，Ω為點(diǎn)x的支持域；x-x′為粒子間距離；h為SPH粒子的光滑長度。從(1)式看出函數(shù)可以通過其核函數(shù)光滑求得，核函數(shù)相當(dāng)于權(quán)函數(shù)。函數(shù)f(x)導(dǎo)數(shù)的核估計(jì)表達(dá)式為：

<▽f(x)>=∫Ωf(x')▽W(xué)(x-x',h)dx'

(2)

圖1 粒子近似示意圖

從(2)式可以看出函數(shù)f(x)的導(dǎo)數(shù)可以轉(zhuǎn)換成其核函數(shù)的導(dǎo)數(shù)，即函數(shù)f(x)的導(dǎo)數(shù)可以由其核函數(shù)導(dǎo)數(shù)光滑求得。

SPH法的粒子近似是將函數(shù)的解域離散成任意分布的粒子，在粒子上進(jìn)行物理量計(jì)算，如圖1所示。

函數(shù)及其導(dǎo)數(shù)的核估計(jì)的積分表達(dá)式轉(zhuǎn)化為粒子的求和的離散形式的表達(dá)式為：

(3)

(4)

式中:mj為粒子j的質(zhì)量，ρj為粒子j的密度；N為在粒子i影響域內(nèi)所有粒子總數(shù)。從(3)和(4)式，可以看出函數(shù)及其導(dǎo)數(shù)可以由離散點(diǎn)上的函數(shù)值經(jīng)過加權(quán)平均而獲得，其權(quán)因子就是核函數(shù)或其導(dǎo)數(shù)。

核函數(shù)決定了SPH法的精度和效率，廣泛使用的核函數(shù)是B樣條函數(shù)，表達(dá)式為：

(5)

SPH法的計(jì)算過程是首先將連續(xù)介質(zhì)離散為SPH粒子，進(jìn)行近鄰粒子搜索，確定每個(gè)粒子影響域中的粒子，根據(jù)動(dòng)量守恒計(jì)算出每個(gè)粒子的速度及坐標(biāo)。接著計(jì)算出粒子的應(yīng)變率，根據(jù)質(zhì)量守恒計(jì)算出粒子的密度，根據(jù)能量守恒計(jì)算出粒子的比內(nèi)能。最后根據(jù)狀態(tài)方程得到靜水壓力，根據(jù)彈塑性本構(gòu)關(guān)系得到應(yīng)力偏張量，依據(jù)應(yīng)力張量是由靜水壓力和應(yīng)力偏張量組成，從而得到某時(shí)刻的應(yīng)力。

2 SPH-FEM的耦合處理

圖2 SPH-FEM耦合算法流程圖

本文采用的SPH-FEM耦合算法是自適應(yīng)耦合算法，即計(jì)算之前都是FEM網(wǎng)格，計(jì)算中根據(jù)網(wǎng)格單元的變形情況決定是否轉(zhuǎn)變?yōu)镾PH粒子。自適應(yīng)耦合算法包括SPH粒子的轉(zhuǎn)換、接觸滑移的計(jì)算、交界面SPH粒子的計(jì)算等內(nèi)容。SPH粒子的轉(zhuǎn)換可以把應(yīng)變作為轉(zhuǎn)換標(biāo)準(zhǔn)，當(dāng)FEM網(wǎng)格單元應(yīng)變超過設(shè)定的應(yīng)變閥值時(shí)，則將其轉(zhuǎn)為SPH粒子。轉(zhuǎn)換過程是將超過應(yīng)變閥值的網(wǎng)格單元?jiǎng)h除，調(diào)整有限元的界面并計(jì)算相關(guān)的有限元節(jié)點(diǎn)質(zhì)量，增加新的SPH粒子并計(jì)算其物理量。接觸滑移的計(jì)算包括不同材料FEM單元間、FEM單元與SPH粒子間、不同材料SPH粒子間的接觸滑移計(jì)算。接觸滑移計(jì)算要求在接觸面上的一方不能侵入另一方的內(nèi)部。如果發(fā)生侵入現(xiàn)象，則按動(dòng)量守恒、動(dòng)量矩守恒及速度協(xié)調(diào)原則調(diào)節(jié)單元或粒子的速度和位置，以保障不發(fā)生侵入現(xiàn)象。交界面SPH粒子的計(jì)算是把交界面一側(cè)處于粒子的影響域的FEM單元作為虛擬粒子參與SPH計(jì)算。但虛擬粒子不作為SPH粒子，而依然保留FEM特征，參與FEM計(jì)算。這樣可以避免FEM一側(cè)由于沒有粒子而導(dǎo)致粒子分布不對(duì)稱引起的計(jì)算誤差。SPH-FEM耦合算法的流程圖如圖2所示。

3 計(jì)算模型

3.1 彎頸刨刀對(duì)45號(hào)鋼切槽的SPH-FEM模型

圖3 彎頸刨刀切槽的FEM模型

在初始時(shí)刻，彎頸刨刀對(duì)45號(hào)鋼切槽全部用FEM建模，如圖3所示。45號(hào)鋼為尺寸為300mm×400mm×200mm的矩形塊，在ABAQUS環(huán)境中繪制。彎頸刨刀在CREO軟件繪制并導(dǎo)入到ABAQUS環(huán)境中。45號(hào)鋼采用C3D8R網(wǎng)格單元?jiǎng)澐?，?duì)切削部位采用比基體更密的網(wǎng)格。FEM網(wǎng)格轉(zhuǎn)換SPH粒子的應(yīng)變閥值設(shè)為0.7，即切削部位的網(wǎng)格應(yīng)變超過0.7時(shí)，自動(dòng)轉(zhuǎn)換為SPH粒子。彎頸刨刀采用C3D4網(wǎng)格單元?jiǎng)澐帧?/p>

采用ABAQUS的EXPLICIT模塊的顯示算法進(jìn)行切削模擬，分析時(shí)長設(shè)置為1s。定義彎頸刨刀與45號(hào)鋼的接觸類型為“通用接觸”。彎頸刨刀沿Z軸直線運(yùn)動(dòng)，速度為400mm/s，約束其X、Y方向自由度及其全部轉(zhuǎn)動(dòng)自由度。45號(hào)鋼的底部及左右兩側(cè)的全部自由度均被約束。

3.2 材料模型

由于在切削時(shí)，45號(hào)鋼的塑性變形、溫度和應(yīng)變率均發(fā)生較大的變化，可采取Johnson- Cook本構(gòu)模型作為材料模型。Johnson-Cook本構(gòu)模型的表達(dá)式為：

(6)

表1 45號(hào)鋼的本構(gòu)參數(shù)

45號(hào)鋼的損傷失效模型采取Johnson- Cook損傷模型，表達(dá)式為：

(7)

式中D為損傷參數(shù)，D∈[0、1]，初始時(shí)D=0。當(dāng)D=1時(shí)材料破壞，在ABAQUS中標(biāo)記相應(yīng)單元的狀態(tài)變量STATUS為1,在隨后的切削變形中該單元即被刪除。 為一個(gè)時(shí)間步的塑性應(yīng)變?cè)隽浚?為當(dāng)前時(shí)間步的破壞應(yīng)變，表達(dá)式為：

(8)

彎頸刨刀的材料為W18Cr4V，材料參數(shù)見表2。

表2 W18Cr4V參數(shù)

4 結(jié)果與分析

當(dāng)彎頸刨刀與45號(hào)鋼發(fā)生接觸進(jìn)行切削時(shí)，由圖4可以看出，刀尖部位的應(yīng)力最大，并且離刀尖的距離越遠(yuǎn)應(yīng)力越小。由于此時(shí)刨削處于彈性變形階段，F(xiàn)EM單元沒有太大的畸變，因此沒有SPH粒子。當(dāng)切削進(jìn)行到中間階段時(shí)，由圖5可以看出，此時(shí)應(yīng)力最大處是緊鄰刀尖處的材料，刀尖處應(yīng)力次之，而已切削的材料則沿刨刀面向上運(yùn)動(dòng)形成切屑。這是因?yàn)樵诘都獾膽?yīng)力作用下，刀尖前部材料達(dá)到強(qiáng)度極限產(chǎn)生裂紋，所以在裂紋尖端處材料應(yīng)力最大。由于切削會(huì)產(chǎn)生殘余內(nèi)應(yīng)力，因此切屑和已切好的溝槽應(yīng)力雖然減少但不為零。由于切削產(chǎn)生較大的變形，促使一些變形大的FEM單元轉(zhuǎn)變?yōu)镾PH粒子分布在切屑和溝槽中，從而避免了網(wǎng)格畸變使計(jì)算出錯(cuò)，以及刪除網(wǎng)格單元所造成的誤差，更加符合切削的真實(shí)過程。切削結(jié)束時(shí)，由圖6可以看出，溝槽應(yīng)力均衡并較小，這是由于溝槽此時(shí)均為殘余內(nèi)應(yīng)力。從圖6還可看出在溝槽的開始處有圓角，在結(jié)束處有毛刺。這是由于開始時(shí)材料產(chǎn)生塑性變形，造成圓角，結(jié)束時(shí)裂紋拉長造成毛刺。從圖7可以看出主切削應(yīng)力在剛開始切削時(shí)，由零迅速增大，然后在一定范圍內(nèi)進(jìn)行上下震蕩，在切削結(jié)束時(shí)，又減少至零。這是因?yàn)樵谇邢鞒跏紩r(shí)，刨刀與工件接觸發(fā)生彈性變形，應(yīng)力迅速增大。隨著切削進(jìn)行，材料達(dá)到屈服極限，發(fā)生塑性變形，然后再達(dá)到強(qiáng)度極限，材料產(chǎn)生裂紋并擴(kuò)展，導(dǎo)致刀尖不與材料直接接觸，切削應(yīng)力有所下降。刀具繼續(xù)向前，刀尖又與材料接觸，又發(fā)生彈性變形，切削應(yīng)力又有所上升。如此循環(huán)，切削應(yīng)力在一定范圍上下震蕩，直至切削完成，切削應(yīng)力降為零。

圖4 0.07秒時(shí)切削應(yīng)力圖

圖5 0.3秒時(shí)切削應(yīng)力圖

圖6 0.8秒時(shí)切削應(yīng)力圖

圖7 主切削應(yīng)力變化曲線

通過彎頸刨刀對(duì)45號(hào)鋼切槽的數(shù)值模擬，證明了采用SPH- FEM耦合方法進(jìn)行切削加工的仿真是可行的。SPH- FEM耦合的方法，吸收了FEM法計(jì)算效率高和SPH法模擬大變形能力強(qiáng)的優(yōu)點(diǎn)，為大變形的材料切削加工提供一種高效、準(zhǔn)確的途徑。

[1] 卞梁,王肖鈞,章杰. SPH/FEM 耦合算法在陶瓷復(fù)合靶抗侵徹?cái)?shù)值模擬中的應(yīng)用[J].高壓物理學(xué)報(bào),2010,24(3): 161- 167.

[2] 于峻宏. 土壤切削大變形的三維數(shù)值仿真[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào),2007,38(4):118- 121.

[3] 宿崇,唐亮,侯俊銘，等. 基于FEM與SPH耦合方法的金屬切削仿真研究[J].系統(tǒng)仿真學(xué)報(bào),2009，21(16):5002- 5005.

[4] 王建明,余豐,劉飛宏，等. SPH和FEM耦合法模擬磨料水射流中單磨粒加速過程[J].山東大學(xué)學(xué)報(bào)(工學(xué)版), 2011,41(5):114- 119.

[5] 郭曉光,魏延軍,張小翼，等.基于FEM- SPH耦合方法的AISI4340鋼超精密切削過程仿真[J].大連理工大學(xué)學(xué)報(bào), 2013,53(4):526- 531．

[6] 林曉東,盧義玉,湯積仁，等.基于SPH- FEM 耦合算法的磨料水射流破巖數(shù)值模擬[J].振動(dòng)與沖擊,2014,33(18):170- 176.

[責(zé)任編輯：張永軍]

Simulation Research of Cutting Based on FEM and SPH

WANG Rui

(School of Engineering, Anhui Agricultural University, Hefei 230036, China)

According to the characteristics of SPH and FEM, the process of cutting is simulated by using the coupled method of SPH and FEM. FEM is employed in the lower distortion regions in order to improve the accuracy and efficiency of calculation. SPH is employed in the higher distortion regions in order to avoid the problem of calculation for mesh distortion. In ABAQUS software environment, the coupled model of SPH and FEM for 45 steel grooving by using bent neck plane knife is established. The simulation of cutting show the stress of material and the cutting stress，and accordingly prove the effectiveness of the coupled method of SPH and FEM.

cutting; simulation; SPH- FEM coupling

2016-09-28

2017-02-20

王 睿(1973— )，男，安徽合肥人，安徽農(nóng)業(yè)大學(xué)工學(xué)院講師；研究方向：模具CAD及數(shù)值仿真。

TP391.9

A

2096-2371(2017)02-0084-05