A simulation: What's the difference between linear MPC and LQG?

Question

I thought it would be good to do a cooperation between basic MPC and basic LQR. So i hook up a state space model for simulation.

$$\dot x = Ax(t) + Bu(t) \\ y(t) = Cx(t) + Du(t)$$

Code:

delay = 0;
A = [0 1; -30 -3];
B = [0; 1];
C = [1 0];

% Create state space
sys = ss(delay, A, B, C); % D = Automatic

Which I convert to a discrete model.

$$x(k+1) = Ax(k) + Bu(k) \\ y(k) = Cx(k) + Du(k)$$

Code:

% Convert to discrete
h = 0.03;
sysd = c2d(sys, h)

To find the LQR gain matrix $L$ and kalman gain matrix $K$ i use the discrete Riccati equation, where $Q$ and $R$ are tuning matrices.

$$0 = A^TXA - X - A^TXB(R + B^TXB)^{-1}B^TXA + Q$$ $$L = (B^TXB + R)^{-1}(B^TXA)$$ $$K = (CXC^T + R)^{-1}(CXA^T)$$ Code:

% Find LQR control law
Q = [2 0; 0 1];
R = [3];
L = lqr(sysd, Q, R)

% Find kalman gain matrix K
Q = [1 0; 0 0.5];
R = [0.3];
K = lqe(sysd, Q, R)

Then I create signals for my LQG system. In this case, I will only write out MATLAB code because I don't know how to write the definition of this:

Code:

% Create signals
r = linspace(5, 5, 800); % Reference signal
d = linspace(0, 0, 800); % Disturbance signal - All zeros
n = 0.2*randn(1, 800); % White noise signal
t = linspace(0, 10, 800); % Time vector

Now I create my LQG state space model with state feedback:

$$\begin{bmatrix} x(k+1)\\ \tilde x(k+1) \end{bmatrix} = \begin{bmatrix} (A-BL) & BL\\ 0 & (A-KC) \end{bmatrix}\begin{bmatrix} x(k)\\ \tilde x(k) \end{bmatrix}+\begin{bmatrix} BK_r & B_d & 0\\ 0& B_d & -KB_n \end{bmatrix}\begin{bmatrix} r\\ d\\ n \end{bmatrix}$$

$$ \begin{bmatrix} y\\ y_m\\ u \end{bmatrix} = \begin{bmatrix} C &0 \\ C & 0\\ -L &L \end{bmatrix}\begin{bmatrix} x(k)\\ \tilde x(k) \end{bmatrix} + \begin{bmatrix} 0 & 0 &0 \\ 0 & 0& B_n\\ K_r & 0 &0 \end{bmatrix}\begin{bmatrix} r\\ d\\ n \end{bmatrix}$$

Where $K_r$ is the precompensator feed forward factor for better tracking.

$$K_r = \frac{1}{C(I-(A-BL))^{-1}B}$$

Code:

% Create LQG state space model
Bd = [0; 0]; % Disturbance matrix
Bn = [1]; % Noise matrix for mesurement
LQG = lqgreg(sysd, L, K, Bd, Bn);

And now I will simulate the system:

Code:

% Simulate
lsim(LQG, [r; d; n], t); % y1 = Filtered, y2 = Noisy, y3 = input signal

The result is:

And we will do the same with MPC. To create the contol law $L$ for MPC, i need to find matrix prediction matrices $F$ and $\Phi$ first.

$$F = \begin{bmatrix} CA\\ CA^2\\ CA^3\\ \vdots \\ CA^{N_p} \end{bmatrix} , \Phi = \begin{bmatrix} CB &0 &0 &\cdots & 0\\ CAB & CB & 0 & \cdots & 0\\ CA^2B& CAB & 0 &\cdots &0 \\ \vdots & \vdots & \vdots & \vdots &\vdots \\ CA^{N_p-1}B & CA^{N_p-2}B & CA^{N_p-3}B & \cdots & CA^{N_p-N_c}B \end{bmatrix}$$

I say that prediction horizon is $N_p = 10$ and control horizon is $N_c = 5$. To get the MPC control law $L$ and precompensator factor $K_r$. Source litterature: "Model Predictive Control System Design and Implementation Using MATLAB" by Liuping Wang, page 16.

