張君
(1.中煤科工集團 太原研究院有限公司, 山西 太原 030006;2.太原重型機械設備協同創新中心, 山西 太原 030024;3.太原科技大學 機械工程學院, 山西 太原 030024)
分數微積分具有悠久的歷史,廣泛應用于工程、科學、應用數學和經濟學等領域[1-5]。許多現實世界的問題,如物理、化學、流體力學、控制和數學生物學,可以通過建立分數本構的模型來建模[6-10]。為在分數階控制系統求解中要求較高的計算精度,就切比雪夫多項式算法應用于分數階控制系統數值模擬,并對切比雪夫多項式算法做了進一步的前期調研。文獻[11]中為了快速精準地建立機床集合誤差項數學模型,提出了一種基于切比雪夫多項式的參數化建模方法,該建模過程簡單且易程序化,切比雪夫多項式的高精度使得建立的模型精度高,同時為機床設計和誤差補償提供了理論依據。文獻[12]中應用了切比雪夫多項式插值代替了傳統CS算法中的泰勒級數展開,得到了新的二維頻譜信號的近似且推導出了完整的成像流程。該算法精度高、誤差小,誤差空變性弱等特點,提升了成像質量,減少了距離向分塊數,驗證了切比雪夫多項式算法在逼近誤差絕對值和空變性上有很大的優勢。文獻[13]中介紹了用切比雪夫逼近多項式分析非線性電路的方法以及用計算機實現的例子,該方法不僅適用于非線性電路的分析和計算,也可用于解決其它工程上的非線性問題。文獻[14]中將多項式函數最佳逼近的代數式運用在大地測量常用計算公式的最佳逼近問題,結果表明其計算速度,子午線孤長正反算為例可比相應的同精度的經典方法提高5倍至數十倍。文獻[15]中采用切比雪夫多項式擬合 GPS 精密星歷,結果表明最終星歷和快速星歷可以達到 mm 級,證明運用切比雪夫點擬合精密衛星星歷方法可行可靠。文獻[16]中利用切比雪夫擬合方法,在GPS衛星信號失鎖情況下,對DR數據進行擬合處理,并把處理結果與未失鎖情況下的GPS導航數據進行比較,驗證切比雪夫擬合效果,且能獲得良好的精度。文獻[17]基于計算全息法檢測圓柱面系統的輔助裝調,提出選用切比雪夫多項式替代傳統的澤尼克多項式進行系統波像差擬合,建立計算機輔助裝調模型指導被測柱面裝調。文獻[18]中基于切比雪夫多項式具有良好的正交性質,展開公式中的冪函數數項更容易進行變分數階微分的計算,具備了進行函數逼近處理的基礎,具有形成算子矩陣的條件。該文獻還進一步探討變分數階微分方程的數值計算方法奠定了理論基礎,具有一定的工程實用價值。
由此可見,上述文獻都闡述了切比雪夫多項式在求解精度上獲得了較高的精度,適用的應用場合較廣,尤其是針對本單位設計、生產的礦用錨桿鉆機機械臂在井下打孔、上錨桿都需要較高精度的定位,該設備在井下變載荷工況下都需要高精度的定位,且均要采用分數階控制算法設計的反饋控制系統。該系統包括受控對象模型、與其串聯鏈接的控制器模型和負反饋連接的反饋模型,同時切比雪夫多項式擬合方法在分數控制系統數值模擬研究中有較強的理論基礎和實際應用場合。因此,對于設備控制系統要解決人工生產效率低、勞動強度大的問題,非常有必要運用切比雪夫多項式算法來研究分數控制系統數值模擬計算。
典型的分數反饋控制系統框架如圖1所示,其中Gc(t)為分數控制器;G0(t)為分數控制器系統的傳遞函數;Gf(t)為分數系統的反饋環傳遞函數;U(t)和Y(t)a為系統的輸入和輸出。

圖1 分數離散控制系統框架圖
該分數控制系統是開關始終處于閉合狀態下的連續系統,其時域模型的建立,即為:
anDαny(t)+an-1Dαn-1y(t)+…+a0Dα0y(t)=
bmDβmu(t)+bm-1Dβm-1u(t)+…+b0Dβ0u(t)
(1)


