Matrices, Rank and Systems of Linear Equations

A system of linear equations is one matrix equation Ax = b. Whether it has a solution, and how many, is decided by a single number: the rank of A.

Matrix operations

For an m × n matrix A (m rows, n columns):

• Product: if A is m × n and B is n × p, then AB is m × p with (AB)ᵢⱼ = Σₖ aᵢₖbₖⱼ. Matrix multiplication is associative and distributive but not commutative: AB ≠ BA in general.
• Transpose: (Aᵀ)ᵢⱼ = aⱼᵢ, with (AB)ᵀ = BᵀAᵀ.
• Inverse: a square A is invertible if there is A⁻¹ with AA⁻¹ = A⁻¹A = I. Then (AB)⁻¹ = B⁻¹A⁻¹ and (Aᵀ)⁻¹ = (A⁻¹)ᵀ. For a 2 × 2 matrix, [[a, b], [c, d]]⁻¹ = (1/(ad − bc)) [[d, −b], [−c, a]].
• Trace: tr(A) = sum of the diagonal entries; tr(AB) = tr(BA).

The product Ax is a linear combination of the columns of A, with the entries of x as coefficients. This view explains most of what follows.

Row echelon form and Gaussian elimination

Elementary row operations do not change the solutions of a system:

1. swap two rows;
2. multiply a row by a non-zero number;
3. add a multiple of one row to another.

Gaussian elimination uses them to reach row echelon form: all zero rows at the bottom, and each row's first non-zero entry (its pivot) to the right of the pivot above it. Continuing until every pivot is 1 and is the only non-zero entry in its column gives the reduced row echelon form (RREF), which is unique.

Rank

The rank of A is the number of pivots in its row echelon form. Equivalently, it is:

• the maximum number of linearly independent rows, which always equals
• the maximum number of linearly independent columns ("row rank = column rank");
• the size of the largest square submatrix with non-zero determinant.

Properties:

• rank(A) ≤ min(m, n); A is full rank if equality holds.
• rank(A) = rank(Aᵀ) = rank(AᵀA).
• rank(AB) ≤ min(rank A, rank B).
• An n × n matrix is invertible ⟺ rank = n ⟺ det ≠ 0.
• A non-zero matrix uvᵀ (an outer product) has rank 1.

Example. A = [[1, 2, 1], [2, 4, 3], [3, 6, 4]]. Row 3 = row 1 + row 2, and rows 1 and 2 are not multiples, so rank(A) = 2. Elimination: R₂ − 2R₁ = (0, 0, 1) and R₃ − 3R₁ = (0, 0, 1), then R₃ − R₂ = (0, 0, 0). Two pivots.

The four subspaces and rank-nullity

• Column space C(A) ⊆ ℝᵐ: all vectors Ax, the span of the columns. dim = rank(A).
• Null space N(A) ⊆ ℝⁿ: all x with Ax = 0. Its dimension is the nullity.
• The row space C(Aᵀ) and the left null space N(Aᵀ) complete the picture. The null space is orthogonal to the row space.

Rank-nullity theorem: for an m × n matrix A,

rank(A) + nullity(A) = n (the number of columns)

In the example, rank 2 + nullity 1 = 3 columns. The null space is spanned by (−2, 1, 0), because −2 × column 1 + column 2 = 0.

Solving Ax = b

Form the augmented matrix [A | b] and compare ranks (the Rouché-Capelli theorem). With n unknowns:

• A homogeneous system Ax = 0 always has the solution x = 0. It has non-zero solutions exactly when rank(A) < n; for a square matrix, when det(A) = 0.
• Every solution of Ax = b is xₚ + xₙ: one particular solution plus anything from the null space.
• A system with more unknowns than equations (n > m) can never have exactly one solution.

Example. For which k and m does this system have no solution, one, or infinitely many?

x + y + z = 6, x + 2y + 3z = 14, x + 3y + kz = m.

Elimination gives y + 2z = 8 and then (k − 5)z = m − 22.

• k ≠ 5: exactly one solution.
• k = 5 and m = 22: infinitely many (the last row is 0 = 0).
• k = 5 and m ≠ 22: no solution (the last row is 0 = non-zero).

from fractions import Fraction as F

def rref(rows):
    M = [[F(x) for x in r] for r in rows]
    pivots, r = [], 0
    for c in range(len(M[0])):
        p = next((i for i in range(r, len(M)) if M[i][c] != 0), None)
        if p is None:
            continue
        M[r], M[p] = M[p], M[r]
        M[r] = [x / M[r][c] for x in M[r]]
        for i in range(len(M)):
            if i != r:
                M[i] = [a - M[i][c] * b for a, b in zip(M[i], M[r])]
        pivots.append(c)
        r += 1
    return M, pivots

A = [[1, 2, 1], [2, 4, 3], [3, 6, 4]]
R, piv = rref(A)
print("RREF:", [[str(x) for x in row] for row in R])
print("rank =", len(piv), " nullity =", len(A[0]) - len(piv))

def classify(k, m):
    aug = [[1, 1, 1, 6], [1, 2, 3, 14], [1, 3, k, m]]
    r_aug = len(rref(aug)[1])
    r_a = len(rref([row[:3] for row in aug])[1])
    if r_a < r_aug:
        return "no solution"
    return "unique solution" if r_a == 3 else "infinitely many"

for k, m in [(4, 20), (5, 22), (5, 20)]:
    print(f"k = {k}, m = {m}: rank(A) and rank([A|b]) -> {classify(k, m)}")
print("solution for k = 4, m = 20:", [str(row[3]) for row in rref([[1, 1, 1, 6], [1, 2, 3, 14], [1, 3, 4, 20]])[0]])

RREF: [['1', '2', '0'], ['0', '0', '1'], ['0', '0', '0']]
rank = 2  nullity = 1
k = 4, m = 20: rank(A) and rank([A|b]) -> unique solution
k = 5, m = 22: rank(A) and rank([A|b]) -> infinitely many
k = 5, m = 20: rank(A) and rank([A|b]) -> no solution
solution for k = 4, m = 20: ['0', '4', '2']

Practice questions