基于CART-LSSVM的球磨機(jī)料位軟測量方法研究*

張　興，李　偉，閻高偉*，龐宇松

（1.太原理工大學(xué)信息工程學(xué)院，太原030024 2.荷蘭代爾夫特理工大學(xué)機(jī)械海運與材料工程學(xué)院，荷蘭）

基于CART-LSSVM的球磨機(jī)料位軟測量方法研究*

張興1，李偉1，閻高偉1*，龐宇松2

（1.太原理工大學(xué)信息工程學(xué)院，太原030024 2.荷蘭代爾夫特理工大學(xué)機(jī)械海運與材料工程學(xué)院，荷蘭）

球磨機(jī)是用于電力、磨礦和冶金等行業(yè)的高能耗設(shè)備，準(zhǔn)確測量其滾筒料位能夠提高運行效率和安全性能。針對其滾筒內(nèi)料位難以實時檢測，球磨機(jī)的軸承振動信號中存在較多的冗余特征，提出了一種基于分類回歸樹和最小二乘支持向量機(jī)的軟測量方法，首先用Welch法獲得振動信號的功率譜密度，并分割得到基本特征，然后建立分類回歸樹模型，根據(jù)最優(yōu)樹模型的分支節(jié)點進(jìn)行特征選擇，最后利用最小二乘支持向量機(jī)實現(xiàn)特征變量與料位間的非線性映射。通過實驗結(jié)果的對比分析，驗證了該模型的有效性和實用性，以及良好的預(yù)測精度。

球磨機(jī)料位；軟測量；特征選擇；分類回歸樹；最小二乘支持向量機(jī)；振動信號

EEACC:7230doi:10.3969/j.issn.1004-1699.2015.09.016

滾筒式球磨機(jī)是一種用于磨礦及制粉工業(yè)等過程的主要設(shè)備，其工作過程是涉及破碎力學(xué)、機(jī)械振動與噪聲學(xué)等多個學(xué)科的復(fù)雜過程，筒體振動信號存在冗余性和互補(bǔ)性［1］，且由于其封閉旋轉(zhuǎn)運行的工作特性，球磨機(jī)料位參數(shù)與筒體振動信號的映射關(guān)系難以直接有效建立，因此經(jīng)常會發(fā)生空磨、飽磨、堵磨等現(xiàn)象，從而影響生產(chǎn)效率或帶來安全問題。因此需要準(zhǔn)確可靠檢測料位，使球磨機(jī)運行在最佳狀態(tài)。

近年來，科研技術(shù)人員針對該問題提出了很多解決方法，其中以軟測量方法為主要代表，其模型的基本結(jié)構(gòu)為特征降維方法結(jié)合機(jī)器學(xué)習(xí)方法。文獻(xiàn)［2］采用主成分分析［3］方法（Principal Component Analysis，PCA）提取振動功率譜特征，再采用支持向量機(jī)（Support Vector Machines，SVM）建立料位和頻譜特征之間的關(guān)系模型。文獻(xiàn)［4］采用主成分分析（Principle Component Analysis，PCA）方法進(jìn)行振動信號的特征提取，然后采用偏最小二乘回歸方法（Partial Least Square Regression，PLSR）建立料位軟測量模型。文獻(xiàn)［5］采用偏最小二乘法（Partial Least Square，PLS）進(jìn)行振動信號的特征提取，然后采用極限學(xué)習(xí)機(jī)（Extreme Learning Machine，ELM）建立料位軟測量模型。以上三種方法利用振動信號，均可以實現(xiàn)球磨機(jī)料位的軟測量，但是存在降維后特征的維數(shù)難以有效準(zhǔn)確確定的問題，通常需要反復(fù)試驗確定一個較好的特征作為模型下一階段輸入。文獻(xiàn)［6］采用自適應(yīng)遺傳算法—偏最小二乘法（Adaptive Genetic Algorithm-Partial Least Squares，AGA-PLS）選擇頻譜特征，然后融合時域電流信號，基于PLS算法建立融合時頻數(shù)據(jù)特征的球磨機(jī)料位軟測量模型，其特征選擇過程結(jié)合了自適應(yīng)遺傳算法（AGA），時間復(fù)雜度有所上升。

