嬰幼兒奶粉壓縮模型參數標定及受力分析

劉少莉 蔣學峰 李歸浦 邢曉偉 熊麗娜 何光華，*

（1浙江科技學院生物與化學工程學院，浙江 杭州 310023；2貝因美（杭州）食品研究院有限公司，浙江 杭州 310053）

現市售嬰幼兒配方奶粉以粉末形式為主，存在易飛散、吸潮、粘勺以及沖調計量不準等弊端［1-5］，大大增加了育兒難度。若能將定質定量的粉狀奶粉通過壓縮成型技術制成統一形狀和重量的塊狀奶粉，不僅可以提高使用便攜性，降低育兒難度，還可以更精準喂養。制備塊狀奶粉時，由于奶粉顆粒會受壓力作用造成顆粒破碎，從而影響了塊狀奶粉抗壓強度及溶解性，因此，需要研究壓縮工藝條件對奶粉顆粒狀態的影響，而奶粉屬于不連續的離散介質，可以通過離散元仿真進行奶粉顆粒的受力分析。

目前離散元仿真涉及的領域范圍非常廣泛，包括食品、制藥、資源、巖土、化工等。Hanley 等［6］以嬰幼兒配方奶粉為對象，使用準靜態載荷試驗校準奶粉團聚體的離散元模型，并利用模型分析了粒子響應力學。Olaleye 等［7］以粘性奶粉為對象，模擬奶粉在管道中輸送，觀察其在管道中的沉積狀況。Afrassiabian 等［8］以乳糖粉為對象，基于離散元對樣品的機械強度進行數值模擬，用于描述機械約束下的樣品行為。離散元模擬仿真參數輸入包括本征參數與模型接觸參數。本征參數是基本特性參數，一般不受仿真場景影響；接觸參數是非常模型化的參數，易受仿真場景影響，而文獻［6-8］中的仿真場景并非壓縮仿真，因此只可采用文獻中奶粉本征參數作為離散元壓縮仿真的參考依據。而非奶粉壓縮模型應用研究表明，不同的仿真對象所采用的離散元接觸模型各不相同［9-12］。因此，需要根據實際需求和奶粉特性來選擇奶粉壓縮接觸模型。

本研究以市售嬰幼兒配方奶粉為研究對象，依次采用Plackett-Burman、最陡爬坡和Box-Behnken 3 個試驗設計方法對奶粉接觸參數進行標定，將標定參數用于離散元壓縮仿真，通過觀察塊狀奶粉的壓縮量、表面力鏈、斷面力鏈，分析塊狀奶粉表層顆粒和內部顆粒受力情況，旨在為塊狀奶粉的制備工藝參數選擇提供參考。

1 材料與方法

1.1 試驗材料

嬰幼兒配方奶粉：市售不同品牌嬰兒配方奶粉，標示為奶粉1～27（以奶粉1 為離散元模型參數標定及分析對象，由于粉體的復雜性使得測定的接觸參數與實際靜態參數存在差異，無法直接用于模擬，因此奶粉2～27 主要目的是用于確定靜摩擦系數的合理范圍，以便水平參數設計以及后續標定）。

1.2 儀器與設備

TA-XT plusC Stable Micro Systems 質構儀，英國Stable Micro Systems 公司；JA1003 電子天平，上海良單儀器儀表有限公司；CMT5000 微機控制電子萬能測試儀，深圳市新三思材料檢測有限公司；FT4 粉體流變儀，英國Freemantechnology 公司；SU1510 電子掃描顯微鏡，日本HITACHI 公司；BT-1001 智能粉體特性測試儀，中國丹東百特儀器有限公司。

1.3 試驗方法

1.3.1 塊狀奶粉及奶粉檢測

1.3.1.1 顆粒間摩擦系數、顆粒與壁面間摩擦系數及休止角的測定 采用文獻［13］的方法測定。用θ1表征內摩擦角，tanθ1表征顆粒間摩擦系數［7］；θ2表征壁面摩擦角，tanθ2表征顆粒與壁面間摩擦系數［7］；奶粉1休止角檢測共測量9次。

1.3.1.2 塊狀奶粉壓縮制備 將奶粉1 分為5 組，每組3 份，每份4.5 g，通過萬能測試儀將每份奶粉在加壓速率為100 N·s-1，壓力分別為1 000、1 500、2 000、2 500、3 000 N，保壓時間為2 s 的條件下進行壓制，獲取塊狀奶粉，然后進行抗壓強度檢測。

