基于最小二乘支持向量機(jī)的結(jié)構(gòu)地震響應(yīng)時(shí)滯控制算法

趙德奇，李春祥，藍(lán)聲寧

(1.上海大學(xué) 土木工程系，上海 200072;2.湖南科技大學(xué) 土木工程學(xué)院，湖南湘潭 411201)

主動(dòng)控制系統(tǒng)受時(shí)滯影響控制效果大大減弱。時(shí)滯影響因素［1］主要有:① 液壓系統(tǒng)或電機(jī)系統(tǒng)執(zhí)行控制動(dòng)作所需時(shí)間;② 結(jié)構(gòu)反應(yīng)從傳感器到控制器傳遞時(shí)間及計(jì)算控制力所耗時(shí)間。因此，如何進(jìn)行時(shí)滯補(bǔ)償以減小時(shí)滯對(duì)主動(dòng)控制效果影響為需探究問(wèn)題。文獻(xiàn)［2-7］提出眾多時(shí)滯補(bǔ)償算法，大致可分三種類(lèi)型:理論型補(bǔ)償、相空間補(bǔ)償及時(shí)域補(bǔ)償。理論型時(shí)滯補(bǔ)償一般在時(shí)滯較小情況下效果較好;時(shí)滯較大時(shí)，需提供較大控制力。實(shí)際控制系統(tǒng)可能出現(xiàn)較大時(shí)滯時(shí)，為獲得與無(wú)時(shí)滯相同控制效果，需大幅提高控制力，既不經(jīng)濟(jì)又不實(shí)用。相空間補(bǔ)償方法前提為準(zhǔn)確獲得延遲狀態(tài)參量與結(jié)構(gòu)主頻率。運(yùn)用不正確的主頻率會(huì)導(dǎo)致結(jié)構(gòu)不穩(wěn)定［6］。針對(duì)地震作用，通過(guò)預(yù)測(cè)系統(tǒng)狀態(tài)向量或反饋增益實(shí)現(xiàn)時(shí)域補(bǔ)償更具有實(shí)用性。因此，本文將線(xiàn)性二次型控制(Linear Quadratic Regulator，LQR)與最小二乘支持向量機(jī)(Least Squares Support Vector Machines，LSSVM)進(jìn)行集成，建立時(shí)滯LSSVMLQR算法。

支持向量機(jī)(Support Vector Machines，SVM)是由Vapnik［8］提出的一種學(xué)習(xí)算法。SVM在分類(lèi)與回歸方面具有優(yōu)勢(shì)，為解決時(shí)滯問(wèn)題提供了新思路。SVM專(zhuān)門(mén)針對(duì)有限樣本情況的學(xué)習(xí)機(jī)器，可實(shí)現(xiàn)結(jié)構(gòu)風(fēng)險(xiǎn)最小化，將實(shí)際問(wèn)題通過(guò)非線(xiàn)性變換轉(zhuǎn)換到高維特征空間，在高維空間中構(gòu)造線(xiàn)性決策從而實(shí)現(xiàn)原空間中的非線(xiàn)性決策函數(shù)，巧妙解決了維數(shù)問(wèn)題［9］;同時(shí)也克服了神經(jīng)網(wǎng)絡(luò)存在的過(guò)學(xué)習(xí)與局部極值點(diǎn)的缺陷，具有小樣本學(xué)習(xí)、全局尋優(yōu)、泛化能力強(qiáng)等特點(diǎn)［10］。然而，在解決大樣本問(wèn)題時(shí)，SVM由于迭代誤差積累而無(wú)法滿(mǎn)足精度要求。為此，Suykens等［11］提出新型的支持向量機(jī)—最小二乘支持向量機(jī) (LSSVM)。LSSVM將支持向量機(jī)中目標(biāo)函數(shù)不等式約束改為等式約束，將求解二次規(guī)劃問(wèn)題轉(zhuǎn)化成求解線(xiàn)性方程組，并將經(jīng)驗(yàn)風(fēng)險(xiǎn)由偏差一次方改為二次方，避免了不敏感損失函數(shù)，大大降低了復(fù)雜度，并解決了其中存在的魯棒性、稀疏性及大規(guī)模運(yùn)算問(wèn)題，運(yùn)算速度高于一般的SVM，在非線(xiàn)性預(yù)測(cè)控制中更具優(yōu)勢(shì)。在本文時(shí)滯LSSVMLQR算法中，運(yùn)用LQR算法計(jì)算結(jié)構(gòu)所需最優(yōu)控制力，并設(shè)為目標(biāo)控制力;采用LSSVM對(duì)其進(jìn)行訓(xùn)練學(xué)習(xí)，對(duì)目標(biāo)控制力進(jìn)行回歸預(yù)測(cè)。取t為結(jié)構(gòu)滯后時(shí)間，采用LSSVM回歸預(yù)測(cè)t時(shí)刻后的目標(biāo)控制力，并實(shí)時(shí)控制驅(qū)動(dòng)力裝置產(chǎn)生的控制力對(duì)結(jié)構(gòu)進(jìn)行及時(shí)控制。

