高速鐵路箱梁?無砟軌道結(jié)構(gòu)溫度場模型的級數(shù)解

戴公連，張強強，劉文碩，黃志斌

(1.中南大學土木工程學院，湖南長沙，410083；2.東南沿海鐵路福建有限責任公司，福建福州，350013)

橋上無砟軌道是我國高速鐵路主要的軌道結(jié)構(gòu)形式，橋梁體系的箱梁軌道系統(tǒng)的溫度作用模式是高速鐵路橋梁軌道結(jié)構(gòu)設計中必須考慮的因素之一。對于橋上軌道結(jié)構(gòu)而言，中國高速鐵路規(guī)范[1]規(guī)定無砟軌道結(jié)構(gòu)的豎向正溫度梯度為90 ℃/m，通過對無砟軌道溫度梯度的實測及分析發(fā)現(xiàn)[2?3]，我國中部地區(qū)的無砟軌道溫度梯度普遍高于80 ℃/m，實測溫度梯度甚至達到100 ℃/m。由于軌道結(jié)構(gòu)的遮擋效應，橋梁頂面的溫度場呈現(xiàn)不均勻性，對處于施工或服役狀態(tài)的橋梁溫度分布及溫度效應的影響也均不同[4?5]。現(xiàn)有的中國規(guī)范、歐洲規(guī)范均沒有給出相應的箱梁?無砟軌道結(jié)構(gòu)的溫度計算模式，因此，箱梁軌道系統(tǒng)的溫度分布特征仍亟待進一步研究。目前，很多學者主要采用數(shù)值模擬、統(tǒng)計分析等方法對軌道板和箱梁的溫度梯度分布模式與溫度取值進行了大量研究。閆斌等[6]利用無砟軌道溫度場的有限元模型，分析了地理位置、熱物性參數(shù)、環(huán)境等因素對無砟軌道溫度梯度的影響。朱勁松等[7]建立了橋?軌溫度場精細化有限元模型，研究了遮擋效應對結(jié)構(gòu)溫度場的影響。戴公連等[8]利用CRTS Ⅱ型軌道板的實測結(jié)果建立了箱梁?無砟軌道結(jié)構(gòu)的溫差GPD 極值模型，給出了具有一定概率的軌道結(jié)構(gòu)溫度取值。申建康等[9]應用機器學習的方法對軌道板溫度場進行預測，具有很高的預測精度。但在實際工程中，研究者大多采用半理論半經(jīng)驗公式或全經(jīng)驗公式預測橋梁結(jié)構(gòu)的溫度。ABID等[10]利用逐步回歸法擬合得到由環(huán)境參數(shù)(太陽輻射量、大氣溫度及風速)預測試驗梁溫差的經(jīng)驗公式，擬合度R2=0.97。劉學毅等[2]給出了一種快速計算高速鐵路道床板溫度的公式，簡化了實際中道床板溫度的預測方法，但沒有涉及箱梁溫度場的計算，無法形成一個整體的箱梁?無砟軌道結(jié)構(gòu)溫度預測方法。在既有的高速鐵路橋梁溫度場研究中，解析法具有較高精度[11]，易于分析環(huán)境變量對結(jié)構(gòu)溫度的影響。為此，本文作者結(jié)合實測溫度建立第三類邊界的箱梁?無砟軌道板溫度場的數(shù)學模型，采用積分變換法推導結(jié)構(gòu)導熱微分方程的級數(shù)解，并與結(jié)構(gòu)的實測結(jié)果和數(shù)值結(jié)果進行對比驗證。同時，參考相關(guān)規(guī)范提出箱梁?無砟軌道結(jié)構(gòu)溫度場的四折線及指數(shù)分布曲線，并給出箱梁?無砟軌道結(jié)構(gòu)升溫模式的計算方法。

1 結(jié)構(gòu)熱交換分析

1.1 結(jié)構(gòu)與環(huán)境的熱交換

橋梁結(jié)構(gòu)與外界環(huán)境的熱交換可分為太陽輻射、對流換熱及輻射換熱，具體熱交換方式如圖1所示。

