固體表面玻璃化冷凍方法降溫結晶特性數值研究

朱文欣 趙芷慧 黃永華 程錦生 陳 威 李 錚

(1 上海交通大學機械與動力工程學院 上海 200240;2 上海東富龍醫療裝備有限公司 上海 201108;3 上海交通大學附屬第一人民醫院泌尿中心男科 上海 200080)

細胞療法被認為是未來醫學的三大支柱之一,生物樣本將成為重要戰略資源[1]。由于生物的活性在低溫環境中受到抑制,代謝等化學反應減慢近乎停滯[2-3],使細胞等生物長期凍存成為可能[4]。然而,冷凍過程溫跨大,造成細胞等生物在冷凍過程中產生不可逆損傷,如氧化性損傷、滲透性損傷和機械損傷等[5-6]。

為減少冷凍損傷,多種冷凍方法在細胞等生物樣本凍存領域獲得應用[7],如直接接觸液氮冷卻法[8]、環網狀冷凍法[9]、微流控冷凍法[10]等。然而上述方法均易造成樣品受熱不均、回收困難并存在樣品污染的風險。而固體表面玻璃化冷凍法(solid surface vitrification,SSV)由于降溫速率快、降溫速率控制性好、樣品不直接與液氮接觸、易回收、無污染等優點,被越來越多地應用于低溫生物凍存[11-12]。目前,針對SSV方法的研究多數為實驗或經驗性報道[13-14],主要聚焦于降溫速率的大小,對于SSV方法在凍存過程中造成的冷凍樣品結晶機理與結晶度變化規律鮮有涉及。此外,由于凍存過程時間短,通過實驗手段以高時間分辨率獲得全冷凍過程樣品全域溫度、結晶度變化仍較為困難,如采用熱電偶測溫不僅侵入性大、易結露、漏熱嚴重,而且無法測量結晶度;低溫顯微系統價格昂貴且無法定量分析冷凍過程的結晶變化[15]等。

鑒于此,本文采用數值計算方法探究SSV方法實現快速降溫的核心要素,并通過建立凍存過程的非等溫結晶預測模型,研究預冷溫度、樣品體積及形狀對冷凍結晶效果的影響,為冷凍方法優化、冷凍對象體積形狀等設計提供理論指導。

1 基于非等溫結晶模型的SSV方法建模

SSV方法采用一定質量的高導熱金屬(銅、銀或鋁合金)塊或片浸入液氮或其他低溫介質預冷,將樣品滴在預冷后的低溫表面,通過瞬間接觸導熱實現快速降溫冷凍,物理過程如圖1所示。

圖1 SSV方法工作原理

聚焦液滴碰觸金屬表面之后的傳熱過程,建立如圖2所示的數學物理模型。鑒于軸對稱性,該模型采用二維軸對稱形式,取一半,其中左上角凸出部分為液滴,其余長方形區域為金屬基底載冷臺。以5 μL液體冷凍為例,真實滴落液滴的幾何尺寸可以近似為半徑1.82 mm,中心高度為0.89 mm的球冠,液滴的底部與金屬塊的上表面直接接觸。此處考慮無氧銅材質的基底載冷臺,其導熱系數在低溫下達400 W/(m·K)。為了消除邊緣效應,取銅塊的寬度為50 mm,高度為20 mm,此時相對液滴而言,可以認為銅塊無限大。

圖2 數學物理模型和邊界條件

基于傅里葉導熱定律,金屬塊內的導熱可表示為:

(1)

式中:T為溫度,K;ρ為密度,kg/m3;λ為導熱系數,W/(m·K);cp為比定壓熱容,J/(kg·K);t為時間,s。

(2)

(3)

式中:L為熔融潛熱,kJ/mol;X為結晶度。忽略結晶界面處由于分子擴散對冰晶前沿生長速率的影響,根據半經驗非等溫結晶模型[16-20]:

(4)