1 最小二乘支持向量機(jī)(LSSVM)原理

對(duì)給定的學(xué)習(xí)樣本集(xi，yi)，xi∈Rn，yi∈R，(i=1，…，n)，xi為n維輸入向量，yi為一維輸出。在高維特征空間中構(gòu)造最優(yōu)線(xiàn)性函數(shù):

式中:φ為輸入空間Rn到高維特征空間的非線(xiàn)性映射;w，b為待求回歸參數(shù)，w為權(quán)向量，b為閥值。

據(jù)結(jié)構(gòu)風(fēng)險(xiǎn)最小化原理，并考慮函數(shù)復(fù)雜度及擬合誤差，回歸問(wèn)題可表示為約束優(yōu)化問(wèn)題，即:

約束條件:

用Lagrange方法求解式(2):

式中:C為懲罰因子，α =(α1，α2，…，αn)為 Lagrange乘子。據(jù)KKT最優(yōu)條件進(jìn)行偏微分處理后所確定的決策函數(shù)為:

式中:α，b為求解線(xiàn)性方程組所得模型參數(shù)。定義核函數(shù)K(x，xi)=〈φ(x)φ(xi)〉，得:

上述算法中核函數(shù)K(x，xi)為高維特征空間內(nèi)積。引入核函數(shù)目的即由原始空間中抽取特征，將原始空間中樣本映射為高維特征空間中的向量［12］，如圖1所示。據(jù)泛函理論，只要滿(mǎn)足Mercer條件的函數(shù)均可作為核函數(shù)，不同核函數(shù)構(gòu)造出不同的支持向量機(jī)。目前，研究較多的有線(xiàn)性核函數(shù)、多項(xiàng)式核函數(shù)、高斯核函數(shù)及sigmoid核函數(shù)等。本文采用徑向基(RBF)核函數(shù):

2 線(xiàn)性二次型最優(yōu)控制(LQR)

2.1 連續(xù)狀態(tài)模型

設(shè)各層安裝主動(dòng)控制裝置的n自由度結(jié)構(gòu)如圖2所示，其運(yùn)動(dòng)方程為:

圖1 輸入空間映射到高維特征空間圖Fig.1 Nonlinear mapping of input data into a high dimensional feature space

圖2 結(jié)構(gòu)主動(dòng)控制系統(tǒng)Fig.2 Active control system for structures

式中:{-1}為元素-1的列向量;為地震地面加速度。式(8)即為狀態(tài)方程:

設(shè)作用在結(jié)構(gòu)上的干擾僅有地震作用，即:

假設(shè)結(jié)構(gòu)全部狀態(tài)均能觀(guān)測(cè)到，定義系統(tǒng)二次型泛函為:

求解式(15)可得狀態(tài)向量Z，代入式(14)可求得最優(yōu)驅(qū)動(dòng)力U。

2.2 離散狀態(tài)模型

干擾荷載僅考慮地震作用時(shí)，結(jié)構(gòu)控制系統(tǒng)離散狀態(tài)方程可表示為:

式中:Z(k)，Z(k+1)分別為結(jié)構(gòu)控制系統(tǒng)第k步(對(duì)應(yīng)時(shí)間為t=k×Δt)與第k+1時(shí)刻的狀態(tài)反應(yīng);U(k)為第k步主動(dòng)控制力(k)為第k步地震輸入;Ad，Bd，Dd分別為離散狀態(tài)系統(tǒng)矩陣、控制力位置矩陣與干擾作用位置矩陣，可通過(guò)A，B，Dg求得。

