基于COMSOL的線(xiàn)性壓縮機(jī)動(dòng)力學(xué)仿真與實(shí)驗(yàn)驗(yàn)證

孔家煊，趙欽宇，晁春陽(yáng)，周文杰，王 博，程 路，甘智華*，徐 軍

（1.浙江大學(xué)制冷與低溫研究所 浙江省制冷與低溫技術(shù)重點(diǎn)實(shí)驗(yàn)室，杭州 310027；2.杭州電子科技大學(xué)能量利用與自動(dòng)化研究所，杭州 310018；3.日照華斯特林科技有限公司，山東 日照 276800）

0 引言

近年來(lái)，由于紅外技術(shù)和超導(dǎo)領(lǐng)域的持續(xù)發(fā)展，低溫制冷機(jī)的需求日益增長(zhǎng)[1?2]。斯特林制冷機(jī)作為低溫制冷機(jī)的重要分支，具有本征效率高、穩(wěn)定性好、結(jié)構(gòu)緊湊和使用壽命長(zhǎng)等優(yōu)點(diǎn)。斯特林制冷機(jī)的驅(qū)動(dòng)源是高頻的交變壓力，壓縮活塞和氣缸間的側(cè)向力不利于延長(zhǎng)壽命和減小噪聲，因此使用直線(xiàn)電機(jī)驅(qū)動(dòng)壓縮機(jī)，不但可以省去運(yùn)動(dòng)轉(zhuǎn)換裝置；同時(shí)能通過(guò)調(diào)節(jié)電壓來(lái)調(diào)節(jié)活塞行程[3]。

建立合理的數(shù)值模型對(duì)斯特林制冷機(jī)以及線(xiàn)性壓縮機(jī)的設(shè)計(jì)非常重要，許多研究者已經(jīng)開(kāi)展了相關(guān)工作。Clark將線(xiàn)性壓縮機(jī)的機(jī)械常量、電參數(shù)、動(dòng)子質(zhì)量等效到一個(gè)電路等效模型中，以匹配一臺(tái)斯特林制冷機(jī)為例，通過(guò)非線(xiàn)性靜磁場(chǎng)的有限元方法，以最小線(xiàn)圈電阻損耗和最小設(shè)備體積為目標(biāo)，對(duì)線(xiàn)性壓縮機(jī)進(jìn)行了設(shè)計(jì)[4]。Almajri等[5]用COMSOL建立了斯特林制冷機(jī)的計(jì)算流體力學(xué)（CFD）模型，通過(guò)實(shí)驗(yàn)驗(yàn)證后使用該模型研究了相角、孔隙率、回?zé)崞鏖L(zhǎng)度對(duì)制冷量的影響。上海技術(shù)物理研究所的鄧偉峰等[6]采用Maxwell建立了自由活塞式線(xiàn)性壓縮機(jī)電磁?機(jī)械?聲耦合模型，模擬了不同負(fù)載條件下壓縮機(jī)電磁參數(shù)和動(dòng)力特性的瞬態(tài)變化規(guī)律。

目前，在回?zé)崾降蜏刂评錂C(jī)設(shè)計(jì)領(lǐng)域，一維數(shù)值計(jì)算軟件SAGE在制冷機(jī)端較為完善，而在壓縮機(jī)端僅提供了幾個(gè)簡(jiǎn)單的接口，簡(jiǎn)化了電磁場(chǎng)的計(jì)算[7]，無(wú)法對(duì)壓縮機(jī)的設(shè)計(jì)起到良好的指導(dǎo)作用。為了彌補(bǔ)SAGE在壓縮機(jī)建模上的不足，本文采用COMSOL對(duì)線(xiàn)性壓縮機(jī)進(jìn)行建模仿真，并通過(guò)實(shí)驗(yàn)對(duì)仿真模型進(jìn)行驗(yàn)證，以期為線(xiàn)性壓縮機(jī)的設(shè)計(jì)和優(yōu)化改進(jìn)提供一種簡(jiǎn)易有效的方案。

1 建模

高頻回?zé)崾降蜏叵到y(tǒng)中各狀態(tài)量在時(shí)域中都可看作交變的正弦波[8]，基本的表示形式為：

利用歐拉公式改寫(xiě)式中動(dòng)態(tài)項(xiàng)Xd：

其中：表示交變的物理量，可根據(jù)正弦波的定義將其分解為Xˉ與動(dòng)態(tài)Xd。對(duì)于制冷系統(tǒng)，運(yùn)行過(guò)程中存在多個(gè)交變量，在時(shí)域中無(wú)法直觀(guān)的表現(xiàn)各物理信號(hào)之間的關(guān)系。一種更常用且直觀(guān)的做法是將時(shí)域上的交變信號(hào)映射到頻域中。在復(fù)平面中，動(dòng)態(tài)項(xiàng)映射的結(jié)果是以幅值Xa為長(zhǎng)度，ωt+φ為方向的向量。

