# Preface

This auxiliary section reminds the reader when a system of linear algebraic equations has a solution. Such systems are called consistent.

Introduction to Linear Algebra with Mathematica

# Mathematica solves linear systems

A system of m equations in n unknowns $$x_1 , \ x_2 , \ \ldots , \ x_n$$ is a set of m equations of the form

$$\label{EqSolve.1} \begin{cases} 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 , \\ \qquad \vdots & \qquad \vdots \\ a_{m1} \,x_1 + a_{m2}\, x_2 + \cdots + a_{mn}\, x_n &= b_m , \end{cases}$$
The numbers $$a_{ij}$$ are known as the coefficients of the system. The m×n matrix $${\bf A} = [\,a_{ij}\,] ,$$ whose $$(i,\, j)$$ entry is the coefficient $$a_{ij}$$ of the system of linear equations is called the coefficient matrix, and is denoted by
${\bf A} = \left[ \begin{array}{cccc} a_{1,1} & a_{1,2} & \cdots & a_{1,n} \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n} \end{array} \right] \qquad\mbox{or} \qquad {\bf A} = \left( \begin{array}{cccc} a_{1,1} & a_{1,2} & \cdots & a_{1,n} \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n} \end{array} \right) .$
Let $${\bf x} =\left[ x_1 , x_2 , \ldots x_n \right]^{\mathrm T}$$ be the vector of unknowns, where "T" stands for transposition. Then the product $${\bf A}\,{\bf x}$$ of the $$m \times n$$ coefficient matrix A and the $$n \times 1$$ column vector x is an $$m \times 1$$ matrix (column-vector)
$\left[ \begin{array}{cccc} a_{1,1} & a_{1,2} & \cdots & a_{1,n} \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n} \end{array} \right] \, \begin{bmatrix} x_1 \\ x_2 \\ \vdots \\ x_n \end{bmatrix} = \begin{bmatrix} a_{1,1} x_1 + a_{1,2} x_2 + \cdots + a_{1,n} x_n \\ a_{2,1} x_1 + a_{2,2} x_2 + \cdots + a_{2,n} x_n \\ \vdots \\ a_{m,1} x_1 + a_{m,2} x_2 + \cdots + a_{m,n} x_n \end{bmatrix} ,$
whose entries are the right-hand sides of our system of linear equations.

If we define another column vector $${\bf b} = \left[ b_1 , b_2 , \ldots b_m \right]^{\mathrm T}$$ whose components are the right-hand sides $$b_{i} ,$$ the given system is equivalent to the vector equation

$$\label{EqSolve.2} {\bf A} \, {\bf x} = {\bf b} .$$

We say that $$s_1 , \ s_2 , \ \ldots , \ s_n$$ is a solution of the system if all $$m$$ equations hold true when

$x_1 = s_1 , \ x_2 = s_2 , \ \ldots , \ x_n = s_n .$
Sometimes a system of linear equations is known as a set of simultaneous equations; such terminology emphasizes that a solution is an assignment of values to each of the $$n$$ unknowns such that each and every equation holds with this assignment. It is also referred to simply as a linear system.

Mathematica has a built-in command LinearSolve[A, b] that solves a linear system of equations, given a coefficient matrix A and a vector of constants b. We need to learn some more theory before we can entirely understand this command, but we can begin to explore its use. For now, consider these methods experimental and do not let it replace row-reducing augmented matrices:

${\bf A}\, {\bf x} = {\bf b} ,$
where A is a given matrix and b is specified column vector. The command LinearSolve[A, b] will provide the vector x of solutions to the linear system $${\bf A}\,{\bf x} = {\bf b}$$ of equations with coefficient matrix A and the vector of constants b.

The system of algebraic equations Ax = b may have a unique solution x, may have infinite many solutions, or may have no solution at all.

Example 1: Let us take a look at the system of algebraic equations

$\begin{cases} \phantom{2\,}x+2\,y &=3 , \\ 2\,x + 4\,y &=6 . \end{cases} \tag{1.1}$
The two equations define the same line on the x-y plane, meaning that we have infinitely many solutions x = 3−2y. The corresponding singular (its determinant is zero) matrix A and nonhomogeneous term b are
${\bf A} = \begin{bmatrix} 1 & 2 \\ 2 & 4 \end{bmatrix} , \qquad {\bf b} = \begin{bmatrix} 3 \\ 6 \end{bmatrix} .$
A = {{1, 2}, {2, 4}}
b = {3, 6}                                 (* define the right-hand side vector *)
Let's apply the LinearSolve[A, b] command to a system with no solutions:
x1 = LinearSolve[A, b]
Out[8]= {3, 0}
So we see that the vector v1 = [ 3, 0 ] gives a solution to system (1.1) of equations. Is this solution unique? To answer this question, we need to find the kernel (or null space) of matrix A:
NullSpace[A]
Out[10]= {{-2, 1}}
Hence, the null space is spanned on the vector [ -2, 1 ]T and the general solution to the given system of linear equations Ax = b is
$x = 3\,t , \quad y = t ,$
where t is an arbitrary real number. We can also find the general solution with the following commands:
A = {{1, 2}, {2, 4}}
b = {3, 6}
Solve[A.{x, y} == b]
{{y -> 3/2 - x/2}}

