基于一維非穩(wěn)態(tài)傳熱模型的爐溫曲線優(yōu)化研究

張曉果，路 月，林宗明，李政君,杜亞冰，蘭奇遜

(1.河南城建學(xué)院 數(shù)理學(xué)院，河南 平頂山 467036； 2.河南城建學(xué)院 能源與建筑環(huán)境工程學(xué)院，河南 平頂山 467036; 3.河南城建學(xué)院 電氣與控制工程學(xué)院，河南 平頂山 467036)

回流焊技術(shù)作為電子產(chǎn)品生產(chǎn)過程中的重要環(huán)節(jié)，通過控制回焊爐各部分的溫度，在保持工藝要求的同時，來提高產(chǎn)品的質(zhì)量一直是現(xiàn)代電子工業(yè)亟待突破的課題[1]。李巖等[2]對某型回流焊爐進行分析，闡述了其爐溫曲線設(shè)定的方法，并引入加熱因子等相關(guān)概念對原溫度參數(shù)進行優(yōu)化；馮志剛等[3]分析了回焊爐中工藝參數(shù)對爐溫曲線的影響，得到爐溫和傳送帶速度兩個關(guān)鍵影響指標(biāo)，并對指標(biāo)進行深層次研究；龔雨兵[4]以最大熱應(yīng)力為優(yōu)化目標(biāo)，對回流焊溫度場進行數(shù)值建模與仿真，得到最優(yōu)溫度分布及熱應(yīng)力分布；朱桂兵等[5]以助焊劑活性導(dǎo)致焊接缺陷的角度反推溫度曲線設(shè)置的不合理性，結(jié)合熱效能理念尋求解決辦法。以上研究大多都是直接從工藝參數(shù)的機理入手，結(jié)合不同加工工藝要求，分析影響溫度曲線的主導(dǎo)因素，以該主導(dǎo)因素為變量進行相應(yīng)的靈敏度分析。但實際上,盲目試驗不僅增加工作量，而且還需通過測試數(shù)據(jù)來設(shè)置工藝參數(shù)，以實現(xiàn)要求的爐溫曲線效果。汪學(xué)軍[6]從實際應(yīng)用出發(fā)，以能量守恒定律與傅里葉定律的熱平衡機理為依據(jù)對焊接過程進行分析，建立了一種全新的多溫區(qū)自動整定策略，通過自學(xué)力能力對控制參數(shù)進行了優(yōu)化，找到了較好的控制參數(shù)。但是，通過數(shù)值模擬對機理模型進一步分析，直接獲取回流焊各控制參數(shù)以及對爐溫曲線進行數(shù)值優(yōu)化的工作目前仍比較欠缺。

本文首先對某回焊爐各溫區(qū)的溫度與時間變化進行線性回歸，并對焊接過程溫度場進行有限元分析，發(fā)現(xiàn)對流傳熱系數(shù)為關(guān)鍵性指標(biāo)，根據(jù)爐溫設(shè)置的特定要求，繪制爐溫曲線。然后通過有限差分法更改初始設(shè)定的爐溫以及傳送帶速度，對原有爐溫曲線進行優(yōu)化，使得傳熱系數(shù)達到最優(yōu)，獲得最佳指標(biāo)值，并通過MATLAB進行蒙特卡洛數(shù)值模擬檢驗。研究發(fā)現(xiàn)，該方法能有效地減少物理試驗次數(shù)以獲得較優(yōu)的爐溫曲線設(shè)置，從而使得爐溫設(shè)置與焊接機理能更好地關(guān)聯(lián)起來。

1 回焊爐介紹

1.1 回焊爐的組成

選取的回焊爐其內(nèi)部由爐前區(qū)域、爐后區(qū)域及11個小溫區(qū)構(gòu)成，其中各小溫區(qū)長度均為30.5 cm，相鄰間隙為5.0 cm,在加熱區(qū)外有爐前、爐后兩長度為25 cm的區(qū)域，如圖1所示。從其截面角度觀察，回焊爐內(nèi)部又分為4個功能區(qū)(大溫區(qū))：預(yù)熱區(qū)、恒溫區(qū)、回流區(qū)、冷卻區(qū)。

圖1 回焊爐截面示意圖

1.2 爐溫分布及設(shè)置

