一種新的土壤凍融特征曲線模型

付子騰， 吳青柏， Miles DYCK， 何海龍

（1.西北農林科技大學資源環境學院，陜西 楊凌 712100； 2.中國科學院西北生態環境資源研究院凍土工程國家重點實驗室，甘肅 蘭州 730000； 3.Department of Renewable Resources，University of Alberta，AB，T6G 2H1，Canada）

0 引言

凍土是指溫度保持在或低于水的冰點，并含有冰的巖石和土壤，分為短時凍土（凍結時間數小時至半月）、季節凍土（凍結時間在半月至數月）和多年凍土（凍結時間在兩年以上）［1］。在凍土中，冰和液態水（未凍水）共存，其含量隨凍融過程中溫度的變化而變化。未凍水含量和含冰量的變化顯著影響凍土的物理和機械性能，從而決定寒區工程建設（如：房屋、道路、石油管線等）的成敗［2-3］。同時，凍融循環過程也會帶來如土壤凍融侵蝕、洪水、污染物的遷移、地下水上升引起的土壤鹽漬化等一系列環境問題［4-7］。因此，掌握凍融過程中土壤水熱耦合的變化機理對生產生活有至關重要的作用。

通常，純水的凍結發生在某一確定的溫度下［8］，具有確定的凝固點，但在濕潤的土壤中，即使當溫度低于自由水的冰點時，土壤中所含水分也不會立刻轉化為固態，而是隨著溫度的降低表現出相對平滑的轉變，以曲線的形式降低［9-12］。由于土顆粒表面能的作用，土體中始終會保持一定數量的液態水被稱作未凍水［1］，當溫度繼續降低，更多的未凍水會轉化成冰，此時毛管力和吸附力會被放大，剩余的未凍水會以薄膜水的形式存在［13］。當溫度低于一定值時，隨著溫度的進一步降低，未凍水含量保持不變，該溫度下的未凍水含量被稱為殘余水含量θres［14-16］。凍融過程中未凍水含量θl隨溫度T（℃）的變化曲線被定義為土壤凍融特征曲線（SFTC），代表了凍土中液態水含量和能量狀態之間的關系［17］，其中未凍水含量θl是指土壤孔隙中沒有凍結的那部分含水量（用體積含水量表示，單位為cm3·cm-3）。一般來說，決定凍土中未凍水含量大小的主要因素是溫度，其次是孔隙水壓力和土質（包括土顆粒的礦物化學成分、分散度、含水量、密度、水溶液的成分和濃度等）［18］。溫度作為主要影響因子，直接決定凍土中未凍水含量的高低，其他因素雖不是主要影響因子，但也顯著影響著凍融特征曲線的變化。在一定的溫度范圍內，未凍水含量隨孔隙水壓力和溶質濃度的增大而增大；相同條件下溶質濃度越大的土壤殘余水含量也越高，此外，土體液態水中溶質的存在還會引起水分凝固點的降低［9］，從而影響土壤凍融特征曲線。土壤的質地和結構包括比表面積、礦物組分、土壤顆粒級配等，這些性質會影響毛管力和吸附力的大小，進而決定凍土中未凍水含量及殘余水含量［19］。

有關凍融特征曲線數學表達的研究已有很多，許多學者也各自通過實驗數據和理論推導建立了相關的模型［15，20-34］。這些模型一般包括兩類：一類是基于土壤基本物理性質提出的半經驗公式［14-15，23］，一類則是通過類比水分特征曲線而提出計算未凍水含量的公式［7，22，24-25］。這些模型雖然在各自文章中能較好地擬合實驗數據，但在實際應用中仍然有其局限性：（1）很多模型只適用于飽和砂土，而且忽視了溶質濃度所引起的凝固點降低現象［14-15］；（2）認為未凍水含量的變化只依賴于溫度和土壤屬性，忽視了土壤中初始含水量的影響［34］，而Tice和Anderson［35］的研究表明，在給定的溫度下，土壤未凍水含量會隨著土壤初始含水量的增加而增加；（3）部分模型形式復雜，參數繁多且只能應用于某種單一的土壤中，如果推廣到其他土壤中，則預測結果偏離實際值［36］。