圖1 箱梁-無砟軌道結(jié)構(gòu)的熱交換示意圖Fig.1 Schematic diagram of heat exchange of box girder-ballastless track structure

以箱梁?軌道結(jié)構(gòu)的水平面熱交換為例，太陽輻射到結(jié)構(gòu)表面上，以輻射平衡等式得到結(jié)構(gòu)表面的第二類邊界條件：

式中：γ為混凝土表面的太陽輻射吸收系數(shù)，本文取0.5[12]；qs為太陽輻射強度，W/m2；v為風速，m/s；ε為混凝土表面的黑度系數(shù)，取0.88[12]；C0為Stefan-Boltzmann常數(shù)，其值為5.67×10?8W/(m2?K4)；Tu和Ta分別為混凝土表面溫度、大氣溫度(℃)；hc為對流換熱系數(shù)，W/(m2?K)，其取值與風速有關(guān)[12?13]；hr為輻射換熱系數(shù)，W/(m2?K)，可根據(jù)結(jié)構(gòu)表面的輻射換熱公式(式(3))計算。

令h=hc+hr，代入式(1)可得

將大氣和太陽輻射與結(jié)構(gòu)的熱交換作用看作等效輻射氣溫[13?14]，用參數(shù)Te表示。

可得到結(jié)構(gòu)表面的第三類邊界條件：

另外，當無太陽輻射情形或結(jié)構(gòu)處于遮陰面時，qs取0 W/m2，則Te=Ta。

1.2 結(jié)構(gòu)熱阻網(wǎng)絡圖

為更好地理解輻射氣溫的換熱作用及物理含義，采用熱阻網(wǎng)絡圖反映結(jié)構(gòu)的綜合對流換熱過程，如圖2所示。根據(jù)熱阻的定義，可得到結(jié)構(gòu)外表面的熱阻Qin表達式(見式(7))。當結(jié)構(gòu)表面的面積A取單元面積，熱源Qst取0 W/m2時，式(7)與對流換熱式(6)相一致。

圖2 結(jié)構(gòu)熱交換的熱阻網(wǎng)絡圖Fig.2 Thermal resistance network of structural heat transfer

式中：Rconv和Rrad分別為結(jié)構(gòu)表面的對流熱阻和輻射熱阻，K/W，可分別用1/(Ahc)和1/(Ahr)計算得到。

同理，結(jié)構(gòu)內(nèi)表面熱交換量Qout為

式中：Td為箱梁內(nèi)表面的待求溫度，℃。對箱梁?無砟軌道結(jié)構(gòu)進行離散，單元熱節(jié)點的導熱量表達式為

式中：Rcond為單元熱節(jié)點Tu1和Tu的導熱熱阻，表達式為δ/(Ak)，與相鄰的熱節(jié)點間的距離δ呈正比。

建立結(jié)構(gòu)熱節(jié)點的能量平衡方程組進行差分求解，求得結(jié)構(gòu)體系各構(gòu)件溫度。

2 結(jié)構(gòu)溫度場模型及求解

根據(jù)熱流邊界條件轉(zhuǎn)換方法，建立對流邊界條件的結(jié)構(gòu)熱傳導微分方程。首先，假定混凝土結(jié)構(gòu)的熱物性參數(shù)不隨時間變化，將結(jié)構(gòu)溫度場簡化為一維瞬態(tài)、無內(nèi)熱源的邊值問題(BVP)。其次，確定流體溫度的函數(shù)表達式。太陽輻射下結(jié)構(gòu)表面的流體溫度呈周期性變化，可用簡諧波描述：

式中：i為結(jié)構(gòu)邊界的個數(shù)；Ti,a和Ti,b分別為第i個邊界處太陽輻射時間內(nèi)流體介質(zhì)的平均溫度和幅值；wi和δi分別為第i個邊界處函數(shù)的頻率及相位；t1和t2分別為日出和日落時間。

結(jié)構(gòu)邊界條件確定后，需要指定結(jié)構(gòu)的初始條件。根據(jù)本文的實測結(jié)果及文獻[15]的研究結(jié)果得到：夏季每天6:00 左右，結(jié)構(gòu)溫度場分布較均勻，因此，以t1時結(jié)構(gòu)均勻溫度為溫度場模型的初始溫度T0。