假設(shè)回流焊爐內(nèi)空氣溫度短時間能達到穩(wěn)定，處于正常工作狀態(tài),各溫區(qū)設(shè)定爐溫即為爐內(nèi)熱風(fēng)溫度。在回焊爐工藝焊接過程中，爐溫曲線需滿足不同工藝的要求，相關(guān)變量有相應(yīng)的制程界限，具體見表1。將傳送帶速度初始設(shè)定為70 cm/min。各溫區(qū)溫度設(shè)定值分別為175 ℃(小溫區(qū)1～5)、195 ℃(小溫區(qū)6)、235 ℃(小溫區(qū)7)、255 ℃(小溫區(qū)8～9)及25 ℃(小溫區(qū)10～11)。

表1 爐溫曲線的制程界限

1.3 實際問題的數(shù)學(xué)化處理

如圖2所示，電路板在回焊爐內(nèi)焊接時，其受熱情況實際上是三維立體溫度擴散，假設(shè)熱輻射的影響極小，因此可忽略平行板的高度和寬度，將立體平行板簡化為只與厚度有關(guān)的一維空間?；睾笭t各功能區(qū)溫度曲線分布如圖3所示。本文需要解決的問題就是通過建立合適的數(shù)學(xué)模型，在設(shè)定傳送帶速度和各溫區(qū)爐溫的前提下，求得爐溫曲線的溫度分布，并在建?；A(chǔ)上進一步優(yōu)化爐溫曲線，在符合制程界限條件下，求得最大過爐速度，從而在保證產(chǎn)品質(zhì)量的同時又提高其加工效率。

圖2 回焊爐中電路板的受熱情況

圖3 爐溫曲線各溫區(qū)分布

2 爐溫曲線建模及機理分析

2.1 基于熱傳導(dǎo)方程的爐溫建模

在回流焊接過程中，其熱量Q的傳遞主要由對流換熱Qc、輻射換熱Qr和熱傳導(dǎo)Qk三種基本傳熱方式組成，滿足能量守恒定律[4]：

Q=Qc+Qr+Qk

(1)

單位時間內(nèi)通過單位面積的流體質(zhì)量由溫度T1流至T2處所傳遞的熱量為：

Qc=CpM(T1-T2)

(2)

其中，Cp是定壓比熱容，M是流體質(zhì)量。

對流傳熱受多種因素影響，主要考慮流體與物體接觸面的熱交換,即電路板表面與熱風(fēng)的對流換熱量。牛頓冷卻公式[6]為：

Qc=hΔTA

(3)

其中，h是對流傳熱系數(shù)，ΔT表示平壁兩側(cè)溫差，A是發(fā)生對流的面積。

為了推導(dǎo)物體熱交換的微分方程,應(yīng)根據(jù)質(zhì)量守恒定律和傅里葉定律來建立物體中的導(dǎo)熱溫度。導(dǎo)熱物體各方向熱流量可分解為x,y,z三個坐標(biāo)軸方向的分熱流量，通過y-z、x-z、x-y三個微元表面導(dǎo)入微元體熱流量分別用qx、qy、qz來表示，根據(jù)傅里葉定律[7]得出：

(4)

假設(shè)各分層內(nèi)不含熱源，故可忽略介質(zhì)比熱和密度的影響，在三維等方向均勻介質(zhì)中的熱傳導(dǎo)方程需滿足：

(5)

為了進一步簡化模型，可將各分層視為無限大，僅考慮在厚度方向即x方向上的溫度變化，便可將式(5)化簡為：

(6)

其中，λ是熱導(dǎo)率，α是介質(zhì)熱擴散率，T表示電子元件的中心溫度。

2.2 分段線性回歸分析

假設(shè)x是時間，為自變量；y是溫度，為因變量。

yi=β0+β1x1+μ

(7)

其中，β0和β1為回歸系數(shù)，μ為無法預(yù)測的且滿足一定條件的擾動項。

原假設(shè)H0：聯(lián)合顯著性檢驗β1=0。

根據(jù)各溫區(qū)溫度隨時間的變化關(guān)系進行分段線性回歸，采用OLS最小二乘估計法，確定各溫區(qū)的回歸系數(shù)，如表2所示。

表2 回歸系數(shù)表

由表2可知:擬合優(yōu)度R2均接近于1，且P小于0.05，說明可以拒絕原假設(shè)，回歸性顯著。按照上述方法將各段位置溫度的不同影響代入不同的控制溫度點，進一步建立線性關(guān)系，用MATLAB繪制爐溫曲線如圖4所示。由圖4可知:在升溫過程中，得到模型結(jié)果的溫度隨時間變化率的線性擬合效果較好；而在降溫過程中，擬合數(shù)據(jù)的溫度下降速度略慢于實際數(shù)據(jù)的下降速度，使得曲線后半部分存在略微差異?？紤]到上文中提到的制程界限對溫度下降速率存在一定影響，因此使得模型中得到的降溫過程中的爐溫曲線受到了相應(yīng)限制，這也是模型的合理之處。

