張 利，陳俊超，陳國(guó)達(dá)，傅昱斐

(浙江工業(yè)大學(xué) 機(jī)械工程學(xué)院，浙江 杭州 310032)

0 引 言

鈦合金人工關(guān)節(jié)是一種植入人體的醫(yī)療器械，預(yù)期使用壽命為10～20年，甚至更長(zhǎng)[1]。人們對(duì)它的表面質(zhì)量要求很高。然而鈦合金材料具有導(dǎo)熱系數(shù)低、塑性弱、硬度高、彈性模量低等特點(diǎn)，導(dǎo)致切削加工性能差[2-3]，人工關(guān)節(jié)表面結(jié)構(gòu)又普遍比較復(fù)雜，所以傳統(tǒng)的機(jī)械加工方法很難達(dá)到理想的效果。

針對(duì)這種表面質(zhì)量需求較高，但構(gòu)造卻比較復(fù)雜的鈦合金工件，學(xué)者們提出了許多新型的加工方法，比如振動(dòng)拋光[4]、電化學(xué)拋光[5]、磁性研磨拋光[6]等。這些加工方法各有優(yōu)缺點(diǎn)：比如振動(dòng)拋光方法所需的裝置簡(jiǎn)單，適用性強(qiáng)，但拋光質(zhì)量相對(duì)較低，且容易產(chǎn)生二次損傷；電化學(xué)拋光為非接觸式拋光，它能有效避免工件表面產(chǎn)生殘余應(yīng)力及二次損傷，但該方法難以保證工件表面材料去除的均勻性，整體拋光效果較差；磁性研磨拋光根據(jù)介質(zhì)不同可以細(xì)分為磁性磨料研磨和磁性流體研磨，實(shí)際拋光中磁性流體研磨效果相對(duì)更好，但它主要適宜加工一些凸形零件的表面，加工范圍存在局限性。

針對(duì)復(fù)雜結(jié)構(gòu)化表面的超精密加工，計(jì)時(shí)鳴等[7]提出了一種軟性磨粒流加工方法。該方法首先仿照復(fù)雜表面的幾何結(jié)構(gòu)設(shè)計(jì)相對(duì)應(yīng)的仿形約束流道，再用弱粘性的流體與磨粒混合制備軟性磨粒流，最后把軟性磨粒流循環(huán)通入仿形約束流道，對(duì)待加工表面進(jìn)行加工[8]。該方法較其他加工方法普適性強(qiáng)，針對(duì)不同的表面結(jié)構(gòu)易于調(diào)整加工參數(shù)，加工效果良好，適合加工各類復(fù)雜的鈦合金表面。

軟性磨粒流是一種固-液兩相流，在兩相流領(lǐng)域，對(duì)于固相和液相的耦合方式有兩種比較典型的描述，分別是歐拉-歐拉描述和歐拉-拉格朗日描述[9]。兩種描述的區(qū)別主要在于對(duì)固相的描述不同。軟性磨粒流的顆粒濃度一般較低，可以忽略固相對(duì)液相的影響，所以一般采用歐拉-拉格朗日法對(duì)其進(jìn)行數(shù)學(xué)建模。

在對(duì)仿真過(guò)程中的湍流運(yùn)動(dòng)進(jìn)行描述時(shí)，雷諾時(shí)均法是一種比較常見(jiàn)的建模方法[10]。該方法的優(yōu)點(diǎn)是計(jì)算量小，得到的結(jié)果已經(jīng)能夠滿足大多數(shù)實(shí)際工程的需要；缺點(diǎn)是對(duì)湍流物理量進(jìn)行時(shí)間平均后抹平了所有的脈動(dòng)值，無(wú)法體現(xiàn)湍流特有的無(wú)序脈動(dòng)這一特點(diǎn)。比如李瓊等[11]在模擬液壓油箱中固體顆粒的運(yùn)動(dòng)和李紹武等[12]在模擬粉塵運(yùn)動(dòng)時(shí)，都采用了基于雷諾時(shí)均法的湍流模型，他們模擬得出的顆粒運(yùn)動(dòng)軌跡雖然無(wú)序，但軌跡本身并不會(huì)隨時(shí)間而改變；從觀測(cè)顆粒運(yùn)動(dòng)的角度出發(fā)，這不是希望得到的結(jié)果。

