瓦斯爆炸在封閉管道內(nèi)沖擊振蕩特征的數(shù)值模擬

朱傳杰，林柏泉，江丙友，劉 謙

(1.中國(guó)礦業(yè)大學(xué) 安全工程學(xué)院，徐州 221116;2.煤炭資源與安全開采國(guó)家重點(diǎn)實(shí)驗(yàn)室，徐州 221008)

瓦斯爆炸在封閉管道內(nèi)沖擊振蕩特征的數(shù)值模擬

朱傳杰1，2，林柏泉1，2，江丙友1，2，劉 謙1，2

(1.中國(guó)礦業(yè)大學(xué) 安全工程學(xué)院，徐州 221116;2.煤炭資源與安全開采國(guó)家重點(diǎn)實(shí)驗(yàn)室，徐州 221008)

為了研究爆炸波在封閉型系統(tǒng)的沖擊和振蕩特征及其特征參數(shù)變化規(guī)律，采用數(shù)值模擬的方法研究了封閉型管道內(nèi)瓦斯爆炸的傳播特征。研究結(jié)果表明：閉口型系統(tǒng)內(nèi)的瓦斯爆炸呈明顯的振蕩特征，對(duì)于爆燃波，反射波有2道，即前驅(qū)沖擊波的振蕩和壓縮波的振蕩。由于沖擊波的振蕩疊加，使其最大爆炸超壓和瞬態(tài)流速峰值與開口型系統(tǒng)相比較高，而且在反射波及稀疏波的影響下，爆炸波超壓分為三個(gè)區(qū)。爆炸溫度和動(dòng)壓同樣呈明顯的振蕩特征，使得爆炸高溫環(huán)境維持較長(zhǎng)時(shí)間，爆炸動(dòng)壓在爆炸傳播方向的動(dòng)壓與其他方向相比明顯較大。研究結(jié)果解釋了受限空間內(nèi)爆炸破壞比開放型系統(tǒng)強(qiáng)烈的原因，為今后受限空間內(nèi)爆炸的預(yù)防與控制提供了基礎(chǔ)理論參考。

爆炸;沖擊;振蕩;瓦斯

煤礦受限巷道系統(tǒng)為井下氣流流動(dòng)的開口循環(huán)系統(tǒng)，因此，以往對(duì)瓦斯爆炸傳播特性的研究多采用開口型管道系統(tǒng)［1－5］。實(shí)際上，采空區(qū)瓦斯爆炸近年來也是煤礦井下爆炸事故的一種重要形式。以美國(guó)為例，由于采用房柱式采煤方式，會(huì)形成良好的爆炸發(fā)生和發(fā)展空間，采空區(qū)瓦斯爆炸事故多發(fā)。從1986年至今，共發(fā)生了12起已知的采空區(qū)瓦斯爆炸事故，最早的見于1986年的Roadfork 1號(hào)煤礦(事故原因?yàn)榉忾]采空區(qū)初期，巷道頂板巖石掉落產(chǎn)生火花引燃瓦斯)［6］。上世紀(jì)90年代，又發(fā)生了7起采空區(qū)爆炸事故，主要是由于密封漏風(fēng)引起。2002年Big Ridge Mine Portal No.2 煤礦、2006 年 Sago、Darby、Jones Fork E－3煤礦等4起爆炸事故都發(fā)生在采空區(qū)［7］。

另外，以往從上述各類事故的分析來看，封閉型系統(tǒng)與開放型系統(tǒng)相比，表現(xiàn)出爆炸超壓導(dǎo)致的破壞明顯加劇現(xiàn)象［7］。因此，對(duì)類似于采空區(qū)的密閉條件下的瓦斯爆炸特征參數(shù)的變化規(guī)律進(jìn)行研究，對(duì)了解破壞加劇原因也是非常必要的，特別對(duì)我們國(guó)家來講，這方面的研究至今未能引起重視，可以引用的研究成果并不多。本文以封閉型管道為研究對(duì)象，采用數(shù)值模擬的方法研究閉口型系統(tǒng)內(nèi)爆炸參數(shù)的變化特征，以期解釋這種現(xiàn)象的原因。

1 數(shù)值模型

