基于COMSOL Multiphysics對(duì)含頂隙液態(tài)罐頭食品熱殺菌的模擬研究

陳舒航，周建偉，劉東紅，4*，劉 懿

（1 浙江大學(xué)生物系統(tǒng)工程與食品科學(xué)學(xué)院 杭州310058 2 浙大寧波理工學(xué)院 浙江寧波315100 3 江蘇省食品先進(jìn)制造裝備技術(shù)重點(diǎn)實(shí)驗(yàn)室 江蘇無錫214122 4 浙江大學(xué)馥莉食品研究院 杭州310058）

罐藏是一種常見的食品貯藏方法，廣泛運(yùn)用于食品工業(yè)中。熱力殺菌因其具有高效、方便的特點(diǎn)，是食品加工中最常用的殺菌手段[1]。罐藏食品的熱殺菌涉及流體傳熱，用傳統(tǒng)方法難以得到罐頭內(nèi)的任意時(shí)刻每個(gè)位置的溫度分布情況，對(duì)罐頭內(nèi)冷點(diǎn)和最低殺菌強(qiáng)度判斷的準(zhǔn)確度不足[2]。借助計(jì)算機(jī)對(duì)熱殺菌過程進(jìn)行仿真模擬，能有效預(yù)測(cè)不同類型食品在該過程中溫度、速度、壓力、最慢加熱區(qū)（Slowing heating zone，SHZ）及殺菌強(qiáng)度的變化[3-8]。將仿真結(jié)果應(yīng)用于調(diào)整幾何尺寸、熱處理溫度、處理時(shí)間、旋轉(zhuǎn)速度等參數(shù)[9-12]，在保證殺菌強(qiáng)度的同時(shí)，減少殺菌時(shí)間，提高產(chǎn)品風(fēng)味和營(yíng)養(yǎng)水平，節(jié)約能耗。精準(zhǔn)的仿真模型可提高試驗(yàn)的精準(zhǔn)性和有效性，為設(shè)備研發(fā)及確定最優(yōu)控制方案提供依據(jù)。

罐頭頂隙是指罐藏食品的表面或液面與罐頭容器上蓋間所留空隙，可維持一定的真空度，平衡罐頭內(nèi)外的壓力差，避免在熱殺菌過程中出現(xiàn)凸角、癟聽等外觀變形的情況[13]。目前對(duì)罐頭熱殺菌過程的模擬大多忽略了頂隙對(duì)傳熱的影響，且未考慮頂隙處存在的液體蒸發(fā)和冷凝現(xiàn)象，因而無法對(duì)殺菌過程中的壓力變化進(jìn)行準(zhǔn)確預(yù)測(cè)[14-15]。忽略頂隙也會(huì)導(dǎo)致對(duì)溫度的預(yù)測(cè)產(chǎn)生偏差，使得最小殺菌值的計(jì)算不夠準(zhǔn)確。

羧甲基纖維素鈉（Carboxymethylcellulose，CMC）具有與一些液態(tài)食品相似的流變學(xué)特性，常被用于液態(tài)罐頭食品熱殺菌的模擬過程[16]。本文選取1.0% CMC 溶液，借助COMSOL Multiphysics軟件，構(gòu)建了包含及不包含頂隙的2 種模型模擬罐頭食品的熱殺菌過程。通過熱穿透試驗(yàn)，由模擬值與真實(shí)值之間的平均相對(duì)誤差驗(yàn)證模型可靠性，以期獲得一種較準(zhǔn)確的罐頭熱殺菌過程的預(yù)測(cè)模型，更精準(zhǔn)地分析熱殺菌過程中速度場(chǎng)、溫度場(chǎng)、致死率場(chǎng)及壓力的變化情況。在含頂隙模型的基礎(chǔ)上，探究不同頂隙高度及保溫溫度對(duì)1.0%CMC 罐頭內(nèi)相對(duì)壓力的影響，為工廠生產(chǎn)罐頭時(shí)所預(yù)留的頂隙高度及不同工作溫度下殺菌釜的壓力控制提供參考。

1 材料與方法

1.1 材料試劑

金屬罐（罐號(hào)7113），山東環(huán)球印鐵制罐有限公司；純凈水，娃哈哈公司；羧甲基纖維素（化學(xué)純級(jí)），國(guó)藥集團(tuán)化學(xué)試劑有限公司。

