Solving Linear Equations

"A batter hit 12 boundaries worth 56 runs. How many fours and how many sixes?" In school you solved this by substitution. With NumPy you write the equations as a matrix and let np.linalg.solve do it — and the same code works for 3, 10 or 1,000 unknowns.

Definition

A system of linear equations is a set of equations where each unknown appears only multiplied by a number (no squares, no products of unknowns). It can always be written as the matrix equation A x = b:

• A — the coefficient matrix: one row per equation, one column per unknown;
• x — the vector of unknowns;
• b — the vector of right-hand-side constants.

np.linalg.solve(A, b) returns the vector x that satisfies A @ x = b, provided A is square and non-singular.

• Square — as many equations as unknowns.
• Non-singular — det(A) ≠ 0, i.e. the equations are independent, so there is exactly one solution.
• Mathematically x = A⁻¹ b, but solve never builds the inverse — it uses LU decomposition, which is faster and more accurate.

From words to matrices

f = number of fours, s = number of sixes

  f +  s = 12      ->   [[1, 1],    [f,    [12,
 4f + 6s = 56            [4, 6]]  @  s]  =   56]
                           A          x        b

Each row of A is one equation's coefficients, in the same column order as the unknowns in x.

Syntax

x = np.linalg.solve(A, b)

Example — fours and sixes

import numpy as np
A = np.array([[1, 1],
              [4, 6]])
b = np.array([12, 56])
x = np.linalg.solve(A, b)
print(x)
print(A @ x)
print(np.allclose(A @ x, b))

[8. 4.]
[12. 56.]
True

How it works: x = [8, 4] → 8 fours and 4 sixes. Check: 8 + 4 = 12 boundaries and 32 + 24 = 56 runs ✔. Always verify with A @ x and np.allclose. The answer is float (8.) even though the inputs were integers.

Example — the inverse method gives the same answer

print(np.linalg.inv(A) @ b)

[8. 4.]

x = A⁻¹ b is the textbook formula and it works, but prefer solve: computing an inverse is slower and amplifies rounding errors on large or nearly-singular systems.

Example — three unknowns

A batter scored 20 scoring shots (singles, twos, fours) for 42 runs, and hit 4 more singles than fours.

# s = singles, d = twos, f = fours
A = np.array([[1, 1, 1],     # s +  d +  f = 20  (scoring shots)
              [1, 2, 4],     # s + 2d + 4f = 42  (total runs)
              [1, 0, -1]])   # s       -  f = 4   (4 more singles than fours)
b = np.array([20, 42, 4])
x = np.linalg.solve(A, b)
print(x)
print(np.allclose(A @ x, b))
s, d, f = x.round().astype(int)
print(f"singles={s}, twos={d}, fours={f}")

[10.  4.  6.]
True
singles=10, twos=4, fours=6

A missing unknown in an equation gets coefficient 0 (no d in the third equation), and "minus" becomes a negative coefficient.

When solve fails

try:
    np.linalg.solve(np.array([[1, 2], [2, 4]]), np.array([3, 6]))
except np.linalg.LinAlgError as e:
    print("LinAlgError:", e)
try:
    np.linalg.solve(np.array([[1, 2, 3], [4, 5, 6]]), np.array([1, 2]))
except np.linalg.LinAlgError as e:
    print("LinAlgError:", e)

LinAlgError: Singular matrix
LinAlgError: Last 2 dimensions of the array must be square

Common mistakes

• ❌ Rows of A in a different unknown order than x. ✅ Fix the order ([f, s]) and use it in every row.
• ❌ Forgetting a 0 coefficient for a missing unknown — the row becomes too short. ✅ Every row has one entry per unknown.
• ❌ np.linalg.inv(A) @ b everywhere. ✅ np.linalg.solve(A, b) is faster and more accurate.
• ❌ Trusting the answer blindly. ✅ Verify with np.allclose(A @ x, b).

Interview tip

Q: Why is np.linalg.solve(A, b) preferred over np.linalg.inv(A) @ b? A: solve factorises A (LU decomposition) and solves directly — about half the work of computing the full inverse, and numerically more stable, so rounding errors are smaller.