Matrices and systems

Why this matters

Real problems rarely have one unknown. A truss has a force in every member, an electrical network has a current in every branch, a chemical equation has a coefficient for every species, and a regression model has a weight for every feature. Each gives a system of linear equationsSeveral linear equations in the same unknowns that must all hold at once., and matrices are how engineers, scientists and software solve them, often with thousands of unknowns at once.

What a matrix is

A matrixA rectangular grid of numbers with m rows and n columns. is a rectangular grid of numbers with m rows and n columns (an m × n matrix):

A = | 2   1 |
    | 1   3 |

A vector can be thought of as a matrix with one column.

Matrix times vector

To multiply A by a vector x, take the dot product of each row of A with x:

| 2  1 | | x |   | 2x + y  |
| 1  3 | | y | = | x + 3y  |

So the whole system

2x +  y = 5
 x + 3y = 10

is one equation: A x = b, with b = (5, 10). The matrix holds the coefficients, x holds the unknowns, and b holds the right-hand sides.

Matrix times matrix

To multiply A (m × n) by B (n × p), entry (i, j) of the result is row i of A dotted with column j of B. The inner sizes must match, and the result is m × p. Order matters: in general AB ≠ BA.

| 1  2 | | 5  6 |   | 1·5 + 2·7   1·6 + 2·8 |   | 19  22 |
| 3  4 | | 7  8 | = | 3·5 + 4·7   3·6 + 4·8 | = | 43  50 |

Solving A x = b

By hand, elimination works: subtract 2 × (second equation) from the first to get −5y = −15, so y = 3 and x = 1. Check: 2·1 + 3 = 5 and 1 + 3·3 = 10.

Conceptually, if A has an inverseThe matrix A⁻¹ that undoes A, so A⁻¹A is the identity. A⁻¹ (the matrix that undoes A), then x = A⁻¹b. For a 2 × 2 matrix:

A = | a  b |      A⁻¹ = 1/(ad − bc) · |  d  −b |
    | c  d |                          | −c   a |

The number ad − bc is the determinantA number computed from a square matrix (ad − bc for 2 × 2); if it is 0 the matrix has no inverse.. If it is 0 the matrix has no inverse: the equations are either contradictory or redundant, so there is no unique solution.

In practice you never compute an inverse to solve a system. You call a solver that uses elimination (LU decompositionThe elimination-based factorisation that solvers such as np.linalg.solve use instead of computing an inverse.), which is faster and more accurate:

import numpy as np

A = np.array([[2.0, 1.0],
              [1.0, 3.0]])
b = np.array([5.0, 10.0])

x = np.linalg.solve(A, b)
print(x)                       # [1. 3.]
print(np.allclose(A @ x, b))   # True: always check the residual

@ is Python's matrix multiplication operator. np.linalg.det(A) gives the determinant (here 5).

A three-unknown example

A lab mixes three stock solutions to hit target amounts of three nutrients. Each column of A is one solution's nutrient content per litre, x is the litres of each, and b is the target. One call to np.linalg.solve(A, b) gives the recipe. If the answer contains a negative volume, the target is not reachable with those stocks, and the maths is telling you something real.

Common mistakes

• Mismatched sizes. An m × n matrix times a vector needs n components. NumPy will raise an error, which is helpful.
• Using * instead of @. In NumPy, A * B multiplies element by element; A @ B is matrix multiplication.
• Assuming AB = BA. Matrix multiplication is not commutative.
• Inverting to solve. Use np.linalg.solve, not np.linalg.inv(A) @ b.
• Ignoring a near-zero determinant. A nearly singular system amplifies measurement errors; the answer may be meaningless.

Quick recap

• A system of linear equations is A x = b: coefficients, unknowns, right-hand side.
• Matrix-vector and matrix-matrix products are built from dot products of rows and columns.
• A unique solution exists when det(A) ≠ 0.
• Solve with np.linalg.solve and check that A @ x is close to b.

Try it

In the Do tab, turn a three-equation system into a NumPy matrix and solve it.