matlab求x(k+1)=Ax(k)+Bu(k),y(t)=Cx(k)+Du(k-t)的增广系统矩阵

假设系统的状态向量为$x(t)$，输入向量为$u(t)$，输出向量为$y(t)$，则有：

$$
\begin{aligned}
\begin{bmatrix}
x(k+1)\\
y(k+1)
\end{bmatrix}
&=
\begin{bmatrix}
A & B\\
C & D
\end{bmatrix}
\begin{bmatrix}
x(k)\\
u(k-t)
\end{bmatrix}\\
&=
\begin{bmatrix}
A & B\\
C & D
\end{bmatrix}
\begin{bmatrix}
x(k)\\
u(k)-u(k-1)+u(k-1)-u(k-2)+\cdots+u(k-t+1)-u(k-t)
\end{bmatrix}\\
&=
\begin{bmatrix}
A & B\\
C & D
\end{bmatrix}
\begin{bmatrix}
x(k)\\
u(k)-u(k-1)+u(k-1)-u(k-2)+\cdots+u(k-t+1)-u(k-t)\\
\end{bmatrix}
\end{aligned}
$$

因此，增广系统矩阵为：

$$
\begin{bmatrix}
A & B\\
C & D
\end{bmatrix}
$$

相关问题

matlab求x(k+1)=Ax(k)+Bu(k),y(t)=Cx(k)+Du(k-t)的增广系统矩阵并求增广系统矩阵的谱半径

对于线性系统$x(k+1)=Ax(k)+Bu(k),y(t)=Cx(k)+Du(k-t)$，其增广系统可以表示为：

$$
\begin{bmatrix}
x(k+1)\\
y(k)
\end{bmatrix}
=
\begin{bmatrix}
A & 0\\
C & D
\end{bmatrix}
\begin{bmatrix}
x(k)\\
u(k-t)
\end{bmatrix}
$$

其中，$0$表示一个全零矩阵，其维度与$A$相同。

增广系统矩阵为$\begin{bmatrix}A & 0\\C & D\end{bmatrix}$，其谱半径为增广系统矩阵的最大特征值的模。可以使用Matlab中的`eig`函数求解增广系统矩阵的特征值，然后取其模的最大值即可。具体的实现代码如下：

% 系统参数
A = [1 1; 0 1];
B = [1; 0];
C = [1 0];
D = 0;
t = 1; % 时滞

% 增广系统矩阵
M = [A zeros(size(A)); C D];

% 计算谱半径
rho = max(abs(eig(M)));
disp(['增广系统的谱半径为：' num2str(rho)]);

其中，`size(A)`表示矩阵$A$的维度，`zeros(size(A))`表示一个与$A$维度相同的全零矩阵。`abs`函数用于取特征值的模，`max`函数用于取模的最大值。

一个稳定的连续时间 LTI系统,其输入信号为x(t)= cos(10t)+cos(20t)+cos(30t) 这个系统的微分方程如下 y(4)(t)+26.1y(3)(t) +341.4y"(t) +2 613.1y'(t)+10000y(t)= 10000x(t) ; 编写MATLAB程序,画出上述系统的幅度响应和在输入信号为x(t)情况下的零状态响应；

首先，我们可以将微分方程转换为状态空间模型：

$$
\begin{cases}
\dot{x}(t) = Ax(t) + Bu(t) \\
y(t) = Cx(t) + Du(t)
\end{cases}
$$

其中，$x(t)$是状态向量，$u(t)$是输入向量，$y(t)$是输出向量，$A,B,C,D$是系统矩阵。

对于本题中的系统，可以得到：

$$
\begin{cases}
\dot{x}(t) = \begin{bmatrix}
0 & 1 & 0 & 0 \\
0 & 0 & 1 & 0 \\
0 & 0 & 0 & 1 \\
-10000 & -2613.1 & -341.4 & -26.1
\end{bmatrix} x(t) + \begin{bmatrix}
0 \\
0 \\
0 \\
10000
\end{bmatrix} x(t) \\
y(t) = \begin{bmatrix}
1 & 0 & 0 & 0
\end{bmatrix} x(t)
\end{cases}
$$

接下来，我们可以使用MATLAB进行计算和绘图。代码如下：

% 定义系统矩阵
A = [0 1 0 0; 0 0 1 0; 0 0 0 1; -10000 -2613.1 -341.4 -26.1];
B = [0; 0; 0; 10000];
C = [1 0 0 0];
D = 0;

% 定义输入信号
t = 0:0.01:10;
x = cos(10*t) + cos(20*t) + cos(30*t);

% 计算幅度响应
sys = ss(A, B, C, D);
freqs = logspace(0, 3, 100);
[mag, ~] = bode(sys, freqs);
figure;
semilogx(freqs, 20*log10(mag));
grid on;
xlabel('Frequency (rad/s)');
ylabel('Magnitude (dB)');
title('Bode Plot of the System');

% 计算零状态响应
initial = [0; 0; 0; 0];
[y, ~, ~] = lsim(sys, x, t, initial);
figure;
plot(t, y);
grid on;
xlabel('Time (s)');
ylabel('Output');
title('Zero-State Response of the System');

运行以上代码，可以得到幅度响应和零状态响应的图像：

注：由于该系统是一个高阶系统，因此需要使用logspace函数生成频率向量。在绘制幅度响应图像时，使用semilogx函数将横坐标设置为对数坐标轴。在绘制零状态响应图像时，使用lsim函数计算系统的响应。