However, if we take another vector b = [−3, 6]T. Obviously, the new system of equations

$\begin{cases} \phantom{2\,}x+2\,y &=-3 , \\ 2\,x + 4\,y &=\phantom{-}6 . \end{cases} \tag{1.2}$
defines two parallel lines on x-y-plane that do not intersect, meaning there is no solut. Therefore, system (1.2) has no solution. We check whether Mathematica has a capability to determine this conclusion.
b2 = {-3, 6}
x2 = LinearSolve[A, b2]
Out[10]=
LinearSolve::nosol: Linear equation encountered that has no solution.
Therefore, the system has no solution for this vector b2 because matrix rank
MatrixRank[A]
Out[12]= 1
is 1, which does not match the dimensions 2×2 of the given matrix. ■

Example 2: Let us consider the singular matrix $${\bf A} = \begin{bmatrix} 2&\phantom{-}1&2 \\ 0&\phantom{-}1&0 \\ 1&-2& 1 \end{bmatrix} .$$ We want to determine when the linear equation $${\bf A} \, {\bf x} = {\bf b}$$ has a solution; namely, we need to find conditions on the vector b so the equation $${\bf A} \, {\bf x} = {\bf b}$$ has a solution. We use Mathematica to answer this question.

A = {{2, 1, 2}, {0, 1, 0}, {1, -2, 1}};
Det[A]
Out[2]= 0
The command MatrixRank[A] provides the rank of a matrix A.
MatrixRank[A]
Out[3]= 2
Atranspose = Transpose[A]
Out[4]= {{2, 0, 1}, {1, 1, -2}, {2, 0, 1}}
Atranspose = ConjugateTranspose[A]
Out[5]= {{2, 0, 1}, {1, 1, -2}, {2, 0, 1}}
We will need the latter command when matrix A would have complex entries.
x = NullSpace[Atranspose]
Out[6]= {{-1, 5, 2}}
So we get the null space of the transpose matrix to be spanned on the vector [-1, 5, 2 ]. The system of linear equations Ax = b has a solution (in this case, we say the the given system is consistent) if and only if the vector b is orthogonal to every vector from the kernel of the transpose (in general adjoint) matrix. In other words, its dot product with each such vector must be zero. So we ask Mathematica to do this job:
b = {b1, b2, b3};
Dot[x, y] == 0
Out[8]= {-b1 + 5 b2 + 2 b3} == 0
x.y == 0

Of course, we can achieve the same result by solving the system of linear equations Ax = b directly using Gaussian elimination method. Mathematica does not provide this algorithmitically fastest way to solve a linear algebraic equation; instead it uses Gauss--Jordan elimination procedure, which is more computationally demanded (and practically is not used). In this case, Mathematica offers a special command to achieve this: RowReduce[A]. ■

A linear system of algebraic equations is said to be consistent if it has at least one solution and inconsistent if it has no solution.
Thus, a consistent linear system has either one solution or infinitely many solutions---there are no other possibilities.

Theorem 1: Let A be an m × n matrix.

• (Overdetermined case) If m > n, then the linear system $${\bf A}\,{\bf x} = {\bf b}$$ is inconsistent for at least one vector b in $$\mathbb{R}^n .$$
• (Undetermined case) If m < n, then for each vector b in $$\mathbb{R}^m$$ the linear system $${\bf A}\,{\bf x} = {\bf b}$$ is either inconsistent or has infinite many solutions.

Theorem 2: If A is an $$m \times n$$ matrix with $$m < n,$$ then $${\bf A}\,{\bf x} = {\bf 0}$$ has infinitely many solutions.

Theorem 3: A linear system of algebraic equations $${\bf A}\,{\bf x} = {\bf b}$$ is consistent if and only if the vector b is orthogonal to every solution y of the adjoint homogeneous equation A*y = 0, where A* is the adjoint matrix to A. If A is real, then the adjoint matrix is just transposed and the statement is equivalent to $${\bf A}^{\mathrm T} \,{\bf y} = {\bf 0}$$ or $${\bf y}^T {\bf A} = {\bf 0} .$$