為此，本文結(jié)合分類回歸樹和最小二乘支持向量機(jī)，并將該方法用于基于球磨機(jī)軸承振動信號的料位軟測量建模。實驗結(jié)果表明，該方法避免了復(fù)雜的特征提取及維數(shù)選擇過程，具有較好的測量精度及穩(wěn)定性，為軟測量建模中數(shù)據(jù)的降維提供了一種新思路。

分類回歸樹［7］（Classification and Regression Tree，CART）為一種基于統(tǒng)計分析的非參數(shù)識別技術(shù)，具有多元參數(shù)、非參數(shù)統(tǒng)計等優(yōu)點，利用其較高的分類正確率［8］，CART能自動進(jìn)行特征選擇，有效降低變量維數(shù)。這種特征選擇的方法有以下幾個優(yōu)點:①能簡單高效地處理高維特征；②為一個非參數(shù)過程，不依賴于特定的函數(shù)形式；③樹的分支節(jié)點上的變量或者變量的組合是以最小的錯分率獲得的比較顯著的特征［9］。CART的這些特性使其能夠快速準(zhǔn)確的提取有效特征，解決了樣本特征維數(shù)高、難以抽取的問題，是一種嵌入到學(xué)習(xí)算法中的特征選擇方法。支持向量機(jī)（SVM）是由Vapnik［10］最初在統(tǒng)計學(xué)習(xí)理論的基礎(chǔ)上建立起來的一種機(jī)器學(xué)習(xí)方法。它是近年來發(fā)展起來的一種人工智能技術(shù)，SVM較好的解決了小樣本、非線性、高維數(shù)、局部極小點等實際問題，具有很強(qiáng)的泛化能力［11］。最小二乘支持向量機(jī)（LSSVM）是標(biāo)準(zhǔn)支持向量機(jī)的一種擴(kuò)展，是支持向量機(jī)在二次損失函數(shù)下的一種形式［12］，它與支持向量機(jī)的主要區(qū)別在于采用了等式約束，并將求解優(yōu)化問題轉(zhuǎn)成求解線性方程，避免了求解耗時的受約束二次規(guī)劃QP問題，求解速度相對加快。

1　理論與算法

1.1分類回歸樹

分類回歸樹（Classification and Regression Tree，CART）采用二分遞歸法，利用訓(xùn)練樣本構(gòu)造一棵結(jié)構(gòu)簡潔的二叉樹，是一種有監(jiān)督的學(xué)習(xí)算法，將當(dāng)前訓(xùn)練樣本集分為兩個子樣本集，使得生成的決策樹的每個非葉子結(jié)點都有兩個分支。具體構(gòu)建分類回歸樹Tmax的過程如圖1所示。

圖1　CART樹模型簡圖

所有的數(shù)據(jù)樣本R=｛（Xi，yi）|i=1，2，…，l｝都屬于樹根節(jié)點t，尋找第一個拆分規(guī)則S1（即選擇整棵樹根節(jié)點的分支條件）時，首先從變量Xi的第一維x1開始掃描，計算并記錄樣本數(shù)據(jù)中該變量的每一個取值作為拆分閾值時節(jié)點的不純度函數(shù)下降值Δi（S11，x1），接著掃描變量第二維x2，同樣計算并記錄該變量的各個不純度函數(shù)下降值Δi（S12，x2），直至掃描完變量最后一維xn，計算并記錄所有的拆分閾值對應(yīng)的不純度下降值。最后找出不純度函數(shù)下降值最大的拆分規(guī)則S1*，滿足:

Δi（S1*，xj） =max Δi（S1j，xj） （j=1，2，...，n）（1）

然后找到該拆分規(guī)則對應(yīng)的拆分變量xi和拆分閾值C1，將其定義為樹根節(jié)點的拆分變量和拆分閾值。此時，已經(jīng)將整個樣本數(shù)據(jù)集分到兩個子樹t1和t2上，其中

t1= ｛X∈R : xi≤C1｝ ， t2= ｛X∈R : xi＞C1｝（2）

對于每一個子樹，重復(fù)以上尋找樹根節(jié)點拆分規(guī)則的掃描過程，尋找每個子樹所屬子樹的根節(jié)點的拆分規(guī)則。整棵樹的建立過程就是一個尋找更小子樹根節(jié)點的拆分規(guī)則的過程。當(dāng)節(jié)點滿足以下三種條件之一，即被視作葉節(jié)點停止拆分。