(2)
迄今為止,給出了求解分數階微分方程的各種數值方法。這些方法包括小波法[19-20],切比雪夫,勒讓德多項式[21-22]和配方法。在文獻[23]中,L. Pezza就是利用一種近似方法,來研究半線性非局部分數演化方程的部分近似能控性。在文獻[24]中,Ali Lotf使用Epsilon penalty和Ritz方法的擴展,來解決一類混合邊界條件下的分數最優控制問題。
本文利用切比雪夫多項式得到了分數控制系統的數值解,其組織架構如下:在第2部分,介紹了分數階微積分的定義;第3部分給出了切比雪夫多項式的一些相關性質;第4部分介紹了數值方法和數值例子。第五部分得出結論。
Riemann-Liouvilleμ階分數階積分函數[25],定義為:
(3)
式中:ζ∈R+和Γ(·)表示gamma函數。
Caputoμ階的分數階微分函數[25],定義為:
m-1<μ≤m
(4)
式中:μ∈R+,m∈N+。
切比雪夫多項式次數的解析形式[26]:
(5)
式中:Ti(0)=(-1)iandTi(1)=1。
正交性式為

(6)
b0=2,bk=1,k≥1。
假設y(t)∈L2[0,1],則由切比雪夫多項式展開多項式展開得:

(7)
其中:ci系數的表達式為:
(8)
考慮式(5)中的截斷級數,其表達式可如下式(9)所示:

(9)
其中
C=[c0,c1,…,cM]T
Φ(t)=[T0(t),T1(t),…,TM(t)]T
(10)
式中:向量Φ(t)的微分函數表達式可為:
(11)
其中P(1)是(M+1)×(M+1)由運行矩陣求導可得,即
(12)
式中,當M為奇數時,k=1,3,5,…,M;當M為偶數時,k=1,3,5,…,M-1。
由此類似,運行矩陣Pn由Φ(t)求n階導數得:
(13)
式中:Pn=(P(1))n。
主要目的是證明切比雪夫多項式分數階導數的下列定理[26]。
引理1 假設Ti(t)是切比雪夫多項式,則

DμTi(t)=0, i=0,1,2,…,
μ

-1,μ>0
(14)
定理1 設Φ(t)為式(12)中定義的Chebyshev向量,并假設μ>0。
(15)
式中:P(μ)由Caputo定義下的(M+1)×(M+1)運行矩陣μ階微分法求得,且定義為:
P(μ)=

000…0???…?000…0Sμ(μ,0)Sμ(μ,1)Sμ(μ,2)…Sμ(μ,M)???…?Sμ(i,0)Sμ(i,1)Sμ(i,2)…Sμ(i,M)???…?Sμ(M,0)Sμ(M,1)Sμ(M,2)…Sμ(M,M)?è?????????????÷÷÷÷÷÷÷÷÷÷÷
式中:

Sμ(i,j)=∑ik=μ·[(-1)i-k2i(i+k-1)!Γk-μ+12()]/[bjΓk+12()(i-k)!Γ(k-μ-j+1)Γ(k-μ+j+1)]

利用Chebyshev多項式對分數控制系統進行數值模擬。其式(1)的每一項都可以用Chebyshev多項式的基表示為:
(16)
(17)
?
(18)
和
(19)
(20)
?
(21)
式中:C和U均由式(11)獲得。
將式(16)~(21)帶入式(1)中,可以得到:
anCTP(αn)Φ(t)+an-1CTP(αn-1)Φ(t)+…+
a0CTP(α1)Φ(t)=bmUTP(βn)Φ(t)+
bm-1UTP(βn-1)Φ(t)+…+b0UTP(β1)Φ(t)
(22)
測試時,考慮分數控制系統,即
D1.8y(t)+D1.5y(t)+y′(t)+y(t)=u(t),
y′(0)=0,y(0)=-1,t∈[0,1]
(23)

該系統的解析解是y(t)=t2-1。當M=4,6,8時,數值和分析結果的絕對誤差如表1所示。由表1可知,數值解與解析解均隨著M的增長有較好的吻合。

表1 數值和分析結果的絕對誤差
本文提出了一種利用切比雪夫多項式求解分數控制系統的數值方法,即利用對運行矩陣求解n階導數,將原問題轉化為易于求解的線性代數方程組。通過數值實驗結果表明,隨著M的增大,數值解收斂于解析解。