式中:Tm為熔融溫度,K;k1為結晶常數,與液體種類、質量濃度相關;Q為活化能,kJ/mol;R為氣體常數8.314 J/(mol·K)。

由于實際采用的拋光銅表面粗糙度Ra為0.8,視作光滑表面條件,且該銅表面與液滴直接接觸,因此在進行建模時忽略了金屬塊與液滴樣品間的接觸熱阻:

TS=TC

(5)

(6)

在邊界條件方面,將金屬塊與液氮接觸邊界設為第一類邊界條件,由于金屬塊與液氮蒸氣對流換熱,因此其他邊界設為第三類邊界條件,對流換熱表面傳熱系數h取為100 W/(m2·K)[21]。

在熱物性方面,銅塊選用常物性參數,如表1所示。

表1 固體材料熱物性參數

需要說明的是,該模型側重于了解被冷凍液滴所經歷的傳熱過程,并實現與系統的傳質過程解耦,其分析尺度大于冰晶生長的微結構,且不包含溶液傳質濃度變化對冰晶微結構的影響。

在COMSOL軟件中,模型劃分自由三角形網格,在樣品與預冷表面接觸處做細化處理,共劃分了1 500、2 500、3 035、4 000四種密度的網格,網格平均質量均優于0.93。經網格無關性驗證,取3 035個網格進行計算在速度和精度方面均滿足要求。對金屬塊用固體傳熱接口求解,對CPA液滴利用導熱接口添加廣義熱源項,通過域常微分和微分代數方程接口定義非等溫結晶源項,并耦合求解。添加溫度判據以保證結晶在熔融溫度下才發生。為保證對極其微小量數值計算正常,將結晶度初始值設為1×10-16(可識別的接近0的最小正值)。采用瞬態求解器,向后差分方法,時間步長不超過0.01 s,相對容差控制為1×10-6。

采用文獻[22]的實測數據來驗證模型的適用性,針對西林瓶(類薄層式SSV結構)中1.9 mL甘油水溶液置于液氮中冷凍過程溫度變化進行對比驗證。需要說明的是,驗證時模型中固體材料選用玻璃,溶液選用甘油水溶液,與文獻[22]的實驗保持一致。預測結果和文獻數據的對比如圖3所示。由圖3可知,耦合非等溫結晶源項的傳熱區域模型可以較好地預測冷凍過程,最大誤差約為10%,平均誤差小于2%。后期溫度略高于實驗值,可能是由于被測樣品溶液的導熱系數在低溫時有所升高。

圖3 模型計算值與文獻[22]實驗值對比

2 計算結果與分析

模型驗證采用的實驗樣本液滴為生物細胞冷凍操作常用的冷凍保護劑(CPA)丙二醇溶液(質量分數為35%),其玻璃化轉變溫度Tg為160 K,冷凍結晶的危險溫區介于Tg與熔化溫度Tm之間(160～255.5 K)。由于CPA物性參數通常缺乏低溫區的實測數據,采用常溫區參數一般會低估物性參數的變化;而若采用水和冰的物性參數替代低溫保護劑物性參數,會高估低溫區物性參數的變化。在沒有準確合適的物性參數數據時,有必要假設物性參數不隨溫度變化來參考計算[23]。一般而言,隨著溫度的降低,CPA導熱系數增加,而比熱容和密度降低,總體上熱擴散系數增加。因此利用常物性參數計算可能會低估冷卻速率,即實際降溫效果優于計算值。樣品液滴(丙二醇溶液)物性參數來源于文獻[24-25],如表2所示。利用上述模型和物性參數,開展被冷凍樣品滴落后的冷卻和結晶特性研究。

表2 冷凍液滴熱物性參數

2.1 非等溫結晶特性