本文回顧了現有的凍融特征曲線模型，基于水分特征曲線與凍融特征曲線的相似性和平行關系，提出一種新的凍融特征曲線表達式。新模型與van Genuchten［37］水分特征曲線（SWCC）相似，考慮到了初始含水量對未凍水含量的影響，包含兩個參數α、β，其中參數α和土壤自身的物理性質有關，參數β則受初始含水量的影響。此外，本文還根據已有數據對模型作了修正，考慮到溶質濃度對未凍水含量的影響，加入參數Tf，使其能適用于不同溶質濃度的土壤中。下面將對模型的確立過程、參數擬合、以及模型的驗證作詳細的說明，同時作模型參數的敏感性分析，以研究模型參數α和β的具體物理意義。

1 模型回顧

凍融特征曲線既能為研究土壤凍融過程和高寒地區水文模型提供理論支持，同時也是進行寒區工程設計及施工建設的理論基礎，能否精確模擬某種土壤中的水熱耦合運移對于寒區生產建設至關重要。許多學者已基于土壤的各種基本參數建立了相關的經驗模型，包括線性、分段、指數和冪函數模型等［36，38］。但由于早期對于凍結過程中恰當的數學表示形式存在分歧，到現在仍未解決，以至于回顧這些模型和經驗公式時，現有模型結構之間由于適用條件和范圍的不同而表現出非常明顯的差異。

1.1 Anderson and Tice（1972）模型

土壤物理特性（如土壤干密度、土壤顆粒大小、土壤孔隙度、比表面積等）都會影響凍土中未凍水含量，進而影響土壤凍融特征曲線的變化［37］。因此，根據土壤的物理特性，尤其是比表面積，可以建立一個簡單的模型來預測土壤中未凍水含量［28］：

式中：θl為體積未凍水含量（cm3·cm-3）；ρd和ρw分別為土壤干密度和標準狀態下未凍水的密度（g·cm-3）；S為土壤比表面積(m2·g-1)；Ac為塑性指數與土壤顆粒小于2 μm 含量所占百分數的比值；T表示溫度，本文中全部使用攝氏溫度（℃）來表示。

Anderson［20，39］通過等溫量熱計和單位表面積歸一化來確定土壤未凍水含量，根據實際測定的未凍水含量數據總結了不同含水量與比表面積之間的關系，提出了一種更為簡單的冪函數模型：

式中：a和b為土壤特征參數與土壤比表面積有關。Kozlowski［14］給出了參數和土壤表面積之間的關系：

該式成立的隱含條件是假定在無窮低的溫度下殘余水含量為零［40］。聯立式（2）和（3）：

式（4）表明，只需測定土壤比表面積和干密度即可描述未凍水含量隨溫度的變化特征。因此，盡管冪函數模型提出已有幾十年歷史，但由于其形式簡單便于應用，各種文獻及資料中仍普遍采用該方法來估算凍土中的未凍水含量［1］。

1.2 Mckenzie（2007）模型

Mckenzie 等［15］提出了一個線性分段函數來近似地描述土壤凍融特征曲線：

式中：m為線性參數，Tres為未凍水含量減少到殘余水含量的溫度，假定為-12 ℃［15］。此外，Mckenzie等［15］還給出了另一個指數函數來預測未凍水含量：

式中：n為擬合參數；T0為標準狀態下自由水的冰點。相比式（5）、式（6）對未凍水含量與溫度之間關系的描述更符合實際的實驗結果［15］。

1.3 Kozlowski（2007）和Zhang（2008）模型

Kozlowski［14］對6 種 黏 土 進 行 差 示 掃 描 量 熱 法實驗后提出了土壤凍融特征曲線的半經驗模型，該模型由三部分組成：

式中：wl、wtotal、wres分別為用重量含水量表示的未凍水含量、土壤總含水量和殘余水含量；Tf代表凝固點的降低量。用體積含水量替換重量含水量則式（7）可變形為：

在這里Kozlowski 通過土壤塑限和含水總量來計算Tf：

式中：θp為土壤塑限（%）。Kozlowski認為殘余水含量與土壤質地和結構等有關，并通過實驗證明，殘余水含量θres和比表面積S之間存在以下關系：