1.2 儀器設(shè)備

Data Trace MP Ⅲ無線溫度傳感器及工作站，美國(guó)Mesa Laboratories 公司；YJ-C20 型手動(dòng)封罐機(jī)，張家港易捷自動(dòng)化設(shè)備；50 L 安全智能型反壓高溫蒸煮鍋，濟(jì)南蘭光機(jī)電技術(shù)有限公司。

1.3 試驗(yàn)方法

采用水和1% CMC 溶液分別對(duì)數(shù)值模型進(jìn)行試驗(yàn)驗(yàn)證，具體試驗(yàn)步驟如下：

1）水的模型驗(yàn)證：用自制支架將溫度傳感器固定在罐頭中心，探針高度位于距底部15 mm 處。在室溫下，將純凈水灌入馬口鐵罐頭中，使液面高度為95 mm，準(zhǔn)備3 罐，用封罐機(jī)將馬口鐵罐封口后放入反壓高溫蒸煮鍋中，同時(shí)在殺菌釜內(nèi)放置1 個(gè)溫度傳感器。隨后殺菌釜開始進(jìn)行保溫溫度121 ℃-保溫時(shí)間30 min 的自動(dòng)殺菌。冷卻完成后，取下4 個(gè)溫度傳感器，讀取殺菌鍋溫度及罐內(nèi)液體預(yù)估冷點(diǎn)溫度。將溫度探針固定于距底部100 mm 的位置，進(jìn)行相同試驗(yàn)，獲取罐頭頂隙處溫度變化曲線。

圖1 熱滲透試驗(yàn)裝置示意圖Fig.1 Schematic diagram of thermal penetration experimental device

2）1.0% CMC 溶液的模型驗(yàn)證：配制質(zhì)量分?jǐn)?shù)為1.0%的CMC 溶液，用同樣方法測(cè)定1.0%CMC 罐頭液體及頂隙處溫度隨時(shí)間變化的曲線。

計(jì)算試驗(yàn)所得實(shí)際溫度與數(shù)值模型的模擬溫度的平均相對(duì)誤差，以評(píng)價(jià)模型的準(zhǔn)確性。

2 數(shù)值模型的搭建

在COMSOL Multiphysics 5.3a 軟件中構(gòu)建了2 種二維軸對(duì)稱模型，對(duì)液體罐頭食品的熱殺菌過程進(jìn)行模擬。模型一忽略了頂隙的影響，認(rèn)為罐內(nèi)被液體充滿，而模型二液面高度為95 mm，頂部為濕空氣，同時(shí)考慮濕空氣與液體表面水分的蒸發(fā)冷凝。利用層流、濕空氣中的傳熱、空氣中的水分輸送3 個(gè)物理場(chǎng)模擬傳熱傳質(zhì)過程，并添加了熱濕、水分輸送、非等溫流動(dòng)節(jié)點(diǎn)來實(shí)現(xiàn)傳熱與蒸汽輸送之間的多物理場(chǎng)耦合。

2.1 模型假設(shè)

對(duì)2 種模型均做以下假設(shè)以簡(jiǎn)化求解過程：罐體內(nèi)液體、氣體均勻?qū)ΨQ；忽略罐壁、探針及固定裝置等對(duì)傳熱的影響；液體在罐內(nèi)壁無滑移；加熱過程中罐體不發(fā)生形變；罐體外壁的溫度分布均勻且等于殺菌釜的溫度。

2.2 控制方程

用Navier-Stokes 方程模擬空氣域中的流體流動(dòng)并用動(dòng)量守恒方程、能量守恒方程模擬相應(yīng)的對(duì)流傳熱[17]。由于涉及水分蒸發(fā)，需將蒸發(fā)的水蒸氣質(zhì)量作為潛熱源添加到輸運(yùn)方程中[18]。蒸發(fā)的水蒸氣質(zhì)量按式（1）計(jì)算：

式中，mvap——蒸發(fā)的水蒸氣質(zhì)量，kg；K——蒸發(fā)率，s-1；aw——水分活度；csat——飽和水分濃度，mol/L；c——水分濃度，mol/L。

式中，psat——飽和壓力，Pa；R——理想氣體常數(shù)，8.314 J/（mol·K）；T——溫度，K。