圖4 爐溫擬合曲線對比圖

2.3 有限元分析結(jié)果

(1)繪圖

繪制單個焊爐的三維幾何模型，見圖5。焊接回流的內(nèi)部近似一個長方體，高溫氣流從焊爐上方的風(fēng)機流出，對流換熱后被下面的吸收器吸收，電子元件從左側(cè)的間隙進入，接受管內(nèi)高溫氣流對流傳熱，隨傳送帶從右側(cè)間隙進入下一個焊爐，因為爐內(nèi)高溫氣體均勻，因此忽略焊爐高度對氣流的影響。單個焊爐的幾何長度為30.5 cm，橫截面為圓形，直徑為15 cm。

圖5 單個焊爐的三維幾何模型

設(shè)置左側(cè)電子元件入口為低溫入口，右側(cè)為高溫出口，高溫氣體從上方吹出，設(shè)置method的流體條件，設(shè)置除元件入口、出口外的幾何外壁的膨脹層，然后劃分網(wǎng)格。通過網(wǎng)格劃分將幾何體分解為多個單元格，讓計算機在不同的單元格中模擬分析，增加結(jié)果的穩(wěn)定性，使結(jié)果收斂，防止模擬數(shù)據(jù)失真。

(3)前處理

設(shè)置氣流的流速、溫度、動力黏度系數(shù)，設(shè)置入口溫度和物體速度。

(4)后處理

模擬爐內(nèi)穩(wěn)定溫度時的溫度云圖，并由ANSYS軟件求解繪制溫度折線如圖6所示。

圖6 溫度折線圖

通過單個焊爐的溫度輪廓如圖7所示，分析爐中溫度變化，得到電子元件加熱焊接的直觀圖像，使溫度變化梯度可視化，并與模擬數(shù)據(jù)對比，校驗?zāi)M數(shù)據(jù)。

她死得很安詳。四年時間足以讓她面對現(xiàn)實。并且，某一天，她突然驚喜地發(fā)現(xiàn)，她再一次迎來了久違的愛情——從未被設(shè)計和植入的愛情——真正的愛情。她愛上秦川，玩偶愛上玩偶，這感覺讓她戰(zhàn)栗。

圖7 溫度輪廓圖

由圖7可以看出，物體進入開口后受到爐內(nèi)高溫氣體的對流換熱，吸收爐內(nèi)高溫氣流熱量，造成爐內(nèi)溫度分布不均勻，但在工作時爐內(nèi)溫度趨于穩(wěn)定。

3 機理模型下的爐溫曲線優(yōu)化

3.1 一維非穩(wěn)態(tài)傳熱模型的建立

在熱傳導(dǎo)方程基礎(chǔ)上，根據(jù)傅里葉定律和能量守恒定律，在任一時間間隔內(nèi)有如下熱平衡關(guān)系：導(dǎo)入微元體的總熱流量=導(dǎo)出微元體的總熱流量+微元體內(nèi)熱力學(xué)能的增量。

其控制方程[8]為：

(8)

其中，ρ是密度，c是比熱容。

對于整個爐內(nèi)傳熱模型，屬于第三類邊界條件，對流換熱對物體提供熱量，假設(shè)在高溫環(huán)境下，爐內(nèi)溫度分布的初始值為當(dāng)前溫區(qū)溫度。邊界條件和初始條件滿足：

(9)

其中，L是導(dǎo)熱長度，Tlu是爐內(nèi)溫度，T0是環(huán)境溫度。

在非穩(wěn)態(tài)一維傳熱模型中，通過建立最小二乘法估計換熱系數(shù)參數(shù)：

(10)

3.2 有限差分法求解

傳熱模型求解的基本思想是把空間、時間上相同且連續(xù)的物理系數(shù)離散化，利用有限差分法進行數(shù)值計算求解。使用顯式差分對傳熱模型模擬離散，假設(shè)物體各個方向傳熱為各向同性，由牛頓冷卻公式聯(lián)立方程計算求解第(n+1)時間軸上的溫度Tn+1，依靠前一層的溫度大小，控制方程離散其中的未知量。

(11)

對一維非穩(wěn)態(tài)傳熱模型進行離散化后，由設(shè)立的邊界條件和初始值，在不同的時間和空間節(jié)點上進行分層離散求解，建立未知值和元件中心溫度值的關(guān)系，進而搜索元件的熱導(dǎo)率，求解擬合溫度數(shù)據(jù)的最優(yōu)組合，具體步驟為：