對(duì)瓦斯爆炸沖擊波傳播過程模擬的數(shù)學(xué)模型，可以模擬理想氣體的受熱膨脹過程，氣體動(dòng)力學(xué)過程可用質(zhì)量、動(dòng)量和能量守恒方程來表示。在笛卡爾坐標(biāo)系下分別為：

其中：ρ為密度;u為速度;p為壓力;e=CvT+mfuHc為比內(nèi)能;k為湍流動(dòng)能;ε為湍流動(dòng)能的耗散率;mfu為燃料質(zhì)量分?jǐn)?shù);μt=Cμρk2/ε為湍流粘性系數(shù);Rc為燃燒速率;Cμ為系數(shù)取0.09;CV為定容比熱;T為溫度;Hc為燃燒熱;Γ*=μt/(σ)*為湍流擴(kuò)散系數(shù);(σ)*為輸運(yùn)特性的湍流普朗特?cái)?shù);δij為克羅內(nèi)克算子;i，j為坐標(biāo)方向;C1，C2為常數(shù)。

湍流強(qiáng)度 u'和特征長(zhǎng)度 Lt可由 k、ε 表示為［10］：

湍流燃燒速度與湍流參數(shù)及混合物特性參數(shù)有以下關(guān)系［11］：

其中：St為湍流燃燒速度，m/s;Sl為層流燃燒速度，m/s;v為運(yùn)動(dòng)粘度。

本次研究采用長(zhǎng)度為20 m的閉口型管道(兩端都封閉)，管道截面為8 cm×8 cm，從0.5 m起每隔1 m布置一個(gè)測(cè)點(diǎn)，如圖1所示。模擬采用由美國(guó)世紀(jì)動(dòng)力公司(Century Dynamics)和荷蘭TNO公司聯(lián)合開發(fā)完成的三維計(jì)算流體分析軟件AutoReaGas，數(shù)值計(jì)算的初始參數(shù)設(shè)置如表1。

圖1 數(shù)值模擬采用的管道模型及測(cè)點(diǎn)布置Fig.1 Experimental configuration and measured points

表1 數(shù)值計(jì)算初始計(jì)算參數(shù)Tab.1 Initial values for numerical simulation

2 結(jié)果與分析

2.1 模型驗(yàn)證

數(shù)值模型或軟件的驗(yàn)證對(duì)計(jì)算結(jié)果的準(zhǔn)確性具有較大影響，AutoReaGas已被證明是模擬氣體爆炸較為有效的數(shù)值手段，如美國(guó)職業(yè)安全與健康局(NIOSH)在Lake Lynn實(shí)驗(yàn)礦井做了大量實(shí)驗(yàn)來驗(yàn)證數(shù)值模擬方法的有效性［12］，模擬得到的爆炸超壓與試驗(yàn)值相差在±47%以內(nèi)。對(duì)于該類氣體爆炸，Lea［13］分析認(rèn)為該試驗(yàn)誤差完全可以滿足實(shí)際工程需求。為了保證數(shù)值計(jì)算結(jié)果的準(zhǔn)確性，我們也進(jìn)行了部分驗(yàn)證工作，由于參數(shù)設(shè)置前面已參照前人研究成果選取，因此主要利用實(shí)驗(yàn)數(shù)據(jù)與數(shù)值計(jì)算結(jié)果的對(duì)比來驗(yàn)證網(wǎng)格劃分的合理性。在對(duì)比參數(shù)的選取上，采用實(shí)驗(yàn)方法易于獲得的爆炸超壓值［14］，實(shí)驗(yàn)管道與數(shù)值模擬采用的管道是一致的，實(shí)驗(yàn)主要采集了2.5 m、10.5 m和19.5 m(分別位于管道的前部、中部和后部)處的超壓值，網(wǎng)格劃分采用邊長(zhǎng)分別為2 cm和4 cm的立方體形單元格，對(duì)應(yīng)的網(wǎng)格數(shù)量分別為4×4×1 000和2×2×500。兩種不同網(wǎng)格得到的數(shù)值結(jié)果見表2，可以看到當(dāng)網(wǎng)格數(shù)量為4×4×1 000時(shí)，數(shù)值結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好，因此，數(shù)值模擬采用該種網(wǎng)格劃分方法。

表2 網(wǎng)格劃分對(duì)數(shù)值結(jié)果的影響Tab.2 Effect of grid size on the numerical results

2.2 爆炸波前流場(chǎng)和超壓波形及振蕩特征