因此，為了盡可能真實(shí)地模擬顆粒的運(yùn)動(dòng)，本文提出一種湍流分散方法，在雷諾時(shí)均法的基礎(chǔ)上還原湍流的速度脈動(dòng)特征，然后在此基礎(chǔ)上進(jìn)行流場(chǎng)分析，最后進(jìn)行加工實(shí)驗(yàn)進(jìn)行驗(yàn)證。

1 流場(chǎng)數(shù)學(xué)模型

1.1 基于歐拉法的液相湍流方程

湍流方程由兩部分組成：第一部分為運(yùn)用雷諾時(shí)間平均思想在經(jīng)典Navier-Stokes方程(以下簡(jiǎn)稱N-S方程)基礎(chǔ)上[13-14]推導(dǎo)出的湍流基本方程，即雷諾時(shí)均方程；第二部分是應(yīng)用最廣泛的標(biāo)準(zhǔn)k-ε湍流雙方程，該方程是為了使雷諾時(shí)均方程封閉而選用的湍流補(bǔ)充方程。

湍流區(qū)別于層流的關(guān)鍵在于湍流內(nèi)部存在尺度不一的渦旋運(yùn)動(dòng)，這些渦旋運(yùn)動(dòng)一般被認(rèn)為仍然符合N-S方程的描述[15]。在不考慮溫度變化的情況下，理論上對(duì)湍流的每一個(gè)瞬時(shí)狀態(tài)用N-S方程中的連續(xù)性方程和動(dòng)量方程就可以準(zhǔn)確地描述湍流運(yùn)動(dòng)，但這種直接數(shù)值模擬方法涉及的計(jì)算量很大[16]，不能滿足實(shí)際工程的需要。

因此，本文依據(jù)雷諾時(shí)間平均思想在N-S方程基礎(chǔ)上構(gòu)建雷諾時(shí)均方程(RANS方程)。

連續(xù)性方程如下：

(1)

動(dòng)量方程(以x方向?yàn)槔?如下：

(2)

(3)

(4)

(5)

用同樣的方法可推導(dǎo)出湍流在y、z方向上的動(dòng)量方程，然后結(jié)合湍流連續(xù)性方程，整理得出RANS方程如下：

(6)

(7)

式中：ρ—流體密度；▽—梯度算子；V—速度矢量；u，v，w—速度分量；p—壓強(qiáng)；μ—流體粘性系數(shù)；fx，fy，fz—作用在單位質(zhì)量上的體積力分量。另外，各個(gè)變量上方的橫線表示該變量是運(yùn)用雷諾時(shí)均法后得到的時(shí)均值。

(8)

針對(duì)新增的k與μt兩項(xiàng)，本文選用LAUNDER和SPALDING[18]提出的標(biāo)準(zhǔn)k-ε模型(Standardk-εModel)建立其內(nèi)在聯(lián)系：

(9)

(10)

其中，方程中各模型常數(shù)一般取值為：

C1ε=1.44，C2ε=1.92，Cμ=0.09，σk=1，σε=1.3。

1.2 基于拉格朗日法的顆粒運(yùn)動(dòng)模型

根據(jù)牛頓第二運(yùn)動(dòng)定律可以得出單顆顆粒的運(yùn)動(dòng)方程：

(11)

式中：mp—顆粒質(zhì)量；Vp—顆粒速度；t—時(shí)間；∑Fi—作用在顆粒上的各種力之和。

顆粒在流體中運(yùn)動(dòng)時(shí)會(huì)受到很多力的作用，其中曳力是在任何情況下都不可以忽略的力。本文只考慮顆粒受到曳力的情形。

曳力的大小與流體密度ρ、迎流面積Sp、速度差(v-vp)有關(guān)。參考空氣阻力的經(jīng)驗(yàn)公式有:

(12)

可以得到曳力計(jì)算公式：

(13)

假設(shè)顆粒為球形并考慮曳力的作用方向，上式變?yōu)椋?/p>

(14)

式中：ρ—密度；v—速度；μ—流體粘度；r—顆粒半徑。下標(biāo)p—該物理量為顆粒物理量；未標(biāo)下標(biāo)p—除了某些無(wú)歧義的物理量外皆為流體物理量；Cd—曳力系數(shù)，其值由顆粒雷諾數(shù)Rep決定:

(15)

(16)