考慮控制系統(tǒng)存在時(shí)滯時(shí)，需采用狀態(tài)方程求解結(jié)構(gòu)控制系統(tǒng)第k+1的狀態(tài)反應(yīng)，即:

式中:i為時(shí)滯步長(zhǎng)，對(duì)應(yīng)iΔt的時(shí)滯。

將U(k-i)=-GZ(k-i)代入式(17)得:

據(jù)結(jié)構(gòu)初始位移及初始速度，用式(18)逐步迭代可求解結(jié)構(gòu)控制系統(tǒng)的地震反應(yīng)。

式中:Q，R分別為半正定及正定權(quán)矩陣，即結(jié)構(gòu)控制系統(tǒng)中的控制參數(shù)。

在閉環(huán)控制情況下，即tf=∞時(shí)的所謂無(wú)限時(shí)間最優(yōu)控制全狀態(tài)反饋，據(jù)Riccati矩陣代數(shù)方程:

求出矩陣后，可得所需增益矩陣G:

進(jìn)而得最優(yōu)驅(qū)動(dòng)力為:

得結(jié)構(gòu)控制系統(tǒng)狀態(tài)方程為:

3 時(shí)滯LSSVM-LQR控制算法

時(shí)間滯后導(dǎo)致結(jié)構(gòu)主動(dòng)控制不能及時(shí)對(duì)結(jié)構(gòu)響應(yīng)做出有效控制，大大降低了LQR最優(yōu)控制效果，結(jié)構(gòu)穩(wěn)定性變差，甚至產(chǎn)生負(fù)反應(yīng)。為此，將LSSVM引入LQR中，先對(duì)采集的最優(yōu)控制力進(jìn)行LSSVM訓(xùn)練，并回歸預(yù)測(cè)第k步時(shí)最優(yōu)控制力U^(k)，則第k+1步狀態(tài)向量計(jì)算式為:

再由式(14)得:

因此，對(duì)LQR最優(yōu)控制力的學(xué)習(xí)訓(xùn)練及回歸預(yù)測(cè)可轉(zhuǎn)化為對(duì)結(jié)構(gòu)系統(tǒng)狀態(tài)向量的學(xué)習(xí)訓(xùn)練及回歸預(yù)測(cè)。形成一個(gè)閉環(huán)控制，既可保證信息采集方便，亦可保證系統(tǒng)的循環(huán)性。時(shí)滯LQR控制系統(tǒng)與時(shí)滯LSSVM-LQR控制系統(tǒng)模塊見(jiàn)圖3。

圖3 結(jié)構(gòu)時(shí)滯LQR與時(shí)滯LSSVM-LQR控制系統(tǒng)模塊Fig.3 Block diagram of LQR and SVM-LQR control systems with time-delay

LSSVM訓(xùn)練時(shí)，選取訓(xùn)練樣本數(shù)量至關(guān)重要。訓(xùn)練樣本量越大，訓(xùn)練速率越慢。故本文取定量數(shù)據(jù)作為訓(xùn)練樣本，在保持訓(xùn)練效果的同時(shí)提高程序執(zhí)行效率及速度。電腦配置:主板/ASUS P8P67、CPU/I7-2600 K 3.4 GHz、內(nèi)存/3 GB，為考察算法訓(xùn)練效率，分別取50、100、150、200個(gè)數(shù)據(jù)作為樣本進(jìn)行訓(xùn)練，產(chǎn)生的時(shí)滯量見(jiàn)表1。

表1 算法訓(xùn)練效率Tab.1 Training efficiency of the algorithm

表1中數(shù)據(jù)顯示，樣本數(shù)量愈大，算法耗時(shí)越長(zhǎng)。樣本量達(dá)到某一值時(shí)算法預(yù)測(cè)的準(zhǔn)確性已有保證。本文在綜合考察預(yù)測(cè)效果與訓(xùn)練效率基礎(chǔ)上選取100個(gè)數(shù)據(jù)作為訓(xùn)練樣本。由于在初采集數(shù)據(jù)的一段時(shí)間里樣本量較少，LSSVM回歸預(yù)測(cè)效果會(huì)有較大誤差，因此可約定樣本量下限為10。樣本量小于10時(shí)，不對(duì)驅(qū)動(dòng)裝置實(shí)施時(shí)滯補(bǔ)償。LSSVM-LQR控制算法流程見(jiàn)圖4。