閉口型管道內(nèi)的爆炸流場(chǎng)和超壓波形變化要比開口型管道復(fù)雜得多，這一點(diǎn)可以從圖2顯示的波形特征看出，圖中顯示為爆炸超壓和流場(chǎng)的波形呈明顯振蕩特征(“鋸齒”狀)。這種振蕩現(xiàn)象對(duì)爆炸超壓和流場(chǎng)的發(fā)展影響非常明顯，以1.5m處的情況為例，爆炸超壓波形在初始階段與開口型系統(tǒng)是完全一樣的［15］，即沖擊波后緊跟壓縮波，爆炸流場(chǎng)也呈突躍上升后出現(xiàn)慣性回流現(xiàn)象。但是在稀疏波還未形成時(shí)，來自封閉端的反射波使爆炸超壓上升到很高的水平，同時(shí)產(chǎn)生反方向的氣流流動(dòng)。此后，這道反射波在點(diǎn)火端迅速反射，進(jìn)一步提高爆炸超壓，這使得爆炸超壓達(dá)到很高的破壞水平。

圖2 閉口型管道內(nèi)的爆炸超壓和流速演化特征Fig.2 Evolution of explosion overpressure and gas velocity in an closed end pipe

上面提到的波形呈“鋸齒”狀，實(shí)際上是由沖擊波在管道內(nèi)的往復(fù)反射引起的。從爆炸波形來看，除前驅(qū)沖擊波外，還存在若干次帶有周期性的沖擊波峰值，見圖3(a)(以10.5 m測(cè)點(diǎn)處獲得的數(shù)據(jù)為例)。對(duì)于爆燃波來講，通常存在2個(gè)峰值，即前驅(qū)沖擊波峰值和火焰熱輻射形成的壓縮波峰值。這2個(gè)峰值在管道末端分別形成2道反射波，圖3(a)中虛線箭頭代表前驅(qū)沖擊波的反射波，實(shí)線箭頭代表火焰熱輻射形成的壓縮波的反射波。其中箭頭向下代表從管道末端向點(diǎn)火端傳播的反射波，向上代表點(diǎn)火端向管道末端傳播。由于反射波從性質(zhì)上來講屬于P波，無論方向如何都會(huì)使到達(dá)點(diǎn)的超壓呈正值增大。正是這兩道反射波使原本很弱的爆炸波上升到很高的破壞水平。例如，圖3(a)中前驅(qū)沖擊波的超壓為0.25 MPa，在前驅(qū)沖擊波的第1道反射波作用下上升到0.65 MPa，然后反射波在點(diǎn)火端反射使超壓上升到0.84 MPa，隨后火焰輻射造成的壓縮波向點(diǎn)火端反射傳播使超壓上升到1.4 MPa，此時(shí)，達(dá)到了該點(diǎn)的最大爆炸超壓，其值是前驅(qū)沖擊波的1.4/0.25=5.6倍，可見在反射波的作用下，爆炸破壞性大大加強(qiáng)。

此外，反射波的到達(dá)時(shí)間存在一定的周期性，從而造成了爆炸波的往復(fù)反射振蕩，本研究中兩種類型的反射波的周期大約都為40 ms，如圖2(b)圖3(b)所示。在理想情況下(不考慮壁面散熱和粗糙度等因素)，這種振蕩會(huì)持續(xù)較長(zhǎng)時(shí)間，但實(shí)際實(shí)驗(yàn)過程中這種振蕩可能會(huì)很快衰減。

圖3 閉口型管道內(nèi)振蕩頻譜圖Fig.3 Oscillation pectrogram plot of reflection wave in closed end pipe

2.3 爆炸超壓和波前流速峰值變化特征

圖4給出了沿管道內(nèi)各測(cè)點(diǎn)獲得的最大正向和方向氣流速度以及最大爆炸超壓，氣流速度從整體上來講，基本呈增大趨勢(shì)，而爆炸超壓的變化情況要稍微復(fù)雜。盡管前面分析了振蕩對(duì)爆炸波的影響非常大，從最大爆炸超壓角度來講，圖1爆炸超壓和流速演化特征表明：影響其值大小的為前驅(qū)沖擊波在管道末端和點(diǎn)火端的首次反射波、火焰輻射壓縮波的首次反射波，以及反射波與波后稀疏波的發(fā)展情況。在3次反射波及稀疏波的影響下，爆炸波超壓分為三個(gè)區(qū)：