從曳力的計(jì)算公式中可以看出，流體速度是影響曳力大小的主要因素。前文的雷諾時(shí)均法抹平了湍流中所有尺度的脈動(dòng)情況，只體現(xiàn)了所有物理量的平均結(jié)果。從研究顆粒運(yùn)動(dòng)的角度來(lái)說(shuō)，這種情況體現(xiàn)不出湍流對(duì)于顆粒的無(wú)序驅(qū)動(dòng)作用。因此，本文在模擬顆粒運(yùn)動(dòng)時(shí)，擬在一定程度上還原湍流運(yùn)動(dòng)特有的無(wú)序狀態(tài)，即湍流分散方法。

在標(biāo)準(zhǔn)k-ε湍流模型中湍流動(dòng)能k的表達(dá)式為：

(17)

(18)

（3）依據(jù)前期設(shè)計(jì)要求，在確定爆破孔及預(yù)裂孔位置無(wú)誤后，可進(jìn)行鉆孔設(shè)備合理選擇。為保證鉆孔精度，本次工程選擇110B型潛孔鉆機(jī)作為預(yù)鉆機(jī)，其具有鉆孔基面小、精度高、操作便捷等優(yōu)點(diǎn)；而在緩沖孔、主爆破孔鉆孔設(shè)備選擇時(shí)，則主要采用CM353型鉆機(jī)。結(jié)合角度尺的應(yīng)用，可以保證緩沖孔坡度的有效控制。

(19)

記速度脈動(dòng)矢量為U′，那么有：

(20)

式中：I—隨機(jī)單位向量，用于表示速度脈動(dòng)的方向。

對(duì)于上式中由3個(gè)隨機(jī)參數(shù)λ1，λ2，λ3組成的系數(shù)，采用一個(gè)新的隨機(jī)參數(shù)λ代替，即：

(21)

上式的關(guān)鍵是隨機(jī)參數(shù)λ的取值，理論上來(lái)說(shuō)λ的取值有無(wú)限的可能，但實(shí)際應(yīng)用中卻不能在無(wú)限范圍內(nèi)對(duì)λ取值，具體的取值需要結(jié)合所模擬的流動(dòng)情況而定。本文中λ的取值是一個(gè)均值為0、標(biāo)準(zhǔn)差為1的隨機(jī)分布。

因此，曳力的計(jì)算公式中流體的速度用重新構(gòu)建的湍流仿真速度U代替，即：

(22)

2 數(shù)值模擬

2.1 模型及參數(shù)設(shè)定

實(shí)際加工中的仿型約束流道如圖1所示。

圖1 仿型約束流道示意圖1-蓋板；2-人工膝關(guān)節(jié)；3-流道；4-流道出口；5-底板；6-螺釘7-固定件；8-密封膠填充槽；9流道入口

圖1中，蓋板1、底板5以及固定件7作為仿型約束配件，三者組成仿型約束流道，人工膝關(guān)節(jié)三維模型與仿真流道模型，如圖2所示。

圖2 物理模型和流道網(wǎng)格劃分1-仿型面A；2-仿型面B；3-過(guò)度面

人工膝關(guān)節(jié)模型是參考實(shí)際膝關(guān)節(jié)結(jié)構(gòu)而設(shè)計(jì)，其中1和2分別是包含6個(gè)曲率參數(shù)的仿型面，3是仿型面之間的過(guò)度面；仿真流道模型是根據(jù)圖1的仿型約束流道設(shè)計(jì)，用COMSOL軟件對(duì)其進(jìn)行自動(dòng)網(wǎng)格劃分，網(wǎng)格形狀為正四面體。

為了便于觀測(cè)磨粒流流道內(nèi)的流場(chǎng)分布情況，本研究對(duì)圖2(b)所示的流道模型進(jìn)行截面分割，得到的結(jié)果如圖3所示。

圖3 截面示意圖

仿真中的一些關(guān)鍵參數(shù)如表1所示。

表1 關(guān)鍵仿真參數(shù)

2.2 仿真結(jié)果

0.1 s時(shí)刻流道截面速度分布如圖4所示。

圖4 截面速度分布

從整體來(lái)看，流道中間部分的速度基本高于初始速度5 m/s，靠近出口部分的速度高于靠近入口部分的速度。從單個(gè)仿型面來(lái)看，仿型面A的截面1、截面2、截面3的速度分布依次呈現(xiàn)高-低-高的趨勢(shì)，結(jié)合物理模型分析，從相對(duì)位置較低的截面1和截面3到相對(duì)位置較高的截面2，途中損失了部分速度，導(dǎo)致速度分布呈現(xiàn)高-低-高這一結(jié)果。仿型面B的情況類似。

0.1 s時(shí)刻流道截面壓力分布如圖5所示。

圖5 截面壓力分布