實驗結果表明，當土壤溫度在0~-1 ℃之間時，土壤中的未凍水含量急劇下降，近似呈線性趨勢。因此，當負溫接近于0 ℃時，未凍水含量又可以表示為土壤溫度的線性函數［41］：

1.4 Bai and Lai（2018）模型

Bai 等［16］考慮到初始含水量和殘余水含量，提出一個分段函數模型來預測未凍水含量：

式中：T0為標準狀態下自由水的冰點，σ為模型的擬合參數。該模型可用于不同初始含水量的土壤。

2 研究方法和材料

2.1 新的凍融特征曲線模型

在融土中，隨著土壤含水量的變化，土壤水吸力也會發生變化［42］。與凍融特征曲線的定義類似，土壤含水量與土壤水能量狀態即土壤水的基質勢之間的曲線關系被稱為土壤水分特征曲線（SWCC）［36］，是土壤持水特性的表征。很多學者已比較了SFTC和SWCC之間的相似性［14-15，23，43-45］。早在20世紀60年代，Koopmans和Miller［17］就證實當土壤水分完全由毛細管力或吸附力控制時，SWCC 和SFTC 之間存在平行關系；Black 和Tice［40］的研究表明，在飽和土體下，相同密度和體積的土壤的SWCC和SFTC可進行定量比較，有學者結合Clapeyron與SWCC 導出了非飽和土壤中的SFTC［36］。基于歸一化的無量綱土壤含水量與土壤基質勢之間的關系可表示為［18］：

式中：a、n、m為未確定的參數；m= 1-1/n；h為土壤基質勢；用壓力水頭表示；Θ可由下式計算：

式中：θ、θr、θs分別為實際含水量、飽和度和殘余水含量。聯立式（13）和式（14），即可得到土壤水分特征曲線的Van Genuchten［37］表達式：

注意到土壤凍融特征曲線和水分特征曲線之間的相似性：（1）都是描述土壤液態含水量隨土壤某個物理量變化的函數；（2）長期的實驗數據表明，土壤SFTC 的變化特征同SWCC 相似（圖1）；（3）凍結-融化過程的SFTC 和脫濕-吸濕過程中的SWCC都會產生滯后現象［42-43］；（4）SWCC 中土壤含水量居于殘余水和飽和水含量之間，類似地SFTC 中未凍水含量居于殘余水和凍土中初始含水量（或總含水量）之間。因此，歸一化理論同樣可以用于描述SFTC：

式中：θinit、θres分別為初始含水量和殘余水含量（cm3·cm-3）；θl為未凍水含量。將式（15）中的自變量h替換為溫度T，即得：

由式（16）和（17）即可導出形式類似于Van Genuchten［18］的歸一化公式：

