滑模變結(jié)構(gòu)控制:Simulink仿真設(shè)計與抖振抑制實踐)
滑模變結(jié)構(gòu)控制這名字聽著唬人但在機器人運動學(xué)仿真里它其實是個特別“直接”的控制思路——說白了就是設(shè)計一條誤差收斂的“滑模面”然后用不連續(xù)的切換控制把系統(tǒng)狀態(tài)“按”到這條面上再沿著面滑到零點。我當年第一次在Simulink里把二連桿機械臂的運動學(xué)模型和滑??刂破鹘悠饋頃r看到示波器里那條誤差曲線穩(wěn)穩(wěn)收斂才真正理解課本上那句“對參數(shù)攝動和外部擾動具有魯棒性”是什么意思。這篇就圍繞“機器人運動學(xué)控制 Simulink仿真模型 滑模變結(jié)構(gòu)控制”三件事展開把建模、控制器設(shè)計、仿真實操和踩坑經(jīng)驗一次講透。適合剛接觸滑??刂频难芯可?、做畢業(yè)設(shè)計的本科生以及想快速驗證控制算法的工程師。1. 先搞清楚機器人運動學(xué)模型1.1 為什么從運動學(xué)入手很多人一上來就懟動力學(xué)方程結(jié)果被慣性矩陣、科氏力、重力項搞得焦頭爛額?;?刂齐m然是魯棒控制里的常客但學(xué)習(xí)路徑上完全可以先從運動學(xué)層面切入。運動學(xué)只研究關(guān)節(jié)角度與末端位置之間的映射關(guān)系不涉及力/力矩與加速度的因果模型形式簡單控制目標也直觀——就是讓末端或關(guān)節(jié)角度精確跟蹤期望軌跡。在這個層面上把滑模設(shè)計的邏輯跑通再去碰動力學(xué)會輕松得多。從工程角度看機器人運動學(xué)控制場景也很常見視覺引導(dǎo)的末端定位、關(guān)節(jié)空間軌跡跟蹤、冗余機械臂的運動規(guī)劃等本質(zhì)上都需要先解決“讓關(guān)節(jié)角度按期望曲線走”的問題。所以用二連桿機械臂作為仿真對象既保留了多軸耦合的特征又不會像六軸機械臂那樣把雅可比矩陣推導(dǎo)變成災(zāi)難。等你把二連桿的運動學(xué)滑模控制吃透換成其他構(gòu)型只是換一套正解公式而已。1.2 二連桿機械臂的運動學(xué)方程我采用最常見的平面二連桿模型兩個連桿長度分別為 (L_1)、(L_2)關(guān)節(jié)角度為 (\theta_1)、(\theta_2)。末端位置 ((x, y)) 由正運動學(xué)給出[ x L_1 \cos(\theta_1) L_2 \cos(\theta_1 \theta_2) ] [ y L_1 \sin(\theta_1) L_2 \sin(\theta_1 \theta_2) ]對時間求導(dǎo)得到末端速度與關(guān)節(jié)角速度的關(guān)系[ \dot{x} -L_1 \sin(\theta_1)\dot{\theta}_1 - L_2 \sin(\theta_1 \theta_2)(\dot{\theta}_1 \dot{\theta}_2) ] [ \dot{y} L_1 \cos(\theta_1)\dot{\theta}_1 L_2 \cos(\theta_1 \theta_2)(\dot{\theta}_1 \dot{\theta}_2) ]寫成矩陣形式就是 (\dot{\mathbf{p}} \mathbf{J}(\boldsymbol{\theta}) \dot{\boldsymbol{\theta}})其中 (\mathbf{J}) 是雅可比矩陣。運動學(xué)控制的目標通常是給定期望末端軌跡 (\mathbf{p}_d(t))求關(guān)節(jié)角速度指令 (\dot{\boldsymbol{\theta}})。如果只用逆雅可比 (\dot{\boldsymbol{\theta}} \mathbf{J}^{-1}(\dot{\mathbf{p}}_d \mathbf{K}\mathbf{e}))那是經(jīng)典的“運動學(xué)PID”但面對模型誤差和擾動時跟蹤性能會打折。滑??刂圃谶@里的用法就是針對跟蹤誤差設(shè)計魯棒控制器不依賴精確的逆雅可比。在Simulink里搭模型時不少教程喜歡直接用積分模塊把 (\dot{\boldsymbol{\theta}}) 變成 (\boldsymbol{\theta})再用角度計算末端位置。我建議把這些運算封裝進一個“正運動學(xué)子系統(tǒng)”輸入關(guān)節(jié)角度輸出末端位置和雅可比矩陣這樣后期換三連桿、換球坐標都很方便。模塊層面用Fcn模塊或MATLAB Function都可以但MATLAB Function寫雅可比更清晰還能順便輸出條件數(shù)以觀察奇異性。2. 滑模變結(jié)構(gòu)控制的核心思路2.1 滑模面的設(shè)計滑模控制的設(shè)計分兩步選滑模面定控制律。先說滑模面。對于關(guān)節(jié)角度跟蹤問題定義誤差 (e \theta_d - \theta)常見的線性滑模面為[ s \dot{e} \lambda e ]其中 (\lambda 0)。這個面的物理意義很直觀如果系統(tǒng)狀態(tài)能保持在 (s 0) 上那么 (\dot{e} -\lambda e)誤差指數(shù)收斂到零收斂速度由 (\lambda) 決定。你可以把 (s) 理解成“誤差空間的綜合指標”它同時包含位置誤差和速度誤差。當 (s) 不為零時控制器要驅(qū)動系統(tǒng)讓 (s) 趨向零這個過程叫“到達階段”一旦 (s 0)就進入“滑模階段”。之所以說滑??刂啤皩ζヅ鋽_動不敏感”是因為在滑模面上系統(tǒng)運動由滑模面方程決定原有的動力學(xué)參數(shù)被“替換”了。當然前提是控制增益能克服擾動上界否則系統(tǒng)會被擾動推出滑模面。這個“上界”的概念在設(shè)計符號函數(shù)增益時特別關(guān)鍵后面細說。2.2 趨近律與抖振抑制光讓 (s 0) 還不夠得規(guī)定 (s) 怎么走向零。最粗暴的做法是用符號函數(shù)[ \dot{s} -\eta ,\text{sgn}(s) ]對應(yīng)控制律里會出現(xiàn) (\eta ,\text{sgn}(s)) 這樣的切換項。這種“指數(shù)趨近律”簡單但會導(dǎo)致一個經(jīng)典問題——抖振。因為符號函數(shù)在零附近高速切換理想狀態(tài)下是無限頻率實際離散仿真里就會表現(xiàn)為高頻振蕩反映在關(guān)節(jié)角度曲線上就是毛刺或極限環(huán)。很多初學(xué)者第一次跑仿真看見角度曲線抖得像心電圖第一反應(yīng)是模型錯了其實只是符號函數(shù)增益太大。抑制抖振的常規(guī)方法有幾種用飽和函數(shù) (\text{sat}(s/\phi)) 代替 (\text{sgn}(s))在邊界層內(nèi)做線性過渡代價是犧牲一點收斂精度。用超螺旋算法等高階滑模本質(zhì)是讓切換項作用在 (s) 的高階導(dǎo)數(shù)上抖振大幅減弱但參數(shù)調(diào)節(jié)復(fù)雜些。適當調(diào)小切換增益 (\eta)只要能覆蓋模型誤差和擾動的上界就行沒必要給得過大。我個人在做運動學(xué)仿真時先用飽和函數(shù)把邏輯跑通再換成符號函數(shù)對比抖振差異這個對比過程本身就是學(xué)習(xí)滑??刂谱詈玫慕滩?。你可以把飽和函數(shù)邊界層厚度 (\phi) 設(shè)成 0.01~0.05增益 (\eta) 先設(shè)個保守值比如 0.5然后在仿真里一點點加。2.3 為什么要用“趨近律”而不是直接設(shè)計控制輸入運動學(xué)模型的輸入是關(guān)節(jié)角速度 (\dot{\boldsymbol{\theta}})不是力矩。所以控制律的形式不能照搬動力學(xué)里的 (u \tau)。這里需要用“運動學(xué)級”的虛擬控制量 (v \dot{\boldsymbol{\theta}}_c)把它當作關(guān)節(jié)角速度指令。我們設(shè)計[ \dot{\boldsymbol{\theta}}_c \mathbf{J}^{-1}\left(\dot{\mathbf{p}}_d \lambda \mathbf{e} \eta ,\text{sgn}(s)\right) ]其中 (\mathbf{e} \mathbf{p}_d - \mathbf{p}) 是末端位置誤差(s \dot{\mathbf{e}} \lambda \mathbf{e})。這塊設(shè)計邏輯是如果用 (\dot{\boldsymbol{\theta}}_c) 驅(qū)動真實機械臂那么末端速度會趨近期望速度誤差沿滑模面收斂。在Simulink里這個 (\dot{\boldsymbol{\theta}}_c) 通常不能直接作為物理關(guān)節(jié)速度輸入需要再串一個底層速度環(huán)或者直接把 (\dot{\boldsymbol{\theta}}_c) 當作指令給理想速度源。做運動學(xué)級仿真時我用的是“積分器理想速度驅(qū)動”把 (\dot{\boldsymbol{\theta}}_c) 積分得到實際角度再反饋給正運動學(xué)模塊這樣簡單且能看清控制器核心性能不摻動力學(xué)干擾。當然這種理想化模型忽略了很多執(zhí)行器特性但它作為學(xué)習(xí)載體非常合適——你可以單獨評估滑??刂频聂敯粜员热缃o雅可比矩陣加5%的參數(shù)偏差看誤差是否還能收斂這就比純粹的PID更能體現(xiàn)滑模優(yōu)勢。3. Simulink仿真模型搭建全流程3.1 模型架構(gòu)與模塊選型整個Simulink模型分四塊軌跡生成、滑??刂破?、被控對象正運動學(xué)雅可比、信號觀測。我從一個大框架說起你按這個結(jié)構(gòu)搭不容易亂。軌跡生成模塊我用MATLAB Function生成圓形軌跡圓心 ((0.6, 0.4))半徑 0.1角頻率 0.5 rad/s仿真時長 10 秒。輸出期望末端位置 (\mathbf{p}_d) 和期望速度 (\dot{\mathbf{p}}_d)。也可以在Constant模塊里用Clock和表達式實現(xiàn)但MATLAB Function最好改參數(shù)?;?刂破髂K輸入期望位置、期望速度、實際位置、實際角度輸出 (\dot{\boldsymbol{\theta}}_c)。這個模塊內(nèi)部用MATLAB Function實現(xiàn)雅可比計算、誤差計算、趨近律和控制律。注意雅可比矩陣可能奇異二連桿在某些位形下 (\det(\mathbf{J}) 0)需要在函數(shù)里加個判斷如果行列式絕對值小于閾值就用偽逆或阻尼最小二乘。常見做法是加個參數(shù) (\epsilon)[ \mathbf{J}^ \mathbf{J}^T(\mathbf{J}\mathbf{J}^T \epsilon \mathbf{I})^{-1} ]這個 (\epsilon) 設(shè)成 0.01 左右能避免仿真中途因奇異而報NAN。被控對象模塊理想速度源就是兩個積分器輸入 (\dot{\boldsymbol{\theta}}_c)輸出實際關(guān)節(jié)角度 (\boldsymbol{\theta})反饋給控制器和正運動學(xué)。如果想更貼近真實可以在這個積分器前加一階慣性環(huán)節(jié)模擬執(zhí)行器延遲但學(xué)習(xí)階段別加先看純運動學(xué)下的控制效果。信號觀測模塊用Scope觀察末端位置跟蹤、關(guān)節(jié)角度曲線、誤差和滑模面 (s)。建議把誤差和 (s) 單獨放一個Scope方便看收斂速度。3.2 機器人運動學(xué)子系統(tǒng)實現(xiàn)我習(xí)慣把正運動學(xué)封裝成一個子系統(tǒng)內(nèi)部用MATLAB Function寫function [pos, J] fk_2link(theta, L1, L2) % 正運動學(xué) q1 theta(1); q2 theta(2); x L1*cos(q1) L2*cos(q1q2); y L1*sin(q1) L2*sin(q1q2); pos [x; y]; % 雅可比矩陣 J [-L1*sin(q1)-L2*sin(q1q2), -L2*sin(q1q2); L1*cos(q1)L2*cos(q1q2), L2*cos(q1q2)]; end輸入是 (\theta_1)、(\theta_2) 合成的一個向量信號用Vector Concatenate模塊或直接兩路輸入到MATLAB Function都行。輸出pos和J。注意MATLAB Function里不要用全局變量參數(shù)通過對話框傳遞。運行時長度參數(shù) (L_1 0.5)(L_2 0.4)單位米。這里有個容易踩的坑MATLAB Function默認輸出類型會被推斷如果初始化時沒有給輸出變量賦值Simulink可能報“輸出未定義”。所以在函數(shù)開頭最好先寫一行 (\text{pos} zeros(2,1); J zeros(2,2);)或者用coder.extrinsic調(diào)用不了就直接在代碼里給默認值。3.3 滑模控制器子系統(tǒng)實現(xiàn)控制器的MATLAB Function寫法大概是這樣function [theta_dot_cmd] smc_controller(pd, pd_dot, p_current, theta, L1, L2, lambda, eta, phi) % pd: 期望位置2x1, pd_dot: 期望速度2x1, p_current: 實際位置2x1 % theta: 關(guān)節(jié)角2x1 [~, J] fk_2link(theta, L1, L2); e pd - p_current; edot pd_dot - J * theta_dot_current; % 注意這里需要當前關(guān)節(jié)角速度 s edot lambda * e; % 滑模面 % 飽和函數(shù)替代符號函數(shù) sat_s min(1, max(-1, s / phi)); % 阻尼最小二乘逆 Jt J; JtJ J * Jt 0.01 * eye(2); Jinv Jt * (JtJ \ eye(2)); % 控制律 theta_dot_cmd Jinv * (pd_dot lambda * e eta * sat_s); end這里隱含一個依賴需要當前關(guān)節(jié)角速度 (\dot{\boldsymbol{\theta}})。我們的被控對象是理想積分器所以 (\dot{\boldsymbol{\theta}}) 正好等于 (\dot{\boldsymbol{\theta}}_c)也就是控制器的輸出。這就形成一個代數(shù)環(huán)控制器輸出決定實際速度實際速度又反饋到控制器計算。Simulink會報警告“代數(shù)環(huán)”雖然在小步長下能跑但最好顯式消除。消除代數(shù)環(huán)的辦法是把積分器輸出的狀態(tài) (\boldsymbol{\theta}) 反饋給控制器但在控制器里用差分近似關(guān)節(jié)速度比如 (\dot{\theta}[k] (\theta[k] - \theta[k-1]) / h)通過Memory模塊或Unit Delay實現(xiàn)。這樣可以打破代數(shù)環(huán)代價是速度估計有一點延遲和噪聲。運動學(xué)仿真中這個噪聲可接受。更優(yōu)雅的方案是在被控對象里加一個一階慣性濾波但學(xué)習(xí)階段我建議用Unit Delay近似簡單實用。另一個更直接的辦法既然被控對象是積分器控制器輸出就是速度那就別把“當前速度”當反饋直接用期望速度的誤差做滑模面。末端位置誤差 (e) 是位置量其導(dǎo)數(shù)可以從期望速度減去“當前速度”但我們可以把滑模面定義為[ s \dot{\mathbf{p}}_d - \dot{\mathbf{p}} \lambda e ]其中 (\dot{\mathbf{p}} \mathbf{J} \dot{\boldsymbol{\theta}})。如果控制器輸出 (\dot{\boldsymbol{\theta}}_c) 直接等于被控對象速度那么 (\dot{\mathbf{p}} \mathbf{J}\dot{\boldsymbol{\theta}}_c)。代入控制律會得到關(guān)于 (\dot{\boldsymbol{\theta}}_c) 的隱式方程。所以還是得打破代數(shù)環(huán)。我在實踐中用如上Unit Delay估計速度方案效果穩(wěn)定。3.4 參數(shù)設(shè)置與仿真配置仿真步長選擇很關(guān)鍵?;?刂茙Х柡瘮?shù)切換時如果步長太大控制器輸出在幾個步長之間反復(fù)跳變誤差曲線會呈現(xiàn)鋸齒狀。我推薦用變步長ode45最大步長設(shè) (0.001) 秒相對誤差 (1e-4)。如果用的是飽和函數(shù)可以放寬到 (0.005) 秒。控制器參數(shù)我試用過一組不錯的初始值(\lambda 3)(\eta 0.8)(\phi 0.05)。仿真10秒期望軌跡是半徑0.1米的圓。初始關(guān)節(jié)角度可以從末端位置反解也可以直接設(shè) (\theta_1 0.5)(\theta_2 0.8) 開始讓控制器自己收。注意期望軌跡的起始點最好和機械臂實際末端位置一致否則初始誤差很大符號函數(shù)增益可能讓速度指令瞬間飽和出現(xiàn)超調(diào)。我一般用MATLAB腳本先算初始關(guān)節(jié)角對應(yīng)末端位置把期望軌跡的起點移到那里比如自定義軌跡函數(shù)[ x_d 0.6 0.1\cos(0.5t - \phi_0) ] [ y_d 0.4 0.1\sin(0.5t - \phi_0) ]其中 (\phi_0) 由初始末端位置的極角決定這樣初始誤差接近零。仿真配置里還要注意解法器是否支持信號代數(shù)環(huán)。用Unit Delay后應(yīng)該沒有代數(shù)環(huán)警告如果有可以再插入一個Memory模塊在反饋路徑上。模型里還建議把所有Scope數(shù)據(jù)記錄到工作區(qū)用logsout方便后續(xù)分析。4. 仿真結(jié)果分析與調(diào)試4.1 跟蹤效果怎么看仿真跑完先看末端位置跟蹤曲線。把期望圓和實際軌跡畫在一張圖里如果軌跡重疊得比較好說明跟蹤精度高。然后看關(guān)節(jié)角度曲線理想情況下是平滑曲線不應(yīng)有高頻分量。最后看誤差曲線 (e_x)、(e_y)穩(wěn)態(tài)誤差應(yīng)該在 (10^{-3}) 量級。還有一個關(guān)鍵指標是滑模面 (s) 的收斂過程。如果設(shè)計正確且增益合適(s) 會在很短時間比如0.1秒內(nèi)從初始值衰減到零隨后一直在零附近小幅波動。這個波動幅度取決于飽和函數(shù)邊界層 (\phi) 和擾動如果波動太大說明增益設(shè)置不合理或者速度估計噪聲大。4.2 抖振問題與參數(shù)調(diào)節(jié)經(jīng)驗最典型的問題是用符號函數(shù)時關(guān)節(jié)速度指令呈高頻切換角度曲線有毛刺。這時候不要急著加大 (\lambda)先檢查 (\eta) 是否過大。(\eta) 的物理意義是“對抗擾動的強度”如果模型精確、無擾動(\eta) 只需要大于系統(tǒng)名義部分的誤差上界。在運動學(xué)級仿真里唯一的擾動是數(shù)值誤差和速度估計誤差(\eta) 設(shè)成 0.2~0.5 通常足夠。把 (\eta) 調(diào)到 2 以上抖振會非常明顯除非你故意要觀察抖振。第二個常見問題是初始誤差大導(dǎo)致速度指令飽和。我在控制器中加過飽和模塊限制 (\dot{\theta}_c) 幅值在 3 rad/s 以內(nèi)否則仿真早期速度會飆到幾十數(shù)值發(fā)散。這個限幅不影響穩(wěn)態(tài)性能但能顯著提高仿真穩(wěn)定性。第三個問題是速度估計用Unit Delay引入的相位滯后。滯后會導(dǎo)致滑模面計算不準確誤差穩(wěn)態(tài)值可能達不到理論精度。解決方法是在控制器里對估計的角速度做一階低通濾波或者用Kalman濾波的簡化版。不過對于學(xué)習(xí)項目Unit Delay足夠讓你理解問題所在。5. 常見問題速查與避坑清單5.1 高頻報錯與解決代數(shù)環(huán)警告反饋路徑插入Unit Delay或Memory模塊注意采樣時間要與信號一致。MATL AB Function輸出未定義函數(shù)開頭給輸出賦零值。仿真速度極慢積分步長太小或控制器里有高頻切換可改用飽和函數(shù)并適當增大最大步長。末端軌跡發(fā)散多半是雅可比逆接近奇異或初始誤差太大改用阻尼最小二乘逆并做速度限幅。Scope顯示沒有信號檢查信號線是否連錯以及MATLAB Function是否設(shè)置成了“每步更新”而不是“每幀更新”如果是變步長仿真最好把Scope采樣時間設(shè)為“繼承”。5.2 獨家實戰(zhàn)心得我做完這個項目后最大的體會是學(xué)習(xí)滑??刂撇荒苤槐弛吔晒揭欢ㄒH手在Simulink里改變 (\lambda)、(\eta)、(\phi)觀察誤差和滑模面的變化。建議你做一個對比實驗其他參數(shù)不變把 (\eta) 從 0.1 調(diào)到 2看誤差收斂速度變快但抖振變強然后把符號函數(shù)換成飽和函數(shù)再看抖振被抑制但穩(wěn)態(tài)誤差略增大。這個過程10分鐘就能完成但比讀十遍論文都有用。另一個心得是運動學(xué)級滑模控制只是入門真正的工業(yè)應(yīng)用中滑??刂聘嘤迷趧恿W(xué)層或電機電流環(huán)。但是你把運動學(xué)模型下的滑模設(shè)計邏輯理清后再去看“滑模變結(jié)構(gòu)控制”的經(jīng)典文獻很多符號就不會那么嚇人了。比如文獻里的 (u u_{eq} u_{sw})等效控制項在運動學(xué)級就對應(yīng)把 (s) 的導(dǎo)數(shù)置零求出的控制量切換項對應(yīng) (\eta \text{sgn}(s))。這套對應(yīng)關(guān)系搞明白滑??刂凭退阏嬲腴T了。最后再分享一個小技巧仿真結(jié)束后用MATLAB命令窗口運行plot(logsout)可以快速查看所有記錄信號。如果你準備寫報告或論文建議把誤差均方根值也計算一下用rms(e_x)就能得到這個指標比肉眼觀察更有說服力。按照上面的步驟搭好模型參數(shù)可以先按我給的初始值跑一遍然后再改動各項參數(shù)觀察效果。相信我親手調(diào)過 (s) 曲線之后“滑?!边@兩個字就再也不會讓你犯怵了。