圖5中顯示的結(jié)果與截面速度分布的結(jié)果相反：磨粒流從入口到出口這一過(guò)程所受到的壓力減小，從截面1到截面6這一過(guò)程中，物理位置相對(duì)較高的地方壓力相對(duì)較大。綜合上述分析可以得出：磨粒流速度分布和壓力分布符合實(shí)際加工經(jīng)驗(yàn)，兩者之間存在一定程度的負(fù)相關(guān)聯(lián)系，該情況與伯努利方程的內(nèi)涵相符。

為了直觀反映前文提出的湍流分散方法的應(yīng)用效果，設(shè)置如下簡(jiǎn)單情形：在10 mm的直管一端以5 m/s的速度通入水，在該端部中間以0.001 s的時(shí)間間隔釋放單顆粒子，粒子在水的曳力作用下運(yùn)動(dòng)，分別模擬未應(yīng)用湍流分散方法時(shí)，顆粒的運(yùn)動(dòng)過(guò)程與運(yùn)用湍流分散方法時(shí)顆粒的運(yùn)動(dòng)過(guò)程，得到的結(jié)果如圖6所示。

圖6 不同模擬方法下顆粒的運(yùn)動(dòng)軌跡

從圖6(a)中可以看到：在曳力作用下相同位置釋放的顆粒都沿著相同的軌跡運(yùn)動(dòng)；但圖6(b)中相同位置釋放的顆粒在曳力作用下的運(yùn)動(dòng)軌跡卻不相同，軌跡之間也相對(duì)無(wú)序。兩圖中顆粒運(yùn)動(dòng)軌跡的差別證實(shí)了湍流分散方法的有效性。

本研究結(jié)合湍流分散方法，對(duì)磨粒在人工膝關(guān)節(jié)加工流道模型內(nèi)的運(yùn)動(dòng)進(jìn)行模擬仿真，結(jié)果如圖7所示。

圖7 不同時(shí)刻的顆粒運(yùn)動(dòng)軌跡

從圖7中可以看到：t=0時(shí)刻磨粒在入口處就緒，流道內(nèi)并沒(méi)有磨粒；t=0.01 s時(shí)刻有兩波磨粒進(jìn)入了流道，最開(kāi)始釋放的磨粒已經(jīng)到了流道中部；t=0.02 s時(shí)刻已有4波磨粒進(jìn)入了流道，最初釋放的磨粒已到達(dá)流道出口；t=0.02 s時(shí)刻到t=0.1 s時(shí)刻之內(nèi)磨粒不斷進(jìn)入流道，同時(shí)先前進(jìn)入流道的磨粒也不斷從出口離開(kāi)，磨粒運(yùn)動(dòng)的軌跡逐漸增多，直到20波磨粒都進(jìn)入了流道。圖中每條軌跡不同部分的顏色表示粒子運(yùn)動(dòng)到該位置時(shí)所具有的速度大小。

0.1 s時(shí)刻流道內(nèi)湍流動(dòng)能分布如圖8所示。

圖8 0.1 s時(shí)刻的湍流動(dòng)能分布

圖8中粒子軌跡的顏色表示0.1 s時(shí)刻不同位置的流體所具有的湍流動(dòng)能，從圖中可以看到：絕大部分區(qū)域的流體所具有的湍流動(dòng)能均低于0.4 m2/s2，只有出口部分的湍流動(dòng)能高于0.4 m2/s2。從湍流動(dòng)能的分布可以看出幾乎整個(gè)流道內(nèi)的流體都已處于湍流狀態(tài)，但絕大部分區(qū)域的湍流脈動(dòng)程度相對(duì)較小。

綜合上述分析可以得到：磨粒在進(jìn)入流道后受到流體曳力的作用而做加速運(yùn)動(dòng)，運(yùn)動(dòng)過(guò)程中與壁面碰撞發(fā)生反彈，反彈后的速度方向也因?yàn)榱黧w曳力的影響而逐漸變?yōu)榕c流體流動(dòng)方向一致。但在湍流分散方法的影響下，以相同條件釋放的磨粒的仿真運(yùn)動(dòng)軌跡不再重合。

考慮到實(shí)際加工中造成磨粒無(wú)序運(yùn)動(dòng)的因素有很多，比如泵的間歇性推動(dòng)作用、流道入口之前的管路狀況等，因此磨粒的仿真運(yùn)動(dòng)結(jié)果雖然與磨粒的實(shí)際運(yùn)動(dòng)結(jié)果仍存在一定差距，但對(duì)于前文所建立的數(shù)學(xué)模型來(lái)說(shuō)，該仿真結(jié)果是符合預(yù)期的。