圖4 LSSVM-LQR控制算法流程圖Fig.4 Flowchart of the LSSVM-LQR control algorithm

4 LSSVM-LQR控制算法數(shù)值驗(yàn)證與比較

考慮3層框架結(jié)構(gòu)，質(zhì)量mi=2×106(kg)，層間剛度ki=1 ×109(N/m)(i=1，2，3)，其質(zhì)量、剛度矩陣為:

由表2～表5看出，時(shí)滯對(duì)結(jié)構(gòu)控制效果產(chǎn)生不利影響，且隨滯后時(shí)間的增長(zhǎng)，控制效果越差。Taft波與上海人工波1在時(shí)滯達(dá)0.18 s時(shí)，因時(shí)滯而增大的結(jié)構(gòu)響應(yīng)約為60%，若不對(duì)時(shí)滯進(jìn)行控制，必導(dǎo)致最終控制結(jié)果發(fā)散。采用時(shí)滯LSSVM-LQR算法對(duì)結(jié)構(gòu)進(jìn)行控制后，因時(shí)滯產(chǎn)生的不利影響得到明顯改善。另外，由表中不難發(fā)現(xiàn)其規(guī)律:① 時(shí)滯較短時(shí)，時(shí)滯LSSVMLQR算法補(bǔ)償效果較明顯，不僅幾乎完全彌補(bǔ)因時(shí)滯產(chǎn)生增大部分響應(yīng)，甚至可減小到較無(wú)時(shí)滯時(shí)效果更好，因在循環(huán)迭代過(guò)程中好的補(bǔ)償效果得以積累與加強(qiáng);而隨著時(shí)滯加長(zhǎng)，時(shí)滯LSSVM-LQR算法補(bǔ)償效果逐漸下降，因隨著時(shí)滯的加長(zhǎng)，預(yù)測(cè)效果下降，導(dǎo)致補(bǔ)償效果下降，但依然保持較好的補(bǔ)償效果。② 時(shí)滯LSSVM-LQR算法對(duì)峰值位移補(bǔ)償效果較好，而對(duì)峰值加速度及最大控制力補(bǔ)償效果稍差。因時(shí)滯LSSVM-LQR算法基于位移與速度進(jìn)行反饋預(yù)測(cè)控制，即算法重點(diǎn)補(bǔ)償了位移及速度因時(shí)滯產(chǎn)生的不利影響，因此對(duì)位移的補(bǔ)償效果較好，而對(duì)加速度與控制力的補(bǔ)償效果稍差。③ 時(shí)滯較短時(shí)，時(shí)滯LSSVM-LQR算法補(bǔ)償效果的離散性較大。因時(shí)滯較短時(shí)，時(shí)滯對(duì)結(jié)構(gòu)影響較小(即計(jì)算時(shí)滯補(bǔ)償效果時(shí)分母較小)，時(shí)滯補(bǔ)償效果對(duì)LSSVM預(yù)測(cè)數(shù)據(jù)敏感性較大，會(huì)出現(xiàn)最大時(shí)滯補(bǔ)償達(dá)700%而最差的時(shí)滯補(bǔ)償達(dá)-210%情形。④ 時(shí)滯LSSVM-LQR算法對(duì)四種波的補(bǔ)償效果各不相同。主要原因?yàn)?與核函數(shù)形式的選取有一定相關(guān)性，不同波對(duì)同一種核函數(shù)敏感程度不同;核函數(shù)參數(shù)通過(guò)試算方式選取，因時(shí)間限制，所選參數(shù)不一定最優(yōu)化。

圖5 地震波Fig.5 Seismic waves

表2 El Centro波作用下結(jié)構(gòu)控制比較Tab.2 Comparison of structural control under the El Centro wave

表3 Loma Prieta波作用下結(jié)構(gòu)控制比較Tab.3 Comparison of structural control under the Loma Prieta wave

表4 Taft波作用下結(jié)構(gòu)控制比較Tab.4 Comparison of structural control under the Taft wave

表5 上海人工波1作用下結(jié)構(gòu)控制比較Tab.5 Comparison of structural control under the Shanghai artificial wave 1

圖6 四種地震波作用下時(shí)滯一定時(shí)三種控制算法下結(jié)構(gòu)的峰值位移Fig.6 Peak displacements of the structure with three control algorithms with given time-delay under four kinds of waves