①節(jié)點很小:分支后的葉節(jié)點的樣本數(shù)小于給定值Nmin（一般取默認(rèn)值Nmin=10）；

②純節(jié)點:分支后的葉節(jié)點中的樣本為同一類；

③空屬性向量集:無屬性向量用作分支選擇。

停止拆分后，根據(jù)節(jié)點ti中來自每個類的樣本比例，可以對每個子葉節(jié)點wi（i=1，2，…，n）指定類別標(biāo)簽yi，即本文中的料位。此時就完成了分類回歸樹Tmax的構(gòu)建。樹完全生長過程中會出現(xiàn)過度擬合現(xiàn)象，在得到最大分類回歸樹的基礎(chǔ)上，采用最小代價—復(fù)雜性剪枝算法對生成的樹進(jìn)行剪枝。經(jīng)過剪枝之后，得到最優(yōu)分類回歸樹。

最后選擇出所有分支節(jié)點上的拆分變量，得到模型選擇出來的有效特征為Xcart

Xcart= xi（t）?xj（t1）?…?xm（tn）（3）

xi（t）表示在節(jié)點t選擇出來的第i維變量。

本文的不純度采用基尼不純度（Gini Impurity），表示一個隨機(jī)選中的樣本在子集中被分錯的可能性。當(dāng)一個節(jié)點中所有樣本都是一個類時，基尼不純度為零。假設(shè)i的可能取值為｛1，2，…，m｝，令fi是樣本被賦予i的概率，則基尼不純度可通過如下公式計算:

1.2最小二乘支持向量機(jī)

SVM是由Vapnik提出的一種用于解決分類與函數(shù)逼近問題的機(jī)器學(xué)習(xí)方法，與神經(jīng)網(wǎng)絡(luò)的經(jīng)驗風(fēng)險最小化相比具有更好的泛化能力。LSSVM是在標(biāo)準(zhǔn)SVM上的一種擴(kuò)展，其需要預(yù)先確定的參數(shù)也比標(biāo)準(zhǔn)SVM算法少。

設(shè)訓(xùn)練樣本集為D=｛（xi，yi）|i=1，2，…，l｝，在原始空間中的優(yōu)化問題可以描述為

約束條件為

yi=wTΦ（xi）+b+ei（6）

式中:xi為輸入變量，yi為輸出變量，w為權(quán)值，b為閾值，γ為正則化參數(shù)，ei為誤差，優(yōu)化問題的拉格朗日形式為

式中:λi為拉格朗日系數(shù)，優(yōu)化問題轉(zhuǎn)化為解線性方程

優(yōu)化問題的決策函數(shù)為

在選定核函數(shù)之后，LSSVM需要進(jìn)一步確定其超參數(shù):正則化參數(shù)γ和核參數(shù)σ，目前還沒有確定γ和σ的統(tǒng)一方法。本文采用耦合模擬退火法（Coupled Simulated Annealing，CSA）與單純形法（Simplex Method，SM）相結(jié)合的方法優(yōu)化LSSVM的超參數(shù)，其主要優(yōu)化過程分為兩步，首先，參數(shù)初始化:對于每個超參數(shù)，利用耦合模擬退火法（CSA）計算合適的初始值，其搜索范圍為［exp（-10），exp （10）］。其次，參數(shù)尋優(yōu):根據(jù)得到的初始值，應(yīng)用單純形法進(jìn)一步對超參數(shù)進(jìn)行優(yōu)化。

2　建模及實驗過程

2.1 基于CART-LSSVM軟測量建模

軟測量技術(shù)主要依據(jù)對可測易測過程變量（稱為輔助變量）與難以直接檢測的待測變量（稱為主導(dǎo)變量）的數(shù)學(xué)關(guān)系認(rèn)識，采用各種計算方法實現(xiàn)待測變量的測量或估計。本文以球磨機(jī)軸承振動信號為輔助變量，結(jié)合基于分類回歸樹（CART）的特征選擇方法和最小二乘支持向量機(jī)（LSSVM）建立球磨機(jī)料位的軟測量模型。該模型進(jìn)行軟測量的過程主要分為兩個階段，離線建模階段和在線測量階段，其流程圖如圖2所示。

圖2　CART-LSSVM系統(tǒng)流程圖