(1)在靠近點(diǎn)火端的位置，在前驅(qū)沖擊波的反射波及其在點(diǎn)火端的反射作用下，最大爆炸超壓較高，但隨著距離的增長(zhǎng)，由于前驅(qū)沖擊波后的稀疏波的發(fā)展，從點(diǎn)火端傳來的反射波越來越晚，因此，爆炸超壓呈降低趨勢(shì)。

(2)隨著距離向管道末端靠近，火焰輻射造成的壓縮波的反射波(FCW)與點(diǎn)火端的反射波相向運(yùn)動(dòng)的距離越來越短，這樣FCW在點(diǎn)火端反射波的未充分發(fā)展稀疏波的基礎(chǔ)上發(fā)生疊加，從而使爆炸超壓越來越大。

(3)在靠近管道末端的位置，首先在管道FSW的影響下出現(xiàn)第一次峰值增大，而后在FCW的影響下產(chǎn)生二次增大，此時(shí)從點(diǎn)火端的FSW的反射波傳回的較慢，稀疏波得以發(fā)展較長(zhǎng)時(shí)間，因此，其最大超壓值反而呈下降趨勢(shì)。但在管道末端由于FSW和FCW相距極近，稀疏波的作用可以忽略，所以其最大爆炸超壓較大。

圖4 最大流速和爆炸超壓Fig.4 Peak gas velocity and peak overpressure

2.4 燃燒速率和溫度

在閉口型系統(tǒng)內(nèi)，火焰反應(yīng)區(qū)內(nèi)的預(yù)混氣燃燒速率(Rc)與開口型［16］稍有不同，在 17.5 m 前 Rc，max(各測(cè)點(diǎn)最大燃燒速率)呈線性增長(zhǎng)。但在靠近管道末端位置，即18.5 m和19.5 m處，Rc，max呈明顯增大趨勢(shì)，Rc，max(18.5 m)=524.31 kg/s，Rc，max(18.5 m)=889.84 kg/s，而在 17.5 m 處 Rc，max(17.5 m)=26.86 kg/s。如圖5所示。

爆炸溫度與開口型系統(tǒng)［16］相比，存在明顯不同，在開口型管道爆炸溫度在上升到最大值后，很快趨于穩(wěn)定，而圖6顯示的爆炸溫度呈明顯的振蕩特征，即存在持續(xù)高溫特征，說明閉口型巷道內(nèi)溫度的破壞性極大。

圖5 閉口性系統(tǒng)內(nèi)燃燒速率Fig.5 Combustion rate in closed end pipe

圖6 閉口型系統(tǒng)內(nèi)爆炸溫度Fig.6 Explosion temperature in closed end pipe

2.5 爆炸動(dòng)壓

動(dòng)壓完全是由波前氣體的瞬時(shí)流動(dòng)產(chǎn)生的，而且可以用公式pV=1/2ρu2gas進(jìn)行計(jì)算，因此，動(dòng)壓的大小和來源完全受波前氣體的運(yùn)動(dòng)情況控制，其值大小由瞬態(tài)流速?zèng)Q定，兩者呈正比關(guān)系，而且在時(shí)間上存在對(duì)應(yīng)關(guān)系，圖7(a)為X方向的動(dòng)壓和流速對(duì)應(yīng)關(guān)系圖，在首次流速出現(xiàn)峰值的同時(shí)，X方向動(dòng)壓對(duì)應(yīng)出現(xiàn)峰值，這一點(diǎn)與開口型管道一致［16］。與開口型管道內(nèi)動(dòng)壓的不同之處體現(xiàn)在第2次動(dòng)壓峰值的出現(xiàn)，從方向上看，流速出現(xiàn)負(fù)值，說明管道內(nèi)氣體出現(xiàn)回流現(xiàn)象，而且流速呈“沖擊”性上升，說明回流完全由管道末端的反射波引起。因此，動(dòng)壓峰值同時(shí)出現(xiàn)第2個(gè)峰值，其值為負(fù)說明由反向氣體流動(dòng)引起。其后，回流的氣體在點(diǎn)火端產(chǎn)生2次正向反射，流速和動(dòng)壓出現(xiàn)2次正向峰值，此后正負(fù)峰值交替出現(xiàn)，即呈現(xiàn)明顯的“振蕩現(xiàn)象”。