(1)代入爐溫的初始值，通過模型離散方程逐層求解，得到元件中心溫度的h值；

(2)更新溫度值，迭代離散方程進行求解，得到優(yōu)化后的h值；

(3)遍歷新的溫度值，全局搜索擬合程度最佳且滿足制程界限的對流傳熱系數(shù)，并確定最優(yōu)指標(biāo)值；

(4)根據(jù)最佳對流傳熱系數(shù)組合，得到元件中心溫度分布。

3.3 結(jié)果分析

根據(jù)上述方法求得各小溫區(qū)最佳對流傳熱系數(shù)下的最優(yōu)溫度及最快過爐速度，見表3。

表3 最優(yōu)爐溫曲線的指標(biāo)值

優(yōu)化目標(biāo)函數(shù) (最優(yōu)傳熱系數(shù)) 的優(yōu)化歷程如圖8所示。

圖8 最優(yōu)傳熱系數(shù)下的優(yōu)化曲線

由圖8可看出,各個小溫區(qū)溫度(除冷卻區(qū)之外)相較于最初設(shè)定溫度均有相應(yīng)的改變，原因是冷卻區(qū)溫度確定之后，對應(yīng)的溫度時間便可求得，利用有限差分法和遞推公式對其他溫區(qū)進行刻畫，使得優(yōu)化曲線更為平滑，連續(xù)性更好。與基于熱傳導(dǎo)方程的模型相比，該模型保證了各功能區(qū)溫度過渡的連貫性。

3.4 蒙特卡洛模擬檢驗

考慮到傳送帶的速度限制會對爐溫曲線的變化造成影響，為了提高模型數(shù)據(jù)的準(zhǔn)確性還可以進行蒙特卡洛模擬。先將速度區(qū)間離散化，然后基于傳熱模型，遍歷傳熱過程中每個傳送帶速度，得到相應(yīng)的溫度分布，代入約束條件進行檢驗，最后求得決策后的最優(yōu)傳熱系數(shù)下的最大速度。

最優(yōu)傳熱系數(shù)下的約束條件：

(12)

其中，T(tj)為在第j個時間微元下的溫度(j=1,2,3…)。

根據(jù)上述約束條件，提出一種搜索算法，其過程如下：

Step1：確定速度初始區(qū)間v∈[65,100]，令v=100；

Step2：繪制當(dāng)前傳送帶速度下的爐溫曲線；

Step3：驗證約束條件(式(12))，直到當(dāng)前速度為允許最大過爐速度時結(jié)束，否則Step4；

Step4：設(shè)置步長為1，即v=v-1，更新搜索區(qū)間，回到Step2繼續(xù)循環(huán)。

根據(jù)以上搜索算法，依次迭代求得滿足加工條件的最優(yōu)爐溫曲線，對升、降速度進行求解，利用MATLAB軟件編程得到?jīng)Q策后的最優(yōu)傳熱系數(shù)下的最大過爐速度為76 cm/min，并繪制出當(dāng)前速度下溫度分布曲線，如圖9所示。結(jié)果發(fā)現(xiàn)與有限差分法所求解的結(jié)果基本相符，從而通過了一致性檢驗，同時也說明基于有限差分法的傳熱模型符合優(yōu)化目標(biāo)參數(shù)的求解。

圖9 v=76 cm/min時的爐溫曲線圖

4 結(jié)論

通過建立一維非穩(wěn)態(tài)傳熱模型，以最優(yōu)傳熱系數(shù)為優(yōu)化目標(biāo)，利用有限差分法和遞推公式對已有的爐溫曲線進一步優(yōu)化，再通過蒙特卡洛模擬證明模型的有效性，主要結(jié)論如下：

(1)利用有限元分析溫度輪廓圖可知，物體的進入使得溫度場紊亂，在工作情況下溫度趨于穩(wěn)定，且溫度變化趨勢很小，可以忽略溫度場不均勻的導(dǎo)熱對物體溫度變化的影響；

(2)非穩(wěn)態(tài)傳熱模型相較于基于熱傳導(dǎo)方程下的分段線性回歸模型，其溫度隨時間的變化率在各功能區(qū)內(nèi)不定，隨環(huán)境變化的自適應(yīng)性較強，優(yōu)化后的曲線更為平滑；

(3)基本上解決了回歸過程中存在的內(nèi)生性問題，保證了各功能區(qū)溫度過渡的連貫性，借用該模型可大大減少物理試驗次數(shù)而獲得可行的爐溫曲線，并能通過優(yōu)化獲得特定制程界限下的最優(yōu)爐溫曲線設(shè)置。