Theorem 3 provides us an algorithm of determination whether a system of equations is consistent or not.
Example 3: We consider a singular 3×3 matrix and 3-column vector b
${\bf A} = \begin{bmatrix} 1&\phantom{-}1&1 \\ 2&\phantom{-}3&1 \\ 0&-1&1 \end{bmatrix} \qquad \mbox{and} \qquad {\bf b} = \begin{bmatrix} 3\\ 4 \\ 2 \end{bmatrix} .$
A3 = {{1, 1, 1}, {2, 3, 1}, {0, -1, 1}}
Det[A3]
0
Next, we calculate its rank (number of linearly independent either rows or columns)
MatrixRank[A3]
2
To determine whether the system of linear equations Ax = b is consistent, we need to find the basis of null space of the adjoint matrix A*. Mathematica easily provides this answer:
NullSpace[Transpose[A3]]
{{-2, 1, 1}}
So we see that the kernel of AT is spanned on the vector
${\bf y} = \left[ -2, 1, 1 \right]^{\mathrm T} .$
Finally, we check the consistency of the system by calculating the dot product:
y = NullSpace[Transpose[A3]] // Flatten
b = {1, 2, 3}
y.b
3
Hence, we conclude that the given system Ax = [1, 2, 3]T is inconsistent. Mathematica confirms
LinearSolve[A3, b]
However, if we take another input vector b = [3, 4, 2]T, the corresponding system becomes consistent.    ■
Example 4: Let us consider the matrix A and vector b:
${\bf A} = \begin{bmatrix} \phantom{-}1&-2&3&-4 \\ -1&\phantom{-}1&1&\phantom{-}2 \\ \phantom{-}3&-5&5&-10 \\ \phantom{-}1 &-4&11&-8 \end{bmatrix} \qquad\mbox{amd} \qquad {\bf b} = \begin{bmatrix} 1 \\ 2 \\ 3 \\ 4 \end{bmatrix} .$
The matrix A has rank 2.
A = {{1, -2, 3, -4}, {-1, 1, 1, 2}, {3, -5, 5, -10}, {1, -4, 11, -8}}
MatrixRank[A]
2
In order to determine whether the given linear system Ax = b is consistent, we need to find the basis of null space of the adjoint matrix A*. Mathematica can help
A4 = {{1, -2, 3, -4}, {-1, 1, 1, 2}, {3, -5, 5, -10}, {1, -4, 11, -8}}
NullSpace[Transpose[A4]]
{{-3, -2, 0, 1}, {-2, 1, 1, 0}}
So we see that the kernel of AT is spanned on the vectors
${\bf y}_1 = \left[ -3, -2, 0, 1 \right]^{\mathrm T} \qquad\mbox{and}\qquad {\bf y}_2 = \left[ -2, 1, 1, 0 \right]^{\mathrm T} .$
Therefore, for solvability of the given system, we need to check two dot products.
y1 = { -3, -2, 0, 1 }
y2 = { -2, 1, 1, 0 }
b = {1, 2, 3, 4}
y1.b
-3
y2.b
3
Therefore, for the given column vector [1, 2, 3, 4]T, the given system is inconsistent. However, if we take another vector b4 = [2, 1, 3, 8]T, the situation will be different.
b4 = {2, 1, 3, 8}
y1.b4
0
y2.b4
0
So we conclude that the system with b4 = [2, 1, 3, 8]T is consistent. Finally, we check with Mathematica
LinearSolve[A4, b4]
{-4, -3, 0, 0}
■

Theorem 4: The general solution x of a consistent linear system of algebraic equations $${\bf A}\,{\bf x} = {\bf b}$$ can be expressed as the sum of a particular solution x0 of that system and the general solution of the corresponding homogeneous system

${\bf x} = {\bf x}_0 + c_1 {\bf v}_1 + c_2 {\bf v}_2 + \cdots + c_k {\bf v}_k .$
Here vector $$\{ {\bf v}_1 , \ldots , {\bf v}_k \}$$ is a basis for the kernel (= null space) of A and $$c_1 , \ldots , c_k$$ are arbitrary constants.

Theorem 5: Let A be an $$m \times n$$ matrix. A linear system of equations $${\bf A}\,{\bf x} = {\bf b}$$ is consistent if and only if the rank of the coefficient matrix A is the same as the rank of the augmented matrix $$\left[ {\bf A}\,\vert \, {\bf b} \right] .$$

The actual use of the term augmented matrix appears to have been introduced by the American mathematician Maxime Bôcher (1867--1918) in 1907.

Example 5: The linear system

