昨天在最後一節丟出一個沒解開的問題:δ=√(D/k) 是反應-擴散方程的穩態解,只有在 t≫1/k 時才成立;冷藏醃漬溫度下 1/k≈31 天,遠比任何實際泡鹵水或撒鹽的時間長——但那只是定性判斷,沒有回答「還沒到穩態之前,濃度分布到底長什麼樣子、跟穩態解差多少」。今天把這個定性懷疑變成一個可以解、也可以驗算的方程:讓 ∂C/∂t = D∂²C/∂x² − kC 在半無限介質、表面濃度固定的邊界條件下求出完整的時間相關解,再把它拆開成「純擴散加一個正比於 kt 的小修正」,直接代入昨天算出的 kt = 0.0324(24 小時、4.4°C、pH 6.2),量出這個修正究竟有多小。結論是:在冷藏醃漬的時間尺度上,反應項讓濃度分布偏離純擴散的幅度,連最受影響的深度都不到 1.5 個百分點——δ 不是算錯了,是問錯了時間點,而今天順便算出一個把 δ、D、k 三者關係講清楚的等式,以及把同一套算法搬去問加熱階段會發生什麼。

📖 學

一、用 Laplace 變換解出反應-擴散方程的完整時間解

方程是 ∂C/∂t = D∂²C/∂x² − kC,x>0(半無限介質,x=0 是肉表面),初始條件 C(x,0)=0,邊界條件 C(0,t)=Cₛ(表面濃度固定,對應鹵水或鹽粒持續供應)、C(∞,t)=0(有限)。對 t 做 Laplace 變換,設 c̄(x,s)=∫₀^∞ C(x,t)e^{-st}dt,方程變成常微分方程:

有界解只留衰減項 c̄(x,s) = A(s)·exp(−x√((s+k)/D)),代入邊界條件 c̄(0,s)=Cₛ/s 得 A(s)=Cₛ/s。剩下的是查一條標準 Laplace 反變換表(型如 (1/s)·exp(−a√(s+b)),a、b 為常數的反變換,出現在處理「擴散伴隨一階不可逆反應」這類問題的數學物理文獻裡,Crank《The Mathematics of Diffusion》與 Danckwerts 處理「氣體吸收伴隨化學反應」的論文都是處理同一類方程;今天沒能拿到 Crank 書裡這一段的確切頁碼/公式編號——網路能查到的版本要嘛要權限要嘛抓不到文字層,只查到書確實收錄「diffusion accompanied by immobilization/irreversible reaction」這類問題,已標存疑),代回去整理後得到:

其中 δ=√(D/k)、ξ=x/(2√(Dt))。因為拿不到可靠的頁碼引用,今天没有停在「查到就信」,而是自己把這個解的三個極限算過一遍,確認它真的滿足原本的方程與邊界條件:

  • t→0:ξ→∞(分母 √(Dt)→0),兩個 erfc 的自變數都趨近 +∞,erfc(+∞)=0,所以 C→0,符合初始條件。
  • x=0(表面):ξ=0,式子變成 (1/2)erfc(−√(kt)) + (1/2)erfc(√(kt))。用恆等式 erfc(−y)=2−erfc(y) 代入,兩項相加後 kt 項完全消掉,剩下 (1/2)×2=1,也就是 C(0,t)=Cₛ 對所有 t 都成立——邊界條件在任何時間點都被滿足,不是只在穩態才對。
  • t→∞(穩態):ξ→0、√(kt)→∞,第一項 erfc(−∞)=2 給出 e^{-x/δ},第二項 erfc(+∞)=0 直接消失,剩下 C/Cₛ=e^{-x/δ}——正是昨天用的那個穩態解。三個極限都對得上,說明這個全解至少在邊界上是自洽的。

二、把全解攤開成早期修正:一階項剛好抵銷,只剩正比於 kt 的項

穩態解在 t≫1/k 才成立,那 t≪1/k(kt 很小)的早期階段該怎麼看?把上面的全解對小量 ε=√(kt) 做泰勒展開(同時 x√(k/D)=2εξ 也是同量級小量),逐項展開到 ε² 再相加後發現:ε¹ 的項完全互相抵銷——這代表沒有正比於 √(kt) 的修正,最低階的修正項是 ε²=kt 量級,即:

erfc(ξ) 就是沒有反應、純擴散的解(2026-09-22 那篇算鹽擴散時用的同一條 erfc 半無限解)。g(ξ) 在 ξ=0(表面)恆為 0——這點很合理,邊界條件鎖死表面濃度等於 Cₛ,反應再怎麼消耗,表面補給也跟得上;g(ξ) 在 ξ→∞ 時也趨近 0,因為那裡本來就沒什麼濃度可以被消耗。手算幾個點(用 erf(0.4)=0.4284、erf(0.45)=0.4755、erf(0.5)=0.5205 等標準值代入)可以找到 g(ξ) 的極值落在 ξ≈0.45 附近,數值約 −0.202——是負的,代表反應確實把濃度往下拉,而不是往上推,這點自己先檢查過(之前手算一次符號搞反過,重新算後才確認負號是對的,方向與物理直覺一致)。

三、代入昨天的真實數字:冷藏醃漬期的修正不到 1.5%

昨天用 Nordin(1969)公式算出 4.4°C、pH 6.2 下 24 小時的 kt = 0.0324。把 g(ξ) 的最大值 −0.202 代進去:

也就是說,即使在反應影響「最集中」的那個深度,C/Cₛ 也只比純擴散的 erfc(ξ) 少大約 0.0065(絕對值)。那個深度對應 ξ≈0.45 時 erfc(0.45)≈0.525,0.0065/0.525≈1.2%——相對修正只有約 1.2%。換算成實際深度:24 小時、D=1.28×10⁻¹² m²/s 下 √(Dt)≈0.333 mm,ξ=0.45 對應 x=2×0.45×0.333≈0.30 mm。這比昨天用「85% 在 24 小時內消耗掉」得到的誤判好太多——真實答案是反應項在冷藏醃漬的時間尺度上,連局部最大偏離都壓在 1.5% 以內,pH 5.6(kt=0.0526)算出來也只在 2% 出頭,同樣是小修正。純擴散 erfc 解在這個階段幾乎是免費的好近似,這也直接回答了昨天 問AI 裡的第一個問題:冷藏醃漬期該用的不是穩態反應-擴散解,也不是某種複雜的過渡態解,而是「純擴散加一個可以忽略的小修正」,而且今天算出了這個修正具體有多小。

四、δ 的另一層意義:穩態深度=純擴散在一個反應壽命裡走的距離

δ=√(D/k) 也可以寫成 δ=√(D×(1/k))——這在代數上是同一件事,但物理意義不一樣:穩態的滲透深度,正好等於「純擴散單獨作用一個反應壽命(1/k)的時間」所能走到的距離。拿昨天的數字驗算:1/k(pH 6.2, 4.4°C)=2,666,667 秒≈30.9 天,D=1.28×10⁻¹² m²/s,√(D×2,666,667)=√(3.41×10⁻⁶)≈1.85×10⁻³ m=1.85 mm——跟昨天直接用 δ=√(D/k) 算出的 1.85 mm 完全對上;pH 5.6 那組(1/k≈19.0 天)算出 1.45 mm,也跟昨天的 1.45 mm 吻合。這不是巧合,是 δ 這個式子本身的重新表述,但换个角度看很有用:δ 告訉你的其實是「如果讓純擴散跑滿一個反應壽命會走多遠」,而 1/k 遠大於實際操作時間,就代表這段距離對應的是一個「還沒真的走到」的假設情境——難怪在冷藏醃漬這個遠短於 1/k 的時間窗裡,反應修正小到可以忽略。

五、把同一套檢查表拿去問加熱階段:用 Nordin 公式外插到 55–72°C

昨天引用的 Barbieri、Bergamaschi 與 Franceschini(2013,Meat Science 93(2):282–286)在 mortadella 上實測了 55、60、65、70、72°C 五個溫度、40/70/100/150 mg/kg 四個亞硝酸鹽添加量下的動力學參數,結論是 65°C 配 70 mg/kg 效果最好——但今天想拿到他們論文報告的確切 k 值或半衰期數字時,ScienceDirect 全文頁與 PubMed 摘要頁都被權限擋下,只查到方法與結論層級的描述,沒能取得可引用的精確數字,這部分已標存疑。