高速鐵路橋梁結(jié)構(gòu)為有限截面，則一維混凝土結(jié)構(gòu)的導熱微分方程為

式中：α為混凝土材料的熱擴散系數(shù)；d為一維混凝土結(jié)構(gòu)的截面高度。

為使結(jié)構(gòu)溫度場的初始擾動為零，引入溫度變量θ(x,t)：

將式(10)代入導熱方程(11)，得到結(jié)構(gòu)溫度的控制方程為

則結(jié)構(gòu)上表面(x=0)的對流邊界條件為

結(jié)構(gòu)下表面(x=d)的對流邊界條件為

結(jié)構(gòu)的初始條件為

式中：k為混凝土導熱系數(shù)；h1和h2分別為結(jié)構(gòu)外表面、內(nèi)表面的綜合對流換熱系數(shù)；fi′(t)為結(jié)構(gòu)表面流體溫度式(式(10))的變形式。

為簡化求解過程，對式(10)進行變換：

式中：ΔTf,i=Ti,a-T0；Tf,i(t)=-Ti,bcos(wi(t-δi))。

通過引用溫度變量對原導熱方程進行轉(zhuǎn)換，得到的控制方程可采用積分變換法進行求解。由雙對流邊界條件式(14)和(15)可以確定本征函數(shù)式X(βn,x)和本征值βn的超越方程為[16]：

采用分離變量法對控制方程進行變量分離，可得到其解的正、逆函數(shù)變換式[16]：

式(20)為逆變公式，式(21)為積分變換式。

將式(21)代入逆變式(20)，可得到控制方程的一般解[16]：

將式(12)代入式(23)，可得原方程的級數(shù)解為

式(24)中，T(x,t)為一維混凝土溫度場模型的疊加表達式，由于大體積或截面過高的結(jié)構(gòu)初始溫度分布不均勻，會產(chǎn)生一定的計算誤差；當初始溫度取值確定后，結(jié)構(gòu)溫度場的計算精度與級數(shù)展開項數(shù)n有關(guān)，一般項數(shù)越多，精度越高[16]。僅需要確定結(jié)構(gòu)邊界的流體溫度表達式，即可計算箱梁?無砟軌道結(jié)構(gòu)中任一構(gòu)件的溫度。

3 模型驗證

在江西境內(nèi)的1座32 m高速鐵路箱梁?無砟軌道板中預埋多個溫度傳感器，并連續(xù)進行4 a 觀測。采用預埋BGK-3700 溫度計作為溫度測量元件，利用分布式網(wǎng)絡測量系統(tǒng)，每隔30 min 進行實時采集，溫度測點布置如圖3所示。

圖3 結(jié)構(gòu)截面的溫度傳感器分布Fig.3 Distribution of temperature sensors of structural section

3.1 定解條件

根據(jù)現(xiàn)場實測結(jié)果統(tǒng)計得到無砟軌道板的豎向溫差在夏至日(7月23日)前后達到最大值，約17 ℃，因此，本文僅對夏至前后幾天的實測溫度進行理論建模及分析。利用2014?07?19 至2014?07?25 連續(xù)升溫的實測溫度建立箱梁?無砟軌道結(jié)構(gòu)和箱梁底板的溫度場模型。

為避免結(jié)構(gòu)的初始溫度對溫度場計算產(chǎn)生誤差，選取連續(xù)升溫的第3天(2014?07?22)的實測溫度確定溫度場模型的邊界條件及初始條件，如圖4所示。2014?07?22，太陽輻射較有規(guī)律地發(fā)生變化，最大輻照度為800 W/m2左右。

圖4 2014年輻射氣溫和太陽輻照度的時程曲線Fig.4 Time history curves of radiation temperature and solar irradiance in 2014

太陽升起時刻為6:00，據(jù)此計算得到箱梁?軌道板和箱梁底板的均勻溫度分別為34 ℃和25 ℃。忽略溫度場分布的微小差異(偏差約2 ℃)，令初始溫度T0分別為34 ℃和25 ℃。