$\begin{cases} x_1 + x_2 + x_3 + x_4 + x_5 &= 3 , \\ 2\,x_1 + x_2 + x_3 + x_4 + 2\, x_5 &= 4 , \\ x_1 - x_2 - x_3 + x_4 + x_5 &= 5 , \\ x_1 + x_4 + x_5 &= 4 , \end{cases} \tag{3.1}$
is an example of a system of four equations in five unknowns, $$x_1 , \ x_2 , \ x_3 , \ x_4 , \ x_5 .$$ One solution of this system is
$x_1 = -1 , \ x_2 = -2 , \ x_3 =1, \ x_4 = 3 , \ x_5 = 2 ,$
as you can easily verify by substituting these values into the equations. Every equation is satisfied for these values of $$x_1 , \ x_2 , \ x_3 , \ x_4 , \ x_5 .$$ However, there is not the only solution to this system of equations. There are many others.

On the other hand, the system of linear equations

$\begin{cases} x_1 + x_2 + x_3 + x_4 + x_5 &= 3 , \\ 2\,x_1 + x_2 + x_3 + x_4 + 2\, x_5 &= 4 , \\ x_1 - x_2 - x_3 + x_4 + x_5 &= 5 , \\ x_1 + x_4 + x_5 &= 6 , \end{cases} \tag{3.2}$
has no solution. There are no numbers we can assign to the unknowns $$x_1 , \ x_2 , \ x_3 , \ x_4 , \ x_5$$ so that all four equations are satisfied. We can solve these two equations simultaneously by considering the augmented matrix. An efficient way to do this is to form the partitioned matrix $$\left[ {\bf A} \, \vert \, {\bf b}_1 \, \vert {\bf b}_2 \right]$$ in which the coefficient matrix A is "augmented" by all 2 vectors. Then this augmented matrix is reduced using Gauss--Jordan elimination to
$\left( \begin{array}{ccccc|cc} 1&1&1&1&1&3&3 \\ 2&1&1&1&2&4&4 \\ 1&-1&-1&1&1&5&5 \\ 1&0&0&1&1&4&6 \end{array} \right) \quad \Longrightarrow \quad \left( \begin{array}{ccccc|cc} 1&0&0&0&1&1&0 \\ 0&1&1&0&0&-1&0 \\ 0&0&0&1&0&3&0 \\ 0&0&0&0&0&0&1 \end{array} \right) .$
The sixth column gives the required solution for the former system of equations. The last row in reduced form shows that the latter system of equations has no solution. ■

Equations in blocks

Suppose that we have the linear system $${\bf M}\,{\bf z} = {\bf r} ,$$ in which the given coefficient matrix M is partitioned into 2 by 2 blocks: $${\bf M} = \begin{bmatrix} {\bf A} & {\bf B} \\ {\bf C} & {\bf D} \end{bmatrix} ,$$ where A and D are square invertible matrices. Suppose that the right-hand side given vector is also broken as r = [ p q ]T, while column vector z = [ x y ]T is unknown. The system is equivalent to the pair of vector equations:

$$\label{EqSolve.3} \begin{split} {\bf A}\, {\bf x} + {\bf B}\, {\bf y} = {\bf p} , \\ {\bf C}\, {\bf x} + {\bf D}\, {\bf y} = {\bf q} . \end{split}$$
If D is nonsingular, eliminating y and solving for x yields
$$\label{EqSolve.4} {\bf x} = \left( {\bf A} - {\bf B}\,{\bf D}^{-1} {\bf C} \right)^{-1} \left( {\bf p} - {\bf B} \,{\bf D}^{-1} {\bf q} \right) = \left( {\bf M} / {\bf D} \right)^{-1} \left( {\bf p} - {\bf B} \,{\bf D}^{-1} {\bf q} \right) ,$$
where $${\bf M} / {\bf D} = {\bf A} - {\bf B}\,{\bf D}^{-1} {\bf C}$$ is the Schur complement. Similarly, if A is nonsingular, eliminating x and solving for y yields
Solve[{A x + B y == p, CC x + DD y == q}, {x, y}]
{{x -> -((DD p - B q)/(B CC - A DD)), y -> -((-CC p + A q)/(B CC - A DD))}}
$$\label{EqSolve.5} {\bf y} = \left( {\bf D} - {\bf C}\,{\bf A}^{-1} {\bf B} \right)^{-1} \left( {\bf q} - {\bf C} \,{\bf A}^{-1} {\bf p} \right) = \left( {\bf M} / {\bf A} \right)^{-1} \left( {\bf q} - {\bf C} \,{\bf A}^{-1} {\bf p} \right) ,$$
where $${\bf M} / {\bf A} = {\bf D} - {\bf C}\,{\bf A}^{-1} {\bf B}$$ is called the Schur complement. An important use of Schur complements is the partitioned solution of linear systems. In fact this is the most important application in Finite Element Method (FEM for short), in connection with superelement analysis.