1.3.1.3 抗壓強度檢測 抗壓強度采用質構儀檢測，具體操作如下：1）設置質構分析（texture profile analysis，TPA）測試條件，采用A/MORS 刀片探頭，探頭測試前速度2 mm·s-1，測試中速度0.1 mm·s-1，測試后速度2 mm·s-1，加載距離3 mm，停留時間5 s，測試壓力5 g［14-16］。2）在時間-壓力曲線中，以頂峰值作為塊狀奶粉抗壓強度［17-19］。3）取2塊塊狀奶粉，每塊檢測3個不同部位的抗壓強度，取平均值。

1.3.1.4 塊狀奶粉表層電鏡掃描 采用電子掃描顯微鏡分別對2 000、2 500 N 壓力下制取的塊狀奶粉表層拍攝200倍的電鏡掃描圖片。

1.3.2 離散元接觸模型確定及顆粒模型建立

1.3.2.1 離散元接觸模型選擇與分析 （1）離散元接觸模型選擇。根據粉末能否壓縮和是否具有黏性兩方面來選擇接觸模型［20］。本研究采用的原料奶粉具有可壓縮和一定黏性的特征，所以奶粉離散元壓縮仿真采用非線性粘結彈性塑形（edinburgh elasto-plastic adhesion，EEPA）模型。

（2）EEPA 模型分析。在奶粉壓縮過程中，主要考慮奶粉顆粒在法向上的受力情況。如圖1 所示，EEPA的法向力與位移的關系主要受6 個因素綜合影響，分別是恒定脫拉力F0、加載分支剛度k1、加載-卸載分支剛度k2、最小作用力Fmin、粘性分支剛度ka和塑性重疊量δp。

圖1 EEPA模型法向力與位移關系Fig.1 Normal force and displacement relationship of the EEPA model

法向的總力（Fn）為：

Fn=FEEPAn+Fdn

式中，FEEPAn代表EEPA 力的法向分量；Fdn代表粘性阻尼分量。而FEEPAn由上述6個因素定義：

根據EEPA 模型分析可知，在奶粉離散元壓縮仿真前，需向仿真軟件輸入表面能、脫拉力、塑性變形比、拉伸指數、切向剛度因子5 個參數，對這5 個參數合理選擇，可以獲取塊狀奶粉在壓縮過程中法向上受力分布狀態，為塊狀奶粉法向上的受力分析提供依據。

1.3.2.2 顆粒模型建立 將奶粉1在200倍下進行電鏡掃描，觀察形態，利用Image-j 軟件對框中的四棱臺形態顆粒進行特征尺寸標量（圖2-A），獲取尺寸比。正面尺寸比值：上底：下底：高≈1∶3∶1.5；側面尺寸比值：上底：下底：高≈2∶1∶1.5。

圖2 顆粒真實形態與虛擬形態Fig.2 Real form and virtual form of particles

根據特征尺寸比值，采用八個等半徑球體內切于四棱臺的方式進行填充［21］，填充結構示意圖如圖3所示。根據填充結構幾何關系、坐標系對稱性以及顆粒縮放原則［22-23］，確定8 個填充球體的球心坐標（如下），構建出體積與半徑1 mm 球體體積相近的顆粒模型（圖2-B）。

O1：(0.656 2，0.161 9，0.58)；O2：(-0.656 2，0.161 9，0.58)；

O3：(0.656 2，-0.161 9，0.58)；O4：(-0.656 2，-0.161 9，0.58)；

O5：(0.269 5，0.355 2，1.16)；O6：(-0.269 5，0.355 2，1.16)；

O7：(0.269 5，-0.355 2，1.16)；O8：(-0.269 5，-0.355 2，1.16)；

1.3.3 休止角模擬測定 構建BT-1001 漏斗三維模型，其構造尺寸如圖4-A 所示。在EDEM 軟件中導入圖2-B 的顆粒模型，顆粒生產方式為Dynamic，顆粒工廠生產時間為7 s，生成總量為70 g，待顆粒全部落到平臺后，采用Matlab軟件對仿真圖像進行采集處理，之后剔除干擾數據，再線性擬合獲得一次項系數λ，arctanλ即為休止角角度。圖像從+x軸方向和-y軸方向采集，測量粉堆兩端角度，共測4 次角度，取平均值代表休止角。

圖4 三維模型Fig.4 Three-dimensional model

1.3.4 接觸參數標定 采用Plackett-Burman 試驗從多個接觸參數中篩選出對休止角模擬影響最大的參數，之后通過最陡爬坡試驗確認休止角模擬最佳參數區間，最后通過Box-Behnken 試驗確定休止角模擬的最佳接觸參數組合［23-25］。