3 加工實(shí)驗(yàn)

3.1 平臺(tái)搭建及參數(shù)設(shè)定

搭建的實(shí)驗(yàn)平臺(tái)如圖9所示。

圖9 實(shí)驗(yàn)平臺(tái)示意圖1-渣漿泵；2-磨?；厥粘兀?-加工流道；4-輸送管

渣漿泵是整個(gè)磨粒流加工過(guò)程的動(dòng)力源，開(kāi)式磨粒回收池是整個(gè)磨粒流加工循環(huán)的起點(diǎn)和終點(diǎn)。磨粒流的加工流道，內(nèi)部安裝著表面待加工的人工膝關(guān)節(jié)實(shí)物。兩條輸送管，其中與泵相連的是加工流道的磨粒流入口管，與磨?；厥粘叵噙B的是加工流道的磨粒流出口管。

實(shí)驗(yàn)分為粗加工和精加工兩個(gè)過(guò)程，實(shí)驗(yàn)參數(shù)如表2所示。

表2 加工實(shí)驗(yàn)主要參數(shù)

實(shí)驗(yàn)加工時(shí)，磨粒流循環(huán)流動(dòng)加工這一過(guò)程會(huì)產(chǎn)生大量熱能，包括渣漿泵、輸送管、磨粒流、加工流道、磨?；厥粘卦趦?nèi)的所有實(shí)驗(yàn)裝置和實(shí)驗(yàn)材料的溫度都會(huì)緩慢升高?？紤]到離心泵的使用壽命，粗加工和精加工過(guò)程都以1 h為加工間隔，到達(dá)時(shí)間后暫停加工，更換基液，等待0.5 h后繼續(xù)加工。同時(shí)，因?yàn)槟チＡ髦形醇臃稚枰〞r(shí)攪拌磨?；厥粘?，防止磨粒沉積。

3.2 實(shí)驗(yàn)結(jié)果及分析

加工前后表面效果對(duì)比如圖10所示。

圖10 加工前后表面效果對(duì)比

從圖10中可以明顯看到：加工前工件表面很粗糙，劃痕、凹坑等缺陷很多；而加工后工件表面比較平整，沒(méi)有明顯的表面缺陷。

在圖10(e～f)中，采用“倒映文字”的方法對(duì)工件表面質(zhì)量情況進(jìn)行宏觀對(duì)比分析，從圖中可以清楚地看到，加工前零件的下表面無(wú)法映照出紙上的“Mechanical”字樣，而加工后零件的下表面已經(jīng)能清晰地映照出各個(gè)字母。此外，憑借肉眼也能清楚地判斷出人工膝關(guān)節(jié)的表面質(zhì)量在經(jīng)過(guò)實(shí)驗(yàn)加工后得到了很好的改善。

4 結(jié)束語(yǔ)

本文在經(jīng)典N-S方程基礎(chǔ)上對(duì)RANS方程進(jìn)行了推導(dǎo)，然后根據(jù)湍流動(dòng)能重新構(gòu)建湍流速度脈動(dòng)，結(jié)合標(biāo)準(zhǔn)k-ε方程、顆粒受力方程共同建立磨粒流加工流場(chǎng)的數(shù)學(xué)模型，以此為基礎(chǔ)借助Comsol軟件對(duì)磨粒流加工流場(chǎng)進(jìn)行了數(shù)值模擬，參考仿真結(jié)果搭建了實(shí)驗(yàn)平臺(tái)，進(jìn)行了加工實(shí)驗(yàn)。得到的結(jié)論如下：

(1)磨粒流加工流場(chǎng)的速度分布與壓力分布各自符合實(shí)際加工經(jīng)驗(yàn)，兩者之間存在一定程度的負(fù)相關(guān)聯(lián)系，與伯努利方程的內(nèi)涵相符；

(2)經(jīng)仿真結(jié)果驗(yàn)證，湍流分散方法能在一定程度上還原湍流速度脈動(dòng)這一特點(diǎn)，從觀測(cè)顆粒運(yùn)動(dòng)的角度出發(fā)，該方法可以有效提高磨粒流仿真分析的可靠性；

(3)分別從微觀角度和宏觀角度對(duì)實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比分析，得出磨粒流加工方法對(duì)人工膝關(guān)節(jié)表面具有良好的加工效果。