MATLAB中 控制工程会用到的
- 基础
- 传递函数表达
- 传递函数 零极点式 状态空间表达式 相互转化
- 画响应图线
- 根轨迹
- Nyquist图和bode图
- 现控部分
- 求约旦
- 判能控能观
- 极点配置和状态观测
基础
传递函数表达
% 拉普拉斯变换
syms t s a
f = exp(a*t) %e的a次方
l = laplace(f) %拉氏变换
f2 = ilaplace(f) %拉氏反变换
%传递函数
%多项式形式
num=[1,2] %表示传递函数分子系数 从右往左表示0次 1次 0次为2 1次为1
den=[1,2,3] %表示分母上的系数
sys = tf(num,den) %连续时间传递函数
%零极点形式
z = [1 2 3] %零点 用逗号空格都可以
p = [4 5 6] %极点
k = 15 %增益
sys = zpk(z,p,k) %传递函数
传递函数 零极点式 状态空间表达式 相互转化
%状态空间表达式
A = [1,2,3;4,5,6;7,8,9]
B = [0;1;1]
C = [1,2,3]
D = 3
sys = ss(A,B,C,D)
%状态空间表达式到多项式形式传递函数--- ss2tf
A = [1,2,3;4,5,6;7,8,9]
B = [0;1;1]
C = [1,2,3]
D = 3
[num,den] = ss2tf(A,B,C,D)
sys = tf(num,den)
%状态空间表达式生成零极点--- ss2zp
A = [1,2,3;4,5,6;7,8,9]
B = [0;1;1]
C = [1,2,3]
D = 3
[z,p,k] = ss2zp(A,B,C,D)
sys = zpk(z,p,k)
%多项式形式到零极点形式--- tf2zp
num = [1,2]
den = [1,2,3]
[z,p,k] = tf2zp(num,den)
%多项式形式到状态空间表达式-- tf2ss
num = [1,2]
den = [1,2,3]
[A,B,C,D] = tf2ss(num,den)
sys = ss(A,B,C,D)
%零极点形式得到状态空间表达式-- zp2ss
z = [1 2 3] %零点
p = [4 5 6] %极点
k = 15
[A,B,C,D] = zp2ss(z,p,k)
sys = ss(A,B,C,D)
%零极点形式得到多项式形式-- zp2tf
z = [1;2;3] %零点
p = [4 5 6] %极点
k = 15
sys = zpk(z,p,k)
[num,den]= zp2tf(z,p,k)
sys = tf(num,den)
画响应图线
%时域响应部分 求单位阶跃响应和单位脉冲响应
num = [1];
den = [1,0.4];
t = [0:0.1:20];
% step 函数用于计算由分子 num 和分母 den 定义的系统在时间向量 t 上的阶跃响应
y = step(num,den,t);
%y = impulse(num,den,t); % 单位脉冲响应
plot(t,y) % 以时间 t 为横轴，阶跃响应 y 为纵轴绘制曲线
grid on % grid on 用于在绘图中显示网格线，方便观察数据
xlabel('时间/s')
ylabel('y')
title('单位阶跃响应曲线')
legend('单位阶跃响应曲线') %在坐标图例添加标签
根轨迹
%根轨迹
num = [1,2]
den = [1,2,3]
sys = tf(num,den)
pzmap(sys) %画出零极点位置
rlocus(sys) %画出根轨迹
num = [1,2]
den = conv([1,0],conv([0.5,1],[4,1]))
%卷积 实际上就是算S*[(0.5s+1)(4s+1)]
sys = tf(num,den)
k = [0,0.1,10] %就是给k一个范围
rlocus(sys,k)
Nyquist图和bode图
%频率特性奈奎斯特图
num = 5
den = [3,1]
G = tf(num,den)
figure % figure 函数用于打开一个新的图形窗口，后续的绘图操作将在这个窗口中进行
% 在图形窗口中创建子图布局
% subplot(2,1,1) 表示将图形窗口划分为 2 行 1 列的子图网格，并激活第 1 个子图（即第一行的子图）
subplot(2,1,1)
nyquist(G)
title('nyquist图')
grid on % 在当前子图中显示网格线
%绘制伯德图
subplot(2,1,2)
bode(G)
title('bode图')
grid on
%频域性能指标：赋值裕度（没有Db单位得） 相角裕度，穿越频率，截止频率。
num = 5
den = [3,1]
G = tf(num,den)
[Gm,pm,Wcg,wcp] = margin(G) %一串直接求了
现控部分
求约旦
判能控能观
%判断能控能观性
A = [1,2,3;4,5,6;7,8,9]
B = [0;1;1]
C = [1,2,3]
D = 1
M = [B,A*B,A*A*B]
M = ctrb(A,B) %计算可控性矩阵
a = rank(M)
if a == 3
disp('可控');
else
disp('no')
end
%判断能观性
N = obsv(A,C) %计算能观判别矩阵
b = rank(N)
极点配置和状态观测
%状态反馈 极点配置
A = [1,2,3;4,5,6;7,8,9]
B = [0;1;1]
C = [1,2,3]
D = 1
%期望闭环极点
p = [1,2,4]
%适用place函数进行任意极点配置
%当然也要先判能控
K = place(A,B,p)
%状态观测器
A = [1,2,3;4,5,6;7,8,9]
B = [7;2;8]
C = [3,1,3]
D = 1
N = obsv(A,C)
N = [C;C*A;C*A*A]
b = rank(N)
p = [-5,-6,-7]
G = place(A',C',p)