采用Fourier函數(shù)擬合日照時間(6:00—18:00)范圍內(nèi)軌道板表面的輻射氣溫，并轉(zhuǎn)化為余弦形式：

箱梁底板無太陽輻射，日氣溫變化可用以下經(jīng)驗公式[17]表示：

3.2 綜合換熱系數(shù)取值

采用式(29)計算箱梁中不同位置的綜合換熱系數(shù)：

箱梁外表面的風速為1 d 中實測風速的平均值；箱梁內(nèi)的風速取v=0 m/s；hr利用式(3)計算。

文獻[13]研究了綜合換熱系數(shù)對結(jié)構(gòu)預估溫度場的影響，發(fā)現(xiàn)1 d中綜合換熱系數(shù)為平均值時引起的溫度計算誤差不超過6%。2014?07?22的實測平均風速為2 m/s，按照上述方法進行計算可得到日照時長內(nèi)換熱系數(shù)平均值分別如下：軌道板表面的換熱系數(shù)為20.7 W/(m2·K)，箱梁內(nèi)表面的換熱系數(shù)為20.2 W/(m2·K)，箱梁底板處的換熱系數(shù)為11.7 W/(m2·K)。

3.3 方法驗證

將定解條件代入式(24)，得到箱梁?無砟軌道結(jié)構(gòu)中不同位置(d1～d6)和箱梁底板不同位置(d7～d10)的溫度級數(shù)解，并與有限差分法的數(shù)值解進行驗證。結(jié)構(gòu)溫度的級數(shù)解與數(shù)值解如圖5所示，其中，溫度級數(shù)解用d表示，數(shù)值解用d′表示。

圖5(a)所示為2014?07?22 太陽輻射時長內(nèi)結(jié)構(gòu)溫度的級數(shù)解(級數(shù)展開式前8項求和的結(jié)果)和數(shù)值解的時程曲線，兩者溫度時空變化規(guī)律相一致。對比兩者的溫度波動特征發(fā)現(xiàn)：1)在太陽升起時刻(6:00)，不同深度處(d1～d6)的溫度計算結(jié)果均收斂于結(jié)構(gòu)的初始溫度；2)導熱方程的級數(shù)解在非線性邊界條件出現(xiàn)不收斂的情形[14]，結(jié)構(gòu)外表面(d1=0 m)溫度的級數(shù)解與數(shù)值解的誤差較大，最大絕對誤差達到2.3 ℃。

圖5 結(jié)構(gòu)溫度的級數(shù)解與數(shù)值解對比Fig.5 Comparison between series solutions and numerical solutions of structural temperature

圖5(b)所示為2014?07?22 太陽輻射時間內(nèi)箱梁底板溫度的級數(shù)解和數(shù)值解的時程曲線。從圖5(b)可見：由于底板處無太陽輻射，邊界的熱流擾動很小，級數(shù)解和數(shù)值解的最大誤差小于0.5 ℃，兩者的時程曲線幾乎重合，表明太陽輻射是影響混凝土溫度場的主要因素。

經(jīng)對比可知，箱梁?軌道結(jié)構(gòu)溫度場的級數(shù)解與數(shù)值解的最大相對誤差小于5%，表明采用積分變換法求解箱梁?軌道結(jié)構(gòu)溫度場模型是準確、可行的。

3.4 計算結(jié)果與實測結(jié)果對比

圖6所示為2014?07?22 太陽輻射時長內(nèi)軌道板溫度的實測值(T1～T4)和計算值(T′1～T′4)時程曲線。從圖6(a)可見，由于溫度預估模型的初始值高于實測值約2 ℃，致使特征點T′1和T′2的計算升溫速率均比實測溫度高；道床板內(nèi)特征點T′1距軌道板表面0.05 m，其最高溫度發(fā)生的時間保持一致，峰值溫度的誤差在1 ℃左右。從圖6(b)可見：底座板距結(jié)構(gòu)表面0.25 m，其溫度時程曲線具有升溫晚、溫度波動小的特點，底座板內(nèi)特征點(T′3和T′4)的溫度計算值與實測特征點(T3和T4)的溫度平均誤差均為1 ℃左右。