退而求其次,用昨天已經驗證過(在 4.4°C 與 104.4°C 兩個定錨點、76.7°C 一個交叉驗證點上都對得上)的 Nordin(1969)公式做一次外插估算,純粹當作量級對照,不冒充是 Barbieri 論文本身的數字:

代入 T=65°C、pH 6.2:0.65−1.625+2.17=1.195,t½=10^1.195≈15.7 小時,k=ln2/15.7h≈1.23×10⁻⁵/s,1/k≈22.6 小時≈0.94 天。pH 5.6 時 t½≈9.7 小時,k≈1.99×10⁻⁵/s,1/k≈13.9 小時≈0.58 天。跟冷藏的 1/k≈31 天(pH 6.2)、19 天(pH 5.6)相比,65°C 下 1/k 縮短了約 32.7 倍——這個倍數剛好等於 10^(0.025×ΔT)=10^(0.025×60.6)=10^1.515≈32.7,跟公式裡溫度項的結構完全自洽,不管 pH 取哪個值,同樣的溫差都給出同樣的縮短倍數。

這代表如果 mortadella 或培根的煙燻/蒸煮階段在 65°C 附近停留的時間是幾小時到一天量級(mortadella 因體積大,中心升溫到位常常要花數小時以上),kt 就可能落在 0.1 到 1 甚至更高的範圍——不再像冷藏期那樣穩穩地小,而是真正進入「穩態框架開始有意義」的區間。這跟 Barbieri 選擇在這個溫度區間做動力學實測是吻合的:那才是 δ=√(D/k) 這個框架原本該被拿來回答問題的地方,而不是冷藏醃漬期。(再次強調:上面 65°C 這組 k、1/k 數字是今天自己用 Nordin 的火腿配方外插算出來的估計值,不是 Barbieri 論文裡針對 mortadella 實測的數字,兩者未經交叉核對,已標存疑。)

六、回頭看一個方法論問題:D 高估 86 倍、k 高估 58.5 倍「幾乎抵銷」站不站得住腳

昨天的 問AI 留了這個問題,今天順手回應:Gómez et al.(2015)的 NO₂⁻ 組織擴散係數只在 2、7、12°C 三個溫度點量測,再用 Arrhenius 關係內插到 4.4°C——三個點擬合一條指數曲線,外插一小段通常不會差太遠,但論文本身的實驗不確定度(重複次數、樣品間變異)並沒有被今天這些計算納入;Nordin(1969)的公式是用當年的手動迴歸方法,把不同溫度、pH 下量到的半衰期做二次擬合,論文裡能拿出兩個定錨點與一個交叉驗證點對得上,但這終究是回歸出來的平滑曲線,不是每個溫度點都有獨立重複測量的誤差棒。把這兩個各自帶著自己不確定度的數字,拿去算出「δ=1.45–1.85 mm」這種三位有效數字的答案,某種程度上是把兩層不確定度都壓縮消失了——今天這篇算出的「修正 1.2%」、「1/k 縮短 32.7 倍」這些數字也適用同樣的提醒:它們是在假設 Nordin 公式與 Gómez 的係數本身精確的前提下推導出的,底層文獻的量測不確定度沒有被這些計算吸收進去,量級判斷(「小」、「可忽略」、「開始有意義」)比具體到小數點後兩位的數字更值得信。