式中：θinit、θres分別為初始含水量和殘余水含量(cm3·cm-3）；T為溫度（℃）；α、β為模型的擬合參數。

未凍水含量的變化是多種因素綜合作用的結果，以往的模型大都忽略了初始含水量和溶質濃度的影響。然而研究表明，當土壤凍結時，隨著土壤水分逐漸向冰相的轉變，土壤溶液的濃度會增加（假設冰相中不含溶質），溶質的存在會降低土壤水分的凝固點，從而影響未凍水含量［47］。此外，凍融過程中土壤的凝固點還會隨初始含水量的增加而升高，隨溶質含量的增加而降低［10，21，48-50］。在式（18）中，我們考慮到土體初始含水量θinit的影響，等式中的兩個參數α、β或許能解釋溶質濃度對未凍水含量的影響，下文將進行參數α、β敏感性分析，以探索其具體物理意義。

2.1.1 對參數α敏感性分析

如圖2（a）顯示，若保持初始含水量不變，當固定參數β不變時，α的變動導致了整個模型曲線發生變化。保持參數β等于2.0 不變，令α分別取0.5、1.0、2.0、3.0、5.0。隨著α的增加，曲線整體下凹并向左移動。結合圖2（b）可以發現，新模型曲線隨參數α的變化規律與新模型在不同溶質濃度土壤中的曲線變化規律相似這表明參數α或許和土壤溶質濃度有關，隨著溶質濃度的增大，α逐漸減小。

圖2 式（18）中參數α的敏感性分析（a）；式（18）在不同溶質濃度土壤中的曲線變化（b）Fig.2 Unfrozen water content changing with temperature for various α in the Eq.（18）（a）and for various solute concentration of the Eq.（18）（b）

2.1.2 對參數β敏感性分析

保持α固定不變，曲線隨參數β的變化如［圖3（a）］所示，保持α等于2.0 不變，令β分別取1.2、1.5、2.0、3.0、5.0，從圖3 中可以看出，隨著參數β的增加，土壤凍融特征曲線向下移動。保持初始含水量不變，曲線的最低點，即殘余水含量隨β的增加而減少。結合圖3（b），在相同初始含水量下，隨著黏粒含量的減小，殘余水含量也越來越低，表明β的取值可能同土質有關。受土壤物理性質（如比表面積、干密度、黏粒含量等）的影響，黏粒含量越小，β越大。另一方面，參數β的增加使得凍融特征曲線在0~-2 ℃的溫度范圍內變得越來越陡峭，當β取到5.0 時曲線幾乎垂直。這與圖2 中曲線隨α的變化規律相似，說明參數β同樣與溶質濃度有關。

圖3 式（18）中參數β的敏感性分析（a）；式（18）在不同質地土壤中的曲線變化（b）Fig.3 Unfrozen water content changing with temperature for various β in the Eq.（18）（a）and unfrozen water content changing with temperature for various soils with different textures of the Eq.（18）（b）

對參數α、β的敏感性分析表明，α受溶質濃度的影響，β則與土質和溶質濃度之間存在聯系，可以據此建立相關的參數轉換函數來計算不同土壤中α、β的值。

2.2 數據的來源與處理

為了驗證新模型的可行性及研究模型中參數α、β的實際物理意義，本文引用了文獻中已有的試驗數據（見文后附表1），進行參數的擬合與模型效果的評價。在進行數據篩選時需要考慮到：1）數據的完整性。必須要能獲取到完整凍結或融化過程的數據資料；2）數據的多樣性。所選擇的數據需要來自于盡可能多的土壤，便于模型驗證；3）優先選擇土壤物理性質比較詳細的數據。如土壤比表面積、黏粒含量、容重、孔隙度等，便于模型參數物理性質的分析。基于以上原則，本文所選土樣數據包括Smith 和Tice［51］文中的6 種粉質土（Hot Springs 粉土、Leda 粉土、Suffield 粉質黏土、Niagara 粉土、蘭州粉土），11 種黏土以及4 種代表性土壤礦物組分（膨潤土、高嶺石、鋰輝石、礫石）的凍結實驗數據；Spaans 和Baker［52］文中分別處于干燥和濕潤兩種狀態下Waukegan 粉砂土的凍結試驗數據；Suzuki［10］中濕地土壤的數據；Wen 等［53］青藏高原粉質黏土的凍結試驗數據；Christ 和Park［54］砂土、粉砂土、粉土的試驗數據；Wang 等［55］中Harbin 黏土、Lanzhou 黃土以及Ma 等［56］粉土和黏土的凍結實驗數據（數據選用的詳細信息見附表2）。

考慮到模型效果評價的合理性，將以上所有數據按照來源和土壤類型整理后分為兩組，其中一組包含70%的數據，用于模型的訓練和參數擬合，另一組包含30%的數據，用于模型的驗證和評價。此外，本文還選取了Patterson 和Smith［57］粉質黏土的凍結實驗數據以研究模型是否適用于不同溶質濃度的土壤中。該數據包含4 個對照組，即分別用0、10、20、35 g·L-1NaCl溶液飽和，然后進行凍結處理，得到不同溶質濃度下的凍結試驗數據。在本文中，同樣將其分為訓練和驗證數據兩部分。

2.3 參數擬合預評價指標

本文重點研究新模型的可行性及適用范圍。通過對實驗數據的分析，研究新模型能否適用于不同的飽和度和不同含水量的土壤中，檢驗新模型在不同類型土壤中的表現，前文中提到的五種模型被用于模型間的比較。此外，文中還討論了模型中參數的變化對擬合結果的影響。通過MathCAD（Prime 3.0，Parametric Technology Corporation），利用非線性最小二乘算法來進行參數α、β的擬合。

采用表1所示三種模型評價指標對模型性能進行評價：（1）均方根誤差（RMSE），用于衡量預測值與真實值間的偏差，反映測量的精密度。均方根誤差的值越小越好；（2）平均偏差（AD），用于來測定預測值與真實值算數平均值間的差異；（3）納什效率系數（NSE），用于模型的效率評價和對擬合結果的評判，NSE 的變化范圍從-∞~1，其值越接近1，則模型的可信度越高，預測值越靠近真實值。此外，文中還給出了各模型在部分土壤中預測值與實測值的1∶1線性圖以便更直觀地展示測試數據和模型預測值之間的差異。

表1 復合介質混合模型的性能評價指標Table 1 Performance evaluation indexes of the composite medium hybrid model

3 模型評價

3.1 新模型在不同類型土壤中的表現

使用文獻中已有的實驗數據擬合參數α、β以研究新模型在不同的類型土壤中的適用情況。圖4展示了新模型在7 種代表性土壤中的擬合結果，其中粉土和粉砂土是Christ 和Park［54］的野外實驗數據；Harbin 黏土來自Wang 等［55］，其余四種土樣均來自文獻Patterson［58］。從圖4可以看出，對于這7種土壤樣品，用新模型所得到的擬合曲線結果與實測值較為吻合，納什效率系數（NSE）都在0.95 以上，這表明新的模型可廣泛適用于不同類型土壤中。

圖4 不同土壤中未凍水含量與隨溫度變化的擬合曲線，其中散點為實測數據，實線為式（18）計算的模擬值Fig.4 Unfrozen water contents of various soils changing with temperature，where the dots are measured data and lines are calculated value by the Eq.（18）

表2 給出了相關的擬合指標，均方根誤差（AD）表明，在大部分土樣中，預測值稍低于實測值，但在統計學意義上誤差均在可接受范圍內。盡管不同土壤中模型參數有所差異，但對于給定的土壤，只要能通過實驗數據得到參數α、β的值，新模型就能準確預測相同條件下該種土的土壤凍融特征曲線。

3.2 新模型不同初始含水量下的表現

為了研究新模型在不同初始含水量土壤中的表現，本文選用了Wang等［53］的實驗數據。土樣是來自青藏高原的粉質黏土，比表面積為8.12 m2·g-1，干密度1.56 cm3·cm-3，保持其他條件不變，設置土壤初始 含 水 量 分 別 為0.254、0.313、0.357、0.398、0.459 cm3·cm-3，分別擬合出每個初始含水量下的參數α、β的值，然后通過該參數去計算所有初始含水量下的凍融特征曲線，觀察其差異性。從圖5 中可以看出，同一土壤分別在五種初始含水量下得到的模型參數對整體的擬合結果并顯示出沒有出現太大的差異，但相對來說，在初始含水量較低的條件下得到的模型參數在擬合整體時會更準確一些。表3 顯示不同初始含水量下得到的參數α、β值不同但差異較小，而且隨著初始含水量的變化β和α并沒有表現出明顯的線性的變化規律。總體來看，預測結果較好地吻合了實測值。表3給出的相關擬合指標顯示，模型的NSE都在0.96以上，也就是說，在任意土壤類型和初始含水量下得到的模型參數都是有效的。

圖5 不同初始含水量下新模型的表現，散點為實測數據，實線為新模型的模擬結果Fig.5 Unfrozen water content changing with temperature for various θinit（Note：the scattered points are observed and the solid lines are simulated of the new model）

表2 新模型在7種土壤樣本中的參數α、β 及擬合指標Table 2 A summary of the parameters α、β for the seven soils and performance measures of the new model to fit the 7 soils samples

3.3 新模型在不同溶質濃度下的表現

本文選擇Patterson 和Smith［57］粉質黏土來進行模型在不同溶質濃度條件下的檢驗，土樣干密度為1.19 g·cm-3，保持初始含水量0.523 cm3·cm-3不變，取溶質濃度為10、20、35 g·L-1NaCl 溶液，設置0 g·L-1NaCl 溶液作為對照。圖6 為式（18）在不同濃度溶液中的表現，從圖中可以看出，模型在低濃度時的表現依舊能很好地同實測值吻合，但隨著溶質濃度的增加，擬合曲線逐漸偏離實測值，且濃度越大，預測值與實測值的差異越大，可能是因為該等式同樣沒有包含溶質濃度對土壤水分凝固點降低的影響，考慮到這一點，引入一個新的參數Tf，式（18）可變為：

圖6 式（18）在不同濃度溶液中的擬合情況，散點為實測值，實線為新模型模擬結果Fig.6 Unfrozen water content changing with temperature of different solute concentrations of the Eq.（18）（Note：the scattered points are observed，the solid lines are simulated by the new model）

Tf表示土壤水分凝固點的降低，在模型中有兩種含義：（1）類似于α、β，看作一種擬合參數，沒有實際的物理意義；（2）表示土壤水分凝固點隨溶質濃度的降低，具有實際的物理意義，并且有相應的計算公式。本文考慮到兩種情況，對比作為參數擬合得到的Tf和作為變量計算得到的Tf，選出最適合的模型形式。

表3 五個初始含水量下分別得到的參數α、β值及其在不同初始含水量情況下的統計指標（加黑數值表示該含水量下對應的指標值）Table 3 The values of parameter α and β corresponding to the five initial water content θinit and their statistical indexes（Note：the corresponding index under initial moisture content θinit is indicated in bold type）

計算Tf的經驗公式有很多，Banin 和Anderson［49］描述了土壤中溶質的存在會降低冰點溫度，隨著土壤水分的凍結，溶質逐漸被滯留在越來越小的空間中，其中，凝固點的降低Tf與溶質濃度之間存在以下關系［59］：

式中：Tref為參考溫度，與溶質種類有關，NaCl 可取62 ℃［31］。Sn為溶質濃度（g·L-1NaCl）。Kozlowski［9］修改了式（9），給出了另一種計算Tf的方法：

式中：w為重量含水量（%）；wL為土壤液限（%）。此后，Kozlowski［60］又提出可根據土壤干密度和初始含水量來求得Tf：

式中：α和β為給定土壤的系數，是一個常數。Bodnar［61］根據試驗數據，回歸出了一個三次多項式來描述Tf與溶質濃度的關系：

用上述四種方法計算Tf，并分別與新模型耦合。表4顯示四種經驗公式的擬合結果在統計學意義上相差不大，考慮到經驗公式的不確定性，新模型采用第一種方法［式（21）］來計算Tf值。圖7顯示在加入Tf后，低濃度條件下，擬合曲線并未有太大變動，但在溶質濃度較高時，式（19）的擬合結果明顯要優于式（18），而且相比于用經驗公式計算的Tf值，參數擬合得到的結果要更加精確（表5），而濃度為35 g·L-1NaCl溶液的樣本在三種方法的計算下NSE都沒有超過0.9，這可能是由于在高濃度條件下，實測數據的不準確所導致。

圖7 加入Tf后的新模型在不同濃度溶液中的擬合情況，Tf用經驗公式和參數擬合兩種方法獲取。其中實線未考慮Tf影響，短劃線Tf通過經驗公式計算，虛線的Tf通過擬合獲得Fig.7 Unfrozen water content changing with temperature while considering effects of Tf under different solute concentrations（Note：solid line is not considered the influence of Tf，short line Tf is calculated by the empirical formula and the dashed line Tf is obtained by fitting）

綜合來看，使用式（18）和式（19），均能準確地描述不同溶質濃度下的未凍水含量隨溫度變化的關系，但式（19）能更好地反映出溶質濃度的變化對土壤凝固點降低的影響。因此，在實際應用過程中，當土壤溶質濃度較低時，可直接通過式（18）得到準確的凍融特征曲線，而當土壤溶質濃度較高時，通過給模型加入Tf［式（19）］可有效提高其擬合精度。

表4 四種經驗公式計算的Tf值在不同土壤溶質濃度下的表現，其中參數α、β通過公式（19）擬合得到Table 4 The statistical indexes of the four empirical formulas calculating the Tf at different soil solute concentrations，of which the parameters α，β are fitted by the Eq.（19）

表5 式（18）和式（19）在不同溶質濃度土壤中的擬合指標，式（19）中Tf的值分別用經驗公式（20）計算和參數擬合Table 5 The statistical indexes of Eq.（18）and Eq.（19）in different soil solute concentrations，the values of Tf in Eq.（19）are calculated by the Eq.（20）and by parameter fitting，respectively

3.4 模型間的對比

本文選擇了前文中提到五種較為典型的模型與新模型對比，包括兩種線性模型［15-41］，兩種分段函數模型［14］以及早先提出的指數模型［20］。四種土壤（包括砂土、粉砂土、粉土、哈爾濱黏土）被用于對比，砂土、粉砂土、粉土以及一種哈爾濱黏土。其中粉土、粉砂土及砂土凍結數據均來Christ 和Park［54］，凍結過程未凍水含量的測定方法為時域反射儀（TDR）法；哈爾濱黏土數據來自Wang等［55］，測量方法為核磁共振法（NMR）。

圖8 和表6 顯示了各模型在四種土壤中的擬合結果。在砂土中，除了Zhang等［41］的線性模型外，其余五種包括新模型都能較為準確地描述未凍水含量和溫度之間的關系。在其他三種土壤中，新模型的擬合精度最高，Bai 等［16］的模型次之，Anderson等［20］的預測結果與實測值偏離最大，預測值往往低于測值，且隨著土壤黏粒成分的增加和土壤顆粒半徑的減小，偏離值增大，可能是由于該模型忽視了初始含水量，只考慮到了土壤比表面積的影響所導致的。Zhang 的線性模型考慮到了初始含水量，其曲線的斜率與初始含水量和殘余水含量有關，因此相對于觀測值，該模型在未凍水含量接近殘余水含量之前總是顯示出過高地估計。總體上看，除了新模型和Bai 的模型外，其余模型在應用到砂土以外的土壤中時均會較大程度地偏離真實值，而且隨著土壤顆粒級配的變化這種差異也會越發顯著。此外，部分土壤樣品被用來擬合模型的1∶1散點圖（圖9），研究預測值和實測值之間的差異。從圖9 可以看出，在所有的土壤樣品中，新模型的預測值都非常接近真實數據。其他模型中，Bai 模型的表現良好，預測值與真實值也較為接近；Kozlowski 模型由于可用數據太少，無法得出準確結論，但從圖中擬合結果來看，在砂土和粉砂土中表現較好。Anderson 模型在未凍水含量較低時表現較好，在較高的未凍水含量下則出現明顯的偏離。

圖8 新模型與另外五種模型分別在砂土、粉砂土、粉土及黏土中的表現Fig.8 Unfrozen water content changing with temperature for sand（a），silt-sand（b），silt（c）and Harbin clay（d），measured and simulated by the new model and other five quoted models

圖9 6種模型在部分土壤數據的預測值和實測值1∶1散點圖Fig.9 The 1∶1 scatter plots of unfrozen water content simulated by the new model and the other five quoted models against the measured

4 結論

本文考慮到初始含水量和溶質濃度對凍融特征曲線的影響，類比水分特征曲線模型，提出一種包含兩個參數的新凍融特征曲線模型。利用7種典型土壤樣本對模型適用性進行了驗證，結果表明新模型可適用于任意類型和飽和度的土壤中。通過加入新的參數Tf，可在統計學上顯著提高模型在高溶質濃度條件下的預測精度。此外，明確模型參數α和β在模型中的物理意義有助于模型的完善和應用，敏感性分析表明，參數α與土壤溶質濃度有關；參數β與土壤溶質濃度和土壤自身性質如比表面積、干密度和黏粒含量有關，可根據這些影響因子可建立參數轉換函數，方便模型在實際中的應用。由于收集數據的不完整，本文未能建立有效的參數轉換函數，但給出了用于測試模型的38種土壤樣品的參數值（附表1）以供參考。

表6 新模型和其他五種模型的擬合指標Table 6 The statistical indexes of the new model and the other five quoted models

附表1 38種土壤樣品的參數α和β值及模型擬合指標Attached Table 1 The parameters α and β of the 38 soil samples，together with determination method and index

附表2 土樣數據的來源、測定方法和基本性質Attached Table 2 Sources，determination methods and basic properties of the soil samples