基于Matlab／Simulink非線性系統(tǒng)相軌跡的繪制

張 婧，蓋文東，高宏巖

（山東科技大學(xué)電氣與自動(dòng)化工程學(xué)院，山東青島266590）

0 引言

系統(tǒng)可分為線性系統(tǒng)和非線性系統(tǒng)，但嚴(yán)格地說(shuō)，實(shí)際系統(tǒng)都含有非線性［1］。若系統(tǒng)的非線性程度較弱且運(yùn)行在平衡工作點(diǎn)附近較小范圍內(nèi)時(shí)，可忽略次要的非線性因素或經(jīng)小偏差法將其線性化，從而可按照線性系統(tǒng)理論來(lái)處理；若非線性程度較明顯或系統(tǒng)的非線性因素不能忽略時(shí)，只能采用非線性系統(tǒng)理論來(lái)分析和設(shè)計(jì)［2］。

目前常用的非線性系統(tǒng)分析的工程近似方法有：相平面法和描述函數(shù)法等。其中，描述函數(shù)法是基于頻域分析法和非線性特性諧波線性化的一種圖解分析法，該方法一般只討論運(yùn)動(dòng)的穩(wěn)定性、是否存在自激振蕩及振蕩的頻率和振幅，不涉及具體的輸入信號(hào)形式，也不涉及在給定初始條件下的具體運(yùn)動(dòng)狀態(tài)，但對(duì)系統(tǒng)結(jié)構(gòu)、非線性環(huán)節(jié)的特性和線性部分的性能有一定的要求，因此該方法的應(yīng)用有一定的限制條件；而相平面法是求解二階非線性常微分方程較精確的時(shí)域分析的圖解法，這種方法不僅可以判斷系統(tǒng)的穩(wěn)定性和自激振蕩，還可以討論系統(tǒng)在不同初始條件以及在不同輸入信號(hào)作用下的運(yùn)動(dòng)狀態(tài)［2～3］。但是由于相平面法是一種圖解法，需要繪制大量曲線，在授課過(guò)程中，教師在黑板上繪制曲線的準(zhǔn)確度難以保證，系統(tǒng)的性能便也難以體現(xiàn)；教學(xué)實(shí)驗(yàn)中，著重采用相平面法分析含有典型非線性環(huán)節(jié)的非線性系統(tǒng)，但在實(shí)驗(yàn)中經(jīng)常會(huì)受到元件參數(shù)不穩(wěn)定的影響，使實(shí)驗(yàn)值與理論值不太吻合，容易導(dǎo)致學(xué)生產(chǎn)生錯(cuò)誤的理解［4］。因此如何準(zhǔn)確地繪制相軌跡，是分析和設(shè)計(jì)非線性系統(tǒng)的前提。

Matlab是將數(shù)值分析、信號(hào)處理和仿真等諸多功能集合在一起的編程工具，在各領(lǐng)域中得到廣泛應(yīng)用。應(yīng)用Matlab／Simulink可以精確繪制非線性系統(tǒng)的相軌跡，并可得到相應(yīng)的時(shí)間響應(yīng)曲線，便于對(duì)系統(tǒng)進(jìn)行仿真和分析。

對(duì)于可用非線性微分方程描述的非線性系統(tǒng)，可在Matlab中調(diào)用ode及plot函數(shù)求解系統(tǒng)微分方程高階數(shù)值，再繪制相應(yīng)的相軌跡。

本文主要討論基于Matlab／Simulink繪制由一個(gè)非線性環(huán)節(jié)與一個(gè)線性環(huán)節(jié)組成的分區(qū)非線性系統(tǒng)的相軌跡，系統(tǒng)結(jié)構(gòu)圖如圖1所示。

圖1 非線性系統(tǒng)結(jié)構(gòu)圖

1 非線性系統(tǒng)相軌跡的繪制

相軌跡的繪制方法可采用解析法、圖解法、計(jì)算機(jī)仿真計(jì)算法等。解析法是通過(guò)求解微分方程獲得相軌跡，對(duì)于難以采用解析法的系統(tǒng)可采用圖解法（等傾線法），但圖解法是一種耗時(shí)多的近似作圖方法，且當(dāng)求出的等傾線是復(fù)雜的曲線時(shí)，就不宜采用。而計(jì)算機(jī)具有極大的靈活性，可以編制程序適應(yīng)各種非線性系統(tǒng)，因此對(duì)于較為復(fù)雜的系統(tǒng)，或需要繪制精確相軌跡時(shí)，可采用計(jì)算機(jī)仿真計(jì)算方法，即采用Matlab／Simulink來(lái)對(duì)系統(tǒng)進(jìn)行仿真計(jì)算。

1．1 二階非線性系統(tǒng)相軌跡的繪制

例1．已知非線性系統(tǒng)的結(jié)構(gòu)圖如圖2所示，輸入r（t）＝R·1（t），系統(tǒng)初始狀態(tài)為0，△ ＝1，T＝1s，k·K＝1，試在相平面e- 上繪制系統(tǒng)的相軌跡。