基于以上原理，在復(fù)平面內(nèi)建立壓縮活塞的質(zhì)量?彈簧?阻尼模型，該模型在復(fù)平面中某一時(shí)刻的受力分析如圖 1 所示[9]，其中點(diǎn)折線(xiàn)箭頭表示此時(shí)刻位移、速度和加速度的相位，實(shí)線(xiàn)箭頭表示壓縮活塞在工作狀態(tài)下所受的力，表示電機(jī)推力，表示彈性力，表示慣性力，表示阻尼力，表示氣體力。

圖1 復(fù)平面活塞受力分析圖Fig.1 Force analysis of piston in the complex plane

式中：pa和xa分別為壓力和位移在時(shí)域中的幅值，A代表活塞截面積；ω為線(xiàn)性壓縮機(jī)的角頻率；θ為回?zé)崞鳠岫速|(zhì)量流領(lǐng)先壓力波的相位角，對(duì)于斯特林制冷機(jī)，實(shí)際設(shè)計(jì)中應(yīng)將其控制在40°至45°，從而減小回?zé)崞鲹p失；kp和Rp為合并常數(shù)項(xiàng)后的結(jié)果，分別代表氣體彈性剛度和氣體阻尼系數(shù)。二者表征了制冷機(jī)側(cè)對(duì)線(xiàn)性壓縮機(jī)的影響，通過(guò)改變壓力幅值和相位角可以模擬壓縮機(jī)在連接不同工況制冷機(jī)時(shí)的運(yùn)行情況。

COMSOL是一款基于有限元方法的多物理場(chǎng)仿真軟件，通過(guò)其內(nèi)置的模塊可以便捷地實(shí)現(xiàn)對(duì)電磁學(xué)、流體流動(dòng)、傳熱學(xué)等領(lǐng)域的仿真和對(duì)不同物理場(chǎng)的耦合。本文采用COMSOL對(duì)一款商用Cryo S 100斯特林制冷機(jī)動(dòng)磁式線(xiàn)性壓縮機(jī)工作過(guò)程中的電磁場(chǎng)及壓縮機(jī)在電磁場(chǎng)中的運(yùn)動(dòng)狀況進(jìn)行了動(dòng)態(tài)分析。

線(xiàn)性壓縮機(jī)的裝配模型與仿真模型的對(duì)應(yīng)關(guān)系如圖2所示，運(yùn)動(dòng)活塞由板彈簧支撐，與氣缸內(nèi)壁形成間隙密封；永磁體2是活塞的一部分，位于內(nèi)定子1和齒形磁軛3之間的氣隙內(nèi)。當(dāng)齒形磁軛包裹的線(xiàn)圈4中通入正弦電流時(shí)，永磁體會(huì)受到交變的軸向電磁力從而帶動(dòng)活塞做往復(fù)直線(xiàn)運(yùn)動(dòng)。

圖2 線(xiàn)性壓縮機(jī)裝配模型與仿真模型對(duì)照?qǐng)DFig.2 The assembly model and simulation model of linear compressor

建立仿真模型的參數(shù)如表1所示。為了提高計(jì)算效率，對(duì)裝配模型進(jìn)行以下的假設(shè)和簡(jiǎn)化：

（1）忽略模型在徑向上的運(yùn)動(dòng)；

（2）簡(jiǎn)化活塞的完整形狀，將所有的質(zhì)量和受力集中到永磁體上；

（3）假定板彈簧對(duì)活塞提供的軸向力與位移是線(xiàn)性相關(guān)的，且摩擦因素被簡(jiǎn)化為固定的機(jī)械阻尼系數(shù)；

（4）使用相對(duì)磁導(dǎo)率為1的氣體介質(zhì)填充模型組件以外的區(qū)域。

表1 線(xiàn)性電機(jī)的基本參數(shù)Tab.1 Primary structure parameter of linear motor

使用不同的磁性本構(gòu)關(guān)系描述磁場(chǎng)中軟磁材料和永磁體的磁性特征。對(duì)于內(nèi)定子和磁軛使用的Somaloy軟磁材料，通過(guò)導(dǎo)入圖3所示的外部的材料磁化屬性表，對(duì)其磁性能進(jìn)行了描述。對(duì)于永磁體，指定了剩余磁感應(yīng)強(qiáng)度關(guān)系式，如式（4）所示，其中表示剩余磁感應(yīng)強(qiáng)度，μ0表示真空磁導(dǎo)率，μr表示永磁體的相對(duì)磁導(dǎo)率，表示徑向矯頑力。對(duì)應(yīng)該系統(tǒng)中使用的永磁體，μr=1.02，