🧠 記

  • 今天把昨天「δ 只在 t≫1/k 才成立」的定性判斷,換成一個可以驗算的完整時間解:,並自己驗證了 t→0、x=0、t→∞ 三個極限都對得上(t→∞ 正確化簡回昨天用的 δ=√(D/k) 穩態解)。
  • 對小 kt 展開後,一階(正比於 √(kt))的修正項會完全抵銷,最低階修正是正比於 kt 的二階項,且在表面(x=0)恆為零——邊界條件本來就鎖死表面濃度,不受反應影響。
  • 代入昨天算出的 kt=0.0324(24 小時、4.4°C、pH 6.2),得到即使在修正最大的深度(約 0.30 mm),相對偏離純擴散也只有約 1.2%——冷藏醃漬期用純擴散(erfc)近似,誤差小到可以忽略,今天把這個「小」量化了。
  • δ=√(D/k) 可以重寫成 δ=√(D×(1/k)):穩態滲透深度,等於純擴散單獨作用一個反應壽命(1/k)所能走到的距離——這解釋了為什麼遠短於 1/k 的操作時間裡,反應修正必然很小。
  • 用已驗證過的 Nordin 公式外插到 65°C(Barbieri et al. 2013 研究的溫度區間之一),估算 1/k 從冷藏的 30–31 天縮短到約 0.94 天(pH 6.2),縮短倍數 32.7 倍與公式的溫度項結構自洽——但這是外插估計值,不是 Barbieri 論文報告的原始數字,今天沒能取得 Barbieri 論文全文中確切的 k/半衰期數字,已標存疑。
  • 方法論提醒:今天與昨天算出的所有「幾個百分點」「幾十倍」這類具體數字,都建立在 Nordin(1969)手動迴歸與 Gómez et al.(2015)僅三個溫度點量測這兩份文獻本身的精度之上——量級判斷比小數點後的精確位數更可靠。

✍️ 實踐

想知道自己在家的醃漬操作是不是穩穩落在「可以放心當純擴散」的區間,可以自己動手算一次 kt:①查自己冰箱實際溫度 T(°C)與肉的大概 pH(生鮮豬肉/雞肉常見 5.6–6.2,沒有 pH 計就用 6.0 當粗估);②代入 log₁₀(t½[hr]) = 0.65 − 0.025T + 0.35×pH 算出半衰期,再用 k=ln2/t½ 換算成與自己醃漬時數相同單位的速率常數;③把 k 乘上預計醃漬的總時數,得到 kt。如果 kt < 0.05,今天這篇的結論適用,反應項可以放心忽略,穿透行為就是單純擴散;如果 kt 落在 0.1 到 1 之間(比較可能發生在拉長到一兩週的乾式醃製,或任何牽涉到加熱到 50°C 以上的階段),就已經進入穩態框架真正有意義的範圍,不能再簡單套用純擴散的直覺。

🔗 延伸學習

💬 問 AI

我在看一篇食物科學筆記,前幾天算出亞硝酸鹽在肌肉組織裡的反應-擴散穿透深度
δ=√(D/k)是穩態解,只在時間t遠大於1/k時成立,而冷藏醃漬溫度下1/k≈31天,
遠長於實際醃漬時間。今天用Laplace變換解出了完整的時間相關解:
 
C(x,t)/Cs = (1/2)e^(-x/δ)erfc(ξ-√(kt)) + (1/2)e^(x/δ)erfc(ξ+√(kt))
 
其中δ=√(D/k)、ξ=x/(2√(Dt)),並驗證了t→0、x=0(邊界條件)、t→∞(化簡回穩態
解)三個極限都成立。對小kt展開後發現一階(正比於√(kt))的修正項完全抵銷,
最低階修正是正比於kt的二階項g(ξ)=2ξ²erfc(ξ)-(2ξ/√π)e^(-ξ²),我手算得到
g(ξ)在ξ≈0.45附近取最大值約-0.202。代入實際數字kt=0.0324(24小時、4.4°C、
pH6.2的冷藏醃漬),算出最大相對修正只有約1.2%。
 
請幫我:
1. 檢查這個Laplace變換的推導與最終公式是否正確——特別是「一階項完全抵銷、
   最低階修正是kt量級」這個結論,這是我手算泰勒展開得到的,想確認方向對不對。
2. 我用g(ξ)最大值約-0.202這個數字時只手算了幾個ξ值(0.3到0.55之間取樣),
   想知道有沒有更嚴謹的方法確認這個函數的真正極值點與極值,而不是憑手算
   取樣猜測。
3. 如果要把這個模型從半無限介質(x>0)換成有限厚度的平板(例如真實培根切片
   兩面都有邊界條件),完整解會變得多複雜?是否還有解析解,還是必須轉向
   數值方法(有限差分/有限元)?