禹忠 楊劉偉 柯熙政
(1. 西安理工大學(xué),西安 710048;2. 西安郵電大學(xué),西安 710121)
徑向點(diǎn)插值法局部形參的改進(jìn)校準(zhǔn)方法
禹忠1,2楊劉偉2柯熙政2
(1. 西安理工大學(xué),西安 710048;2. 西安郵電大學(xué),西安 710121)
徑向點(diǎn)插值法(Radial Point Interpolation Method,RPIM)中基函數(shù)形參、支持域大小和平均節(jié)點(diǎn)間距等因素直接影響算法的精度與計(jì)算效率,而過(guò)小形參會(huì)引起插值不穩(wěn)定現(xiàn)象. 針對(duì)此問(wèn)題,提出一種基于局部形參校準(zhǔn)法(Local Shape Factor Calibration Method,LSFCM)改進(jìn)的形參優(yōu)化算法,研究RPIM應(yīng)用于電磁場(chǎng)問(wèn)題中插值精度的影響因素,在不同形參、支持域大小和節(jié)點(diǎn)距離時(shí)的全局均方根插值誤差曲線上,根據(jù)插值精度和計(jì)算效率靈活選擇全局形參,簡(jiǎn)化形參設(shè)置與節(jié)點(diǎn)插值計(jì)算過(guò)程,提高計(jì)算效率. 數(shù)值試驗(yàn)結(jié)果驗(yàn)證了所提方法的有效性.
無(wú)網(wǎng)格法;徑向點(diǎn)插值;徑向基函數(shù);支持域;形參校準(zhǔn)
近年來(lái),徑向點(diǎn)插值法(Radial Point Interpolation Method,RPIM)成為時(shí)域電磁場(chǎng)數(shù)值計(jì)算方法的研究熱點(diǎn)之一. 相對(duì)于目前廣泛應(yīng)用的時(shí)域有限差分[1](Finite Difference Time Domain ,FDTD)法和有限元法[2](Finite Element Method,FEM),其不需要生成復(fù)雜的網(wǎng)格結(jié)構(gòu), 可以靈活設(shè)置節(jié)點(diǎn)位置和密度,網(wǎng)格法中的網(wǎng)格重建亦可以通過(guò)自適應(yīng)的節(jié)點(diǎn)分布來(lái)代替,因此在處理復(fù)雜多尺度結(jié)構(gòu)問(wèn)題時(shí)具有很多優(yōu)勢(shì)[3-6].
RPIM是一種高效的數(shù)值計(jì)算方法,其計(jì)算點(diǎn)場(chǎng)值由支持域中其他節(jié)點(diǎn)場(chǎng)值采用徑向基函數(shù)插值近似獲得,因其優(yōu)異的插值特性被引入于電磁場(chǎng)領(lǐng)域求解偏微分方程問(wèn)題,但RPIM中基函數(shù)形參、支持域大小和平均節(jié)點(diǎn)間距等因素對(duì)插值精度影響較大. 傳統(tǒng)RPIM參數(shù)選擇優(yōu)化算法的研究集中在交叉驗(yàn)證(Leave-One-Out Cross Validation,LOOCV)法,但其多適用于全局插值或分區(qū)域插值[7-10].局部形參校準(zhǔn)法(Local Shape Factor Calibration Method,LSFCM)雖可獲取各支持域最優(yōu)形參[11],但每個(gè)支持域設(shè)置不同形參較為繁瑣.
本文中,采用改進(jìn)的LSFCM設(shè)置校準(zhǔn)函數(shù),RPIM計(jì)算支持域插值,結(jié)合均方根誤差(Root Mean Square Error,RMSE)函數(shù),計(jì)算獲得不同情況下全局插值RMSE曲線,確定合適的最大誤差,即可獲得全局最優(yōu)形參或符合精度的形參范圍. 通過(guò)分析RPIM應(yīng)用于二維交錯(cuò)節(jié)點(diǎn)分布時(shí),不同支持域大小、不同節(jié)點(diǎn)間距對(duì)插值精度的影響,測(cè)試各種情況下獲得的最優(yōu)形參范圍,為RPIM的形參選擇提供參考.
計(jì)算域中位置X處場(chǎng)分量u(X)可由其支持域中N個(gè)節(jié)點(diǎn)插值近似[12]為
=rT(X)a+pT(X)b.
(1)
式中:rn(X)為徑向基函數(shù);pm(X)為單項(xiàng)式基函數(shù);an和bm為相應(yīng)插值系數(shù).
實(shí)際應(yīng)用中存在多種徑向基函數(shù),而求解電磁場(chǎng)問(wèn)題時(shí),高斯徑向基函數(shù)具有更好的性能,其為節(jié)點(diǎn)距離的指數(shù)函數(shù),即
(2)
式中:形參αc控制基函數(shù)平滑度及插值精度;dc表示支持域中節(jié)點(diǎn)的平均距離.
低階多項(xiàng)式基函數(shù)項(xiàng)pm(X),可改善局部基函數(shù)的收斂特性,增加插值近似的精度.
(3)
式中,插值系數(shù)矢量a和b由下式計(jì)算:
(4)
式(4)代入式(1),即可獲得RPIM形函數(shù)Φ(X):
=Φ(X)Us.
(5)
形函數(shù)Φ(X)滿足Kronecker delta函數(shù)特性,保證了對(duì)節(jié)點(diǎn)場(chǎng)值的準(zhǔn)確擬合插值. 形函數(shù)沿κ=x,y方向空間導(dǎo)數(shù)可表示為
=?κΦ(X)Us.
(6)
RPIM應(yīng)用于一階麥克斯韋微分方程時(shí),其形式類似于FDTD,時(shí)域微分由中心差分形式離散,而空間導(dǎo)數(shù)則由式(6)插值近似獲得,二維TM波更新方程形式如下:
(7)
2.1徑向基函數(shù)
目前,應(yīng)用于無(wú)網(wǎng)格RPIM中的徑向基函數(shù)有很多種[13]:復(fù)合二次(Multi-Quadrics, MQ)函數(shù)、高斯(Exponential, EXP)函數(shù)、薄板樣條(Thin Plate Spline, TPS)函數(shù)、對(duì)數(shù)型函數(shù)和緊支徑向基函數(shù)(Compactly Supported Radial Basis Function, CSRBF)等. 其中MQ和EXP應(yīng)用最為廣泛,MQ基函數(shù)多用于固體力學(xué)領(lǐng)域,而EXP基函數(shù)更多應(yīng)用于電磁場(chǎng)領(lǐng)域. 這些基函數(shù)均含有部分形參,且形參對(duì)其插值精度影響較大. EXP基函數(shù)中,形參αc控制徑向基函數(shù)的平滑度,如圖1所示[10].
圖1 形參對(duì)EXP基函數(shù)的影響及支持域選擇
研究表明:αc越小,EXP基函數(shù)越平坦,RPIM插值近似精度越高;但隨著αc趨于0,式(4)中求解形函數(shù)的矩量矩陣G條件數(shù)急劇增大(指數(shù)階),從而導(dǎo)致矩陣求逆困難,徑向點(diǎn)插值會(huì)在αc趨于零時(shí)產(chǎn)生嚴(yán)重的不穩(wěn)定現(xiàn)象. 因此,最優(yōu)形參αc的選擇是對(duì)穩(wěn)定性和精度的綜合考慮.
2.2最優(yōu)形參
經(jīng)過(guò)研究,RPIM全局最優(yōu)形參αc通常選擇接近于0的較小值,一般小于1,如取αc=0.01[12],但由于矩陣求逆困難問(wèn)題,此形參值并非對(duì)每一個(gè)支持域都可以獲得較好的插值精度. 因此,Rippa、Kaufmann等[5-6,14-15]采用了LOOCV法進(jìn)行RMSE分析,獲得了全局或分區(qū)域插值最優(yōu)形參. 而Machdo等[11]通過(guò)實(shí)驗(yàn)獲得最優(yōu)形參,分別測(cè)試了αc=0.1、αc=7.4和αc=8.5的百分比誤差,并提出了LSFCM,獲得各支持域的局部最優(yōu)形參,避免矩陣求逆困難問(wèn)題,減小了插值誤差.
相對(duì)于之前的LOOCV法,本文中改進(jìn)的形參校準(zhǔn)方法選取插值誤差曲線極小值點(diǎn)作為最優(yōu)形參,避免了過(guò)小形參引起的矩陣求逆困難問(wèn)題;又不同于局部形參校準(zhǔn)法,采用了全局插值誤差曲線,獲得計(jì)算域統(tǒng)一的最優(yōu)形參,簡(jiǎn)化了形參設(shè)置與支持域插值計(jì)算過(guò)程.
本文采用改進(jìn)的LSFCM設(shè)置校準(zhǔn)函數(shù),通過(guò)RPIM計(jì)算支持域插值,結(jié)合RMSE函數(shù),計(jì)算不同情況下全局插值誤差,確定合適的最大誤差即可獲得最優(yōu)形參值選擇范圍. 并進(jìn)一步分析二維交錯(cuò)節(jié)點(diǎn)分布時(shí),電磁場(chǎng)計(jì)算中不同支持域大小、不同節(jié)點(diǎn)距離對(duì)插值精度的影響,測(cè)試各種情況下獲得的最優(yōu)形參范圍,為RPIM的應(yīng)用提供參考.
3.1改進(jìn)的形參校準(zhǔn)方法
改進(jìn)的LSFCM首先設(shè)置二維計(jì)算域Ω中校準(zhǔn)函數(shù):
C(ui)=sin(Kxi)cos(Kyi).
(8)
采用RPIM局部方案對(duì)計(jì)算域各節(jié)點(diǎn)插值近似,為確定最優(yōu)形參值αc,定義各節(jié)點(diǎn)插值誤差為
ξi(c)=C(c,ui)-C(ui).
(9)
式中,C(c,ui)為形參αc取c時(shí),對(duì)校準(zhǔn)函數(shù)C(ui)計(jì)算域中第i個(gè)節(jié)點(diǎn)的插值,由式(1)獲得. 取計(jì)算域各節(jié)點(diǎn)插值誤差,計(jì)算全局插值誤差:
(10)
獲得插值誤差曲線,形參優(yōu)化方法通過(guò)檢測(cè)誤差函數(shù)極小值點(diǎn),即可獲得全局的最優(yōu)形參值或符合精度要求的形參選擇范圍.其整個(gè)計(jì)算域只需設(shè)置同一個(gè)最優(yōu)形參,可以簡(jiǎn)化形參設(shè)置與節(jié)點(diǎn)插值計(jì)算過(guò)程,提高計(jì)算效率,優(yōu)化RPIM.
3.2數(shù)值分析
本文采用電磁場(chǎng)問(wèn)題交錯(cuò)節(jié)點(diǎn)分布[16-17],設(shè)置校準(zhǔn)函數(shù)中最高頻率fmax=3 GHz,計(jì)算域大小為600×600節(jié)點(diǎn).
3.2.1 支持域半徑對(duì)插值精度影響
首先考慮支持域半徑對(duì)插值精度的影響,支持域分別選擇4節(jié)點(diǎn)(1.0dc)、12節(jié)點(diǎn)(2.0dc)、32節(jié)點(diǎn)(3.0dc)、80節(jié)點(diǎn)(5.0dc),節(jié)點(diǎn)間距選擇λ/20(5 mm),獲得了不同支持域半徑時(shí)插值誤差曲線,如圖2所示.
圖2 支持域半徑對(duì)插值精度影響
分析數(shù)據(jù)可知:形參αc選擇0.01,支持域選取12節(jié)點(diǎn)時(shí),插值誤差為1.26×10-4,而32節(jié)點(diǎn)和80節(jié)點(diǎn)時(shí)均出現(xiàn)不穩(wěn)定現(xiàn)象;形參αc選取0.1時(shí),支持域半徑越大,節(jié)點(diǎn)數(shù)越多,插值誤差越小,80節(jié)點(diǎn)時(shí)可達(dá)到3.31×10-5;檢測(cè)誤差曲線極小值點(diǎn),支持域取12節(jié)點(diǎn),極小值點(diǎn)αc=2.4時(shí),其插值誤差為5.37×10-4,為三者中最小,32及80節(jié)點(diǎn)時(shí),支持域半徑越大,極小值點(diǎn)越靠近零值,且插值誤差越小,分別達(dá)到3.65×10-3(αc=1.6)和1.35×10-4(αc=1.3);從總體分布來(lái)看,支持域半徑越大,誤差函數(shù)在形參取零值到極小值點(diǎn)間的誤差越小,而取大于極小值點(diǎn)時(shí),誤差函數(shù)正好相反. 支持域選取4節(jié)點(diǎn)時(shí),由于支持域半徑過(guò)小,其插值誤差為1.22×10-2,插值精度極低,不建議采用.
綜合以上分析可獲得以下結(jié)論:1) 形參取接近零值時(shí),αc越小,插值精度越高;2) 增大支持域半徑,可提高形參取接近零值時(shí)的插值精度,但αc過(guò)小會(huì)引起插值不穩(wěn)定現(xiàn)象;3) 形參優(yōu)化方法中,形參取誤差函數(shù)極小值點(diǎn),可有效避免過(guò)小形參引起的插值不穩(wěn)定現(xiàn)象,并可達(dá)到近零點(diǎn)的計(jì)算精度;4) 增大支持域半徑,可減小誤差函數(shù)極小值點(diǎn)插值誤差,但對(duì)插值精度的提升不明顯;5) 綜合考慮計(jì)算精度和計(jì)算效率,12節(jié)點(diǎn)(2.0dc)為支持域半徑最優(yōu)選擇,而交錯(cuò)節(jié)點(diǎn)分布下支持域選擇1.5dc~3.0dc時(shí),計(jì)算效率較高,插值誤差較小.
3.2.2 平均節(jié)點(diǎn)間距對(duì)插值精度的影響
首先選取支持域半徑為2.0dc(12節(jié)點(diǎn)),平均節(jié)點(diǎn)間距分別設(shè)置為λ/10(10 mm)、λ/14(7 mm)、λ/20(5 mm)、λ/30(3 mm)、λ/50(2 mm),應(yīng)用形參優(yōu)化方法,分別獲得不同平均節(jié)點(diǎn)間距時(shí)插值誤差曲線,如圖3所示.
圖3 平均節(jié)點(diǎn)間距對(duì)插值精度影響
分析數(shù)據(jù)可知:形參αc選取0.01時(shí),插值誤差最小達(dá)2.07×10-5,最大為3.93×10-3,均可穩(wěn)定插值;而選取極小值點(diǎn)時(shí),不同平均節(jié)點(diǎn)間距的極小值點(diǎn)集中于2.3~2.4,且插值誤差最小為4.69×10-5,最大為9.17×10-4,均可達(dá)到較高精度.總體來(lái)看,平均節(jié)點(diǎn)間距越小,插值誤差越小.
綜合以上分析可得以下結(jié)論:1) 平均節(jié)點(diǎn)間距越小,插值誤差越小;2) 平均節(jié)點(diǎn)間距對(duì)極值點(diǎn)位置影響較小,即對(duì)最優(yōu)形參值影響較小;3) 采用形參優(yōu)化方法獲得極值點(diǎn)最優(yōu)形參可獲得較高的插值精度.
3.2.3 數(shù)值驗(yàn)證
綜合考慮以上結(jié)論,選取支持域半徑為2.0dc(12節(jié)點(diǎn)),平均節(jié)點(diǎn)間距為λ/20(5 mm),其插值誤差為5.37×10-4,具有較高插值精度,對(duì)校準(zhǔn)函數(shù)應(yīng)用RPIM插值近似,獲得了良好效果,如圖4所示.
(a) 校準(zhǔn)函數(shù)
(b) 徑向點(diǎn)插值圖4 校準(zhǔn)函數(shù)驗(yàn)證
為驗(yàn)證形參優(yōu)化方法,參考文獻(xiàn)[4]設(shè)置二維計(jì)算域中高斯點(diǎn)輻射源:
(11)
式中:τ=5.33e-9s;t0=0.8τ. 采用規(guī)則交錯(cuò)節(jié)點(diǎn)分布,支持域半徑為2.0dc(12節(jié)點(diǎn)),節(jié)點(diǎn)距離d=0.01 m,計(jì)算域設(shè)置足夠大(600×600節(jié)點(diǎn))以避免反射誤差,在距離輻射源一倍節(jié)點(diǎn)距離處設(shè)置觀察點(diǎn). 參考解析解形式如下:
(12)
由RPIM和FDTD方法獲得的時(shí)域結(jié)果如圖5所示.與FDTD對(duì)比,RPIM更為接近解析解,具有較好精度.
圖5 觀察點(diǎn)時(shí)域波形
表1中列出RPIM和FDTD的計(jì)算消耗及節(jié)點(diǎn)插值的均方根誤差.對(duì)比數(shù)據(jù)可知,在相同的時(shí)間步進(jìn)和節(jié)點(diǎn)數(shù)時(shí),RPIM內(nèi)存消耗略微高于FDTD,但計(jì)算時(shí)間遠(yuǎn)小于FDTD,且插值誤差亦小于FDTD.
表1 計(jì)算效率
注:文中所有結(jié)果均使用配置為主頻2.93 GHz、Intel i5 CPU、4GB RAM的計(jì)算機(jī)獲得,采用Intel MKL LAPACK進(jìn)行矩陣運(yùn)算和方程組求解.
本文采用LSFCM改進(jìn)的形參優(yōu)化方法,解決RPIM應(yīng)用于求解電磁場(chǎng)問(wèn)題中的最優(yōu)形參選擇問(wèn)題. 主要分析了采用二維交錯(cuò)節(jié)點(diǎn)分布時(shí),支持域半徑和平均節(jié)點(diǎn)間距等因素對(duì)插值精度的影響,并通過(guò)數(shù)值試驗(yàn)驗(yàn)證了形參優(yōu)化方法的有效性.
相對(duì)于目前主要采用的實(shí)驗(yàn)方法或交叉驗(yàn)證法,改進(jìn)的形參校準(zhǔn)方法由RMSE函數(shù)極小值點(diǎn)獲取最優(yōu)形參,更準(zhǔn)確快捷,避免了過(guò)小形參引起的插值不穩(wěn)定現(xiàn)象;又不同于之前的形參校準(zhǔn)法,其計(jì)算域選取統(tǒng)一的最優(yōu)形參,簡(jiǎn)化了誤差計(jì)算與形參設(shè)置,提高了計(jì)算效率.
數(shù)值試驗(yàn)結(jié)果表明,二維交錯(cuò)節(jié)點(diǎn)分布時(shí),支持域半徑選擇1.5dc~3.0dc可獲得較高精度,而支持域半徑選取2.0dc(12節(jié)點(diǎn))為最優(yōu)選擇. 平均節(jié)點(diǎn)間距越小精度越高,可根據(jù)計(jì)算效率和內(nèi)存消耗靈活選擇,建議平均節(jié)點(diǎn)間距選取λ/20,可減少內(nèi)存消耗,同時(shí)可獲得較高的插值精度. 本方法亦可擴(kuò)展到三維領(lǐng)域,而完全不規(guī)則節(jié)點(diǎn)分布的應(yīng)用需要進(jìn)一步研究.
[1] 蘇卓, 譚峻東, 張俊, 等. 基于高階時(shí)域有限差分算法的電磁波傳播計(jì)算[J]. 電波科學(xué)學(xué)報(bào), 2014, 29(3): 431-436.
SU Z, TAN J D, ZHANG J, et al. An electromagnetic wave propagator based on higer-order FDTD method[J]. Chinese journal of radio science, 2014, 29(3): 431-436. (in Chinese)
[2] 彭文峰, 宛汀, 郁美艷. 復(fù)雜媒質(zhì)目標(biāo)電磁散射問(wèn)題的有限元分析[J]. 電波科學(xué)學(xué)報(bào), 2013, 28(1): 76-81.
PENG W F, WAN T, YU M Y. Applying the finite element method to analyze EM scattering of the complex media targets[J]. Chinese journal of radio science, 2013, 28(1): 76-81. (in Chinese)
[3] ZHANG X, CHEN Z D, YU Y. Numerical simulation of coupled transient electrothermal fields with the meshless RPIM[C]//IEEE International Symposium on Antennas and Propagation &USNC/URSI National Radio Science Meeting, 2016: 625-626.
[4] ZHANG X, CHEN Z, YU Y. An unconditional stable meshless ADI-RPIM for simulation of coupled transient electrothermal problems[J]. IEEE journal on multiscale and multiphysics computational techniques, 2016, 1: 98-106.
[5] TANAKA Y, WATANABE S, OKO T. Study of eddy current analysis by a meshless method using RPIM [J]. IEEE transactions on magnetics, 2015, 51(3): 1-4.
[6] ZHU H, GAO C, CHEN H. An unconditionally stable radial point interpolation method based on Crank-Nicolson scheme [J]. IEEE antennas and wireless propagation letters, 2017, 16: 393-395.
[7] RIPPA S. An algorithm for selecting a good value for the parameter c in radial basis function interpolation[J]. Advances in computational mathematics, 1999, 11(2): 193-210.
[8] FASSHAUER G E, ZHANG J G. On choosing “optimal” shape parameters for RBF approximation [J]. Numerical algorithms, 2007, 45(1): 345-368.
[9] WEI YX, XU L, CHEN X Q. The Radial Basis Function shape parameter chosen and its application in engineering [C]//IEEE International Conference on Intelligent Computing and Intelligent Systems. Shanghai, November 20-22, 2009: 79-83.
[10] KAUFMANN T, ENGSTROM C, FUMEAUX C, et al. Eigenvalue analysis and longtime stability of resonant structures for the meshless radial point interpolation method in time domain[J]. IEEE transactions on microwave theory and techniques, 2010, 58(12): 3399-3408.
[11] MACHADO P L, OLIVEIRA R, SOUZA W C, et al. An automatic methodology for obtaining optimum shape factors for the radial point interpolation method [J]. Journal of microwaves, optoelectronics and electromagnetic applications, 2011, 10: 389-401.
[12] WANG J G, LIU G R. A point interpolation meshless method based on radial basis functions [J]. International journal for numerical methods in engineering, 2002, 54(11): 1623-1648.
[13] LIU G R, GU Y T. An introduction to meshfree methods and their programming [M]. Netherlands: Springer, 2005: 276-283.
[14] FASSHAUER G F. Meshfree approximation methods with MATLAB[M]. World Scientific, 2007: 376-382.
[15] KAUFMANN T, ENGSTROM C, FUMEAUX C. Residual-based adaptive refinement for meshless eigenvalue solvers[C]//International Conference on Electromagnetics in Advanced Applications. Sydney, September 20-24, 2010: 244-247.
[16] LEE J F, LEE R, CANGELLARIS A. Time-domain finite-element methods[J]. IEEE transactions on antennas & propagation, 1997, 45(3): 430-442.
[17] BARBER C B, DOBKIN D P, HUHDANPA A H. The quickhull algorithm for convex hulls[J]. ACM transactions on mathematical software, 1996, 22(4): 469-483.
禹忠(1973—),男,陜西人,副教授,主要研究方向?yàn)橛?jì)算電磁學(xué)、天線與微波器件、無(wú)線通信信道、移動(dòng)互聯(lián)網(wǎng)技術(shù).
楊劉偉(1990—),男,陜西人,碩士研究生,研究方向?yàn)橛?jì)算電磁學(xué).
Theoptimizedcalibrationmethodforthelocalshapefactorsofradialpointinterpolationmethod
YUZhong1, 2YANGLiuwei2KEXizheng2
(1.Xi’anUniversityofTechnology,Xi’an710048,China;2.Xi’anUniversityofPosts&Telecommunication,Xi’an710121,China)
In the radial point interpolation method (RPIM), the shape factors such as the parameters of basis function, the support domain size and the average node distance directly affect the interpolation precision and computational efficiency, and the small shape parameter may lead to interpolation instability. In order to optimize parameter selection in RPIM, the algorithm based on the local shape factor calibration method(LSFCM)is proposed in the electromagnetic field problem. The global root mean square (RMS) interpolation error curve with different support domain size, node distance and different parameters can simplify the process of the shape parameter selection and node interpolation. The numerical results show the algorithm has improved the interpolation precision and computational efficiency.
meshless method;radial point interpolation method; radial basis functions; support domain; shape parameter calibration
禹忠, 楊劉偉, 柯熙政.徑向點(diǎn)插值法局部形參的改進(jìn)校準(zhǔn)方法[J]. 電波科學(xué)學(xué)報(bào),2017,32(4):410-415.
10.13443/j.cjors.2017033001
YU Z, YANG L W, KE X Z. The optimized calibration method for the local shape factors of radial point interpolation method[J]. Chinese journal of radio science,2017,32(4):410-415. (in Chinese).DOI: 10.13443/j.cjors.2017033001
O441.4
A
1005-0388(2017)04-0410-06
DOI10.13443/j.cjors.2017033001
2017-03-30
國(guó)家自然科學(xué)基金(No.61377080)
聯(lián)系人: 禹忠 E-mail: zhongyu@mail.xjtu.edu.cn