饒 宇,趙 根,吳新霞,李 鵬,胡英國,黎衛(wèi)超
(長江科學(xué)院 水利部巖土力學(xué)與工程重點實驗室,武漢 430010 )
應(yīng)力波入射滑移節(jié)理傳播特性的數(shù)值模擬
饒 宇,趙 根,吳新霞,李 鵬,胡英國,黎衛(wèi)超
(長江科學(xué)院 水利部巖土力學(xué)與工程重點實驗室,武漢 430010 )
基于彈性波理論,建立應(yīng)力波入射滑移節(jié)理的數(shù)值模擬模型,并計算應(yīng)力波的透射系數(shù)、反射系數(shù)和吸收系數(shù),并據(jù)此研究應(yīng)力波垂直入射滑移節(jié)理的傳播特性。研究表明:縱波垂直入射滑移節(jié)理的傳播特性主要受節(jié)理面的抗拉強度影響,當(dāng)縱波產(chǎn)生的拉應(yīng)力大于抗拉強度的瞬間,節(jié)理面產(chǎn)生分離而破壞,透射波瞬間消失,反射波增強,吸收能量增加,發(fā)生全反射;剪切波垂直入射滑移節(jié)理的傳播特性主要受節(jié)理面的切向剛度與抗剪強度、等效切向剛度乘積的比值λ的影響,隨著λ增大,反射波逐漸增強,透射波逐漸減弱,吸收能量先增大后減小,當(dāng)λ足夠大時,應(yīng)力波在節(jié)理面處發(fā)生全反射而不產(chǎn)生透射。
應(yīng)力波傳播;滑移節(jié)理;反射系數(shù);透射系數(shù);吸收系數(shù);能量;數(shù)值模擬
天然巖體的不連續(xù)性,諸如包含的裂隙、軟弱夾層以至節(jié)理、斷層帶等,使得天然巖體內(nèi)的應(yīng)力波傳播更為復(fù)雜。一方面巖體不連續(xù)性導(dǎo)致應(yīng)力波傳播過程中發(fā)生復(fù)雜的透射、反射甚至能量損失;另一方面應(yīng)力波傳播也使得巖體內(nèi)部發(fā)生調(diào)整甚至導(dǎo)致巖體參數(shù)的改變。因此,研究應(yīng)力波在巖體節(jié)理處的傳播規(guī)律對于判定巖體穩(wěn)定性具有重要的工程意義。
目前,國內(nèi)外學(xué)者對應(yīng)力波在節(jié)理巖體內(nèi)傳播進行研究時,主要將節(jié)理分為黏結(jié)型和滑移型。對于黏結(jié)型的節(jié)理面,無論應(yīng)力波以何種角度入射節(jié)理面,都不會引起節(jié)理面處巖體的相對滑移和能量損失,即在任何時間段內(nèi),節(jié)理面的入射波能流恒等于反射波能流加上透射波能流。Pyrak-Nolte等[1]、Li等[2]、楊風(fēng)威等[3]、劉立波等[4]、劉婷婷等[5]在這方面都做了較為深入的研究。
自然狀態(tài)下,巖體內(nèi)節(jié)理多為滑移型節(jié)理。應(yīng)力波傾斜入射在節(jié)理面處產(chǎn)生剪應(yīng)力,當(dāng)剪應(yīng)力大于抗剪強度時,兩側(cè)節(jié)理面將產(chǎn)生相對滑動,從而消耗入射波能量。Li 等[6-7]考慮應(yīng)力波在節(jié)理面處的能量損失,研究了節(jié)理巖體的等效黏彈性介質(zhì)模型;宋林等[8]采用位移不連續(xù)理論,研究了節(jié)理面滑移狀態(tài)下的透射系數(shù)、反射系數(shù)及能量消耗系數(shù);王帥等[9]運用頻域分析方法,研究了應(yīng)力波在層狀節(jié)理中的傳播,并運用離散元程序UDEC進行驗證;盧文波[10]以節(jié)理剛度作為判斷節(jié)理黏結(jié)程度的指標(biāo),來區(qū)分完全黏結(jié)型節(jié)理和滑移型節(jié)理。
應(yīng)力波入射滑移節(jié)理時的傳播較完全黏結(jié)節(jié)理更為復(fù)雜,目前的研究側(cè)重于將巖體節(jié)理視為彈性或半彈性材料來計算應(yīng)力波的透射系數(shù)和反射系數(shù),而對應(yīng)力波入射滑移型節(jié)理的能量損失研究偏少。本文基于彈性波動理論,推導(dǎo)了應(yīng)力波傾斜入射滑移節(jié)理時的透射系數(shù)、反射系數(shù)以及吸收系數(shù)的計算公式,并對節(jié)理不同抗拉強度和抗剪強度參數(shù)下,產(chǎn)生滑移的過程及機理采用數(shù)值模擬的方法進行驗證,并探討透射系數(shù)、反射系數(shù)以及吸收系數(shù)與節(jié)理面參數(shù)的關(guān)系,拓展了滑移型節(jié)理面的動態(tài)響應(yīng)規(guī)律的研究。
2.1 透射系數(shù)、反射系數(shù)及吸收系數(shù)定義
應(yīng)力波在斜入射不連續(xù)界面時將產(chǎn)生波型轉(zhuǎn)換,即產(chǎn)生透射同類型波、透射轉(zhuǎn)換波、反射同類型波和反射轉(zhuǎn)換波。本文采用能量法來定義透射系數(shù)、反射系數(shù)及吸收系數(shù),即取節(jié)理處單位周期內(nèi)通過單位面積的透射波、反射波與入射波的能量比值來定義透射系數(shù)、反射系數(shù),在透射系數(shù)、反射系數(shù)基礎(chǔ)上定義吸收系數(shù)。定義如下:
(1)
(2)
(3)
(4)
(5)
(6)
式中:下標(biāo)c=P,S,表示入射波類型;下標(biāo)k=P,S,表示透射或反射波類型(其中:P表示縱波,S表示剪切波);ETk,ERk,EIc分別表示1個周期內(nèi)透射、反射和入射波通過單位面積的能量;Tkc,Rkc分別為透射系數(shù)、反射系數(shù);A為不連續(xù)界面處的能量吸收系數(shù);σ=ρCv,表示應(yīng)力波作用于節(jié)理界面的應(yīng)力(其中:P波作用為法向應(yīng)力,S波則為剪應(yīng)力);ρ為傳播介質(zhì)密度;C為應(yīng)力波傳播速度;v(vIc,vTk,vRk)為入射波、透射波和反射波的質(zhì)點振動速度;T0為振動周期;t1為應(yīng)力波到達節(jié)理界面后某一時刻(作為計算開始時刻)。
圖1 模型示意圖Fig.1 Schematic diagrams of numerical models
2.2 FLAC3D模型的建立
由于應(yīng)力波(P波或SV波)傾斜入射不連續(xù)結(jié)構(gòu)面時將產(chǎn)生波形轉(zhuǎn)換,為了更好地區(qū)分不同類型的透射波、反射波,模型的尺寸需要足夠大,以通過應(yīng)力波的不同到達時間來對波型進行區(qū)分。如圖1所示,取模型的尺寸為400 m×4 m×400 m,在模型中心位置設(shè)置一條橫向貫穿節(jié)理面,其沿y軸旋轉(zhuǎn)到與xy平面成一定角度,即為入射角。動力加載方式采用的是應(yīng)力時程加載,應(yīng)力波加載在模型的底部,沿z方向傳播;其中S波入射時沿x方向加載,P波入射時沿z方向加載。固定模型的y方向,模型上、下邊界設(shè)置為z向和x向的靜態(tài)邊界條件;模型械右邊界S波入射時設(shè)置z方向的黏性邊界,而P波入射時設(shè)置x方向的黏性邊界。
為更好地區(qū)分透射和反射波形,在模型內(nèi)部共布置4個速度監(jiān)測點,監(jiān)測點A 布置在模型左下角,監(jiān)測點B布置在模型左上角,監(jiān)測點C,D布置在模型內(nèi)部,通過監(jiān)測A,B,C,D的x方向和z方向的速度,以此來計算S波和P波的透射系數(shù)、反射系數(shù)。
2.3 模型參數(shù)的確定
為了更好地研究滑移節(jié)理對應(yīng)力波傳播的影響,考慮在巖體結(jié)構(gòu)內(nèi)設(shè)置一組滑移型節(jié)理面,而傳播介質(zhì)為完整的巖塊,視為彈性材料,巖塊的參數(shù)取值見表1。
表1 巖塊物理力學(xué)參數(shù)
Table 1 Physical and mechanical parameters of rock mass
密度ρ/(kg·m-3)彈性模量En/GPa剪切模量Ea/GPa體積模量K/GPa265075.0530.0250.03
則P波和S波的傳播速度如下:
(7)
(8)
節(jié)理面初始參數(shù)取為無滑移彈性節(jié)理的物理力學(xué)參數(shù),如表2所示。
表2 無滑移彈性節(jié)理的物理力學(xué)參數(shù)
Table 2 Physical and mechanical parameters of non-slip elastic joint
法向剛度EV/GPa切向剛度EM/GPa內(nèi)摩擦角φ/(°)黏聚力c/MPa抗拉強度T/MPa3.873.00402.5100
FLAC3D動力分析中,動荷載的施加方法有加速度曲線、速度曲線、位移曲線以及應(yīng)力曲線4種方法。本文采用入射波幅值為0.1 m/s、頻率為50 Hz的正弦函數(shù)波,由于采用Lysmer和 Kuhlemeyer[11]的黏性吸收邊界,動荷載輸入必須為應(yīng)力曲線。應(yīng)力輸入應(yīng)為2倍的黏滯邊界法向黏滯力或切向黏滯力,則對于入射P波和入射S波等價施加的應(yīng)力函數(shù)如下:
(9)
(10)
式中:σP,σS為入射波等價應(yīng)力函數(shù)(MPa);ρ為傳播介質(zhì)密度;CP,CS分別為P波、S波傳播速度;ω為質(zhì)點振動角頻率;vP,vS為P波、S波的質(zhì)點振動速度函數(shù),vP=vS=0.1sin(ωt)。模型網(wǎng)格尺寸4 m,是入射P波波長的1/29,入射S波波長的1/17,符合網(wǎng)格尺寸要求[11]。
P波入射或者S波入射節(jié)理面時,在節(jié)理面法向和切線方向?qū)⒎謩e產(chǎn)生拉應(yīng)力和切應(yīng)力,當(dāng)拉應(yīng)力超過抗拉強度或剪切力超過抗剪強度時,節(jié)理兩側(cè)的巖體介質(zhì)將產(chǎn)生相對運動,此時即存在能量損失。本文選取節(jié)理面不同黏聚力以及抗拉強度參數(shù),對P波和S波垂直入射滑移節(jié)理時進行數(shù)值模擬,獲得相應(yīng)的傳播特性。
3.1 節(jié)理抗拉強度對應(yīng)力波傳播特性的影響
完整巖石的抗壓性能良好,抗拉性能差,這就是動力荷載容易造成巖質(zhì)邊坡的拉裂破壞的原因。而巖體內(nèi)由于存在不良結(jié)構(gòu)面,使得其抗拉強度更差。當(dāng)拉伸應(yīng)力波入射到節(jié)理面時,在節(jié)理面的法向?qū)a(chǎn)生拉應(yīng)力。當(dāng)法向拉應(yīng)力小于節(jié)理面的法向抗拉強度時,接觸面如同膠合模型,不會產(chǎn)生張開或滑移;當(dāng)應(yīng)力波產(chǎn)生的拉應(yīng)力一部分超過節(jié)理面的法向抗拉強度,接觸面分離產(chǎn)生破壞,從而影響到應(yīng)力波的傳播特性?;乒?jié)理的抗拉強度按表3取值。
表3 滑移節(jié)理抗拉強度選取
Table 3 Tensile strength of slip joint
節(jié)理類型工況抗拉強度/MPa無滑移彈性節(jié)理100.012.0滑移節(jié)理20.530.0
圖2 T=2 MPa時入射、反射和透射波形Fig.2 Incidence, reflection and transmission waveforms when T=2 MPa
拉伸P波入射時,在節(jié)理面處產(chǎn)生應(yīng)力方向與P波振動方向一致,即垂直入射時,應(yīng)力方向為節(jié)理面的外(內(nèi))法線方向,根據(jù)式(9),可得其應(yīng)力大小為
(11)
由式(11)可知應(yīng)力峰值σmax=1.545 MPa。
圖2—圖4分別為抗拉強度為2,0.5,0 MPa時,測點C(即入射波形和反射波形),D(即透射波形)速度時程曲線。
圖3 T=0.5 MPa時入射、反射和透射波形Fig.3 Incidence, reflection and transmission waveforms when T=0.5 MPa
圖4 T=0 MPa時入射、反射和透射波形Fig.4 Incidence, reflection and transmission waveforms when T=0 MPa
由圖2—圖4可知,當(dāng)應(yīng)力波沿底面內(nèi)法向加載時,前半周期在節(jié)理面上產(chǎn)生壓應(yīng)力,此時節(jié)理面抗拉強度大小不影響應(yīng)力波傳播,后半個周期應(yīng)力波在節(jié)理面產(chǎn)生拉應(yīng)力;當(dāng)應(yīng)力波沿底面外法向加載時,前半周期在節(jié)理面上產(chǎn)生拉應(yīng)力。
無論應(yīng)力波沿內(nèi)法向加載還是沿外法向加載,當(dāng)抗拉強度T>應(yīng)力波應(yīng)力峰值σmax時,抗拉強度不影響應(yīng)力波的傳播,應(yīng)力波在節(jié)理面處正常透射和反射;當(dāng)抗拉強度T<應(yīng)力峰值σmax時,在應(yīng)力波產(chǎn)生拉應(yīng)力的半周期內(nèi)(內(nèi)法向加載的后半周期、外法向加載的前半周期),抗拉強度將影響應(yīng)力波傳播。
當(dāng)應(yīng)力波在節(jié)理面處產(chǎn)生拉應(yīng)力σ
3.2 節(jié)理剪切參數(shù)對應(yīng)力波傳播的影響
垂直入射時,剪切波在節(jié)理面上產(chǎn)生切向應(yīng)力。因此,節(jié)理面的抗剪強度控制著剪切波垂直入射節(jié)理面的傳播。本節(jié)主要針對剪切波垂直入射滑移節(jié)理的傳播特性進行探討。定義剪切應(yīng)力波入射時節(jié)理的等效切向剛度Ks為
Ks= ks/(ωzs) 。
(12)
式中:ks為節(jié)理切向剛度;zs為剪切波的波阻抗,zs=ρCS。則入射剪切波在節(jié)理面上產(chǎn)生的最大剪應(yīng)力為τ=ksU,U為應(yīng)力波入射產(chǎn)生的最大位移值。
定義無量綱常數(shù)λ為
(13)
式中τs為節(jié)理面的抗剪強度。
采用摩爾-庫倫強度理論的抗剪強度公式,即
τs=σtanφ+c 。
(14)
式中:σ為節(jié)理面的法向應(yīng)力;φ為內(nèi)摩擦角;c為節(jié)理面的黏結(jié)強度。
為了更好地了解透射系數(shù)、反射系數(shù)和吸收系數(shù)與無量綱常數(shù)λ的關(guān)系,節(jié)理的初始物理力學(xué)參數(shù)可按照表2取值,考慮內(nèi)摩擦角φ=0°,節(jié)理的抗剪強度全部由節(jié)理面的黏結(jié)強度c提供。采用式(10)所示的剪切應(yīng)力波入射,則剪切S波入射產(chǎn)生的剪應(yīng)力為
(15)
由式(15)可知剪應(yīng)力峰值τmax=0.892 MPa。
為了研究節(jié)理面抗剪強度對應(yīng)力波傳播特性的影響,節(jié)理面抗剪強度取值如表4所示,圖5為抗剪強度分別取表4所示值時的入射、透射和反射波形。
表4 滑移節(jié)理抗剪強度選取
Table 4 Shear strength of slip joint
節(jié)理類型工況抗剪強度/MPa無滑移彈性節(jié)理2.5011.00 滑移節(jié)理20.2030.0540.00
圖5 不同抗剪強度下入射、反射和透射波形Fig.5 Incidence, reflection and transmission waveforms at different shear strength
圖5表明,當(dāng)入射波剪切應(yīng)力峰值τmax<節(jié)理面抗剪強度τs時,應(yīng)力波能很好地通過節(jié)理面,產(chǎn)生透射和反射,而無能量的損失,此時節(jié)理面的抗剪強度不影響應(yīng)力波的透射和反射;當(dāng)τmax>τs時,節(jié)理面產(chǎn)生滑移和能量的損失,隨著抗剪強度τs減小,反射波逐漸增強,透射波逐漸減弱。當(dāng)抗剪強度τs減小到0時,節(jié)理面被破壞而如同自由面,應(yīng)力波在節(jié)理面處發(fā)生全反射而不產(chǎn)生透射。
如圖1所示,通過監(jiān)測點C,D的速度,采用式(1)—式(6)來計算剪切應(yīng)力波入射滑移節(jié)理面時透射系數(shù)、反射系數(shù)及吸收系數(shù)。根據(jù)設(shè)置節(jié)理面的不同抗剪強度和切向剛度獲得不同的無量綱常數(shù)λ,研究無量綱常數(shù)λ與剪切應(yīng)力波的透射系數(shù)、反射系數(shù)和吸收系數(shù)的關(guān)系。圖6分別為節(jié)理等效切向剛度Ks=0.5和Ks=5時無量綱常數(shù)λ與剪切應(yīng)力波的透射系數(shù)、反射系數(shù)和吸收系數(shù)的關(guān)系。
圖6 不同節(jié)理等效切向剛度下反射系數(shù)、透射系數(shù)和吸收系數(shù)Fig.6 Reflection, transmission and absorption coefficients at different normalized tangential joint stiffness
由圖6可知,主要影響剪切波在節(jié)理處透反射的因素是無量綱常量λ,垂直入射剪切應(yīng)力波在節(jié)理處的傳播主要與剪切波的波阻抗zs、入射剪切波的角頻率ω、應(yīng)力波入射產(chǎn)生的最大位移值U以及節(jié)理面的抗剪強度τs有關(guān)。
當(dāng)無量綱常量λ≤1時,節(jié)理面上不會產(chǎn)生滑移,從而沒有能量損失,隨無量綱常量λ增大,節(jié)理面的剪應(yīng)力大于抗剪強度,能量吸收系數(shù)先增大后減小,反射系數(shù)逐漸增大至1,透射系數(shù)逐漸減小到零。即節(jié)理切向剛度ks→∞或節(jié)理抗剪強度τs→∞時,應(yīng)力波能很好地透射節(jié)理而反射波很小,節(jié)理近似于完整巖體而對應(yīng)力波傳播無影響;節(jié)理切向剛度ks→0或節(jié)理抗剪強度τs→0時(節(jié)理面被破壞),應(yīng)力波幾乎全部被反射而透射波很小,節(jié)理面如同自由界面而發(fā)生全反射。
本文建立應(yīng)力波入射滑移節(jié)理的FLAC3D模型,研究表明縱波垂直入射滑移節(jié)理的傳播特性主要受節(jié)理面的抗拉強度影響,剪切波垂直入射滑移節(jié)理的傳播特性主要受節(jié)理面的抗剪強度和切向剛度的影響,并獲得如下結(jié)論:
(1) 縱波垂直入射,當(dāng)產(chǎn)生的峰值應(yīng)力小于抗拉強度,應(yīng)力波很好地通過節(jié)理面產(chǎn)生透射和反射,而無能量損失;當(dāng)縱波產(chǎn)生的拉應(yīng)力大于抗拉強度的瞬間,節(jié)理面產(chǎn)生分離而破壞,透射波瞬間消失,反射波增強,吸收能量增加,節(jié)理面如同自由界面而發(fā)生全反射。
(2) 剪切波垂直入射時,當(dāng)剪切波在節(jié)理面產(chǎn)生的剪應(yīng)力小于抗剪強度,應(yīng)力波能很好地通過節(jié)理面,產(chǎn)生透射和反射,而無能量的損失;當(dāng)剪應(yīng)力大于抗剪強度時,節(jié)理面產(chǎn)生滑移,隨著抗剪強度τs減小,反射波逐漸增強,透射波逐漸減弱,當(dāng)抗剪強度為0時,應(yīng)力波被全部反射而無透射波。
(3) 節(jié)理面的透射系數(shù)、反射系數(shù)及吸收系數(shù)主要受無量綱常量λ影響,隨λ增大,反射系數(shù)逐漸增大至1,透射系數(shù)逐漸減小到0,吸收系數(shù)先增大后減小;當(dāng)節(jié)理切向剛度或節(jié)理抗剪強度足夠大時,應(yīng)力波能很好地透射節(jié)理而反射波很小,節(jié)理近似于完整巖體而對應(yīng)力波傳播無影響;節(jié)理切向剛度或節(jié)理抗剪強度足夠小時,應(yīng)力波幾乎全部被反射而透射波很小,節(jié)理面如同自由界面而發(fā)生全反射。
[1] PYRAK-NOLTE L J, MYER L R, COOK N G W. Transmission of Seismic Waves Across Single Natural Fractures[J]. Journal of Geophysical Research, 1990, 95(B6):8617-8638.
[2] LI J, MA G. Analysis of Blast Wave Interaction with a Rock Joint[J]. Rock Mechanics and Rock Engineering, 2010, 43(6): 777-787.
[3] 楊風(fēng)威, 李海波, 李建春, 等. 斜入射線彈性節(jié)理應(yīng)力波傳播特征的數(shù)值模擬[J]. 巖土力學(xué), 2013, 34(3):901-907.
[4] 劉立波, 李建春, 李海波, 等. 應(yīng)力波斜入射黏彈性節(jié)理的傳播規(guī)律[J]. 巖石力學(xué)與工程學(xué)報, 2012, 31(增2): 3593-3598.
[5] 劉婷婷, 李建春, 李海波, 等. 應(yīng)力波通過非線性平行節(jié)理的能量分析[J]. 巖石力學(xué)與工程學(xué)報, 2013, 32(8): 1610-1617.
[6] LI J C, MA G W, ZHAO J. An Equivalent Viscoelastic Model for Rock Mass with Parallel Joints[J]. Journal of Geophysical Research: Solid Earth, 2010, 115(B3).DOI:10.1029/2008JB006241.
[7] 李建春, 李海波, MA G W, 等. 節(jié)理巖體的一維動態(tài)等效連續(xù)介質(zhì)模型的研究[J]. 巖石力學(xué)與工程學(xué)報, 2010, 29(增2): 4063-4067.
[8] 宋 林, 邵珠山, 吳敏哲,等. 考慮剪切滑移效應(yīng)時節(jié)理巖體中剪切S波的傳播特性探析[J]. 西安建筑科技大學(xué)學(xué)報(自然科學(xué)版), 2011, 43(4):494-500.
[9] 王 帥, 盛 謙, 朱澤奇, 等. 頻域分析方法在層狀節(jié)理巖體波動問題中的應(yīng)用[J]. 長江科學(xué)院院報, 2013, 30(4): 62-66.
[10]盧文波. 應(yīng)力波與可滑移巖石界面間的相互作用研究[J]. 巖土力學(xué), 1996, 17(3):70-75.
[11]LYSMER J, KUHLEMEYER R L. Finite Dynamic Model for Infinite Media[J]. Journal of the Engineering Mechanics Division, 1969,95(4):859-877.
(編輯:羅 娟)
Numerical Simulation of Propagation Characteristics ofStress Waves across Slip Joint
RAO Yu, ZHAO Gen, WU Xin-xia, LI Peng, HU Ying-guo, LI Wei-chao
(Key Laboratory of Geotechnical Mechanics and Engineering of Ministry of Water Resources,Yangtze River Scientific Research Institute, Wuhan 430010, China)
Based on the elastic wave theory, the numerical simulation model of stress waves across slipping joint has been built, and transmission, reflection and absorption coefficients are calculated to analyze the propagation characteristics of stress waves. Results show that propagation characteristics of normal incident longitudinal wave are affected by joint’s tensile strength, and if the tensile stress produced on the joint by normal incident longitudinal wave exceeds the tensile strength, the interface will separate, transmitted wave disappears instantaneously, and then total reflection occurs, with the reflected wave and absorbed energy increasing. The propagation characteristics of normal incident shear wave are affected by the ratio of tangential joint stiffness and the product of shear strength and normalized tangential joint stiffness. With the increasing of the ratio, the reflection coefficient increases, transmission coefficient decreases, and absorption coefficient increases at first and then decreases, then ends with total reflection.
stress waves propagation; slip joint; reflection coefficient; transmission coefficient; absorption coefficient; energy; numerical simulation
2015-09-17
國家自然科學(xué)基金青年科學(xué)基金項目(51309026)
饒 宇(1990-),男,湖南長沙人,助理工程師,碩士,主要從事巖土動力學(xué)方面的研究工作,(電話)027-82829793(電子信箱)raoyuray@163.com。
10.11988/ckyyb.20150795
2016,33(12):94-98
TU45
A
1001-5485(2016)12-0094-05