圖2 非線性系統(tǒng)的結(jié)構(gòu)圖

若要繪制此二階非線性系統(tǒng)的相軌跡，可根據(jù)非線性環(huán)節(jié)的分段特點(diǎn)，先將相平面e- 分成若干區(qū)域，再應(yīng)用線性系統(tǒng)的相平面分析法在各個(gè)區(qū)域繪制相軌跡，具體分析繪制過(guò)程可參見(jiàn)文獻(xiàn)［1］。除此之外，還可借助Matlab／Simulink來(lái)繪制此系統(tǒng)的相軌跡。

1）利用Matlab文本語(yǔ)言繪制相軌跡

由死區(qū)非線性特性將相平面分成三個(gè)區(qū)域，根據(jù)已知條件可得系統(tǒng)誤差的分區(qū)線性微分方程為

則可編寫(xiě)分區(qū)線性微分方程的M函數(shù)“fundeadzone”。

%含有死區(qū)特性的非線性系統(tǒng)的M函數(shù)

然后調(diào)用ode45函數(shù)求解微分方程組／狀態(tài)方程，設(shè)階躍函數(shù)幅度R分別為2和10，t0＝0 s，tf＝16 s，即可繪制非線性系統(tǒng)的相軌跡。為突出初始條件對(duì)系統(tǒng)性能的影響，將兩條相軌跡在同一圖中繪制出來(lái)，如圖3所示。

2）利用Simulink繪制相軌跡

按照?qǐng)D2所示非線性系統(tǒng)的結(jié)構(gòu)圖，在Simulink中構(gòu)建其仿真模型，由于題目要求在相平面e-上繪制相軌跡，故在偏差處引出e和 的值，如圖4所示。取階躍函數(shù)幅度分別為2和10，即可在XY Graph中繪制相應(yīng)的相軌跡，繪制結(jié)果與圖3相同。

由圖3可知，當(dāng)輸入的階躍信號(hào)幅度不同時(shí)，即誤差的初始條件不同時(shí)，非線性系統(tǒng)的運(yùn)動(dòng)區(qū)域及運(yùn)動(dòng)軌跡皆不同。由相軌跡的運(yùn)動(dòng)形式可知，系統(tǒng)存在穩(wěn)態(tài)誤差，穩(wěn)態(tài)誤差的大小不僅取決于系統(tǒng)的結(jié)構(gòu)和參數(shù)，而且與系統(tǒng)的初始條件和輸入信號(hào)的幅度有關(guān)，其值等于相軌跡的終點(diǎn)與坐標(biāo)原點(diǎn)之間的距離。

圖3 階躍輸入時(shí)非線性系統(tǒng)的相軌跡圖

圖4 非線性系統(tǒng)Simulink仿真模型

若想用相平面法分析含有其它非線性環(huán)節(jié)的非線性系統(tǒng)，只需根據(jù)非線性環(huán)節(jié)的分界線（開(kāi)關(guān)線），求出分區(qū)線性微分方程的表達(dá)式，修改M函數(shù)中sys與狀態(tài)變量x的關(guān)系即可。但此方法需要根據(jù)非線性環(huán)節(jié)的分界線將其線性部分的傳遞函數(shù)寫(xiě)出分區(qū)線性微分方程，仍較為復(fù)雜，下面介紹一種更為簡(jiǎn)便的方法。

1．2 高階非線性系統(tǒng)相軌跡的繪制

對(duì)于如圖1所示的高階非線性系統(tǒng)，其相軌跡的繪制比較困難，在文獻(xiàn)［5～6］中闡述了高階系統(tǒng)可采用δ法或等傾線法繪制其相軌跡，但這兩種方法還是比較復(fù)雜，也有文獻(xiàn)［7～8］借助Matlab／simulink繪制高階非線性系統(tǒng)的相平面。下面主要介紹基于Matlab／Simulink來(lái)繪制高階非線性系統(tǒng)相軌跡的方法。

一般將高階線性環(huán)節(jié)部分寫(xiě)成可觀標(biāo)準(zhǔn)型而非可控標(biāo)準(zhǔn)型，是因?yàn)榭捎^標(biāo)準(zhǔn)型線性環(huán)節(jié)的輸出c（t）恰為最后一個(gè)狀態(tài)變量xn，即狀態(tài)變量圖中最后一個(gè)積分器輸出就是線性環(huán)節(jié)輸出c（t），其輸入為線性環(huán)節(jié)輸出的一階導(dǎo)數(shù) （t）。若已知系統(tǒng)輸出的初始條件，或要求根據(jù)相平面法分析高階非線性系統(tǒng)的運(yùn)動(dòng)特性，即在相平面c-上繪制系統(tǒng)相軌跡，相較可控型，可觀型的表示方法更為簡(jiǎn)便。