Example 6: Let us consider a linear algebraic equation Mz = r with two input vectors r and b, where
${\bf M} = \begin{bmatrix} 2&5&-2&\phantom{-}1&2 \\ 3&4&-2& \phantom{-}1&2 \\ 1&3&\phantom{-}1&\phantom{-}3&1 \\ 2&1&-3&-2&1 \\ 4&2&-1&\phantom{-}3&1 \end{bmatrix} , \qquad {\bf r} = \begin{bmatrix} 1 \\ 2 \\ 3 \\ 4 \\ 5 \end{bmatrix} \qquad\mbox{and} \qquad {\bf b} = \begin{bmatrix} 1 \\ 7 \\ 3 \\ 4 \\ 5 \end{bmatrix} . \tag{6.1}$
We check whether the given matrix is singular or not.
M = {{2, 5, -2, 1, 2}, {3, 4, -2, 1, 2}, {1, 3, 1, 3, 1}, {2, 1, -3, -2, 1}, {4, 2, -1, 3, 1}}
Det[M]
0
Next, we check its rank
MatrixRank[M]
4
To determine the solvability of the system, we need to find a solution of the adjoint vector equation ATy = 0, where AT is the transposed matrix (because A is real, we don't need to take complex conjugate). First, we find the transposed matrix.
Transpose[M]
{{2, 3, 1, 2, 4}, {5, 4, 3, 1, 2}, {-2, -2, 1, -3, -1}, {1, 1, 3, -2, 3}, {2, 2, 1, 1, 1}}
${\bf M}^{\mathrm T} {\bf y} = {\bf 0} \qquad\mbox{or} \qquad \begin{bmatrix} \phantom{-}2&\phantom{-}3&1&\phantom{-}2&\phantom{-}4 \\ \phantom{-}5&\phantom{-}4&3&\phantom{-}1&\phantom{-}2 \\ -2&-2&1&-3&-1 \\ \phantom{-}1&\phantom{-}1&3&-2&\phantom{-}3 \\ \phantom{-}2&\phantom{-}2&1&\phantom{-}1&\phantom{-}1 \end{bmatrix} \begin{bmatrix} y_1 \\ y_2 \\ y_3 \\ y_4 \\ y_5 \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ 0 \\ 0 \\ 0 \end{bmatrix} . \tag{6.2}$
We find the null-space of Eq.{6.2} with a one line command:
NullSpace[Transpose[M]]
{{0, -1, 1, 1, 0}}
So the column-vector y is
${\bf y} = \left[ 0, -1, 1, 1, 0 \right]^{\mathrm T} . \tag{6.3}$
Since its dot product with given vector r is not zero, the system Mz = r is not consistent.
y = {0, -1, 1, 1, 0}
r = {1, 2, 3, 4, 5}
y.r
5
On the other hand, the system Mz = b is consistent because the dot product of vectors b and y is zero.

Equation in blocks

Since the system is too big to handle it manually, we break it into blocks:

${\bf M} = \begin{bmatrix} {\bf A} & {\bf B} \\ {\bf C} & {\bf D} \end{bmatrix} , \tag{6.4}$
where the partition of the given matrix M into submatrices A, B,C, and D can be chosen in four basic different ways. Of course, there are many other choices because the order of equations does not matter, but we stick with these four choices to demonstrate the method (due to Issai Schur). We consider every case separately.

Four-dimensional matrix A

We consider the case when matrix A is a 4×4 matrix:
${\bf A} = \begin{bmatrix} 2&5&-2&\phantom{-}1 \\ 3&4&-2& \phantom{-}1 \\ 1&3&\phantom{-}1&\phantom{-}3 \\ 2&1&-3&-2 \end{bmatrix} , \qquad {\bf B} = \begin{bmatrix} 2 \\ 2 \\ 1 \\ 1 \end{bmatrix} , \qquad {\bf C} = \begin{bmatrix} 4 & 2 & -1 & 3 \end{bmatrix} , \qquad {\bf D} = \begin{bmatrix} 1 \end{bmatrix} . \tag{6.4.1}$
Correspondingly, we partition the inbut vectors r and b:
${\bf r} = \begin{bmatrix} {\bf p}_4 \\ {\bf q}_4 \end{bmatrix} , \qquad\mbox{with} \qquad {\bf p}_4 = \begin{bmatrix} 1 \\ 2 \\ 3 \\ 4 \end{bmatrix} , \quad {\bf q}_4 = \begin{bmatrix} 5 \end{bmatrix} , \tag{6.4.2}$
and
${\bf b} = \begin{bmatrix} {\bf p}_4 \\ {\bf q}_4 \end{bmatrix} , \qquad\mbox{with} \qquad {\bf p}_4 = \begin{bmatrix} 1 \\ 7 \\ 3 \\ 4 \end{bmatrix} , \quad {\bf q}_4 = \begin{bmatrix} 5 \end{bmatrix} . \tag{6.4.3}$
The system is partitoned as
$\begin{cases} 2\,x_1 + 5\,x_2 -2\,x_3 + x_4 + 2\,x_5 &= 1 , \\ 3\,x_1 + 4\,x_2 -2\,x_3 + x_4 + 2\,x_5 &= 2 \quad \mbox{or} \quad 7 , \\ \phantom{-\,}x_1 + 3\,x_2 + x_3 + 3\,x_4 + \phantom{2}x_5 &= 3 , \\ 2\,x_1 + x_2 -3\,x_3 -2\,x_4 +\phantom{2}x_5 &= 4 \\ \hline 4\,x_1 + 2\,x_2 -x_3 + 3\,x_4 + \phantom{2}x_5 &= 5 . \end{cases}$
This leads to two equations for x and y:
$\left( {\bf A} - {\bf B}\,{\bf D}^{-1} {\bf C} \right) {\bf x} = {\bf p} - {\bf B}\,{\bf D}^{-1} {\bf q} \tag{6.4.4}$
and
$\left( {\bf D} - {\bf C}\,{\bf A}^{-1} {\bf B} \right) {\bf y} = {\bf q} - {\bf C}\,{\bf A}^{-1} {\bf p} . \tag{6.4.5}$
With our partition, these Schur complement matrices become
${\bf M} / {\bf D} = {\bf A} - {\bf B}\,{\bf D}^{-1} {\bf C} = \begin{bmatrix} -12 & -9 & -16 & -13 \\ -11&-10 &-16& -13 \\ -13& -11& -13& -11 \\ -12 & -13&-17&-16 \end{bmatrix} \qquad \Longrightarrow \qquad \left( {\bf M} / {\bf D} \right)^{-1} = \frac{1}{224} \begin{bmatrix} -84& 105& 35& 7 \\ 140&-119& -35& 7 \\ 84& -193& 11& 81 \\ -140&223& 43&-111 \end{bmatrix}$
A4 = {{2, 5, -2, 1}, {3, 4, -2, 1}, {1, 3, 1, 3}, {2, 1, -3, -2}}
B4 = {2,2,1,1}
C4 = {4,2,-1,3}
MD = A4 - B4.C4
{{-12, -9, -16, -13}, {-11, -10, -16, -13}, {-13, -11, -13, -11}, {-12, -13, -17, -16}}
Its inverse matrix exists because det(M/D) = -224.
Inverse[%]
{{-(3/8), 15/32, -(5/32), 1/32}, {5/8, -(17/32), -(5/32), 1/32}, {3/ 8, -(193/224), 11/224, 81/224}, {-(5/8), 223/224, 43/ 224, -(111/224)}}
but another Schur complement matrix does not exist because matrix A is singular.
Det[A4]
0
Therefore, our partition does not work: it is impossible to separate the problem into two independent systems of sizes 4 and 1.

Three-dimensional matrix A

We consider the case when matrix A is a 3×3 matrix:
${\bf A} = \begin{bmatrix} 2&5&-2 \\ 3&4&-2 \\ 1&3&\phantom{-}1 \end{bmatrix} , \qquad {\bf B} = \begin{bmatrix} 1 &2 \\ 1&2 \\ 3&1 \end{bmatrix} , \qquad {\bf C} = \begin{bmatrix} 2 & 1 & -3 \\ 4&2&-1 \end{bmatrix} , \qquad {\bf D} = \begin{bmatrix} -2&1 \\ \phantom{-}3&1 \end{bmatrix} .$
First, we check whether our square matrices are not singular:
A3 = {{2, 5, -2}, {3, 4, -2}, {1, 3, 1}}; Det[A3]
-15
D3 = {{-2, 1}, {3,1}}; Det[D3]
-5
Therefore, it is possible to separate the given 5×5 system of equations into two separate systems of sizes 3×3 and 2×2. The corresponding Schur complements are
${\bf M} / {\bf D} = {\bf A} - {\bf B}\,{\bf D}^{-1} {\bf C} = \begin{bmatrix} -4 & 2 & 2 \\ -3&1 & 2 \\ -3& 1& 2 \end{bmatrix}$
B3 = {{1, 2}, {1, 2}, {3, 1}}; C3 = {{2,1,-3}, {4,2,-1}};
MD = A3 - B3.Inverse[D3].C3
{{-4, 2, 2}, {-3, 1, 2}, {-3, 1, 2}}
Matrix M/D is singular
Det[MD]
0
and
${\bf M} / {\bf A} = {\bf D} - {\bf C}\,{\bf A}^{-1} {\bf B} = \frac{1}{3} \begin{bmatrix} 0&\phantom{-}0 \\ 4&-2 \end{bmatrix} .$
MA = D3 - C3.Inverse[A3].B3
{{0, 0}, {4/3, -(2/3)}}
Matrix M/A is also singular:
Det[MA]
0
In order to determine solvabilities of each system
$\left( {\bf M} / {\bf D} \right) {\bf x} = {\bf p} - {\bf B}\,{\bf D}^{-1} {\bf q} \qquad \mbox{and} \qquad \left( {\bf M} / {\bf A} \right) {\bf y} = {\bf q} -{\bf C}\,{\bf A}^{-1} {\bf p} \tag{6.3.1}$
we need to find the null-space of each Schur complement.
NullSpace[Transpose[MD]]
{{0, -1, 1}}
and
NullSpace[Transpose[MA]]
{{1, 0}}
So we see that the kernels of Scur complements are spanned on vectors
${\bf y}_3 = \left[ 0, -1, 1 \right]^{\mathrm T} \qquad \mbox{and} \qquad {\bf y}_2 = \left[ 1, 0 \right]^{\mathrm T} , \tag{6.3.2}$
respectively. Next, we calculate the right-hand sides of Eq.(6.3.1). For vector r, we have
${\bf p} - {\bf B}\,{\bf D}^{-1} {\bf q} = \left[ -8, -7, -2 \right]^{\mathrm T} , \qquad {\bf q} -{\bf C}\,{\bf A}^{-1} {\bf p} = \left[ 5, 1 \right]^{\mathrm T} .$
p3 = {1, 2, 3} - B3.Inverse[D3].{4, 5}
{-8, -7, -2}
This vector [-8, -7, -2]T is not orthogonal to y3 from Eq.(6,3,2). For another 2-dimensional system, we have
{4, 5} - C3.Inverse[A3].{1, 2, 3}
{5, 1}
This vector [5, 1]T is not orthogonal to y2 = [1, 0]T from Eq.(6,3,2), and the system Mz = r is inconsistent.

Hoverver, if we choose vector b, situation becomes different.

${\bf p} - {\bf B}\,{\bf D}^{-1} {\bf q} = \left[ -8, -2, -2 \right]^{\mathrm T} , \qquad {\bf q} -{\bf C}\,{\bf A}^{-1} {\bf p} = \left[ 0, -32 \right]^{\mathrm T} .$
p3 = {1, 7, 3} - B3.Inverse[D3].{4, 5}
{-8, -2, -2}
This vector [-8, -2, -2]T is orthogonal to y3 from Eq.(6,3,2). Another two-dimensional vector is
q2 = {4, 5} - C3.Inverse[A3].{1, 7, 3}
{0, -(32/3)}
which is orthogonal to y2. Therefore, the system Mz = b is consistent.

Two-dimensional matrix A

We consider the case when matrix A is chosen as a 2×2 matrix:
${\bf A} = \begin{bmatrix} 2&5 \\ 3&4 \end{bmatrix} , \qquad {\bf B} = \begin{bmatrix} -2&1&2 \\ -2&1&2 \end{bmatrix} , \qquad {\bf C} = \begin{bmatrix} 1 & 3 \\2&1 \\4&2 \end{bmatrix} , \qquad {\bf D} = \begin{bmatrix} \phantom{-}1&\phantom{-}3&1 \\ -3&-2&1 \\ -1& \phantom{-}3&1 \end{bmatrix} .$
First, we check that our square matrices are not singular:
A2 = {{2, 5}, {3, 4}}; D2 = {{1, 3, 1}, {-3, -2, 1}, {-1, 3, 1}};
Det[D2]
-10
Det[A2]
-7
So we can split the given 5-dimensional systems of equations Mz = r and Mz = b into two separate equations
$\left( {\bf M} / {\bf D} \right) {\bf x} = {\bf p}_2 - {\bf B}\,{\bf D}^{-1} {\bf q}_2 \qquad \mbox{and} \qquad \left( {\bf M} / {\bf A} \right) {\bf y} = {\bf q}_2 -{\bf C}\,{\bf A}^{-1} {\bf p}_2 , \tag{6.2.1}$
where
${\bf M} / {\bf D} = {\bf A} - {\bf B}\,{\bf D}^{-1} {\bf C} = \begin{bmatrix} -1 & 1 \\ \phantom{-}0&0 \end{bmatrix} \qquad \mbox{and} \qquad {\bf M} / {\bf A} = {\bf D} - {\bf C}\,{\bf A}^{-1} {\bf B} = \frac{1}{7} \begin{bmatrix} \phantom{-}15 & \phantom{-}17& -1 \\ -15& -17 & \phantom{-}1 \\ \phantom{-1}5& \phantom{-}15 & -5 \end{bmatrix} , \tag{6.2.2}$
B2 = {{-2, 1, 2}, {-2, 1, 2}}; C2 = {{1, 3}, {2, 1}, {4, 2}};
MD = A2 - B2.Inverse[D2].C2
{{-1, 1}, {0, 0}}
MA = D2 - C2.Inverse[A2].B2
{{15/7, 17/7, -(1/7)}, {-(15/7), -(17/7), 1/7}, {5/7, 15/7, -(5/7)}}
and vectors r and b are partitioned as
${\bf r} = \begin{bmatrix} {\bf p} \\ {\bf q} \end{bmatrix} \qquad\mbox{with} \qquad {\bf p}_2 = \begin{bmatrix} 1 \\ 2 \end{bmatrix}, \quad {\bf q}_2 = \begin{bmatrix} 3 \\ 4 \\ 5 \end{bmatrix} \qquad\mbox{and} \qquad {\bf b} = \begin{bmatrix} {\bf p} \\ {\bf q} \end{bmatrix} \qquad\mbox{with} \qquad {\bf p}_2 = \begin{bmatrix} 1 \\ 7 \end{bmatrix}, \quad {\bf q}_2 = \begin{bmatrix} 3 \\ 4 \\ 5 \end{bmatrix} ,$
The right-hand side vectors in Eq.(6.2.1) are
${\bf p}_2 - {\bf B}\,{\bf D}^{-1} {\bf q}_2 = \begin{bmatrix} -6 \\ -5 \end{bmatrix} , \qquad {\bf q}_2 -{\bf C}\,{\bf A}^{-1} {\bf p}_2 = \frac{1}{7} \begin{bmatrix} 18 \\ 17 \\ 13 \end{bmatrix} \tag{6.2.3}$
{1, 2} - B2.Inverse[D2].{3, 4, 5}
{-6, -5}
{3, 4, 5} - C2.Inverse[A2].{1, 2}
{18/7, 17/7, 13/7}
for vector r and
${\bf p}_2 - {\bf B}\,{\bf D}^{-1} {\bf q}_2 = \begin{bmatrix} -6 \\ \phantom{-}0 \end{bmatrix} , \qquad {\bf q}_2 -{\bf C}\,{\bf A}^{-1} {\bf p}_2 = \frac{1}{7} \begin{bmatrix} \phantom{-}23 \\ -23 \\ -67 \end{bmatrix} \tag{6.2.4}$
{1, 7} - B2.Inverse[D2].{3, 4, 5}
{-6, 0
{3, 4, 5} - C2.Inverse[A2].{1, 7}
{23/7, -(23/7), -(67/7)}
for input vector b. Finally, we check whether vectors (6.2.3) and (6.2.4) are orthogonal to vectors that span the null space of Schur complements:
${\bf y}_2 = \left[ 0, 1 \right]^{\mathrm T} \qquad \mbox{and} \qquad {\bf y}_3 = \left[ 1, 1, 0 \right]^{\mathrm T} \tag{6.2.5}$
NullSpace[Transpose[MA]]
{{1, 1, 0}}
NullSpace[Transpose[MD]]
{{0, 1}}
Since vectors (6.2.3), generated by components p and q of the vector r = [p, q]T, are not orthogonal (their dot products are not zero) to vectors (6.2.5), we claim that system Mz = r is inconsistent. On the other hand, vectors (6.2.4), generated by components p = [1, 7]T and q = ]3, 4, 5]T of the vector b = [p, q]T, are orthogonal (their dot products are zeroes) to vectors (6.2.5), and so we claim that system Mz = b is consistent.

One-dimensional matrix A

Finally, let us consider the case when matrix A is choisen as a 1×1 matrix:
${\bf A} = \begin{bmatrix} 2 \end{bmatrix} , \qquad {\bf B} = \begin{bmatrix} 5&-2&1&2 \end{bmatrix} , \qquad {\bf C} = \begin{bmatrix} 3 \\ 1 \\ 2 \\ 4 \end{bmatrix} , \qquad {\bf D} = \begin{bmatrix} 4&-2& \phantom{-}1&2 \\ 3&\phantom{-}1&\phantom{-}3&1 \\ 1&-3&-2&1 \\ 2&-1&\phantom{-}3&1 \end{bmatrix} .$
First, we check whether matrix D is invertible. Mathematica tells us that this is wrong
D1 = {{4, -2, 1, 2}, {3, 1, 3, 1}, {1, -3, -2, 1}, {2, -1, 3, 1}};
Det[D1]
0
So we conclude that such separation into 1-dimensional and 4-dimensional systems does not lead to an equivalent systems of algebraic equations.    ■