具體測量步驟為:①利用振動傳感器采集球磨機(jī)軸承振動信號。②采用Welch法計算振動信號中時域數(shù)據(jù)的功率譜并進(jìn)行頻譜的分割，不同頻率段的功率譜對應(yīng)不同的特征［13］，需進(jìn)一步進(jìn)行頻段的選擇，即特征的選擇。③根據(jù)離線建模階段得到的CART選擇特征變量Xcart，將其做為LSSVM模型輸入。④將Xcart輸入離線建模階段所建立的LSSVM模型，預(yù)測得出實際料位。

2.2實驗過程及數(shù)據(jù)預(yù)處理

采用小型球磨機(jī)進(jìn)行試驗數(shù)據(jù)的采集，振動傳感器安裝在遠(yuǎn)離電動機(jī)一端的軸承座上，分別采集1 L，2 L，…，20 L共20組振動信號，使用Welch法計算振動信號的功率譜密度（Power Spectrum Density，PSD），將每組PSD均分成22個樣本，則樣本總數(shù)為440個，并把每組的22個樣本隨機(jī)的分成15個訓(xùn)練樣本與7個測試樣本，因此得到的訓(xùn)練集有300個樣本，測試集有140個樣本。對PSD進(jìn)行分析后發(fā)現(xiàn)，振動信號的能量主要集中在600 Hz～6 000 Hz之間，因此取這個區(qū)間的頻譜為有效頻譜。綜合考慮頻譜的波動范圍及計算效率，以20 Hz為單位對有效頻譜進(jìn)行分割，并求其均值，所以單個樣本的維數(shù)是270維。

然后，將訓(xùn)練數(shù)據(jù)（300×270）輸入CART建立分類回歸樹模型，根據(jù)建立的最優(yōu)樹模型，抽取出樹的所有分支節(jié)點上的拆分變量Xcart，這些變量就是CART選擇出的特征對應(yīng)的維度。本文中，將270維的數(shù)據(jù)輸入CART后，得到最優(yōu)的回歸樹模型，其在訓(xùn)練集上的RMSE為0.085 8，選擇樹模型中的拆分變量得到基本的26維特征，從這270維數(shù)據(jù)中抽取出這26維，將300×26的訓(xùn)練數(shù)據(jù)輸入LSSVM建立特征變量與料位間的非線性映射。最終，采用與建模數(shù)據(jù)在同一數(shù)據(jù)空間的測試數(shù)據(jù)集（140×270）對CART-LSSVM模型進(jìn)行測量效果驗證。

2.3結(jié)果分析與討論

為了證實本文所建立的球磨機(jī)料位測量模型的有效性，在實驗的同時采用主成分回歸分析（PCR）、主成分分析-最小二乘支持向量機(jī)（PCALSSVM）、核主成分分析-最小二乘支持向量機(jī)（KPCA-LSSVM）、最小冗余最大相關(guān)-最小二乘支持向量機(jī)（MRMR-LSSVM）和RELIEFF-最小二乘支持向量機(jī)（RELIEFF-LSSVM）5種方法作為對照。采用均方根誤差（Root Mean Square Error，RMSE）、平均絕對誤差（Average Absolute Error，AAE）和最大絕對誤差（Maximal Absolute Error，MAE）作為模型性能評價指標(biāo)，其中RMSE為主要的評價指標(biāo)。

式中:yi為第i個樣本的真實值，y^i為第i個樣本的估計值，N為測試集的樣本個數(shù)。實驗中，LSSVM算法的實現(xiàn)采用J.A.K.Suykens等人開發(fā)的LS-SVM-lab1.8Toolbox工具箱。

最小冗余最大相關(guān)（Min-Redundancy Max-Relevance，MRMR）是一種經(jīng)典的基于互信息（Mutual Information，MI）進(jìn)行特征選擇的方法，能選擇出單個表示能力較好的特征，并使這些單個較好的特征之間的冗余度最?。?4］。RELIEFF是一種啟發(fā)式的特征選擇方法，可以很好的處理有噪音、多標(biāo)簽及有缺失的數(shù)據(jù)［15］。

表1　不同測試方法測試誤差比較

