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

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

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

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

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

1 材料與方法

1.1 材料試劑

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

1.2 儀器設備

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

1.3 試驗方法

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

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

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

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

計算試驗所得實際溫度與數(shù)值模型的模擬溫度的平均相對誤差，以評價模型的準確性。

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

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

2.1 模型假設

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

2.2 控制方程

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

式中，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 計算條件

初始溫度為室溫，初始絕對壓力（-40 kPa）作為罐內(nèi)的真空度。純水的熱物理特性參數(shù)由COMSOL Multiphysics 內(nèi)置材料庫給出，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 種模型均利用細化的物理場控制的自由三角形網(wǎng)格進行剖分，模型一與模型二分別將幾何體剖分為5 429 及5 971 個單元進行求解。

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

表2 模擬值與試驗值的平均相對誤差Table 2 Average relative error between simulated and experimental values

3 結(jié)果與討論

3.1 模型驗證

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

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

圖2 實際測量及不同模型仿真得到的罐內(nèi)冷點及罐頭頂隙處溫度變化曲線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 物理場分析

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

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

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

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

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

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

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

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

以肉毒桿菌為參考，式中，F(xiàn)——致死率，min；T——食品冷點的溫度，℃；Tref——殺菌標準溫度，121℃；Z——微生物耐熱性特征值，10 ℃。

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

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

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

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

表3 兩種模型模擬的致死率最大值及最小致死點坐標Table 3 Maximum and minimumt lethality and minimum lethal coordinates of the two models

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

4 結(jié)論

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

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

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

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