首先從液滴對稱軸的高、中、低三個位點進行局部降溫結晶特性分析,如圖4所示。三個位點的冷凍降溫過程均未呈現溫度平臺期,無明顯的一級相變特征。由于低點最接近預冷表面,傳熱路徑短、熱阻小,因此液滴溫度在低點處最接近預冷表面溫度,降溫速率最快,相比于液滴中點及高點的降溫速率有數量級優勢。不同的降溫特性導致液滴起始結晶時間和最終結晶度不同。根據Uhlmann理論,結晶度小于10-6即為玻璃化[26],低點的結晶度幾乎不變,趨向初始設定值10-16量級,實現了玻璃化轉變。中點降溫速度比高點更快,先于高點降溫到熔融溫度,所以結晶比高點更早產生。同時,由于中點比高點的傳熱路徑更小,傳熱速率更大,能較快通過結晶危險溫區,最終結晶度小于高點。只有以較高的溫變速率快速通過結晶危險溫區,才能實現被冷凍樣品的玻璃化轉變。

圖4 樣品內低、中、高點處溫度與結晶度變化規律

對液滴樣品選取切面進行整體降溫結晶特性分析,如圖5所示(注:圖中只針對實現玻璃化冷凍及處于危險溫區的部分顯示云圖)。液滴的降溫最先從底部發生并實現玻璃化轉變。樣品中間結晶度增大,并向上擴散結晶,由于液滴邊緣薄中間厚,結晶度與等溫線均呈上凹形。但最大結晶度未出現在距離冷源最遠、溫度最高的液滴頂端,而出現于距頂端一定距離的中部。需要說明的是,由于本文結晶模型與傳質解耦,因此,液滴中部存在最大結晶度這一現象與冷凍過程中樣品界面處存在液氣-固氣相變無關,今后可進一步分析液滴內部復雜的溫度變化情況對樣品結晶度的影響。

圖5 液滴內部結晶度和溫度云圖

進一步對冷凍過程中的液滴頂點處與結晶度最大點處開展降溫速率和結晶度變化速率特性分析,如圖6所示。頂點的降溫速率和最大結晶處的降溫速率趨勢一致,但被冷凍液滴內部的降溫速度存在不均勻性,不利于大體積樣品的均勻冷凍。最大結晶度變化速率點比最大降溫速率點滯后約0.3 s,表明結晶過程存在滯后性。在整個液滴均出現結晶度后,頂點的降溫速率始終大于最大結晶點處的降溫速率,頂點的液體更快通過結晶危險溫區,表明更大的降溫速率可以抑制結晶,使被冷凍樣品更容易實現玻璃化。

圖6 液滴頂點、結晶度最大處的降溫速率與結晶度變化速率

需要說明的是,上文通過引入冷凍樣品的結晶度X,更直觀地展示了樣品在被冷凍過程中的內部狀態,可預測樣品的被冷凍效果,為設法實現無損凍存的各種嘗試提供便捷高效的驗證方法。

2.2 固體表面玻璃化冷凍影響因素

針對不同的生物樣本,通常需要設定不同的冷卻預冷溫度。因此,有必要探索預冷溫度對樣品在固體表面被玻璃化冷凍的影響。液氮是冷卻金屬載物臺最常見的工質,其預冷極限溫度Tprec,min為78 K,而35%丙二醇溶液的Tm為255.5 K,Tg為160 K。因此,分別選取Tprec為80、100、120、140、160、200、250 K。圖7所示為預冷溫度對液滴結晶度和穿越危險溫區時降溫最慢處的平均降溫速率vcr的影響。預冷溫度越低,降溫速率越慢,結晶度越小,預冷溫度與降溫速率和結晶度均無線性關系。此外,當預冷溫度低于樣品玻璃化轉變溫度時,樣品的降溫速率均較快,達到103K/min量級,其結晶度在10-5量級;而當預冷溫度高于Tg時,降溫速率在102K/min量級,結晶度在10-5～10-4量級;當預冷溫度越接近熔融溫度,降溫速率越慢,極易結晶。因此,預冷溫度越低,傳熱溫差越大,越有利于實現樣品的玻璃化冷凍。

