Matrix augmentation understanding problem I'm currently implementing a Model Predicitve Control (MPC) for my mobile robot and i'm using this paper as a reference. I'm currently having the problem of understanding the augmentation of a matrix mentioned in this paper (it's been a while since I worked with matrices). We have these two matrices

which should be augmented into these matrices

The first one is easy, it's just calculating A with the corresponding parameter k and append them. But for B I don't undertstand how this should be done. Here A is a 3x3 matrix and B is 3x2 matrix. First I cannot multiply them because of dimension mismatch and secondly even if B would be a 1x3 matrix how would I resize B(k|k) to fit in there?
 A: Let's expand it a little bit.
With a prediction horizon $N_p$ and control horizon $N_c$ (in your case $N_c=N_p=N$), the expanded explanation is as follows
The system dynamics:
\begin{align}
&\boldsymbol x(t+1) = \boldsymbol A \boldsymbol x(t) + \boldsymbol B \boldsymbol u (t) \nonumber \\
&\boldsymbol y(t) =  \boldsymbol C \boldsymbol x(t) \end{align}
The future states:
\begin{align}
&\boldsymbol x(t+1|t)=\boldsymbol A \boldsymbol x(t)+\boldsymbol B \boldsymbol u(t) \label{x-future-calc} \\
&\boldsymbol x(t+2|t)=\boldsymbol A \boldsymbol x(t+1|t)+\boldsymbol B \boldsymbol u(t+1)= \nonumber \\
&~~~~~~~~\boldsymbol A^2 \boldsymbol x(t)+\boldsymbol A\boldsymbol B \boldsymbol u(t) +\boldsymbol B\boldsymbol u(t+1) \nonumber\\
&\dots \nonumber\\
&\boldsymbol x(t+N_p|t)=\boldsymbol A^{N_p} \boldsymbol x(t)+\boldsymbol A^{N_p-1} \boldsymbol B \boldsymbol u(t)+ \nonumber \\
&~~~~~~~~\boldsymbol A^{N_p-2} \boldsymbol B \boldsymbol u(t+1)+\dots+ \nonumber \\
&~~~~~~~~\boldsymbol A^{N_p-N_c} \boldsymbol B \boldsymbol u(t+N_c-1) \nonumber
\end{align}
The future output:
\begin{align}
&\boldsymbol y(t+1|t)= \boldsymbol C\boldsymbol A \boldsymbol x(t)+\boldsymbol C\boldsymbol B \boldsymbol u(t) \label{y-future-calc} \\
&\boldsymbol y(t+2|t)= \boldsymbol C\boldsymbol A \boldsymbol x(t+1|t)+\boldsymbol B \boldsymbol u(t+1)=  \nonumber \\
&~~~~~~~~\boldsymbol C\boldsymbol A^2 \boldsymbol x(t)+\boldsymbol C\boldsymbol A\boldsymbol B \boldsymbol u(t) + \boldsymbol C\boldsymbol B\boldsymbol u(t+1) \nonumber\\
&\dots \nonumber\\
&\boldsymbol y(t+N_p|t)=\boldsymbol C\boldsymbol A^{N_p} \boldsymbol x(t)+\boldsymbol C\boldsymbol A^{N_p-1} \boldsymbol B \boldsymbol u(t)+ \nonumber \\
&~~~~~~~~\boldsymbol C\boldsymbol A^{N_p-2} \boldsymbol B \boldsymbol u(t+1)+\dots+\nonumber \\
&~~~~~~~~\boldsymbol C\boldsymbol A^{N_p-N_c} \boldsymbol B \boldsymbol u(t+N_c-1) \nonumber
\end{align}
Vectors of inputs and outputs:
\begin{align}
\boldsymbol Y=\begin{bmatrix}\boldsymbol y(t+1|t)\\\boldsymbol y(t+2|t)\\ \vdots \\ \boldsymbol y(t+N_p|t)\\\end{bmatrix} \in \mathbb{R}^{(N_pN_{out})\times 1} ,\quad
\boldsymbol U=\begin{bmatrix}\boldsymbol u(t) \\ \boldsymbol u(t+1) \\ \vdots \\ \boldsymbol u(t+N_c-1) \end{bmatrix} \in \mathbb{R}^{(N_cN_{in})\times 1}
\end{align}
Future dynamics in a neat form:
\begin{align} 
\boldsymbol Y=\tilde{ \boldsymbol A}\boldsymbol X(t)+\tilde{ \boldsymbol B} \boldsymbol U
\end{align}
where
$$\tilde{\boldsymbol A}=\begin{bmatrix}\boldsymbol CA \\ \boldsymbol C\boldsymbol A^2 \\ \boldsymbol C\boldsymbol A^3 \\ \vdots \\ \boldsymbol C\boldsymbol A^{N_p}\end{bmatrix}, \quad \tilde{ \boldsymbol B}=\begin{bmatrix}
\boldsymbol C\boldsymbol B & \boldsymbol 0 & \dots & \boldsymbol 0 \\
\boldsymbol C\boldsymbol A\boldsymbol B & \boldsymbol C\boldsymbol B & \dots & \boldsymbol 0 \\
\boldsymbol C\boldsymbol A^2\boldsymbol B & \boldsymbol C\boldsymbol A\boldsymbol B & \dots & \boldsymbol 0 \\
\vdots & \vdots & \ddots & \vdots \\
\boldsymbol C\boldsymbol A^{N_p-1}\boldsymbol B & \boldsymbol C\boldsymbol A^{N_p-2}\boldsymbol B & \dots & \boldsymbol C\boldsymbol A^{N_p-N_c}\boldsymbol B
\end{bmatrix}$$

Here is the point: You should not expect the matrix blocks to be square.
BTW, in your case, $\boldsymbol C=\boldsymbol I$.
Reference: Mohammadi, Arash, et al. "Optimizing Model Predictive Control horizons using Genetic Algorithm for Motion Cueing Algorithm." Expert Systems with Applications 92 (2018): 73-81. [sciencedirect] [researchgate]