圖7 四種地震波作用下時(shí)滯一定時(shí)三種控制算法下結(jié)構(gòu)的峰值加速度Fig.7 Peak accelerations of the structure with three control algorithms with given time-delay under four kinds of waves

圖8 四種地震波作用下時(shí)滯一定時(shí)三種控制算法下結(jié)構(gòu)的最大控制力Fig.8 Maximum control forces of the structure with three control algorithms with given time-delay under four kinds of waves

圖9 El Centro波作用下三種控制算法下結(jié)構(gòu)的頂層響應(yīng)曲線(xiàn)(時(shí)滯0.14 s)Fig.9 Displacement of top storey with three control algorithms with time delay under El Centro wave(time-delay:0.14 s)

圖10 Loma prieta波作用下三種控制算法下結(jié)構(gòu)的頂層響應(yīng)曲線(xiàn)(時(shí)滯0.18 s)Fig.10 Displacement of top storey with three control algorithms with time delay under Loma prieta wave(time-delay:0.18 s)

圖11 Taft波作用下三種控制算法下結(jié)構(gòu)的頂層響應(yīng)曲線(xiàn) (時(shí)滯0.18 s)Fig.11 Displacement of top storey with three control algorithms with time delay under Taft wave(time-delay:0.18 s)

圖12 上海人工波1波作用下三種控制算法下結(jié)構(gòu)的頂層響應(yīng)曲線(xiàn) (時(shí)滯0.18 s)Fig.12 Displacement of top storey with three control algorithms with time delay under Shanghai artificial wave 1(time-delay:0.18 s)

圖6～圖8為四種地震波作用下時(shí)滯一定(時(shí)滯最長(zhǎng)情況)時(shí)三種控制算法下結(jié)構(gòu)的峰值位移、峰值加速度、最大控制力。由圖看出，即使時(shí)滯較長(zhǎng)時(shí)，時(shí)滯LSSVM-LQR算法的控制效果依然介于無(wú)時(shí)滯LQR控制與時(shí)滯LQR控制效果之間，且其中對(duì)Loma Prieta與Taft波的控制效果依然保持較好的控制水平。圖9為El Centro波作用下不同控制算法時(shí)(時(shí)滯0.14 s)結(jié)構(gòu)頂層時(shí)程響應(yīng)曲線(xiàn);圖10為L(zhǎng)oma prieta波作用下不同控制算法時(shí)(時(shí)滯0.18 s)結(jié)構(gòu)頂層時(shí)程響應(yīng)曲線(xiàn);圖11為T(mén)aft波作用下不同控制算法時(shí)(時(shí)滯0.18 s)結(jié)構(gòu)頂層時(shí)程響應(yīng)曲線(xiàn);圖12為上海人工波1作用下不同控制算法時(shí)(時(shí)滯0.18 s)結(jié)構(gòu)頂層時(shí)程響應(yīng)曲線(xiàn)。由圖9～圖12知，時(shí)滯LQR算法的控制效果較無(wú)時(shí)滯LQR算法差很多，而時(shí)滯LSSVM-LQR算法控制效果較無(wú)時(shí)滯LQR控制效果稍差些，但較時(shí)滯LQR算法控制效果好很多，尤其在極值點(diǎn)附近，效果更明顯。只有離極值點(diǎn)較遠(yuǎn)的極少點(diǎn)的控制效果較時(shí)滯LQR控制效果差。值得注意的是，LSSVM-LQR使結(jié)構(gòu)響應(yīng)趨于平坦的效果有利于抗震。時(shí)滯較長(zhǎng)時(shí)，由于LSSVM回歸預(yù)測(cè)點(diǎn)離最近觀(guān)測(cè)點(diǎn)時(shí)間間隔太長(zhǎng)，LSSVM預(yù)測(cè)的最優(yōu)控制力與實(shí)際數(shù)據(jù)誤差較大，LSSVM-LQR的控制效果變得不理想，嚴(yán)重時(shí)也會(huì)導(dǎo)致響應(yīng)發(fā)散。

5 結(jié)論