圖3 軟磁材料B-H曲線(xiàn)Fig.3 B-H curve of Somaloy

根據(jù)圖1，使用“全局常微分和微分代數(shù)方程”接口導(dǎo)入帶氣體力的運(yùn)動(dòng)方程描述活塞的機(jī)械運(yùn)動(dòng)。電磁場(chǎng)和運(yùn)動(dòng)方程通過(guò)位移量x和電磁力FI相互耦合，在線(xiàn)圈內(nèi)指定輸入電流，可以求解電磁力、位移量和速度量隨時(shí)間的變化關(guān)系。

為仿真模型添加移動(dòng)網(wǎng)格時(shí)，考慮到運(yùn)動(dòng)部件只有活塞，可以將模型如圖4所示分為內(nèi)定子、永磁體和包裹線(xiàn)圈的磁軛三個(gè)部分。僅在中間的永磁體區(qū)域使用移動(dòng)網(wǎng)格接口。在該區(qū)域內(nèi)，將圖4中紅色線(xiàn)框所包裹的永磁體和其周?chē)囊徊糠謿怏w作為整體，指定其軸向位移為微分方程接口求解得到的位移量x，上下方的氣體域設(shè)置為被動(dòng)變形的自由變形域，并使用了映射網(wǎng)格來(lái)減少由于網(wǎng)格劃分帶來(lái)的計(jì)算誤差，防止發(fā)散。

圖4 模型網(wǎng)格剖分Fig.4 model meshed geometry

電磁場(chǎng)的模擬必須使用麥克斯韋應(yīng)力張量對(duì)活塞受到的電磁力進(jìn)行精確的求解，因此在永磁體和氣隙部分繪制了更為精密的網(wǎng)格。同時(shí)為了更好地保證移動(dòng)域和非移動(dòng)域之間物理量的連續(xù)性，在二者相交的邊界上，移動(dòng)域一側(cè)必須進(jìn)行更精細(xì)的邊界單元?jiǎng)澐帧Ｗ罱K網(wǎng)格的平均質(zhì)量為0.82。經(jīng)無(wú)關(guān)性驗(yàn)證，進(jìn)一步加密網(wǎng)格對(duì)求解結(jié)果的影響不超過(guò)5%，求解器步長(zhǎng)設(shè)置為0.001 s。

2 實(shí)驗(yàn)流程

2.1 空載實(shí)驗(yàn)裝置

圖5為測(cè)試空載工況活塞位移的裝置示意圖，線(xiàn)性壓縮機(jī)由變頻電源驅(qū)動(dòng)，可直接讀出輸入電壓、輸出電流和有效功率；壓縮機(jī)活塞板簧軸向上方裝有激光位移傳感器，能對(duì)活塞板彈簧的運(yùn)動(dòng)進(jìn)行實(shí)時(shí)測(cè)量，活塞板彈簧的中心孔通過(guò)緊固螺母和活塞相連接，二者的運(yùn)動(dòng)始終具有一致性；在空載條件下設(shè)置輸出電壓和頻率為80 Hz、7 V，記錄此時(shí)壓縮機(jī)活塞的位移結(jié)果，與仿真結(jié)果進(jìn)行對(duì)比。

圖5 空載工況下活塞位移測(cè)試裝置實(shí)物圖和示意圖Fig.5 Physical diagram and schematic diagram of piston displacement test device under on-load condition

2.2 比推力實(shí)驗(yàn)裝置

將壓縮機(jī)直線(xiàn)電機(jī)平穩(wěn)放置在拉壓力測(cè)試儀的托盤(pán)中心。直線(xiàn)電機(jī)由直流電源驅(qū)動(dòng)，可以直接讀出輸入電壓和輸入電流；拉壓力測(cè)試儀在壓縮活塞軸向上方連接有壓力觸頭，可以通過(guò)手輪控制的限位機(jī)構(gòu)精確限制壓縮活塞的位置，并對(duì)活塞的推力進(jìn)行測(cè)量。

測(cè)量時(shí)，設(shè)置直流電源輸出電流依次為1 A、2 A、3 A、4 A、5 A；記錄活塞在不同限位條件下推力與電流的變化關(guān)系，從而對(duì)仿真結(jié)果進(jìn)行對(duì)比驗(yàn)證。

3 討論

3.1 空載實(shí)驗(yàn)