圖6 結(jié)構(gòu)溫度的級數(shù)解與數(shù)值解對比Fig.6 Comparison between measured temperatures and calculation temperatures of ballastless track slab

經(jīng)統(tǒng)計，實測溫度與計算溫度的最大相對誤差為8%，可滿足實際工程中結(jié)構(gòu)的溫度場分析及溫度模式預測要求。

4 箱梁?無砟軌道升溫模式的計算

4.1 軌道板表面溫度的計算方法對比

軌道板中的溫度測點距結(jié)構(gòu)外表面的距離越小，對應的溫度波動越大，如圖6所示。日照時間內(nèi)道床板中的溫度波動幅值均在10 ℃以上，而底座板中溫度測點距軌道板表面0.44 m，溫度的波動幅值僅為0.5 ℃，由此假設：在日照時間內(nèi)無砟軌道板表面的熱擾動范圍僅局限在軌道結(jié)構(gòu)內(nèi)，可按半無限大結(jié)構(gòu)的溫度場模型分析箱梁?無砟軌道結(jié)構(gòu)溫度場。

根據(jù)半無限大固體的一維非穩(wěn)態(tài)溫度場的解析解[16]，同時忽略結(jié)構(gòu)溫度場中瞬態(tài)溫度項的影響。當f1(t)>T(0,t)時，可得到結(jié)構(gòu)表面溫度的簡化計算式：

式中：ta和tb分別為結(jié)構(gòu)表面流體溫度的平均值和振幅；β為流體溫度的相位；φ為結(jié)構(gòu)表面溫度與流體溫度的振幅比值；ψ為結(jié)構(gòu)表面溫度相比流體溫度滯后的相位角。