圖7 爆炸動(dòng)壓Fig.7 Dynamic explosion pressure

圖7(b)為在坐標(biāo)系內(nèi)三個(gè)不同方向(X、Y、Z)的動(dòng)壓分布情況。其中，在爆炸傳播方向(X)的動(dòng)壓與其他方向(Y、Z)相比較大，如在X方向的最大動(dòng)壓為1.07 MPa(19.5 m處)，Y和Z方向的最大動(dòng)壓基本相同，數(shù)值同為1.7 kPa。這與開口型管道內(nèi)動(dòng)壓值相比明顯增大。另外，峰值動(dòng)壓的大都出現(xiàn)在沖擊波的振蕩回傳時(shí)，說明了閉口型系統(tǒng)內(nèi)的沖擊波反射振蕩對(duì)爆炸動(dòng)壓的強(qiáng)度有很大影響。

3 結(jié)論

本文采用數(shù)值模擬的方法來模擬CH4爆炸時(shí)各爆炸參數(shù)在封閉管道內(nèi)的振蕩特征，給出了模擬采用的基本數(shù)值模型，通過研究得出如下結(jié)論：

(1)閉口型系統(tǒng)內(nèi)的瓦斯爆炸呈明顯的振蕩特征，對(duì)于爆燃波，振蕩波有2道，即前驅(qū)沖擊波的振蕩和壓縮波的振蕩。由于沖擊波的振蕩疊加，使其超壓和流速較開口型系統(tǒng)內(nèi)的高。

(2)影響爆炸超壓和流速值大小的為前驅(qū)沖擊波在管道末端和點(diǎn)火端的首次反射波、火焰輻射壓縮波的首次反射波，以及反射波與波后稀疏波的發(fā)展情況，在反射波及稀疏波的影響下，爆炸波超壓和流速分布分為三個(gè)區(qū)。

(3)在閉口型系統(tǒng)內(nèi)，火焰反應(yīng)區(qū)內(nèi)的預(yù)混氣燃燒速率(Rc)與開口型稍有不同，在17.5 m前Rc，max(各測(cè)點(diǎn)最大燃燒速率)呈線性增長(zhǎng)。但在靠近管道末端位置，即18.5 m和19.5 m處，Rc，max呈明顯增大趨勢(shì)。爆炸溫度的振蕩使得存在持續(xù)高溫特征，說明閉口型巷道內(nèi)溫度的破壞性極大。

(4)爆炸動(dòng)壓峰值同樣呈現(xiàn)出明顯的振蕩特征。在爆炸傳播方向(X)的動(dòng)壓與其他方向(Y、Z)相比較大。另外，動(dòng)壓的峰值大都出現(xiàn)在沖擊波的振蕩回傳時(shí)，說明了閉口型系統(tǒng)內(nèi)的沖擊波反射振蕩對(duì)爆炸動(dòng)壓的強(qiáng)度有很大影響。

［1］司榮軍.礦井瓦斯煤塵爆炸傳播實(shí)驗(yàn)研究［J］.中國(guó)礦業(yè)，2008，17(12)：81－84.

SI Rong－jun.Experimental research into methane－coal dust explosion propagation［J］.China Mining，2008，17(12)：81－84.

［2］ Weiss E S，Cashdollar K L，Harteis S P，et al.Explosion effects on mine ventilation stoppings［R］.NOSIH，Pittsburgh Research Laboratory，Pittsburgh，PA，2008.

［3］林柏泉，菅從光，張 輝.管道壁面散熱對(duì)瓦斯爆炸傳播特性影響的研究［J］.中國(guó)礦業(yè)大學(xué)學(xué)報(bào)，2009，38(1)：1－4.

LIN Bai－quan， JIAN Cong－guang， ZHANG Hui.The influence of heat dissipated in aconduitwallon the transmission characteristics of gas explosions［J］.Journal of China University of Mining＆ Technology，2009，38(1)：1－4.

［4］Liu Q M，Bai C H，Li X D，et al.Coal dust/air explosions in a large－scale tube［J］.Fuel，2010，89(2)：329 － 335.