$$K_r = \overset{N_c}{\begin{bmatrix} I_{mxm} & o_{mxm} & o_{mxm} & \cdots & o_{mxm} \end{bmatrix}}(\Phi^T\Phi + \bar R)^{-1}\Phi^T \bar R_s$$

$$L = \overset{N_c}{\begin{bmatrix} I_{mxm} & o_{mxm} & o_{mxm} & \cdots & o_{mxm} \end{bmatrix}}(\Phi^T\Phi + \bar R)^{-1}\Phi^T F$$

Where:

$$\bar R = R_{\omega} I_{N_c x N_c}$$

There $ R_{\omega}$ is a tuning factor and $I_{N_c x N_c}$ is a square identity matrix with dimension $N_c$ and $m$ is the number of output signals.

Also this is important too:

$$\bar R_s = \overset{N_p}{\begin{bmatrix} 1 & 1 & 1 & \cdots & 1 \end{bmatrix}^T}$$

Anyway! We will still use the LQG-state space model, with just a change of control law $L$ and precompensator factor $K_r$.

Code:

Np = 10; % Prediction horizon
Nc = 5; % Control horizon
Rw = 0.0001; % Tuning parameter "R omega".
MPC = mpcreg(sysd, Np, Nc, Rw, K, Bd, Bn);

And then do a new smulation:

Code:

figure(2) % new figure to plot on
lsim(MPC, [r; d; n], t); % y1 = Filtered, y2 = Noisy, y3 = input signal

And the result is:

As you can see. MPC model can predict the future states and do better tracking of trajectories.

Questions:

1. Is it suppose that $R_{\omega}$ should be so low as 0.0001?
2. What's the really difference between LQG and MPC in this case? Isn't LQR also a predictable controller?
3. Does my controller be better if I sett the control horizon and predict horizon even further?
4. Why are the result from MPC much more noisy that the LQG result? The observer is the same.
5. Is it better that I compute the $K_r$ for MPC by using $$K_r = \frac{1}{C(I-(A-BL))^{-1}B}$$

Instead?

Here is the code of mpcreg function.

https://github.com/DanielMartensson/matavecontrol/blob/master/sourcecode/mpcreg.m

You can also download my library and install it easy on your MALTAB or GNU Octave. https://github.com/DanielMartensson/matavecontrol

EDIT:

Here is a simulation with saturation of the states. Only available for Matavecontrol.

% Create matricies
delay = 0;
A = [0 1; -30 -3];
B = [0; 1];
C = [1 0];

% Create state space
sys = ss(delay, A, B, C); % D = Automatic

% Convert to discrete
h = 0.03;
sysd = c2d(sys, h)

% Find LQR control law
Q = [1 0; 2 0.5];
R = [2];
L = lqr(sysd, Q, R)

% Find kalman gain matrix K

K = lqe(sysd, Q, R)

% Create signals
r = linspace(1, 1, 800); % Reference signal 
d = linspace(0, 0, 800); % Discurbance signal - All zeros
n = 0.2*randn(1, 800); % White noise signal
t = linspace(0, 10, 800); % Time vector


% Create LQG state space model
Bd = [0; 0]; % Discurbance matrix
Bn = [1]; % Noise matrix for mesurement
[LQG, Kr] = lqgreg(sysd, L, K, Bd, Bn)

% Simulate
nlsim(LQG, [r; d; n], t, [0;0;0;0], [0 2; 0 2; 0 2; 0 2]); % y1 = Filtered, y2 = Noisy, y3 = input signal

% Now do MPC
Np = 10; % Prediction horizon
Nc = 5; % Control horizon
Rw = 0.0001; % Tuning parameter
[MPC, Lmpc, Krmpc] = mpcreg(sysd, Np, Nc, Rw, K, Bd, Bn);


% Simulate
figure(2) % new figure to plot on
nlsim(MPC, [r; d; n], t, [0;0;0;0], [0 2; 0 2; 0 2; 0 2]); % y1 = Filtered, y2 = Noisy, y3 = input signal

Results:

LQG

MPC