前文非等溫結晶模型計算表明,最大結晶處可能在被冷凍樣品的中上部。因此,有必要探究樣品體積對用SSV方法冷凍的樣本的降溫速率和結晶度的影響。基于現有實際應用和文獻調研情況,樣品體積V分別選取1、5、10、20、50、100、500 μL,定義玻璃化比例VR為樣品中實現玻璃化冷凍體積與樣品總體積之比,計算結果如圖8所示。當體積為1 μL時,樣品完全玻璃化冷凍,因此,完全玻璃化轉變的安全樣品體積為1 μL。隨著樣本體積增加,樣本內部導熱路徑變長,降溫速率變慢,結晶顯著增加,VR不斷下降,且呈非線性關系。在一定范圍內,體積越大,結晶度變化越劇烈;而體積過大,結晶度將完全結晶。因此,樣本體積是影響SSV方法的關鍵因素之一。本文提出的耦合熱源項的非等溫結晶模型可以預測不同體積的樣品內部結晶程度和不同條件下玻璃化安全體積,為實際固體表面玻璃化冷凍設計提供一定依據。

圖8 不同體積樣品中結晶度和玻璃化轉變率

此外,根據傅里葉導熱定律和牛頓冷卻公式,有效傳熱面積是研究傳熱問題的關鍵因素之一。因此,引入傳熱面積與總體積比值A/V參數,研究樣品形狀對SSV方法冷凍的影響。樣品體積為5 μL,改變樣品形狀,如圖9所示。正球缺即為圖 2介紹的模型;反球缺與正球缺尺寸相同倒置在挖空的金屬預冷表面;反半球形樣品半徑為1.34 mm;反圓柱形半徑為1.51 mm,厚0.7 mm;反薄層形半徑為2.3 mm,厚0.3 mm。

圖9 樣品形狀

不同形狀液滴結晶度和玻璃化比例的計算結果如圖10所示。由圖10可知,相同體積下薄層形樣品的A/V更大,玻璃化程度更高,而其他形狀樣品出現部分結晶。因此,A/V一定程度上可以反映傳熱拓撲結構的優劣。相同情況下,A/V值越大,結晶度越小,樣品越容易實現玻璃化冷凍。若以實現玻璃化冷凍為目標,樣品形狀應盡量薄層化。說明本文模型可以簡便高效地研究不同形狀樣品在冷凍時內部結晶情況,可為實際冷凍樣品形狀設計提供依據。

圖10 不同形狀液滴結晶度和玻璃化比例

3 結論

本文對基于固體表面玻璃化冷凍方法的冷凍液滴建立了耦合非等溫結晶源項的傳熱區域模型,研究了液滴降溫特性和冷凍液滴內部的結晶情況,得到如下結論:

1)固體表面玻璃化冷凍方法可以實現快速降溫,降溫速率可達到103K/min量級。其可獲得快速降溫能力的核心要素是固體冷臺得到了充分預冷并儲存足夠冷量。

2)降溫速率越快,結晶程度越小。結晶易出現在遠離預冷表面的樣品上部,結晶最嚴重區域不一定位于樣品頂端。

3)預冷溫度越低、樣品體積越小、樣品面積/體積(即薄層化)越大,結晶程度越小,越容易實現玻璃化冷凍。

4)耦合非等溫結晶源項的傳熱區域模型可以通過結晶度X反映冷凍過程中樣品內部結構變化,可作為設計開發新型冷凍設備和低溫保護劑的輔助手段。

本文受上海市科委創新行動計劃(20S31903400);上海交通大學“交大之星”計劃醫工交叉研究基金(YG2021QN86)項目資助。(The project was supported by Innovation Action Plan of Shanghai Science and Technology Commission (No. 20S31903400), and Shanghai Jiao Tong University Med-Engineering Fund (No. YG2021QN86).)