1.3.4.1 Plackett-Burman 試驗 根據接觸參數范圍確定各因素的兩水平數值（-1 和+1），將數值輸入到Design-Expert 軟件的Plackett-Burman 模塊中進行試驗設計并對試驗結果進行方差分析，篩選出影響顯著的參數。

1.3.4.2 最陡爬坡試驗 將Plackett-Burman 試驗獲取的顯著參數按照合理爬坡步長確定參數組合，根據參數組合進行離散元模擬仿真得出模擬休止角角度，再與實際休止角比較計算出相對誤差，結合相對誤差確定顯著參數最佳響應區域。

1.3.4.3 Box-Behnken 試驗 根據顯著參數最佳響應區域設定各參數的三水平數值（-1、0、+1），將數值輸入到Design-Expert 軟件的Box-Behnken 模塊中分別進行試驗設計、對試驗結果進行方差分析并得出回歸模型，最后獲取最佳參數組合。

1.3.5 壓縮模擬測定 構建內腔空間為23 mm×31 mm×40 mm的模具幾何模型，如圖4-B所示。在EDEM軟件中輸入模擬參數，定義上沖模運動方向為Z 軸正方向，加壓速率100 N·s-1，加壓時間25 s，離散元仿真完畢后，再通過EDEM 軟件的后處理平臺分析壓縮量、力鏈。力鏈是通過線段顏色表示顆粒間或顆粒與壁面間的接觸狀態或接觸力，顏色越暖表示接觸力越大［26］。接觸力是切向力和法向力的合力，由于奶粉壓縮主要考慮法向的受力情況，因此力鏈分析主要分析接觸力的法向力。

2 結果與分析

2.1 嬰兒配方奶粉顆粒間摩擦系數、顆粒與不銹鋼間摩擦系數確定

摩擦系數是仿真過程中的重要參數之一，合理范圍的確定才能保證模擬不會失真。電子附表1 中顆粒與不銹鋼間的靜摩擦系數和顆粒間靜摩擦系數測定值范圍分別為0.07～0.21 和0.39～1.00，其范圍跨度較大，且用于標定的奶粉1實測值也在該范圍內，表明參數范圍確定較合理，為仿真初期提供前提保證。其余參數是根據文獻［6-7，27-28］確定（電子附表2），其可靠性需進一步仿真驗證。

2.2 Plackett-Burman試驗分析

Plackett-Burman試驗是通過比較各個因素之間差異以及整體差異來確定各因素顯著性程度的方法。根據電子附表2 參數范圍確定12 個因素水平（表1），Design Expert 軟件共生成20 組試驗，試驗設計方案與結果如表2 所示。結果顯示，仿真休止角角度范圍為37.43°～55.43°，其范圍值跨度較大，表明該12 個因素中具有對休止角影響顯著的因素，并且實際休止角角度47.92°（電子附表2）在該范圍內，說明各個因素的參數范圍選擇較為合理，能夠保證后續仿真參數的正常標定。

表1 因素水平表Table 1 Factor level table

表2 Plackett-Burman試驗設計及結果Table 2 Plackett-Burman test design and results

利用Design-Expert 軟件對試驗結果進行方差分析，結果如表3 所示。A2、A3和A5的P值小于0.01，表明這3 個因素對顆粒休止角的影響極顯著，可作為最陡爬坡試驗設計對象。A10的P值小于0.05，表明該因素對顆粒休止角的影響顯著。剩余因素的P值大于0.05，表示剩余因素對休止角的影響不顯著。

表3 Plackett-Burman試驗參數方差分析Table 3 Plackett-Burman test parameters analysis of variance

2.3 最陡爬坡試驗分析

最陡爬坡試驗是依據合理的步長增加試驗密集度來確定顯著因素最佳響應區域。將A2、A3、A5的爬坡步長分別設為0.1、0.05、0.03（表4），其余仿真參數具體數值設定參考電子附表2。由表4可知，隨著3個參數增大，相對誤差先減小后增大，其中2號參數組合相對誤差最小，所以最優參數組合范圍在2號組合參數附近。

表4 最陡爬坡試驗設計及結果Table 4 Steepest climbing test design and results

2.4 顯著接觸參數的Box-Behnken試驗分析

由于2 號參數組合相對誤差最小，所以最佳參數組合在1號參數組合和3號參數組合的范圍內，且靠近2號參數組合，因此在Box-Behnken試驗設計中將1號、2 號和3 號參數組合分別作為低水平（-1）、中間水平（0）和高水平（+1）。Box-Behnken 試驗設計及結果如表5所示。