表1給出了幾種不同的方法10次預(yù)測誤差的均值。其中，Ncomp表示降維后所得特征的維數(shù)，這里PCR降維后所得預(yù)測精度最好的維數(shù)為28維，PCA-LSSVM為20維，KPCA-LSSVM為10維。從這幾個最好的結(jié)果可以看出，PCR和PCA-LSSVM表現(xiàn)普通，分析其原因，主要是由于PCA是一種“線性”降維技術(shù)，只能處理線性問題，一方面主成分是原始特征的線性組合，另一方面是由于其對原始數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理，使協(xié)方差矩陣變成相關(guān)系數(shù)矩陣，而相關(guān)系數(shù)矩陣僅反映指標(biāo)間的“線性”相關(guān)程度。經(jīng)過KPCA進(jìn)行非線性降維得到的特征優(yōu)于PCA得到的特征，對振動信號特征的處理能力較PCA有所提高，但還是不能有效減少特征中的冗余信息，對此類特征的非線性降維難以達(dá)到理想的結(jié)果。

為了比較MRMR、RELIEFF和CART三種特征選擇方法與前三種方法所建立模型的效果，將文中的六種方法均降到26維做對比試驗，不同測試方法的最優(yōu)預(yù)測曲線如圖3所示，其橫坐標(biāo)表示每個料位下7個測試樣本組成的測試樣本集，縱坐標(biāo)表示球磨機(jī)的實際料位，單位為升。結(jié)合表1可以看出，利用CART得到的特征建模效果較好，該方法有效降低了特征中的冗余信息，自動對特征變量進(jìn)行選擇，并快速將特征降低到較小的維數(shù)，為非線性降維的一種有效方法，優(yōu)于傳統(tǒng)的特征選擇方法MRMR 和RELIEFF。分析其原因，CART利用基尼系數(shù)逐步找出所含信息最大的屬性，并遞歸的建立樹模型，能選擇出特征中最具區(qū)分度的屬性作為樹的分支，這種機(jī)制有效的減少了特征中的冗余信息。從圖中可以看出，CART-LSSVM模型的預(yù)測精度RMSE= 0.538 5，優(yōu)于其他對比方法，該方法能較好的跟蹤球磨機(jī)料位，達(dá)到了比較理想的預(yù)測精度。

通過對比分析可以發(fā)現(xiàn)，CART可以自動將特征變量快速降低到一個較小的維數(shù)，避免了特征提取方法的維數(shù)選擇過程，同時優(yōu)于傳統(tǒng)的特征選擇方法，其選擇出的特征具有一定的穩(wěn)定性，可有效的應(yīng)用于球磨機(jī)料位的軟測量建模。

圖3　不同方法最優(yōu)測試結(jié)果

3　結(jié)論

針對球磨機(jī)振動信號中存在較多的冗余特征，本文從選擇球磨機(jī)軸承振動信號中的有效特征出發(fā)，采用分類回歸樹（CART）建立樹模型，選擇樹干上的頻域特征變量，根據(jù)得到的特征建立最小二乘支持向量回歸模型來進(jìn)行球磨機(jī)料位的軟測量。對比試驗預(yù)測結(jié)果表明，該方法具有的優(yōu)勢有:①降維過程中維數(shù)的確定不需要人為參與，避免了傳統(tǒng)的特征提取及特征選擇過程中繁瑣的維數(shù)調(diào)整過程；②具有良好的預(yù)測精度，且建模過程簡單，易于在實際中操作，算法復(fù)雜度不高。

同時該模型有以下待改進(jìn)的部分，本文中的分類回歸樹，其拆分規(guī)則都只考慮了單個變量，有些情況下需要考慮多變量組合的拆分規(guī)則，因此，基于該特征選擇方法的軟測量建模為需要進(jìn)一步研究的課題。

［1］ 湯健，柴天佑，叢秋梅，等.基于EMD和選擇性集成學(xué)習(xí)算法的磨機(jī)負(fù)荷參數(shù)軟測量［J］.自動化學(xué)報，2014，40（9）:1853-1866.

［2］ Tang Jian，Zhao Lijie，Yu Wen，et al.Soft Sensor Modeling of Ball Mill Load via Principal Component Analysis and Support Vector Machines［J］.Lecture Notes in Electrical Engineering，2010（67）: 803-810.

［3］ 付華，王馨蕊，王志軍，等.基于PCA和PSO-ELM的煤與瓦斯突出軟測量研究［J］.傳感技術(shù)學(xué)報，2014，12:1710-1715.

