Model prediction gives stochastic inputs - Predictive Control

Question

I have a discrete state space model:

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

And I'm trying to compute the predicted inputs. The first thing I do is that I fist create the extended observability matrix $\Phi$

$$\Phi = \begin{bmatrix} CA\\ CA^2\\ CA^3\\ \vdots \\ CA^{n-1} \end{bmatrix}$$

Then create the lower traingular topelitz matrix $$\Gamma = \begin{bmatrix} CB & 0 & 0 & 0 &0 \\ CAB & CB & 0 & 0 &0 \\ CA^2B & CAB & CB & 0 & 0\\ \vdots & \vdots &\vdots & \ddots & 0 \\ CA^{n-2} B& CA^{n-3} B & CA^{n-4} B & \dots & CA^{n-j} B \end{bmatrix}$$

Where $j=n$

For SISO-case, as in this case. I compute the predicted inputs signals by using:

$$U = (\Gamma)^{-1}(R-\Phi x)$$

For MIMO case. I need to use this equation:

$$U = (\Gamma^T \Gamma)^{-1}\Gamma^T(R-\Phi x)$$

Where $R$ is my reference vector and $x$ is the initial state. I have matlab code for that.

function [U] = predict (A, B, C, x, N, r)

  ## Find matrix
  PHI = phiMat(A, C, N);
  GAMMA = gammaMat(A, B, C, N);
  U = inv(GAMMA'*GAMMA)*GAMMA'*(repmat(r, N, 1) -PHI*x);
endfunction

function PHI = phiMat(A, C, N)

  ## Create the special Observabillity matrix
  PHI = [];
  for i = 1:N
    PHI = vertcat(PHI, C*A^i);
  endfor

endfunction

function GAMMA = gammaMat(A, B, C, N)

  ## Create the lower triangular toeplitz matrix
  GAMMA = [];
  for i = 1:N
    GAMMA = horzcat(GAMMA, vertcat(zeros((i-1)*size(C*A*B, 1), size(C*A*B, 2)),cabMat(A, B, C, N-i+1)));
  endfor

endfunction

function CAB = cabMat(A, B, C, N)

  ## Create the column for the GAMMA matrix
  CAB = [];
  for i = 0:N-1
    CAB = vertcat(CAB, C*A^i*B);
  endfor

endfunction

If I have a state space model, discrete with sample time $h=0.5$.

A =

   0.89559   0.37735
  -0.37735   0.51825

B =

   0.10441
   0.37735

C =

   1   0

When I run this code above. I get this result if I choose $N = 40$

   239.4510
  -301.5661
   301.1320
  -208.4868
   222.4277
  -141.9374
   166.1560
   -94.3562
   125.9231
   -60.3368
    97.1576
   -36.0137
    76.5909
   -18.6233
    61.8862
    -6.1896
    51.3728
     2.7002
    43.8559
     9.0562
    38.4815
    13.6006
    34.6389
    16.8497
    31.8916
    19.1727
    29.9273
    20.8336
    28.5229
    22.0211
    27.5188
    22.8702
    26.8009
    23.4772
    26.2876
    23.9113
    25.9206
    24.2216
    25.6582
    24.4434

The correct input signals that tracks the reference is 25. So the algorithm computes correct at last. But when I heard prediction, I was told that I should select the first value of the input signals, not the last one. The first value is 239.4510 and the last value is 24.4434.

Questions:

• I'm I wrong according to prediction that I should use a lower triangular toeplitz matrix?

• Should prediction, by solving linear least squares, looks like this? In this case, I think I need to select the last value.

• If I have a model with a pole equals 1, e.g integration. How would this affect this prediction?

For Model Predictive Control by using a QP-solver, then I can use the first value. But it this case...no.

EDIT:

Here is a simulation when we say that the input signals cannot be $ u < 0$

Answer 1

I have not calculated the formulation by myself. But, I can certainly tell you that you should use the head of $\boldsymbol{U}$ instead of its tail.

Using the tail of $\boldsymbol{U}$ is definitely wrong. Look at the trend of the optimal control inputs. They start shrinking. The system looks stable. But, if you read the inputs upside-down then the system's input is growing (with oscillation). That's scary and looks unstable.

Ignore the the number 25. Simulate the system as it is and then see how the system behaves.