表5 Box-Behnken試驗設計及結果Table 5 Box-Behnken test design and results

利用Design-Expert 軟件對試驗結果進行方差分析（表6），并對參數因素與休止角之間進行全項二次多項式回歸擬合，建立二階回歸模型1如下：

表6 試驗結果方差分析Table 6 ANOVA of test results

方差分析結果表明，模型1 的P值為0.036<0.05，失擬項PL值大于0.05，決定系數R2=0.843 5，表明回歸模型擬合較好。由表6 可知，A2、A3、A32、A52都具有顯著性，其余沒有顯著性，所以將沒有顯著性的A2A3、A2A5、A3A5、A22從回歸模型中剔除，獲得新的回歸模型2：

模型2 的P值為0.001 7<0.05，失擬項PL值大于0.05，決定系數R2=0.792 4，模型尚佳。模型α2的響應范圍為42.507 5°～52.157 5°，真實休止角測定值（47.92°）被包含在內，該模型可以適用于奶粉參數標定。

2.5 仿真最優參數組合的驗證及確定

以奶粉1 的休止角47.92°為標準，利用Design Expert 軟件中參數優化模塊進行參數尋優，從中獲取最優參數組合，A2、A3、A5分別為0.359、0.099、0.125，其余非顯著參數取值如電子附表2 所示。將上述參數輸入到EDEM 中進行休止角仿真，測定后的仿真值平均值為47.94°（電子附表4），與真實值的相對誤差為0.04%，因此上述標定的參數比較可靠，可用于奶粉壓縮仿真分析。

2.6 奶粉壓縮量與力鏈分析

采用2.5中相同的參數設置，之后按1.3.5步驟進行壓縮仿真，壓縮仿真獲取的壓縮曲線如圖5 所示。壓縮仿真曲線與真實曲線的變化趨勢基本一致，表明EEPA 接觸模型能夠模擬奶粉壓縮過程中的壓縮狀態。在壓縮量上，仿真曲線與真實曲線存在一定差異，可能是某參數未調節適當所致。

圖5 壓縮仿真曲線與真實曲線對比Fig.5 Comparison of the simulation curves with the real curves

由圖6 可知，隨著壓力增加，表面顆粒與上沖模接觸的大力鏈數量逐漸增多。在20～25 s時，大力鏈數量明顯增多，幾乎將塊狀奶粉表面占滿，表明隨壓力增加，塊狀奶粉表面分載壓力的顆粒數量逐漸增多。時間為25 s時，表面顆粒幾乎都在承受巨大壓力，無多余的顆粒繼續分載壓力，由此可能導致表層顆粒破碎。

由圖7 可知，隨壓力增加，內部顆粒間接觸的大力鏈數量逐漸增多，但增加數量較少，這是因為內部顆粒數量較多，分載壓力效果較好，所以力鏈顏色主要分布在藍色和綠色。壓縮時間為20～25 s時，綠色力鏈數量較明顯增多，還有少部分綠色力鏈轉化為紅色力鏈，顆粒破碎具有從表層開始逐漸向內部發展的趨勢。上述結果表明，壓縮時間為25 s（壓力為2 500 N）時，塊狀奶粉表層顆粒產生破碎的可能性較高，而內部顆粒雖未產生破碎，但有破碎的趨勢。

圖7 不同壓縮時間斷面接觸力鏈的分布Fig.7 Distribution of contact force chains between particles at different suppression time sections

2.7 奶粉壓縮試驗驗證

質構儀與電子掃描顯微鏡對塊狀奶粉的檢測結果如圖8和圖9 所示。由圖8可知，塊狀奶粉抗壓強度隨壓力增加而增加，當壓力大于2 000 N后，塊狀奶粉抗壓強度增幅加大。可能是壓力達到2 500 N時，顆粒達到承受極限，開始產生破碎，細碎顆粒致使顆粒間作用力劇烈增加所致。由圖9 可知，壓力為2 000 N 時，表層顆粒較完整，無明顯碎片化顆粒，而壓力為2 500 N時，部分表層顆粒受壓出現破裂，產生碎片化顆粒，表明塊狀奶粉表層顆粒在壓力為2 500 N 時無法分載壓力而產生破碎，而在仿真模型受力分析中同樣推斷出壓力在2 500 N時顆粒產生破碎的可能性高。因此認為，離散元壓縮模型可以對奶粉壓縮過程中的受力狀態進行分析。

圖8 壓力對塊狀奶粉抗壓強度的影響Fig.8 Effect of pressure on the compressive strength of bulk milk powder