將流體溫度函數(shù)式(25)的相關(guān)參數(shù)分別代入式(24)和(30)，可求得不同溫度場模型的表面溫度計算值。日照時間內(nèi)軌道結(jié)構(gòu)表面溫度的時程曲線如圖7所示。從圖7可見：由于日出后2～3 h 范圍內(nèi)，結(jié)構(gòu)表面流體溫度小于結(jié)構(gòu)表面溫度(f1(t)

圖7 結(jié)構(gòu)表面溫度的計算結(jié)果對比Fig.7 Comparison of calculation results of structural surface temperature

4.2 箱梁?無砟軌道溫度模式的擬合結(jié)果

新西蘭規(guī)范[18]、美國AASHTO 規(guī)范[19]分別給出了混凝土結(jié)構(gòu)豎向溫度的五次拋物線形式和三折線分布形式。雖然規(guī)范中均考慮了不同鋪裝厚度(5～10 cm)對結(jié)構(gòu)溫差取值的影響，但既有溫度模式對板厚近50 cm的箱梁?無砟軌道結(jié)構(gòu)不適用。文獻[15]通過對箱梁?無砟軌道板溫度場進行長期測試與統(tǒng)計，發(fā)現(xiàn)在夏季高溫時期，混凝土結(jié)構(gòu)出現(xiàn)年最大正溫差，可采用指數(shù)曲線描述結(jié)構(gòu)豎向溫度分布。

參考混凝土結(jié)構(gòu)溫度的三折線分布形式，利用2014年及2015年夏至日(15:00)的箱梁?無砟軌道結(jié)構(gòu)的實測溫度擬合結(jié)構(gòu)最不利正溫差，得到結(jié)構(gòu)豎向溫度的四折線分布形式，如圖8所示。由于箱梁底板無太陽輻射，溫度波動很小，距底板外表面0.15 m范圍內(nèi)的溫差為1 ℃左右，因此，底板溫度分布也采用折線表示。

同時，采用指數(shù)曲線擬合箱梁?無砟軌道結(jié)構(gòu)的溫度計算值，得到指數(shù)溫度分布與實測溫度、折線溫度分布均比較接近，擬合度較高，如圖8所示。利用2組實測數(shù)據(jù)計算得到道床板的溫度梯度分別為80 ℃/m 和70 ℃/m，兩者均小于中國高速鐵路設計規(guī)范值90 ℃/m[1]。

圖8 箱梁?無砟軌道結(jié)構(gòu)豎向溫度的擬合曲線Fig.8 Fitting curves of vertical temperature of box girderballastless track structure

4.3 箱梁?無砟軌道升溫模式的參數(shù)計算

結(jié)合箱梁?無砟軌道結(jié)構(gòu)溫度場的實測溫度與計算溫度，提出2種箱梁軌道結(jié)構(gòu)的最不利溫度分布曲線，如圖9所示(其中，l為結(jié)構(gòu)截面總高度；d為結(jié)構(gòu)溫度作用范圍；ΔT1～ΔT4，ΔTu1，ΔTd1均為結(jié)構(gòu)溫度特征值)。圖9(a)中，四折線分布形式需要4 個參數(shù)(ΔT1，ΔT2，ΔT3，ΔT4)才能確定，每個參數(shù)之間需按指定間隔距離取值；圖9(b)中，指數(shù)分布形式僅需要3個參數(shù)(ΔTu1，ΔTd1和a)，可適用于不同高度的結(jié)構(gòu)截面，其函數(shù)形式為

圖9 箱梁?無砟軌道結(jié)構(gòu)的溫度作用模式Fig.9 Temperature action model of box girder ballastless track structure

式中：ΔTu1和ΔTd1分別為結(jié)構(gòu)截面上表面和下表面的溫度特征值，其值分別為截面上、下表面溫度與截面內(nèi)部最小溫度的差值；a為指數(shù)溫度分布的衰減參數(shù)。

對于箱梁?無砟軌道結(jié)構(gòu)的指數(shù)溫度分布曲線，衰減參數(shù)a可利用結(jié)構(gòu)溫度計算結(jié)果擬合得到。由于溫度級數(shù)解為時間t的變量，可計算出不同時刻的結(jié)構(gòu)溫度場分布；同時，對于不同截面高度，不同工況下的結(jié)構(gòu)溫度模式均可計算。相比于傳統(tǒng)溫度模式的計算方法[2?3]，利用結(jié)構(gòu)溫度級數(shù)解可以得到不同環(huán)境參數(shù)、材料參數(shù)或截面高度的結(jié)構(gòu)溫度場分布。

表1所示為不同溫度分布曲線的參數(shù)計算結(jié)果。由表1可知：1) 不同環(huán)境下結(jié)構(gòu)表面的溫度ΔT1取值不同；2)結(jié)構(gòu)內(nèi)的溫度ΔT2～ΔT4幾乎保持一致。在沒有結(jié)構(gòu)實測溫度時，利用結(jié)構(gòu)表面溫度的級數(shù)解(式(24))或近似解(式(30))計算表面溫度，ΔT2～ΔT4可按照表1取值。

表1 不同溫度模式的計算結(jié)果對比Table 1 Comparison of calculation results of different temperature models

5 結(jié)論

1)結(jié)合實橋溫度與氣象參數(shù)建立了第三類邊界條件的箱梁?無砟軌道結(jié)構(gòu)溫度場的預估模型，采用積分變換法求得了太陽輻射時長內(nèi)結(jié)構(gòu)溫度場模型的級數(shù)解。

2)箱梁?無砟軌道溫度場的級數(shù)解可由級數(shù)展開式前8項求和近似得到，與數(shù)值解的最大相對誤差為5%，與實測數(shù)據(jù)的最大誤差為8%?；炷两Y(jié)構(gòu)溫度場的級數(shù)解可適用于不同截面高度、不同時刻的結(jié)構(gòu)溫度預測。

3)當?shù)刃л椛錃鉁卮笥诨炷两Y(jié)構(gòu)表面溫度時，有限域結(jié)構(gòu)的表面溫度級數(shù)解與半無限域結(jié)構(gòu)的表面溫度近似解的計算誤差在1 ℃之內(nèi)，可采用結(jié)構(gòu)溫度的近似解快速計算結(jié)構(gòu)表面溫度。

4)高速鐵路箱梁?無砟軌道結(jié)構(gòu)的四折線溫度分布曲線和指數(shù)溫度分布曲線與實測溫度分布曲線擬合度高。