［4］ 湯健，鄭秀萍，趙立杰，等.基于頻域特征提取與信息融合的磨機(jī)負(fù)荷軟測量［J］.儀器儀表學(xué)報，2010，31（10）:2161-2167.

［5］ Tang Jian，Wang Dianhui，Chai Tianyou.Predicting Mill Load Using Partial Least Squares and Extreme Learning Machines［J］.Soft Computing，2012，16（9）:1585-1594.

［6］ 湯健，柴天佑，趙立杰，等.融合時/頻信息的磨礦過程磨機(jī)負(fù)荷軟測量［J］.控制理論與應(yīng)用，2012，29（5）:564-570.

［7］ Loh W Y.Classification and Regression Trees［J］.Wiley Interdisciplinary Reviews:Data Mining and Knowledge Discovery，2011，1 （1）:14-23.

［8］ Bittencourt H R，Clarke R T.Feature Selection by Using Classification and Regression Trees（CART）［J］.The International Archives of the Photogrammetry，Remote Sensing and Spatial Information Sciences，2004:66-70.

［9］ Gomez-Chova L，Calpe J，Soria E，et al.CART-Based Feature Selection of Hyperspectral Images for Crop Cover Classification［C］// IEEE International Conference on Image Processing，2003，3（2）: 589-592.

［10］Vapnik V N.The Nature of Statistical Learning Theory［J］.Neural Networks IEEE Transactions，1995，10（5）:988-999.

［11］賽吉爾呼，戴盛芳，董愛華，等.基于SVM和RBFN的汽車主動降噪系統(tǒng)傳感器故障診斷［J］.傳感技術(shù)學(xué)報，2014，27（4）: 512-517.

［12］Suykens J A K，Vandewalle J.Least Squares Support Vector Machine Classifiers［J］.Neural Processing Letters，1999，9（3）:293-300.

［13］湯健，趙立杰，岳恒，等.濕式球磨機(jī)筒體振動信號分析及負(fù)荷軟測量［J］.東北大學(xué)學(xué)報:自然科學(xué)版，2010，31（11）:1521-1524.

［14］Peng H，Long F，Ding C.Feature Selection Based on Mutual Information Criteria of Max-Dependency，Max-Relevance，and Min-Redundancy［J］.Pattern Analysis&Machine Intelligence IEEE Transactions，2005，27（8）:1226-1238.

［15］Kononenko I，?imec E，Robnik-?ikonja M.Overcoming the Myopia of Inductive Learning Algorithms with Relieff［J］.Applied Intelligence，1997，7（1）:39-55.

張興（1990-），男，碩士研究生，研究方向為機(jī)器學(xué)習(xí)，智能信息處理，zhangxing0265@link.tyut.edu.cn；

閻高偉（1970-），男，博士，教授，主要研究方向為智能信息處理和多傳感器信息融合，yangaowei@tyut.edu.cn。

Soft Sensor for Ball Mill Fill Level Based on CART-LSSVM Model*

ZHANGXing1，LI Wei1，YAN Gaowei1*，PANG Yusong2
（1.College of Information Engineering，Taiyuan University of Technology，Taiyuan 030024，China；2.College of Mechanical，Maritime and Materials Engineering，Delft University of Technology，Delft，Holland）

Ball mill is a high energy consumption equipment used in electricity，grinding and metallurgical industries.Accurate measurement of its fill level（FL）can improve operational efficiency and safety performance.However，The real-time measurement of FL is difficult to realize，and the components of bearing vibration of ball mill are complex and redundant.Aiming at these problems，a new soft sensor approach of FL based on Classification and Regression Tree（CART）and Least Squares Support Vector Machine（LSSVM）is proposed.Firstly，the Power Spectrum density（PSD）of bearing vibration is obtained by welch method，essential features are achieved by partition subsequently.Secondly，these features are adopted to build CART，and branch nodes of the best model is selected as features.Finally，the LSSVM are used to implement the non-linear mapping between features and FL.The comparative experiments verifies that this model is feasible and practical with high prediction accuracy.

ball mill fill level；soft sensor；feature selection；classification and regression tree；least squares support vector machine；vibration signal

TP29

A

1004-1699（2015）09-1361-06

項目來源:國家自然科學(xué)基金項目（61450011）；山西省自然科學(xué)基金（2015011052）

2015-05-07修改日期:2015-06-24