將蒸發(fā)熱作為潛熱源插入到傳熱方程中：

式中，Q——蒸發(fā)熱，J；Hvap——蒸發(fā)潛熱，44 172 J/mol。

2.3 計(jì)算條件

初始溫度為室溫，初始絕對(duì)壓力（-40 kPa）作為罐內(nèi)的真空度。純水的熱物理特性參數(shù)由COMSOL Multiphysics 內(nèi)置材料庫(kù)給出，1.0%CMC 溶液的參數(shù)如表1[19]。

表1 1.0% CMC 溶液的熱物理參數(shù)特性Table 1 Thermal physical properties of 1.0% CMC solution

2.4 網(wǎng)格剖分及求解

2 種模型均利用細(xì)化的物理場(chǎng)控制的自由三角形網(wǎng)格進(jìn)行剖分，模型一與模型二分別將幾何體剖分為5 429 及5 971 個(gè)單元進(jìn)行求解。

在模型驗(yàn)證試驗(yàn)中，將實(shí)際測(cè)量得到的保溫溫度121 ℃，保溫時(shí)間為30 min 的殺菌釜內(nèi)溫度作為內(nèi)插函數(shù)，導(dǎo)入軟件中作為罐頭外壁的溫度函數(shù)，以表征罐體在殺菌釜內(nèi)的受熱情況。對(duì)液面高度為95 mm 的水及1.0% CMC 罐頭進(jìn)行仿真求解及模型驗(yàn)證后，通過改變溫度函數(shù)以適應(yīng)不同的殺菌工藝。利用參數(shù)化掃描，將液面高度由75 mm 以5 mm 為間隔增加至100 mm，保溫溫度由85 ℃以5 ℃為間隔遞增至130 ℃進(jìn)行求解。

表2 模擬值與試驗(yàn)值的平均相對(duì)誤差Table 2 Average relative error between simulated and experimental values

3 結(jié)果與討論

3.1 模型驗(yàn)證

圖2比較了水和1.0% CMC 溶液在罐頭中心高度為15 mm（內(nèi)容物）及高度為100 mm 處（頂隙）的試驗(yàn)溫度變化及不同模型在相同位置的溫度模擬值，表2顯示了模擬值與試驗(yàn)值之間的平均相對(duì)誤差，其誤差均在5%以下。由于實(shí)際生產(chǎn)過程中，罐頭內(nèi)容物主要為黏性液體，故本研究主要對(duì)1.0% CMC 的模擬結(jié)果進(jìn)行分析。

由圖2及表2可看出，模型二對(duì)罐頭頂隙處的溫度變化尤其是冷卻過程描述得更為準(zhǔn)確，說明考慮頂隙的模型能更真實(shí)地反映熱殺菌過程中罐頭內(nèi)的傳熱情況，精準(zhǔn)地分析罐頭內(nèi)物理量的變化，準(zhǔn)確地預(yù)測(cè)罐內(nèi)冷點(diǎn)位置及能達(dá)到的殺菌程度。

圖2 實(shí)際測(cè)量及不同模型仿真得到的罐內(nèi)冷點(diǎn)及罐頭頂隙處溫度變化曲線Fig.2 Temperature variation curve of the cold spot and the head gap of the can obtained by different model simulations and experimental measurements

3.2 物理場(chǎng)分析

3.2.1 速度場(chǎng) 選取保溫溫度131 ℃，升溫時(shí)間10 min，保溫時(shí)間15 min，冷卻時(shí)間10 min 的理想殺菌曲線，對(duì)加熱條件下熱殺菌過程中的物理場(chǎng)進(jìn)行分析。圖3～4 為2 種模型在殺菌全過程中不同階段的速度場(chǎng)模擬結(jié)果，用處處都與速度矢量方向相切的流線及面上箭頭來表征速度大小，流線分布越密，箭頭越長(zhǎng)表明該處速度越快。

模型一僅考慮液體的非等溫流動(dòng)引起的速度場(chǎng)變化。在升溫階段，靠外壁處的液體升溫較快，密度降低而上浮，罐內(nèi)溫度相對(duì)較低的液體下沉。由于液體和罐頭底面的碰撞產(chǎn)生擠壓，在靠近底部中心位置出現(xiàn)速度向上的二級(jí)流場(chǎng)。冷卻過程的流場(chǎng)分布與升溫過程類似，二級(jí)流場(chǎng)出現(xiàn)在頂部靠近中心軸線處。保溫階段的液體流動(dòng)速度較低。