圖9 塊狀奶粉表層電鏡掃描圖Fig.9 SEM image of the surface layer of bulk milk powder

3 討論

通常認為在離散元仿真中，顆粒大小是影響DEM模擬效果和效率的關鍵因素之一［23，29］，而在實際應用中，常需要通過顆粒大小或密度的縮放來降低模擬時間，提高仿真效率，從而降低仿真成本［30-32］。Hanley等［6］研究認為嬰幼兒奶粉在單軸壓縮仿真過程中，奶粉的顆粒密度縮放不會影響法向力與應變之間的曲線變化。Thakur 等［33］研究發現，在DEM 壓縮仿真中，顆粒尺寸放大與未放大的離散元模型所產生的機械加載響應相一致。Coetzee 等［34］研究了改變顆粒大小對離散單元法校準和驗證的影響，得出玉米顆粒縮放比例為4時，可以通過休止角仿真試驗準確標定顆粒-顆粒摩擦系數。上述結果表明，顆粒的合理縮放不會使離散元仿真結果偏離實際。本試驗將奶粉顆粒模型放大為與半徑1 mm 球體體積相近的顆粒，經過參數標定后，休止角仿真結果為47.94°，與實際值47.92°的相對誤差為0.04%，表明顆粒縮放在參數標定后可以保證仿真結果的精準性。

離散元仿真接觸參數的選擇和設置非常重要。Hanley 等［6］采用擴散限制凝聚模型模擬奶粉單軸壓縮的力學行為時，將奶粉顆粒間的摩擦系數設置為0.65。Olaleye 等［7］以Hysteretic spring with JKR 模型作為接觸模型模擬管道運輸奶粉時，將顆粒間的摩擦系數設置為0.33，顆粒與壁面間的摩擦系數設為0.43。上述研究表明，同樣以奶粉為模擬對象，但采用不同的接觸模型以及不同模擬應用場景，都會使仿真參數數值設定不一，為了保證仿真結果的準確性，需進行參數標定試驗確定合理的接觸參數。本研究采用3 種試驗設計標定參數，確定奶粉顆粒間的靜摩擦系數、奶粉顆粒間的滾動摩擦系數和奶粉顆粒-不銹鋼靜摩擦系數分別為0.359、0.099 和0.125，通過塊狀奶粉抗壓強度檢測與電鏡掃描結果驗證了以標定參數為基礎建立的離散元壓縮模型可以模擬奶粉壓縮時的受力狀態。

在仿真結果中，奶粉壓縮仿真曲線與真實曲線雖然變化趨勢相同，但貼合度有待提高，可能是某參數未調節適當所致，因此，在后續仿真研究中，需添加新的仿真標定指標（可壓縮度）進行二次參數標定完善其他參數，從而進一步增加離散元壓縮仿真模型的準確性。

奶粉組成成分復雜多樣，壓縮成型后相較原料奶粉在溶解性、保存性、油脂氧化性等方面上可能發生變化。Pugliese 等［35］和Vishwakarma 等［36］采用相同的不溶指數（Isi）方法分別測定脫脂奶粉與脫脂奶粉片的溶解性，結果顯示脫脂奶粉的溶解性最佳，為92%，而1 000 N 壓力下制取的脫脂奶粉片溶解性在91.2%～93.6%之間，二者溶解性相近。柴田滿穗［1］研究發現，在保存3個月后，塊狀奶粉游離脂肪、風味、溶解性、脫脂、硬度等檢測值與制備時的初期值無差異，而過氧化值在乳粉壓縮前后仍處在相同程度，表明塊狀奶粉能夠長期保存，且壓縮對奶粉的油脂氧化性不產生影響。根據文獻［1，35-36］可知，奶粉的合理壓縮不會影響其溶解性、油脂氧化性及貯藏穩定性，但本研究壓縮條件下制備的塊狀奶粉的貯藏穩定性有待進一步深入研究。

4 結論

本研究確定了最優仿真接觸參數，分別為奶粉顆粒間的靜摩擦系數0.359、奶粉顆粒間的滾動摩擦系數0.099、奶粉顆粒-不銹鋼靜摩擦系數0.125，并經休止角試驗驗證相對誤差為0.04%，表明接觸參數標定具有較高的可靠性。以接觸參數為基礎，結合EEPA接觸模型構建的奶粉離散元壓縮模型可以模擬奶粉的壓縮狀態以及顆粒受力分布情況，并分析得出塊狀奶粉在受壓達到2 500N 時，表層顆粒破碎可能性高，且具有由表及里的破碎趨勢。