例2．已知非線性系統(tǒng)的結(jié)構(gòu)圖如圖5所示，h＝1，a＝1，系統(tǒng)的初始條件為［c（0），（0），（0）］＝［0．6，1，0．7］，試在相平面 c- 上繪制系統(tǒng)的相軌跡。

圖5 非線性系統(tǒng)結(jié)構(gòu)圖

將線性環(huán)節(jié)寫(xiě)成可觀標(biāo)準(zhǔn)型：

由已知條件［c（0），（0），（0）］＝［0．6，1，0．7］及式（2），可得：

求得狀態(tài)變量的初始值為：［x1（0），x2（0），x3（0）］＝［3．8，2．5，0．6］。

1）利用Matlab文本語(yǔ)言繪制相軌跡

令k＝2，基于狀態(tài)方程可編寫(xiě)分區(qū)線性微分方程的M函數(shù)“funsatu”。

然后調(diào)用ode45函數(shù)求解這個(gè)微分方程組／狀態(tài)方程，設(shè)t0 ＝0 s，tf＝40 s，分別令k＝1，2 和3，即可繪制非線性系統(tǒng)的相軌跡。為突出系統(tǒng)參數(shù)對(duì)系統(tǒng)性能的影響，將三條相軌跡在同一圖中繪制出來(lái)，如圖6所示。

圖6 k取不同值時(shí)高階非線性系統(tǒng)的相軌跡圖

在Matlab文本語(yǔ)言中，K雖是變量，但沿用為正體。

2）利用Simulink繪制相軌跡

在Simulink中將圖5所示的線性環(huán)節(jié)以狀態(tài)變量圖表示，再結(jié)合非線性環(huán)節(jié)構(gòu)造仿真模型如圖7所示，其中狀態(tài)變量的初始值可以在積分器模塊參數(shù)中設(shè)置，第一個(gè)積分器的初始值為3．8，第二個(gè)積分器的初始值為2．5，第三個(gè)積分器的初始值為0．6；增益“Gain”k可取不同值，則在XY Graph中即可繪制相應(yīng)的相軌跡，繪制結(jié)果與圖6相同。

由圖6可知，當(dāng)k取不同值時(shí)，相軌跡從初始條件（0．6，1）出發(fā)，以不同形式運(yùn)動(dòng)。當(dāng)k＝1和k＝2時(shí)，相軌跡以螺旋線的方式最終趨于穩(wěn)定，且穩(wěn)定在坐標(biāo)原點(diǎn)；當(dāng)k＝3時(shí)，相軌跡最終趨于極限環(huán)，系統(tǒng)輸出呈現(xiàn)自激振蕩，振幅約為1．5。由描述函數(shù)法分析可知，當(dāng)系統(tǒng)出現(xiàn)自激振蕩時(shí)的臨界值k＝2．5，將k＝2．5代入圖7，相平面法分析結(jié)果與描述函數(shù)法分析結(jié)果一致。

若想用相平面法分析含有其它非線性環(huán)節(jié)的系統(tǒng)，只需調(diào)換非線性環(huán)節(jié)，根據(jù)輸入信號(hào)及初始條件即可繪制相應(yīng)的相軌跡。

通過(guò)上述非線性二階、高階系統(tǒng)相軌跡繪制的示例，可以看出：基于Matlab文本語(yǔ)言的程序編寫(xiě)較為復(fù)雜，基于Simulink的方法更為簡(jiǎn)便；并且針對(duì)某些非線性系統(tǒng)要求分析系統(tǒng)在輸入信號(hào)下的運(yùn)動(dòng)特性，即在相平面上c-繪制相軌跡，基于Simulink繪制非線性系統(tǒng)相軌跡的方法亦更為簡(jiǎn)便。

圖7 高階非線性系統(tǒng)的Simulink仿真模型

2 結(jié)語(yǔ)

針對(duì)非線性二階、高階系統(tǒng)相軌跡繪制的復(fù)雜性和不準(zhǔn)確性，本文主要介紹了基于Matlab文本語(yǔ)言及基于Simulink構(gòu)建非線性系統(tǒng)仿真模型的兩種方法。由于相平面法實(shí)際上是狀態(tài)空間法的應(yīng)用，因此在繪制相軌跡時(shí)將線性環(huán)節(jié)表示成可觀標(biāo)準(zhǔn)型的狀態(tài)空間形式，根據(jù)非線性環(huán)節(jié)、輸入信號(hào)形式及初始條件的不同，即可繪制出符合條件的相軌跡。借助Matlab和Simulink，編寫(xiě)語(yǔ)言和構(gòu)建仿真模型兩種方法都簡(jiǎn)單易懂，但基于Simulink的方法則更為簡(jiǎn)便。Matlab的引入擴(kuò)大了相平面法的適用范圍，具有一定的實(shí)用價(jià)值，為非線性系統(tǒng)的分析和設(shè)計(jì)提供了更為精確的、快速的制圖方法。