A = np.array([[1,2,3]])
B = np.array([[1],[2],[3]])2 Matrices and matrix multiplication
In a vector all components are listed sequentially. In a matrix they are listed in two directions: vertically and horizontally. We will use uppercase letters to represent matrices. If \(A\) is a matrix one writes \(a_{ij}\) for the component in row \(i\) and column \(j\), and we use parenthesis to indicate the start and end of the matrix: \[ A=\begin{pmatrix} a_{11} & a_{12} & \cdots & a_{1n} \\ a_{21} & a_{22} & \cdots & a_{2n} \\ \vdots & \vdots & \vdots & \vdots \\ a_{m1} & a_{m2} & \cdots & a_{mn}\end{pmatrix}. \] It is also assumed that all rows have an equal number of components, and similarly for the columns. We say that \(A\) is an \(m\times n\)-matrix (or has dimension \(m\times n\)) if it has \(m\) rows and \(n\) columns. \(A\) is called square if \(m=n\).
A horizontal segment in \(A\) is called a row vector, a vertical segment a column vector. A vector can be interpreted both as a row vector, and as a column vector. In the code
\(A\) is a row vector (or a \(1\times 3\)-matrix), \(B\) a column vector (or a \(3\times 1\)-matrix). We will as before write \(a_i\) for the components in a vector, regardless of whether it is a row vector or a column vector. Rows and columns will be indexed from \(0\) and upwards, just as for vectors:
A[0,1]np.int64(2)
B[1,0]np.int64(2)
In many math books it is common to interpret a vector as a column vector, if nothing else is said. When writing \(\mathbf{a}=(a_1,a_2,...,a_n)\) we shall therefore interpret \(\mathbf{a}\) as a column vector.
The \(3\times 4\)-matrix \(A=\begin{pmatrix} 1 & 2 & 3 & 4 \\ 5 & 6 & 7 & 8 \\ 9 & 10 & 11 & 12 \end{pmatrix}\) is obtained by writing
A = np.array([[1,2,3,4],[5,6,7,8],[9,10,11,12]])To find the dimensions of \(A\) write
m, n = np.shape(A)
m, n(3, 4)
Sometimes we will write \(\mathbf{a}_{i,:}\) for the row vector obtained from row \(i\) in \(A\), and \(\mathbf{a}_{:,j}\) for the column vector obtained from column \(j\) in \(A\).
A[i,:]
A[:,j]Square brackets are thus used to indicate the part of the matrix we want. Colon here stands for “all components”. Note the differences in the following two calls:
A[1,:]array([5, 6, 7, 8])
A[[1],:]array([[5, 6, 7, 8]])
The latter should be used if the result is to be used in computations involving matrices.
There is also a more general syntax for extracting a sub-matrix of a matrix. The code below extracts a \(2\times 3\)-matrix consisting of all components from row 1 and 3, which also are in column 1, 3, and 4.
A[np.ix_([0,2],[0,2,3])]array([[ 1, 3, 4],
[ 9, 11, 12]])
If \(A\) is an \(m\times n\)-matrix we define the transpose of \(A\) as the \(n\times m\)-matrix where component \((i,j)\) is \(a_{j,i}\). We write \(A^T\) for the transpose of \(A\). For the matrix \(A=\begin{pmatrix} 1 & 2 & 3 & 4 \\ 5 & 6 & 7 & 8 \\ 9 & 10 & 11 & 12 \end{pmatrix}\) above we get that \(A^T=\begin{pmatrix} 1 & 5 & 9 \\ 2 & 6 & 10 \\ 3 & 7 & 11 \\ 4 & 8 & 12 \end{pmatrix}\). Matrices are easily transposed as follows:
A.Tarray([[ 1, 5, 9],
[ 2, 6, 10],
[ 3, 7, 11],
[ 4, 8, 12]])
Hvis \(A=A^T\), \(A\) is said to be symmetric. In particular \[ A = \begin{pmatrix} 5 & 2 & 3 \\ 2 & 6 & 4 \\ 3 & 4 & 7\end{pmatrix}. \] is symmetric. Symmetric matrices are always square, and the components are symmetric about the diagonal. The diagonal in this example consists of the numbers \(5,6,7\).
The transpose satisfies simple rules: \[\begin{align*} (A+B)^T &= A^T + B^T & (A-B)^T&=A^T-B^T & \left(A^T\right)^T&=A & (sA)^T&=sA^T \end{align*}\] Some simple matrices will turn out to be useful.
- The zero-matrix, written \(\mathbf{0}\), is the matrix where all components are \(0\).
- The ones-matrix, written \(\mathbf{1}\), is the matrix where all components are \(1\).
- A diagonal matrix is a matrix where there are only zeroes outside the diagonal.
- The identity matrix, written \(I_n\), is the \(n\times n\)-diagonal matrix where all diagonal entries equal \(1\).
These can be obtained as follows:
np.zeros((m,n)) # Creates an mxn matrix of zeros
np.ones((m,n)) # Creates an mxn matrix of ones
np.diag(a) # Creates a diagonal matrix with the components
# of a on the diagonal.
np.eye(n) # Creates an nxn identity matrixMany operations give meaning both for vectors and matrices:
np.max(A,axis=0)array([ 9, 10, 11, 12])
np.min(A,axis=0)array([1, 2, 3, 4])
np.sum(A,axis=0)array([15, 18, 21, 24])
np.prod(A,axis=0)array([ 45, 120, 231, 384])
The difference from the vector version is that the operations are performed column by column.
2.1 Vectorised operations on matrices
As for vectors, matrix addition and subtraction are vectorised: If \(A\) and \(B\) are both \(m\times n\) then
- \(A+B\) is the matrix with component \((i,j)\) being \(a_{ij}+b_{ij}\),
- \(A-B\) is the matrix with component \((i,j)\) being \(a_{ij}-b_{ij}\),
- \(tA\) is the matrix with component \((i,j)\) being \(ta_{i,j}\) (where \(t\) is a scalar).
We thus have \[A + B = \begin{pmatrix} a_{11} + b_{11} & a_{12} + b_{12} & \cdots & a_{1n} + b_{1n} \\ a_{21} + b_{21} & a_{22} + b_{22} & \cdots & a_{2n} + b_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m1} + b_{m1} & a_{m2} + b_{m2} & \cdots & a_{mn} + b_{mn} \end{pmatrix}\] It is a simple exercise to check that the properties in Proposition 1.1 also hold if we replace vectors with matrices.
Example 2.1 We have that \[\begin{align*} 3\begin{pmatrix} 1 & 4 \\ 2 & 1 \end{pmatrix} +2\begin{pmatrix} -1 & 2 \\ -2 & 1 \end{pmatrix} &= \begin{pmatrix} 3\cdot 1 & 3\cdot 4 \\ 3\cdot 2 & 3\cdot 1 \end{pmatrix} + \begin{pmatrix} 2\cdot(-1) & 2\cdot 2 \\ 2\cdot(-2) & 2\cdot 1 \end{pmatrix} \\ &= \begin{pmatrix} 3 & 12 \\ 6 & 3 \end{pmatrix} +\begin{pmatrix} -2 & 4 \\ -4 & 2 \end{pmatrix} = \begin{pmatrix} 3-2 & 12+4 \\ 6-4 & 3+2 \end{pmatrix} \\ &= \begin{pmatrix} 1 & 16 \\ 2 & 5 \end{pmatrix} \end{align*}\] \(\clubsuit\)
Also for matrices Python allows vectorised operations with a scalar: If we write A + b where A is a matrix and b a scalar, \[
\begin{pmatrix}
a_{11} + b & a_{12} + b & \cdots & a_{1n} + b \\
a_{21} + b & a_{22} + b & \cdots & a_{2n} + b \\
\vdots & \vdots & \ddots & \vdots \\
a_{m1} + b & a_{m2} + b & \cdots & a_{mn} + b
\end{pmatrix}
\] is returned. In particular will A - 1 return the matrix where 1 is subtracted in all components, and 1/A will return the matrix with components \(1/a_{ij}\).
Vectorised functions are defined in the same way as for vectors. Let us consider some examples.
Example 2.2 We saw previously that the operations sqrt, abs, exp, and log, were vectorised. This is also the case for matrix input. This is easily seen for trigonometric functions: \[
\sin\begin{pmatrix} 0 & \pi/2 \\ \pi & 3\pi/2 \end{pmatrix}
= \begin{pmatrix} \sin 0 & \sin(\pi/2) \\ \sin\pi & \sin(3\pi/2) \end{pmatrix}
= \begin{pmatrix} 0 & 1 \\ 0 & -1 \end{pmatrix}
\] This is simple to verify:
np.sin(np.array([[0,np.pi/2],[np.pi,3*np.pi/2]]))array([[ 0.0000000e+00, 1.0000000e+00],
[ 1.2246468e-16, -1.0000000e+00]])
\(\clubsuit\)
Example 2.3 max and min are also vectorised for matrix input.
A = np.array([[1,2,3],[4,5,6]])
B = np.array([[4,5,6],[1,2,3]])
np.maximum(A,B)array([[4, 5, 6],
[4, 5, 6]])
np.minimum(A,B)array([[1, 2, 3],
[1, 2, 3]])
The same applies for multiplication and taking powers.
A*Barray([[ 4, 10, 18],
[ 4, 10, 18]])
A**Barray([[ 1, 32, 729],
[ 4, 25, 216]])
\(\clubsuit\)
A special syntax exists when \(A\) is a column vector (a) and \(B\) a row vector (b): a + b will then return \[
\begin{pmatrix}
a_1 + b_1 & a_1 + b_2 & \cdots & a_1 + b_n \\
a_2 + b_1 & a_2 + b_2 & \cdots & a_2 + b_n \\
\vdots & \vdots & \ddots & \vdots \\
a_m + b_1 & a_m + b_2 & \cdots & a_m + b_n
\end{pmatrix},
\] i.e., the matrix with components \((i,j)\) equal to \(a_i+b_j\). It is not really defined mathematically to add a column vector and a row vector, but some programming languages still allow this. This may sound strange, but it can be useful in some cases, since it can give very compact code. Let us give some examples.
Example 2.4 (The Hilbert matrix) The Hilbert matrix \(H_n\) is defined as the \(n\times n\)-matrix with components \[h_{kl}=1/(k+l-1),\] where \(1\leq k,l\leq n\). In particular \[
H_4=\begin{pmatrix} 1 & 1/2 & 1/3 & 1/4 \\ 1/2 & 1/3 & 1/4 & 1/5 \\ 1/3 & 1/4 & 1/5 & 1/6 \\ 1/4 & 1/5 & 1/6 & 1/7 \end{pmatrix}.
\] The matrix \(H_n\) can be obtained with help of the function scp.hilbert(n), but let us instead create it with the syntax above. First we write
y = np.arange(1,5)
yarray([1, 2, 3, 4])
Here a vector with 4 components was created, without interpretation in rows and columns. With the command reshape this can be interpreted as a row vector (a \(1\times 4\)-matrix):
y.reshape((1,-1))array([[1, 2, 3, 4]])
or a column vector (a \(4\times 1\)-matrix):
y.reshape((-1,1))array([[1],
[2],
[3],
[4]])
The number \(-1\) here means that the number of rows/columns is calculated automatically, from the requirement that there should be only one row/column. The Hilbert matrix can now be obtained as
x = y.reshape((1,-1))
1/(x + x.T - 1)array([[1. , 0.5 , 0.33333333, 0.25 ],
[0.5 , 0.33333333, 0.25 , 0.2 ],
[0.33333333, 0.25 , 0.2 , 0.16666667],
[0.25 , 0.2 , 0.16666667, 0.14285714]])
If we, the other way, would like to convert a row- or column vector back to a vector (without interpretation of rows/columns), this can also be done with the help of reshape:
x.reshape(-1)array([1, 2, 3, 4])
We will return to the Hilbert matrix later, as it will turn out to be quite vulnerable to roundoff errors on a computer. As a result it is difficult to use numerically. \(\clubsuit\)
Example 2.5 (The Vandermonde matrix) The matrix \[
V=
\begin{pmatrix}
1 & x_1 & \cdots & x_1^n \\
1 & x_2 & \cdots & x_2^n \\
\vdots & \vdots & \ddots & \vdots \\
1 & x_m & \cdots & x_m^n
\end{pmatrix}
\tag{2.1}\] is called a Vandermonde matrix, and \(x_1,\dots,x_m\) are called its generators. They are thus obtained by assembling the powers of the generators in each row. Such matrices need not be square - the number of generators (\(m\)) can be different from the number of powers (\(n+1\)). Vandermonde matrices will appear later in connection with interpolation, and numerical differentiation and -integration. If we collect the generators of the vandermonde matrix in the column vector x, and the powers \(0\),…,\(n\) in the row vector p, we get the matrix in Equation 2.1 by writing
p = np.arange(0,n+1).reshape((1,-1));
V = x**p\(\clubsuit\)
Example 2.6 (The discrete Fourier transform) Matrices can, as vectors, have complex components. The DFT matrix, \(F_N\), is an example of this. It is defined as the \(N\times N\)-vandermonde matrix with generators \(\{e^{-2\pi ik/N}\}_{k=0}^{N-1}\). \(\text{DFT}_N\) thus has components \(e^{-2\pi ikn/N}\), where \(0\leq k,n<N\), and it can be obtained with the code
p = np.arange(0,N).reshape((1,-1))
F = np.exp( -2*np.pi*1j*p.T*p )Here there is vectorised multiplication of a column vector with a row vector in the exponent, giving us all possible values for \(kn\). The exponential function is also vectorised, i.e., \(e^A\) gives the matrix with components \((k,n)\) equal to \(e^{a_{kn}}\).
The Fourier matrix has a very useful application to sound. In this setting \(N\) can be very large, so large that the matrix may not fit in the memory of the computer! It turns out that the Fourier matrix can be rewritten so that it can be used for computations after all. To be more precise, the function np.fft.fft(x) computes the multiplication \(F_N\mathbf{x}\), without actually writing down the full matrix \(F_N\) (fft here stands for Fast Fourier Transform). Here we should add some comments.
Both \(\mathbf{x}\) and \(F_N\mathbf{x}\) are column vectors. But it turns out that np.fft.fft(x) actually accepts a vector as input, and returns a vector. If one wants to compute the matrix-vector multiplication \(F_N\mathbf{x}\) (i.e., with both input and output assumed to be column vectors), one should instead call np.fft.fft(x, axis=0). These two things return the same, but the shape is either that of a vector, or that of a matrix:
x = np.array([1,2,3,4])
np.fft.fft(x)array([10.+0.j, -2.+2.j, -2.+0.j, -2.-2.j])
x = x.reshape((-1,1))
np.fft.fft(x, axis=0)array([[10.+0.j],
[-2.+2.j],
[-2.+0.j],
[-2.-2.j]])
Many other functions in numpy accept this form where the inputs and output are flattened to vectors. Since we will heavily use linear algebra notation, we will usually reshape the input to a column vector before we multiply with a matrix to get the resulting column vector. \(\clubsuit\)
2.2 Matrix multiplication
Later we will look at how one can solve linear systems with \(m\) equations and \(n\) unknowns. A general such system can be written \[ \begin{array}{lllll} a_{11}x_1 &+ a_{12}x_2 &+ \cdots &+ a_{1n}x_n &= b_1 \\ a_{21}x_1 &+ a_{22}x_2 &+ \cdots &+ a_{2n}x_n &= b_2 \\ \vdots & \vdots & \vdots & \vdots & \vdots\\ a_{m1}x_1 &+ a_{m2}x_2 &+ \cdots &+ a_{mn}x_n &= b_m \end{array} \tag{2.2}\] Let us write \[\begin{align*} A&= \begin{pmatrix} a_{11} & a_{12} & \cdots & a_{1n} \\ a_{21} & a_{22} & \cdots & a_{2n} \\ \vdots & \vdots & \vdots & \vdots\\ a_{m1} & a_{m2} & \cdots & a_{mn} \end{pmatrix} & \mathbf{x}&=\begin{pmatrix} x_1 \\ x_2 \\ \vdots \\ x_n \end{pmatrix} & \mathbf{b}&=\begin{pmatrix} b_1 \\ b_2 \\ \vdots \\ b_m \end{pmatrix} \end{align*}\] \(A\) is called the coefficient matrix of the system. In order to write Equation 2.2 more compactly the following definition is useful.
Definition 2.1 (Matrix-vector multiplication) Let \(A\) be an \(m\times n\)-matrix, and \(\mathbf{x}\) a column vector with \(n\) components. We define \(A\mathbf{x}\) as the vector \[A\mathbf{x} = \begin{pmatrix} a_{11}x_1 + a_{12}x_2 + \cdots + a_{1n}x_n \\ a_{21}x_1 + a_{22}x_2 + \cdots + a_{2n}x_n \\ \vdots \\ a_{m1}x_1 + a_{m2}x_2 + \cdots + a_{mn}x_n\end{pmatrix} = x_1\mathbf{a}_1 + x_2\mathbf{a}_2 + \cdots + x_n\mathbf{a}_n.\]
With this definition Equation 2.2 can be written more compactly as \(A\mathbf{x}=\mathbf{b}\). The rewriting to \(x_1\mathbf{a}_1 + x_2\mathbf{a}_2 + \cdots + x_n\mathbf{a}_n\) above follows from for the properties of vectors. Note that \(A\mathbf{x}\) is defined if and only if the number of columns in \(A\) equals the number of components in \(\mathbf{x}\), and that \(A\mathbf{x}\) then can be interpreted as av \(m\) scalar products (between the rows in \(A\), and \(\mathbf{x}\)).
Example 2.7 Let us compute \(A\mathbf{x}\) when \(A=\begin{pmatrix} 2 & 1 & -2 \\ 3 & 0 & -1 \\ 7 & -4 &2\end{pmatrix}\) and \(\mathbf{x}=\begin{pmatrix}3 \\ -1 \\ -2\end{pmatrix}\). This is well-defined since the number of columns in \(A\)
equals the number of components in \(\mathbf{x}\) (both are 3). We get \[\begin{align*}
\begin{pmatrix} 2 & 1 & -2 \\ 3 & 0 & -1 \\ 7 & -4 &2\end{pmatrix}
\begin{pmatrix}3 \\ -1 \\ -2\end{pmatrix}
&=
\begin{pmatrix}
2\cdot 3 + 1\cdot(-1) + -2\cdot(-2) \\ 3\cdot 3 + 0\cdot(-1) + -1\cdot(-2) \\ 7\cdot 3 + -4\cdot(-1)+ 2\cdot(-2)
\end{pmatrix}
=
\begin{pmatrix} 6 - 1 + 4 \\ 9 + 0 + 2 \\ 21 + 4 - 4 \end{pmatrix}
= \begin{pmatrix} 9 \\ 11 \\ 21 \end{pmatrix}
\end{align*}\] \(\clubsuit\)
Suppose we want to solve Equation 2.2 for different right hand sides \(\mathbf{b}_1,...,\mathbf{b}_k\): \[\begin{align*} A\mathbf{x}_1 &= \mathbf{b}_1 & A\mathbf{x}_2 &= \mathbf{b}_2 & & \cdots & A\mathbf{x}_k &= \mathbf{b}_k \end{align*}\] We seek a compact notation also for the collection of these systems. For this the following definition is useful.
Definition 2.2 (Matrix multiplication) If \(A\) is an \(m\times n\) matrix, and \(X\) an \(n\times k\)-matrix, we define \(AX\) as the \(m\times k\)-matrix \[ AX = \begin{pmatrix} AX_{:,1} & AX_{:,2} & \cdots & AX_{:,k} \end{pmatrix} \tag{2.3}\]
If we now set \[\begin{align*} X & = \begin{pmatrix} \mathbf{x}_1 & \mathbf{x}_2 & \cdots & \mathbf{x}_k\end{pmatrix} & B &= \begin{pmatrix} \mathbf{b}_1 & \mathbf{b}_2 & \cdots & \mathbf{b}_k \end{pmatrix}, \end{align*}\] Column \(j\) in the equation \(AX=B\) is \(A\mathbf{x}_j=\mathbf{b}_j\). \(AX=B\) can also be interpreted as the assembly of \(k\) systems. \(AX\) is only defined when the number of columns in \(A\) equals the number of rows in \(X\).
Let us see how a matrix product can be computed. From the definition of the matrix-vector product we know that \(AX_{:,j}=\sum_{r=1}^n a_{ir}x_{rj}\). Equation 2.3 now says that this equals column \(j\) in \(AX\), and it follows that \[ (AX)_{ij} = \sum_{r=1}^n a_{ir}x_{rj}. \tag{2.4}\] This is also called the inner form of a matrix product.
Example 2.8 Let us compute the matrix product \[\begin{align*} \begin{pmatrix} 2 & 4 & 3 \\ 1 & -1 & 3 \end{pmatrix}\begin{pmatrix} 2 & -1 \\ 4 & 1 \\ 2 & -1\end{pmatrix} \end{align*}\] The product should be a \(2\times 2\)-matrix. If we use Equation 2.4 to compute component \((1,1)\) in this we get \[2 \cdot 2 + 4\cdot 4 + 3\cdot 2 =4+16+6=26\] If we do the same for component \((1,2)\) we get \[2\cdot(-1) + 4\cdot 1 + 3\cdot(-1) = -2+4-3=-1.\] We have thus computed the first row in the matrix product. It is common to write out these calculations as follows. \[\begin{align*} \begin{pmatrix} 2 & 4 & 3 \\ 1 & -1 & 3 \end{pmatrix}\begin{pmatrix} 2 & -1 \\ 4 & 1 \\ 2 & -1\end{pmatrix} &= \begin{pmatrix} 2 \cdot 2 + 4\cdot 4 + 3\cdot 2 & 2\cdot(-1) + 4\cdot 1 + 3\cdot(-1) \\ 1\cdot 2 -1\cdot 4 + 3\cdot 2 & 1\cdot(-1) -1\cdot 1 + 3\cdot(-1) \end{pmatrix}\\ &= \begin{pmatrix} 4+16+6 & -2 + 4 -3 \\ 2-4+6 & -1-1-3 \end{pmatrix} = \begin{pmatrix} 26 & -1 \\ 4 & -5 \end{pmatrix} \end{align*}\] We can verify this with the built-in method for matrix multiplication.
A = np.array([[2,4,3],[1,-1,3]])
B = np.array([[2,-1],[4,1],[2,-1]])
A @ Barray([[26, -1],
[ 4, -5]])
Matrix multiplication can also be interpreted in terms of the scalar product: Above we computed 4 scalar products: If the left hand side matrix has rows \(\mathbf{a}_1\) and \(\mathbf{a}_2\), and the right hand side matrix har columns \(\mathbf{b}_1\) and \(\mathbf{b}_2\), we computed \[\begin{align*} \mathbf{a}_1\cdot\mathbf{b}_1&=26 & \mathbf{a}_1\cdot\mathbf{b}_2&=-1 & \mathbf{a}_2\cdot\mathbf{b}_1&=4 & \mathbf{a}_2\cdot\mathbf{b}_2&=-5. \end{align*}\] \(\clubsuit\)
Example 2.9 (The outer product) Let \(\mathbf{a}\in\mathbb{R}^m\) and \(\mathbf{x}\in\mathbb{R}^n\) be column vectors. We have that \[\begin{align*} \mathbf{a}\mathbf{x}^T&=\begin{pmatrix} a_1 \\ a_2 \\ \vdots \\ a_m \end{pmatrix}\begin{pmatrix}x_1 & x_2 & \cdots & x_n\end{pmatrix} = \begin{pmatrix} a_1x_1 & a_1x_2 & \cdots & a_1x_n \\ a_2x_1 & a_2x_2 & \cdots & a_2x_n \\ \vdots & \vdots & \ddots & \vdots \\ a_mx_1 & a_mx_2 & \cdots & a_mx_n \\ \end{pmatrix} \end{align*}\] \(\mathbf{a}\mathbf{x}^T\) is also called the outer product of \(\mathbf{a}\) and \(\mathbf{x}\). \(\clubsuit\)
Matrix multiplication satisfies several known properties:
Proposition 2.1 (Properties of matrix multiplication) Let \(A\), \(B\), and \(C\) be matrices. The following hold
- \((AB)C=A(BC)\)
- \(A(B+C)=AB+AC\)
- \((B+C)A=BA+CA\)
- For all scalars \(s\) we have that \(s(AB)=(sA)B=A(sB)\).
- If \(A\) is \(m\times n\) then \(AI_n=A\) and \(I_mA=A\) (so that the identity matrix behaves in the same way as the number \(1\) when multiplying real numbers).
It is required that the matrices have matching dimensions, so that the products are defined.
If we multiply many matrices, 1. here implies that it does not matter which ‘’neighbouring matrices’’ we multiply first. It therefore does not matter where we set the parenthesis when performing matrix multiplication. Therefore, the parenthesis are usually omitted, i.e., one writes \(ABC\) instead of \((AB)C=A(BC)\)). This, however, does not mean that we can change the order of the matrices: If \(A\) is \(m\times n\) and \(B\) is \(n\times m\), it is clear that \(AB\) and \(BA\) both are defined, but they have different dimensions when \(m\neq n\). Even if \(m=n\), we have usually that \(AB\neq BA\), as the following examples shows: \[\begin{align*} \begin{pmatrix} 1 & 0\\ 0 & 0 \end{pmatrix} \begin{pmatrix} 0 & 0 \\ 1 & 0\end{pmatrix} &= \begin{pmatrix} 0 & 0\\ 0 & 0 \end{pmatrix} & \begin{pmatrix} 0 & 0 \\ 1 & 0\end{pmatrix} \begin{pmatrix} 1 & 0\\ 0 & 0 \end{pmatrix} &= \begin{pmatrix} 0 & 0 \\ 1 & 0 \end{pmatrix}. \end{align*}\] Thus, while the order does not matter when multiplying scalars, this is not the case for matrices!
Proof.
- is perhaps the most difficult one to prove, since one obtains a ‘’double sum’’ when multiplying three matrices. Let us therefore take a closer look at this. If \(A\) is \(m\times k\), \(B\) \(k\times l\), \(C\) \(l\times n\), then both \((AB)C\) and \(A(BC)\) are well-defined. Furthermore, \(AB\) is \(m\times l\), and \(BC\) is \(k\times n\), and both \((AB)C\) and \(A(BC)\) are \(m\times n\). The definition of matrix multiplication gives \[\begin{align*} ((AB)C)_{ij} &= \sum_{s=1}^l (AB)_{is}c_{sj} = \sum_{s=1}^l \left( \sum_{r=1}^k a_{ir}b_{rs} \right) c_{sj} = \sum_{s=1}^l \sum_{r=1}^k a_{ir}b_{rs}c_{sj}\\ (A(BC))_{ij} &= \sum_{r=1}^k a_{ir}(BC)_{rj} = \sum_{r=1}^k a_{ir} \left( \sum_{s=1}^l b_{rs}c_{sj} \right) = \sum_{r=1}^k \sum_{s=1}^l a_{ir}b_{rs}c_{sj}. \end{align*}\] A couple of comments are in place:
- in the first equality on the first line the definition of the matrix product of \(AB\) and \(C\) was used,
- in the first equality on the second line the definition of the matrix product of \(A\) and \(BC\) was used,
- in the second equality on the first line the definition of the matrix product \(AB\) was used,
- in the second equality on the second line the definition of the matrix product \(BC\) was used.
The order of summation here does not matter, so that we get the same result. Therefore, all the components in \((AB)C\) and \(A(BC)\) are equal, so that \((AB)C=A(BC)\).
2.-4. is somewhat simpler to prove, and is gives as an exercise.
5.: We have that \[(AI_n)_{ij}=\sum_{r=1}^n a_{ir}(I_n)_{rj} = a_{ij}(I_n)_{jj}=a_{ij},\] where we used that \(I_n\) has only one non-zero in every column. Since all components are equal it follows that \(AI_n=A\). In the same way one shows that \(I_mA=A\).
There is also another useful form for a matrix product, which we also will return to. This form uses the notation \(\mathbf{a}_{:r}\) for column \(r\) in \(A\), and \(\mathbf{x}_{r:}\) for row \(r\) in \(X\):
Proposition 2.2 (Outer form for a matrix product) We have that \[\begin{align*} AX&=\sum_{r=1}^n \mathbf{a}_{:r}\mathbf{x}_{r:}. \end{align*}\] This is also called the outer form of a matrix product.
Proof. Let \(\mathbf{e}_i\) be the column vector which has \(1\) in component \(i\), and \(0\) elsewhere. From Example 2.9 it follows that
- \(\mathbf{a}_{:j} \mathbf{e}_j^T\) is the matrix of the same size as \(A\), where column \(j\) equals column \(j\) in \(A\), and where all other columns are \(\mathbf{0}\).
- \(\mathbf{e}_i\mathbf{x}_{i:}\) is the matrix with the same size as \(X\), where row \(i\) equals row \(i\) in \(X\), and where all other rows are \(\mathbf{0}\).
This means that we can write \(A=\sum_{j=1}^n \mathbf{a}_{:j} \mathbf{e}_j^T\), and \(X=\sum_{i=1}^n \mathbf{e}_i \mathbf{x}_{i:}\). If we use this we get \[\begin{align*} AX &= \left(\sum_{j=1}^n \mathbf{a}_{:j} \mathbf{e}_j^T\right)\left(\sum_{i=1}^n \mathbf{e}_i \mathbf{X}_{i:}\right) =\sum_{i=1}^n\sum_{j=1}^n \mathbf{a}_{:j} \mathbf{e}_j^T \mathbf{e}_i \mathbf{x}_{i:} = \sum_{r=1}^n \mathbf{a}_{:r} \mathbf{x}_{r:}. \end{align*}\] where we used that \(\mathbf{e}_j^T \mathbf{e}_i=1\) if \(i=j\), and \(0\) ellers. This completes the proof.
In the proof above several of the properties from Proposition 2.1 were used. Which ones?
2.3 Linear transformations and matrix multiplication
A mapping \(\mathbf{F}\) from \(\mathbb{R}^n\) to \(\mathbb{R}^m\) is called linear if we have that, for all \(\mathbf{a},\mathbf{a}_1,\mathbf{a}_2\in\mathbb{R}^n\), and all scalars \(c\),
- \(\mathbf{F}(\mathbf{a}_1+\mathbf{a}_2)=\mathbf{F}(\mathbf{a}_1) + \mathbf{F}(\mathbf{a}_2)\)
- \(\mathbf{F}(c\mathbf{a})=c\mathbf{F}(\mathbf{a})\)
From Proposition 2.1 it follows that the mapping \(\mathbf{x}\to A\mathbf{x}\), where \(A\) is a matrix, is linear. The other way, we will now show that all linear transformations can be written on this form:
Proposition 2.3 Suppose \(\mathbf{F}\) is a linear transformation from \(\mathbb{R}^n\) to \(\mathbb{R}^m\). There exists an \(m\times n\)-matrix \(A\) so that \(\mathbf{F}(\mathbf{x})=A\mathbf{x}\) for all \(\mathbf{x}\in\mathbb{R}^n\). More precisely, \[ A=\begin{pmatrix} \mathbf{F}(\mathbf{e}_1) & \mathbf{F}(\mathbf{e}_2) & \cdots & \mathbf{F}(\mathbf{e}_n) \end{pmatrix}, \tag{2.5}\] where \(\mathbf{e}_i\) was the column vector with \(1\) in component \(i\), and \(0\) otherwise.
Proof. For the proof we will use that \[(x_1,x_2,...,x_n) = x_1\mathbf{e}_1+x_2\mathbf{e}_2+\cdots+x_n\mathbf{e}_n,\] which is easy to show. This gives \[\begin{align*} \mathbf{F}(x_1,x_2,...,x_n) &= \mathbf{F}(x_1\mathbf{e}_1+x_2\mathbf{e}_2+\cdots+x_n\mathbf{e}_n) \\ &= x_1\mathbf{F}(\mathbf{e}_1)+x_2\mathbf{F}(\mathbf{e}_2)+\cdots+x_n\mathbf{F}(\mathbf{e}_n). \end{align*}\] By the definition of matrix-vector multiplication this can be written \[\begin{pmatrix} \mathbf{F}(\mathbf{e}_1) & \mathbf{F}(\mathbf{e}_2) & \cdots & \mathbf{F}(\mathbf{e}_n) \end{pmatrix}\begin{pmatrix} x_1 \\ x_2 \\ \vdots \\ x_n \end{pmatrix}.\] It follows that all linear transformations can be written on the form \(\mathbf{x}\to A\mathbf{x}\), with \(A\) given by Equation 2.5.
Thus, one only needs compute \(\mathbf{F}(\mathbf{e}_i)\) to find the matrix of a linear transformation. Let us take a look at some useful linear transformations, and find their matrices.
Example 2.10 (Projections) From the properties of the scalar product it follows that the projection operator \(\text{proj}_{\mathbf{b}}(\mathbf{a})=\frac{\mathbf{a}\cdot\mathbf{b}}{|\mathbf{b}|^2}\mathbf{b}\) is linear. The matrix of \(\text{proj}_{\mathbf{b}}\) is \[\begin{align*}
&\begin{pmatrix} \text{proj}_{\mathbf{b}}(\mathbf{e}_1) & \text{proj}_{\mathbf{b}}(\mathbf{e}_2) & \cdots & \text{proj}_{\mathbf{b}}(\mathbf{e}_n) \end{pmatrix} \\
&= \begin{pmatrix} \frac{\mathbf{e}_1\cdot\mathbf{b}}{|\mathbf{b}|^2}\mathbf{b} & \frac{\mathbf{e}_2\cdot\mathbf{b}}{|\mathbf{b}|^2}\mathbf{b} & \cdots & \frac{\mathbf{e}_n\cdot\mathbf{b}}{|\mathbf{b}|^2}\mathbf{b}\end{pmatrix}
= \begin{pmatrix} \frac{b_1}{|\mathbf{b}|}\frac{\mathbf{b}}{|\mathbf{b}|} & \frac{b_2}{|\mathbf{b}|}\frac{\mathbf{b}}{|\mathbf{b}|} & \cdots & \frac{b_n}{|\mathbf{b}|}\frac{\mathbf{b}}{|\mathbf{b}|}\end{pmatrix} \\
&= \begin{pmatrix} \frac{\mathbf{b}}{|\mathbf{b}|}\frac{b_1}{|\mathbf{b}|} & \frac{\mathbf{b}}{|\mathbf{b}|}\frac{b_2}{|\mathbf{b}|} & \cdots & \frac{b_n}{|\mathbf{b}|}\frac{\mathbf{b}}{|\mathbf{b}|}\end{pmatrix}
= \frac{\mathbf{b}}{|\mathbf{b}|} \begin{pmatrix} \frac{b_1}{|\mathbf{b}|} & \frac{b_2}{|\mathbf{b}|} & \cdots & \frac{b_n}{|\mathbf{b}|} \end{pmatrix}
= \frac{\mathbf{b}}{|\mathbf{b}|} \left(\frac{\mathbf{b}}{|\mathbf{b}|}\right)^T.
\end{align*}\] The second last equality used the definition of matrix multiplication, while the last equality used that \(\frac{\mathbf{b}}{|\mathbf{b}|} = \begin{pmatrix} b_1/|\mathbf{b}| \\ b_2/|\mathbf{b}| \\ \vdots \\ b_n/|\mathbf{b}|\end{pmatrix}\), so that
\(\left(\frac{\mathbf{b}}{|\mathbf{b}|}\right)^T = \begin{pmatrix} \frac{b_1}{|\mathbf{b}|} & \frac{b_2}{|\mathbf{b}|} & \cdots & \frac{b_n}{|\mathbf{b}|} \end{pmatrix}\). It follows that, if \(\mathbf{b}\) is a vector with length \(1\), the matrix \(\mathbf{b}\mathbf{b}^T\) computes the projection onto the line through the origin, with direction vector \(\mathbf{b}\). Such matrices are therefore called projection matrices. \(\clubsuit\)
Example 2.11 (Rotations) Consider the mapping \(\mathbf{F}\) which rotates a vector in \(\mathbb{R}^2\) with angle \(\theta\) . If we use polar form we get \[\begin{align*} \mathbf{F}\begin{pmatrix}r\cos\alpha \\ r\sin\alpha\end{pmatrix} &=\begin{pmatrix} r\cos(\alpha+\theta) \\ r\sin(\alpha+\theta)\end{pmatrix} =\begin{pmatrix} \cos\theta r\cos\alpha-\sin\theta r\sin\alpha \\ \sin\theta r\cos\alpha + \cos\theta r\sin\alpha \end{pmatrix} \\ &= \begin{pmatrix} \cos\theta & -\sin\theta \\ \sin\theta & \cos\theta \end{pmatrix} \begin{pmatrix} r\cos\alpha\\r\sin\alpha\end{pmatrix} \end{align*}\] It follows that \(\mathbf{F}\) is linear, with matrix \(\begin{pmatrix} \cos\theta & -\sin\theta \\ \sin\theta & \cos\theta\end{pmatrix}\). Assuming that \(\mathbf{F}\) is linear, it is easy to find the matrix of \(\mathbf{F}\): \[\begin{align*} \mathbf{F}(\mathbf{e}_1) &= \begin{pmatrix}\cos\theta \\ \sin\theta \end{pmatrix} & \mathbf{F}(\mathbf{e}_2) &= \begin{pmatrix}\cos(\theta+\pi/2) \\ \sin(\theta+\pi/2) \end{pmatrix} = \begin{pmatrix}-\sin\theta \\ \cos\theta \end{pmatrix}. \end{align*}\] \(\clubsuit\)
2.4 Plotting polynomials
Matrix multiplication can be used naturally to plot polynomials. The function plt.plot can be used to plot curves in the plane. The common syntax for plt.plot takes two arguments x and y, which contain the \(x\)- and \(y\)-coordinates of the plot points, respectively. The plt.plot function accepts many variants. For instance, if we want to plot the polynomials \[\begin{align*} q_1(x)&=x & q_2(x)&=x^2 & q_3(x)&=x^3 \end{align*}\] together, with \(m\) points distributed over the interval \([0,1]\), this can be achieved as follows:
x = np.linspace(0,1,m).reshape((-1,1))
plt.plot(x, x)
plt.plot(x, x**2)
plt.plot(x, x**3)
plt.legend(['q_1','q_2','q_3'])
The plt.plot commands above can also be assembled into one call as follows.
plt.plot(x, x, x, x**2, x, x**3)where 6 arguments are given, and where each pair represents a plot. The result can be seen in Figure 2.1. If the polynomials are of higher degree, and with many terms, things get more complicated. Let us take a look at how we can do this systematically by assembling the polynomial coefficients in vectors. Let the first polynomial be \[q_1(x)=c_0+c_1x+\cdots+c_nx^n,\] and let the plot points be \(x_1\),…,\(x_m\). We write \[\begin{align*} q_1(x_1) &= c_0+x_1c_1+\cdots+x_1^nc_n \\ q_1(x_2) &= c_0+x_2c_2+\cdots+x_2^nc_n \\ &\vdots \\ q_1(x_m) &= c_0 + x_mc_1 +\cdots + x_m^nc_n. \end{align*}\] By the definition of matrix-vector-multiplication, this can be written \[ \begin{pmatrix} q_1(x_1) \\ q_1(x_2) \\ \vdots \\ q_1(x_m) \end{pmatrix} = \begin{pmatrix} 1 & x_1 & \cdots & x_1^n \\ 1 & x_2 & \cdots & x_2^n \\ \vdots & \vdots & \ddots & \vdots \\ 1 & x_m & \cdots & x_m^n \end{pmatrix} \begin{pmatrix} c_0 \\ c_1 \\ \vdots \\ c_n \end{pmatrix}, \tag{2.6}\] where we recognise the Vandermonde matrix from Example 2.5. If we now have another polynomial, \(q_2(x)=d_0+d_1x+\cdots+d_nx^n\), we can in the same way write \[ \begin{pmatrix} q_2(x_1) \\ q_2(x_2) \\ \vdots \\ q_2(x_m) \end{pmatrix} = \begin{pmatrix} 1 & x_1 & \cdots & x_1^n \\ 1 & x_2 & \cdots & x_2^n \\ \vdots & \vdots & \ddots & \vdots \\ 1 & x_m & \cdots & x_m^n \end{pmatrix} \begin{pmatrix} d_0 \\ d_1 \\ \vdots \\ d_n \end{pmatrix}. \tag{2.7}\] By the definition of matrix multiplication Equation 2.6 and Equation 2.7 can be assembled into \[ \begin{pmatrix} q_1(x_1) & q_2(x_1) \\ q_1(x_2) & q_2(x_2) \\ \vdots & \vdots \\ q_1(x_m) & q_2(x_m) \end{pmatrix} = \begin{pmatrix} 1 & x_1 & \cdots & x_1^n \\ 1 & x_2 & \cdots & x_2^n \\ \vdots & \vdots & \ddots & \vdots \\ 1 & x_m & \cdots & x_m^n \end{pmatrix} \begin{pmatrix} c_0 & d_0 \\ c_1 & d_1 \\ \vdots & \vdots \\ c_n & d_n \end{pmatrix}, \tag{2.8}\] and in the same way if we have even more polynomials. The Vandermonde matrix on the right side can be computed as follows:
p = np.arange(0,n+1).reshape((1,-1))
V = x**pWe collect the coefficients of the three polynomials \(q_1\), \(q_2\), and \(q_3\), in column vectors q1, q2, and q3, and compute Equation 2.8:
q1 = np.array([[0],[1],[0],[0]])
q2 = np.array([[0],[0],[1],[0]])
q3 = np.array([[0],[0],[0],[1]])
y = V @ np.block([q1,q2,q3])This gives a matrix with three columns, where column 1 contains values for \(q_1\), column 2 for \(q_2\), and column 3 for \(q_3\). The polynomials can now be plotted together by writing
plt.plot(x, y)
plt.legend(['q_1','q_2','q_3'])This is a new syntax for the plt.plot-function, where
- the first argument is the column vectors with plot points, and
- the second argument is the matrix with as many rows as in the first argument, and 3 columns.
If x is a column vector and Y a matrix with as many rows as x and 3 columns, then plt.plot(x,Y) is the same as writing
plt.plot(x, V @ q1, x, V @ q2, x, V @ q3)This also applies more generally if we have more columns with polynomials: plt.plot(x,Y) plots all columns of Y against the vector x.
2.5 Block matrices
Matrices can also be used as building blocks for bigger matrices. If we have four matrices \(A\), \(B\), \(C\), \(D\), we write \[ \left(\begin{array}{c|c} A & B \\ \hline C & D \end{array}\right) \tag{2.9}\] for the matrix where the components in these four matrices are placed over, under, to the left, or to the right of each other as indicated. Lines are used to separate the blocks. We call this a block matrix. Block matrices can be created as follows:
np.block([[A,B],[C,D]])Equation 2.9 only gives meaning if
- \(A\) and \(B\) have equally many rows,
- \(C\) and \(D\) have equally many rows,
for then these can be put next to each other. In the same way \(A\) and \(C\) must have equally many columns, \(B\) and \(D\) have equally many columns, for then they can be placed under one-another. It is easy to generalise to block matrices with more than four blocks.
Example 2.12 If \[\begin{align*} A &=\begin{pmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \end{pmatrix} & B &= \begin{pmatrix} 2 & 2 \\ 2 & 2 \end{pmatrix} & C &= \begin{pmatrix} 3 & 3 & 3 \\ 3 & 3 & 3 \\ 3 & 3 & 3 \end{pmatrix} & D &= \begin{pmatrix} 4 & 4 \\ 4 & 4 \\ 4 &4 \end{pmatrix} \end{align*}\] then \[ \left(\begin{array}{c|c} A & B \\ \hline C & D \end{array}\right) = \begin{pmatrix} 1 & 1 &1 & 2 & 2 \\ 1 & 1 & 1 & 2 & 2 \\ 3 & 3 & 3 & 4 & 4 \\ 3 & 3 & 3 & 4 & 4 \\ 3 & 3 & 3 & 4 & 4 \end{pmatrix} \] This is easily verified:
A = np.array([[1,1,1],[1,1,1]])
B = np.array([[2,2],[2,2]])
C = np.array([[3,3,3],[3,3,3],[3,3,3]])
D = np.array([[4,4],[4,4],[4,4]])
np.block([[A,B],[C,D]])array([[1, 1, 1, 2, 2],
[1, 1, 1, 2, 2],
[3, 3, 3, 4, 4],
[3, 3, 3, 4, 4],
[3, 3, 3, 4, 4]])
\(\clubsuit\)
The product \[ \left(\begin{array}{c|c} A & B \\ \hline C & D \end{array}\right) \left(\begin{array}{c|c} E & F \\ \hline G & H \end{array}\right) \tag{2.10}\] of two block matrices give meaning if
- the number of columns in \(A\) equals the number of rows in \(E\),
- the number of columns in \(B\) equals the number of rows in \(G\).
Then the product can be computed as
np.block([[A,B],[C,D]]) @ np.block([[E,F],[G,H]])Equation 2.10 can also be computed in terms of ordinary matrix multiplication: \[ \left(\begin{array}{c|c} AE+BG & AF+BH \\ \hline CE+DG & CF+DH \end{array}\right) \] The code then becomes
np.block([[A@E +B@G,A@F+B@H],[C@E+D@G,C@F+D@H]]) Later we will in particular use the following, where the order of the matrices is switched: \[ \left(\begin{array}{c|c} A & B \end{array}\right) \left(\begin{array}{c} E \\ \hline G \end{array}\right) = \left(\begin{array}{c|c} B & A \end{array}\right) \left(\begin{array}{c} G \\ \hline E \end{array}\right) = AE+BG \tag{2.11}\]
Quiz
Exercises
Exercise 2.1 Compute
- \(\begin{pmatrix} 2 & 4 & 3 \\ -1 & 2 & 4 \\ 4 & 1 & -1\end{pmatrix} - \begin{pmatrix} - 1 & 2 & 8 \\ -1 & -2 & -5 \\ 2 & 4 & -2\end{pmatrix}\)
- \(2\begin{pmatrix} 3 & 2 \\ -1 & -2 \end{pmatrix} +4\begin{pmatrix} -3 & 3 \\ -1 & -1 \end{pmatrix}\).
Verify the result with your own python code:
Exercise 2.2 Transpose the following matrices:
- \(\begin{pmatrix} 2 & 4 \\ 3 & -1 \\ 0 & 1\end{pmatrix}\)
- \(\begin{pmatrix} 1 & 4 & 7 \\ 2 & 4 & 3 \\ -1 & -2 & 1 \end{pmatrix}\).
Verify the result with your own python code:
Exercise 2.3 Compute the following matrix products.
- \(\begin{pmatrix} 2 & 4 & 3 \\ 1 & -1 & -2 \\ -3 & -1 & 1 \end{pmatrix}\begin{pmatrix} 2 \\ 5 \\ 3 \end{pmatrix}\)
- \(\begin{pmatrix} 1 & 3 \\ 2 & 4\end{pmatrix}\begin{pmatrix} -1 & 0 \\ 3 & 2 \end{pmatrix}\)
- \(\begin{pmatrix} -1 & 2 & 4 & 1 \\ 1 & 3 & 2 & 4 \end{pmatrix}\begin{pmatrix} 2 & 1 \\ 2 & 0 \\ 1 & 0 \\ 3 & 4 \end{pmatrix}\)
Verify the result with your own python code:
Exercise 2.4 Let \(A\) be an \(n\times n\)-matrix. Show that \(\mathbf{x}^TA\mathbf{y}=0\) for all \(\mathbf{x},\mathbf{y}\in\mathbb{R}^n\) if and only if \(A=0\).
Exercise 2.5 Let \(A\) be an \(n\times n\)-matrix. Show that \((A\mathbf{x})\cdot\mathbf{y}=\mathbf{x}\cdot(A\mathbf{y})\) for all \(\mathbf{x},\mathbf{y}\in\mathbb{R}^n\) if and only if \(A\) is symmetric.
Exercise 2.6 The function np.random.rand(m,n) returns an \(m\times n\)-matrix with random components between 0 and 1, Consider the code
m = 10
n = 7
A = np.random.rand(m,n)
B = np.random.rand(m,n)
C = np.random.rand(m,n)
print( ((A+B)+C) - (A+(B+C)) )Which matrix property does this code verify?
Exercise 2.7 Show that the properties in Proposition 1.1 also hold if we replace the vectors \(\mathbf{a}\) and \(\mathbf{b}\) with matrices.
Exercise 2.8 Prove property 2.-4. in Proposition 2.1.
Exercise 2.9 How many arithmetic operations are required by matrix multiplication in general?
Exercise 2.10 We write \(A^n=AA\cdots A\), i.e., the matrix where we multiply in \(A\) \(n\) times. Compute \(A^2\), \(A^3\) and \(A^4\) when \(A=\begin{pmatrix} 1 & 4 \\ -2 & -1\end{pmatrix}\).
Verify the result with your own python code:
Exercise 2.11 (Mirroring) Let \(\mathbf{F}\) be the linear transformation from \(\mathbb{R}^2\) to \(\mathbb{R}^2\) which mirrors a point about the line \(y=x\). Find the matrix of \(\mathbf{F}\).
Exercise 2.12 Show that \((AB)^T=B^TA^T\), given that \(A\) has as many columns as \(B\) has rows.
Exercise 2.13 Assume that \(\mathbf{x}=(x_0,x_1,...,x_{N-1})\) is a real vector, and let \(\mathbf{y}=F_N\mathbf{x}\). Show that \(y_{N-n}=\overline{y_n}\) for \(1\leq n\leq N-1\).