模型二模擬結(jié)果的液體部分流線分布情況與模型一的模擬結(jié)果類似。升溫階段，頂隙處的濕空氣與底部液體大致呈逆時(shí)針流動(dòng)，降溫階段與之相反，且存在二級(jí)流場(chǎng)。與底部液體的流動(dòng)場(chǎng)相比，頂隙處流線較密，箭頭較長(zhǎng)，可能是由于氣體密度低，動(dòng)力黏度小，流動(dòng)速度遠(yuǎn)大于底部液體。同時(shí)，受水分蒸發(fā)冷凝的影響，液面處流線較密，速度較大，且存在垂直液面的速度分量。

圖3 模型一模擬的速度場(chǎng)分布圖Fig.3 Velocity field distribution simulated by model 1

圖4 模型二模擬的速度場(chǎng)分布圖Fig.4 Velocity field distribution simulated by model 2

3.2.2 溫度場(chǎng) 液態(tài)罐頭食品的傳熱方式以對(duì)流傳熱為主，流體速度場(chǎng)會(huì)對(duì)溫度分布產(chǎn)生影響。圖5為2 種模型模擬生成的131 ℃-15 min 的殺菌過程中不同階段的溫度場(chǎng)分布云圖。受自然對(duì)流傳熱影響，罐頭頂部升溫較快，保溫后期SHZ 位于離底部約10%罐高，離罐壁10%～15%半徑處而非罐頭中心位置。對(duì)比模型一與模型二得到的溫度場(chǎng)分布情況，發(fā)現(xiàn)由于頂部濕空氣的比熱與密度較小，傳熱速度更快，使得模型二的頂隙部分升溫更快，并由此導(dǎo)致靠液體更快升溫。同時(shí)水分蒸發(fā)吸熱也會(huì)對(duì)液面附近的溫度分布造成影響，在水蒸氣未飽和時(shí)靠近液面處的空氣溫度較低。

圖5 不同模型模擬的溫度場(chǎng)分布圖Fig.5 Temperature distribution simulated by different model

3.2.3 致死率場(chǎng) 工業(yè)生產(chǎn)中通常用累計(jì)致死率評(píng)價(jià)特定殺菌條件下的總體殺菌效果，致死率按式（4）計(jì)算：

以肉毒桿菌為參考，式中，F(xiàn)——致死率，min；T——食品冷點(diǎn)的溫度，℃；Tref——?dú)⒕鷺?biāo)準(zhǔn)溫度，121℃；Z——微生物耐熱性特征值，10 ℃。

2 種模型溫度場(chǎng)的差別最終會(huì)對(duì)致死率場(chǎng)產(chǎn)生影響，圖6為131℃-15 min 的殺菌過程完成后，罐內(nèi)的微生物致死率場(chǎng)分布情況，表3為罐內(nèi)內(nèi)容物的最大、最小致死率值及最小致死率所在的坐標(biāo)。結(jié)果表明，2 種模型模擬得到的致死率最小值不同，且出現(xiàn)的位置有所差別。模型二模擬的致死率值總體較模型一高，分布更為均勻。隨著保溫溫度的升高及罐頭幾何尺寸的增大，不同模型對(duì)致死率分布的模擬差別更為顯著，說明建立模型時(shí)考慮頂隙與否會(huì)影響罐頭致死率場(chǎng)的分布。為更準(zhǔn)確地尋找罐內(nèi)冷點(diǎn)位置，更精準(zhǔn)地確定最佳保溫溫度及最短殺菌時(shí)間，需應(yīng)用含頂隙的模型。