［5］朱傳杰，林柏泉，江丙友，等.瓦斯爆炸火焰和沖擊波在并聯(lián)巷網(wǎng)的傳播特征［J］.中國(guó)礦業(yè)大學(xué)學(xué)報(bào)，2011，40(3)：385－389.

ZHU Chuan－jie，LIN Bai－quan，JIANG Bing－you，et al.Characteristics of flame and blast－wave propagation during a gas explosion in paralle tunnels［J］.Journal of China University of Mining＆ Technology，2011，40(3)：385－389.

［6］South J S.Report of investigation(underground coal mine)，Non－injury methane explosion，Roadfork No.1 mine(I.D.15 10753)［R］.US：U.S.Department of Labor，Mine Safety and Health Administration，1987.

［7］Karl Z R.Explosion pressure design criteria for new seals in U.S.coal mines［R］.US Pittsburgh：Department of Health and Human Services，CentersforDiseaseControland Prevention，National Institute for Occupational Safety and Health，Pittsburgh Research Laboratory，2007.

［8］ Bone W A，F(xiàn)raser R P，Wheeler W H.A photograhpic investigation of flame movements in gaseous explosions.VII.The phenomenon of spin in detonations［J］.Philosophical Transactions of the Royal Society，1936，235(746)：29－68.

［9］ Launder B E，Spalding D B.Mathematical models of turbulence［M］.London：Academic Press，1972.

［10］ Favre A. Problems ofhydrodynamics and continuum mechanics［M］.US Philadelphia：SIAM (Society for Industry and Applied Mathematics)，1969.

［11］ Bray K N C.Studies of the turbulent burning velocity［J］.Proceedingsofthe RoyalSociety ofLondon， 1990，431(1882)：315－335.

［12］ Zipf R K，Sapko M J，Brune J F.Explosion pressure design criteria for new seals in U.S.coal mines［R］.USA：Department of Health and Human Services，Centers for Disease Controland Prevention，NationalInstitute for Occupational Safety and Health，2007.

［13］ Lea C J，Ledin H S.A review of the state of the art in gas explosion modeling［R］.U.K.，Buxton：Health ＆ Safety Laboratory，2002.

［14］ Salzano E，Marra F S，Russo G，et al.Numerical simulation of turbulent gas flames in tubes［J］.Journal of Hazardous Materials，2002，95(3)：233 247.

［15］ Bjerketvedt D，Bakke J R，Wingerden K van.Gas explosion handbook ［J］.Journal of Hazardous Materials，1997，52(1)：1－150.

［16］朱傳杰.爆炸沖擊波前流場(chǎng)揚(yáng)塵特征及其多相破壞效應(yīng)［D］.徐州：中國(guó)礦業(yè)大學(xué)，2011.

Numerical simulation on oscillation and shock of gas explosion in a closed end pipe

ZHU Chuan－jie1，2，LIN Bai－quan1，2，JIANG Bing－you1，2，LIU Qian1，2
(1.Faculty of Safety Engineering，China University of Mining＆ Technology，Xuzhou 221116，China;2.State Key Laboratory of Coal Resources and Safe Mining，Xuzhou 221008，China)

In order to study the effect of the blast wave oscillation on explosion parameters，the gas explosion propagation characteristics in a closed end pipe was studied by numerical simulation.The results show that the blast wave oscillates obviously in the pipe.There are two different reflection waves during deflagration.One is the front shock wave's reflection wave，and the other is the compression wave.The peak overpressure and the gas velocity in closed end pipe are higher than those in open end pipe.In addition，the evolution of explosion overpressure along the pipe can be divided into three parts under the influence of the reflection wave and the rarefaction wave.The temperature and dynamic pressure oscillation were also observed.The high temperature environment may maintain in a rather long time.The dynamic pressure in the direction of explosion propagation is higher than in other two directions.The above results explain why the explosion damage is more violent in the limited spaces than in open systems.The results give guidance for explosion prevention and control in the limited spaces.

explosion;shock;oscillation;methane

TB936

A

國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃(973)項(xiàng)目(2011CB201205);中國(guó)礦業(yè)大學(xué)人才引進(jìn)和青年教師啟航計(jì)劃資助項(xiàng)目(2011RC07);中央高?；究蒲袠I(yè)務(wù)費(fèi)專項(xiàng)資金資助(2012QNB01)

2011－01－09 修改稿收到日期：2012－02－29

朱傳杰 男，講師，1982年生