Computing matrices for Model Predictive Control

user_6014171 0
Viewed 94

In model predictive control, an optimization problem is solved at every time instant and it is very common to write down the matrices in a compact form. Without going into the details of the optimization problem, suppose I have matrices $A \in \mathbb{R}^{n\times n}$ and $B \in \mathbb{R}^{n\times m}$. I need to compute the matrices $\mathcal{A}$ and $\mathcal{B}$ defined as $$\mathcal{A}=\left(\begin{array}{c} I \ A \ A^{2} \ \vdots \ A^{N_p-1} \end{array}\right),: \mathcal{B}=\left(\begin{array}{ccccc} 0 & 0 & 0 & \ldots & 0 \ B & 0 & 0 & \ldots & 0 \ A B & B & 0 & \ldots & 0 \ \vdots & \ddots & \ddots & \ddots & \vdots \ A^{N_p-2} B & \ldots & A B & B & 0 \end{array}\right)$$

Note that N_p is called "prediction horizon" and it is not the order of matrix A. How can I compute these matrices in a fast and efficient way? In Matlab, I have done the following, but maybe there is a more efficient way to compute these matrices:

A_cal = zeros(length(A)*Np, length(A)); %calligraphic A matrix
B_cal = zeros(size(B,1)*Np, size(Bd,2)*Np); %calligraphic B matrix

temp = eye(size(A));

for j = 1:Np
    A_cal(1+(j-1)*length(A):j*length(A),:)= temp;
    if j > 1
        %The current row is obtained as shift of the previous row, and only the block in the first column is computed
        B_cal(1+(j-1)*size(B,1):j*size(B,1),:) = circshift(B_cal(1+(j-2)*size(B,1):(j-1)*size(B,1),:),size(B,2),2);
        B_cal(1+(j-1)*size(B,1):j*size(B,1),1:size(B,2)) = temp_prev*B_cent;
    end
    temp_prev = temp; %this variable contains A^(j-1)
    temp = temp * A_cent; %use temp variable to speed up the matrix power computation
end
matrix matlab
1 Answers

I assume you already solved your problem, but here is the code from my exercises sheet, creating your matrices.

S_x being the first matrix.

function S_x = compute_Sx(A,N)
  S_x = eye(size(A));
  for i=1:N
    S_x = [S_x;A^i];
  end
end

function S_u = compute_Su(A,B,N)
   S_u = zeros(size(A,1)*N,size(B,2)*N);
 for i=1:N
   S_u = S_u + kron(diag(ones(N-i+1,1),-i+1),A^(i-1)*B);
 end
S_u = [zeros(size(A,1),size(B,2)*N);S_u];
end

Best Regards

user_182547340
Related