3.2.4 壓力 為保證罐頭內(nèi)外壓力差小于罐頭可承受的臨界壓力差，避免出現(xiàn)凸罐或癟聽等罐頭變形甚至破損等問題，實(shí)際熱殺菌過程中需對(duì)壓力進(jìn)行控制，同時(shí)根據(jù)罐內(nèi)外壓力差選擇合適的罐壁厚度以確保罐頭的耐壓值大于殺菌過程的最大壓力差。由于模型一不考慮頂隙，無法模擬熱殺菌過程中罐內(nèi)的壓力變化，只能通過理想氣體方程，利用溫度變化推算罐內(nèi)的壓力變化情況。圖7為121 ℃-30 min 滅菌過程中，通過模型二模擬和通過罐內(nèi)溫度推算的罐內(nèi)相對(duì)壓力及殺菌釜壓力表讀數(shù)，結(jié)果表明由溫度計(jì)算罐內(nèi)的相對(duì)壓力不準(zhǔn)確，而模型二模擬得到的罐內(nèi)絕對(duì)壓力與實(shí)際情況較為接近。對(duì)該殺菌過程而言，罐內(nèi)外壓差最大值出現(xiàn)在剛開始進(jìn)行補(bǔ)壓時(shí)，差值為0.025 MPa，需保證罐頭的耐壓值大于該值。

圖6 不同模型模擬殺菌結(jié)束時(shí)微生物致死率場(chǎng)的分布圖Fig.6 Microbial lethality field distribution at the end of sterilization simulated by different models

工業(yè)生產(chǎn)中以純水的飽和蒸氣壓為參考對(duì)殺菌釜進(jìn)行壓力控制，然而由于分壓效果及水分活度的影響，完全用飽和蒸氣壓控制壓力并不夠準(zhǔn)確。利用模型二模擬得到不同保溫溫度及不同頂隙高度下保溫階段罐頭內(nèi)的絕對(duì)壓力及相應(yīng)溫度下的飽和水蒸氣壓。由表4可知，隨著溫度的升高，頂隙高度的減小，罐內(nèi)的相對(duì)壓力逐漸增大，這與理想氣體方程相吻合。頂隙過小，罐內(nèi)壓力大，容易產(chǎn)生凸角、假脹等現(xiàn)象，而頂隙過大則可能在冷卻過程中出現(xiàn)癟罐。生產(chǎn)廠家可通過模型模擬輔助設(shè)定不同頂隙高度下殺菌工藝中各階段的壓力，以避免出現(xiàn)脹罐、癟聽等罐頭外觀變形問題。

表3 兩種模型模擬的致死率最大值及最小致死點(diǎn)坐標(biāo)Table 3 Maximum and minimumt lethality and minimum lethal coordinates of the two models

圖7 罐內(nèi)相對(duì)壓力模擬值及殺菌釜相對(duì)壓力變化曲線Fig.7 Simulated relative pressure curve in cans and relative pressure curve of sterilizer

4 結(jié)論

1）利用COMSOL Multiphysics 軟件構(gòu)建的2種二維軸對(duì)稱模型，都能較為準(zhǔn)確地模擬水、1.0%CMC 溶液的熱殺菌過程。其中模型二對(duì)該過程模擬的溫度變化過程與實(shí)測(cè)值更為接近，是一種能精準(zhǔn)分析物理場(chǎng)變化的仿真模型。

表4 不同保溫溫度及不同頂隙高度下罐內(nèi)壓力（MPa）與相應(yīng)溫度下飽和水蒸氣壓Table 4 Pressure （MPa）in the cans at different holding temperatures and different headspace heights and saturated water vapour pressure at corresponding temperature

2）為更準(zhǔn)確地預(yù)測(cè)罐頭內(nèi)速度場(chǎng)、溫度場(chǎng)及微生物致死率場(chǎng)的分布情況，不可忽略頂隙的存在。本文建立了更準(zhǔn)確地分析殺菌過程中各物理場(chǎng)的變化的模型，適用于不同罐形大小及不同內(nèi)容物的液態(tài)罐頭食品，有助于在熱滲透測(cè)試進(jìn)行前準(zhǔn)確尋找冷點(diǎn)位置，最終實(shí)現(xiàn)殺菌工藝的優(yōu)化。

3）含頂隙模型可對(duì)罐內(nèi)壓力變化情況進(jìn)行分析，本文給出不同頂隙高度及不同保溫溫度下的罐內(nèi)相對(duì)壓力，頂隙高度越小、溫度越高，罐內(nèi)壓力越大。該模型可用于多種罐形及不同熱殺菌工藝液態(tài)罐頭食品物理場(chǎng)的預(yù)測(cè)，以便生產(chǎn)廠家對(duì)殺菌釜壓力變化情況及罐頭最大耐壓值進(jìn)行控制，避免出現(xiàn)脹罐或癟罐等質(zhì)量問題。