在模型中輸入與實(shí)驗(yàn)相同的空載工況，對(duì)活塞的位移情況進(jìn)行了計(jì)算。采用測(cè)試裝置對(duì)多臺(tái)壓縮機(jī)進(jìn)行了空載測(cè)試，測(cè)試結(jié)果如圖6所示。在隨機(jī)抽取的50臺(tái)線(xiàn)性壓縮機(jī)中，有四臺(tái)壓縮機(jī)的活塞位移與計(jì)算值出現(xiàn)較大偏差，其余樣機(jī)的測(cè)試結(jié)果平均分布在計(jì)算值附近。進(jìn)一步比較后發(fā)現(xiàn)，出現(xiàn)較大偏差的四臺(tái)線(xiàn)性壓縮機(jī)的活塞與排出器之間的位移相位差遠(yuǎn)大于其余線(xiàn)性壓縮機(jī)的位移相位差，推測(cè)是在裝配過(guò)程中零件安裝不當(dāng)產(chǎn)生的影響。剔除錯(cuò)誤樣本后，測(cè)試值與計(jì)算值的平均偏差為3.49%。

圖6 空載工況下活塞位移的結(jié)果分析圖Fig.6 The result of displacement in no-load condition

3.2 比推力測(cè)試

在SAGE模型中，壓縮機(jī)端的影響是通過(guò)指定比推力引入的，為驗(yàn)證壓縮機(jī)模型的合理性，使用前文所述測(cè)量方式，依次測(cè)量了活塞在距離平衡位置0 mm、1 mm、2 mm、3 mm、4 mm處的比推力。同時(shí)在模型中將壓縮活塞的永磁體固定在對(duì)應(yīng)的位置，設(shè)定線(xiàn)圈的電流參數(shù)，從穩(wěn)態(tài)的角度對(duì)壓縮活塞在不同位置上的比推力進(jìn)行求解，結(jié)果如圖7所示，實(shí)測(cè)值略大于計(jì)算值，最大偏差出現(xiàn)在離初始位置1 mm處，為7.2%；最小偏差出現(xiàn)在離初始位置4 mm處，為4.37%，無(wú)論是實(shí)測(cè)值還是計(jì)算值，比推力的最大變化不超過(guò)0.4 N∕A。在模型中進(jìn)一步增大與初始點(diǎn)的距離后，發(fā)現(xiàn)當(dāng)壓縮活塞位移超過(guò)6 mm時(shí)，活塞會(huì)由于逐漸脫離線(xiàn)圈生成的磁場(chǎng)，導(dǎo)致比推力迅速衰減。

3.3 負(fù)載條件下與Sage模型比較

該模型中，制冷機(jī)側(cè)的影響被簡(jiǎn)化為氣體力，為了驗(yàn)證這種方法的可行性，使用SAGE模型計(jì)算了該款制冷機(jī)在100 K、193 K、253 K和293 K四個(gè)工況點(diǎn)下的動(dòng)態(tài)壓力幅值、壓縮機(jī)出口的壓力?質(zhì)量流相位、阻抗等參數(shù)，將它們換算為氣體力的形式輸入COMSOL模型進(jìn)行計(jì)算，并對(duì)比了用兩種模型計(jì)算的活塞位移幅值。該制冷機(jī)使用表1所列的線(xiàn)性電機(jī)。計(jì)算結(jié)果如圖8和表2所示，在不同的冷頭極限溫度下，對(duì)壓縮活塞位移的計(jì)算結(jié)果均吻合的較好，最大偏差不超過(guò)3%。

圖7 比推力隨位移的變化關(guān)系比較曲線(xiàn)Fig.7 The measurement and simulated result of specific thrust vs.displacement

圖8 兩種模型計(jì)算結(jié)果對(duì)比曲線(xiàn)Fig.8 Comparison of displacement in two models

表2 基于壓力幅值的位移幅值比較結(jié)果Tab.2 Comparison of displacement base on pressure amplitude

4 結(jié)論

本文基于斯特林制冷機(jī)電?機(jī)械?聲的耦合原理，將聲阻抗拆分為位移和速度方向分量，建立了帶負(fù)載條件的直線(xiàn)電機(jī)COMSOL模型，并進(jìn)行了空載測(cè)試、比推力測(cè)試和與SAGE模型的對(duì)比驗(yàn)證。其中空載測(cè)試的50臺(tái)電機(jī)中，測(cè)試值與計(jì)算值的平均偏差為3.49%；對(duì)位于磁場(chǎng)不同位置處活塞的比推力測(cè)試表明，測(cè)試值和計(jì)算值的最大偏差出現(xiàn)在距離平衡位置1 mm處，為7.2%，最小偏差出現(xiàn)在距離平衡位置4 mm處，為4.37%。計(jì)算結(jié)果表明，當(dāng)壓縮活塞距離平衡位置6 mm以上時(shí)，比推力將迅速衰減；在不同負(fù)載條件下，該模型的活塞位移計(jì)算結(jié)果和SAGE的計(jì)算結(jié)果均吻合良好，最大偏差不超過(guò)3%。