本文通過(guò)建立結(jié)構(gòu)時(shí)滯LSSVM-LQR地震響應(yīng)智能控制算法，基于 MATLAB平臺(tái)，編制時(shí)滯 LSSVMLQR算法的執(zhí)行程序;用3層框架結(jié)構(gòu)進(jìn)行數(shù)值驗(yàn)證與比較。結(jié)論如下:

(1)時(shí)滯對(duì)結(jié)構(gòu)控制效果與穩(wěn)定性具有極不利影響，對(duì)不同地震波敏感程度亦不相同，其敏感度與地震加速度正反變換速率相關(guān)。

(2)時(shí)滯LSSVM-LQR可有效減小時(shí)滯對(duì)結(jié)構(gòu)控制效果影響，但控制效果隨時(shí)滯的過(guò)度加長(zhǎng)而變差。

(3)在時(shí)滯LSSVM-LQR中，徑向基核函數(shù)參數(shù)的選取對(duì)結(jié)構(gòu)地震響應(yīng)控制效果有很大影響。

［1］歐進(jìn)萍.結(jié)構(gòu)振動(dòng)控制-主動(dòng)、半主動(dòng)和智能控制［M］.北京:科學(xué)出版社，2003.

［2］Lin C C，Chang C C，Chen H L.Optimal H∞output feedback control systems with time delay［J］.Journal of Engineering Mechanics，2006，132(10):1096-1105.

［3］蔡國(guó)平，黃金枝.控制存在時(shí)滯的線(xiàn)性系統(tǒng)主動(dòng)控制的滑移模態(tài)方法［J］.力學(xué)季刊，2002，23(2):164-172.

CAI Guo-ping，HUANG Jin-zhi.Sliding-mode control method for linear systems with time-delay in control［J］.Chinese Quarterly of Mechanics，2002，23(2):164-172.

［4］陳龍祥，蔡國(guó)平.地震作用下建筑結(jié)構(gòu)的多時(shí)滯最優(yōu)控制［J］.振動(dòng)與沖擊，2008，27(4):63-65，169.

CHEN Long-xiang，CAI Guo-ping.Optimal control method for scimically excited structures with multiple time delays in control［J］.Journal of Vibration and Shock，2008，27(4):63-65，169.

［5］劉軍龍，代晶輝，張春巍，等.基于位移反饋控制的主動(dòng)增加時(shí)滯補(bǔ)償方法及其試驗(yàn)驗(yàn)證［J］.振動(dòng)與沖擊，2011，30(6):185-191.

LIU Jun-long，DAI Jing-hui，ZHANG Chun-wei，et al.Time delay compensation method based on displacement feedback and active increasing of time delays and its experimental verification［J］.Journal of Vibration and Shock，2011，30(6):185-191.

［6］ McGreevy S，Soong T T，Reinhorn A M.An experimental study of time delay compensation in active structural control［C］.Proc.6th Int.Modal Analysis Conf.and Exhibits，Society for Experimental Mechanics，1988:1733-1739.

［7］徐龍河，李忠獻(xiàn)，錢(qián)稼茹.半主動(dòng)預(yù)測(cè)控制系統(tǒng)的時(shí)滯與補(bǔ)償［J］.工程力學(xué)，2011，28(9):79-83.

XU Long-he，LI Zhong-xian，QIAN Jia-ru.Time delay and compensation for semi-active predictive control system［J］.Engineering Mechanics，2011，28(9):79-83.

［8］Vapnik V.Statistical learning theory［M］.New York:Wiley，1998.

［9］ Cheng J S，Yu D J，Yang Y.Application of support vector regression machines to the processing of end effects of Hilbert-Huang transform［J］.Mechanical Systems and Signal Processing，2007，21(3):1197-1211.

［10］ David M，Bart B，Tony V G.Comprehensible credit scoring models using rule extraction from support vector machines［J］.European Journal of Operational Research，2007，183(3):1466-1476.

［11］ Suykens J A K，Vandewalle J.Least squares support vector machine classifiers［J］.Neural Processing Letters，1999，9:293-300.

［12］顧燕萍，趙文杰，吳占松.最小二乘支持向量機(jī)的算法研究［J］.清華大學(xué)學(xué)報(bào)(自然科學(xué)版)，2010，50(7):1063-1066，1071.

GU Yan-ping， ZHAO Wen-jie， WU Zhan-song. Least squares support vector machine algorithm［J］.J Tsinghua Univ(Sci＆Tech)，2010，50(7):1063-1066，1071.