DEV Community

Michael Lip
Michael Lip

Posted on • Originally published at zovo.one

Solving Linear Equations Programmatically Without a Math Library

#javascript #tutorial #beginners #math

Every developer eventually needs to solve a system of linear equations. Maybe you are building a pricing model, fitting a trendline, balancing a physics simulation, or solving a constraint satisfaction problem. Reaching for a math library is the pragmatic choice, but understanding the underlying algorithm makes you a better debugger when things go wrong.

What a linear equation system looks like

A system of two linear equations in two variables:

2x + 3y = 8
4x - y = 2
Enter fullscreen mode Exit fullscreen mode

In matrix form, this is Ax = b:

| 2  3 | | x |   | 8 |
| 4 -1 | | y | = | 2 |
Enter fullscreen mode Exit fullscreen mode

The goal is to find values of x and y that satisfy both equations simultaneously.

Solving a 2x2 system with Cramer's Rule

For a 2x2 system, Cramer's Rule is the fastest approach:

function solve2x2(a1, b1, c1, a2, b2, c2) {
  // a1*x + b1*y = c1
  // a2*x + b2*y = c2
  const det = a1 * b2 - a2 * b1;

  if (Math.abs(det) < 1e-10) {
    return null; // No unique solution
  }

  const x = (c1 * b2 - c2 * b1) / det;
  const y = (a1 * c2 - a2 * c1) / det;
  return { x, y };
}

// 2x + 3y = 8, 4x - y = 2
const result = solve2x2(2, 3, 8, 4, -1, 2);
// { x: 1, y: 2 }
Enter fullscreen mode Exit fullscreen mode

The determinant check is critical. If the determinant is zero (or very close to zero due to floating-point), the system either has no solution (parallel lines) or infinite solutions (same line). In either case, there is no unique answer.

Gaussian elimination for larger systems

Cramer's Rule becomes computationally expensive for systems larger than 3x3. Gaussian elimination scales better and is the foundation of most numerical linear algebra.

The algorithm has two phases:

Forward elimination converts the matrix to upper triangular form:

function gaussianElimination(A, b) {
  const n = A.length;

  // Augmented matrix
  const aug = A.map((row, i) => [...row, b[i]]);

  // Forward elimination
  for (let col = 0; col < n; col++) {
    // Find pivot
    let maxRow = col;
    for (let row = col + 1; row < n; row++) {
      if (Math.abs(aug[row][col]) > Math.abs(aug[maxRow][col])) {
        maxRow = row;
      }
    }
    [aug[col], aug[maxRow]] = [aug[maxRow], aug[col]];

    if (Math.abs(aug[col][col]) < 1e-10) {
      return null; // Singular matrix
    }

    for (let row = col + 1; row < n; row++) {
      const factor = aug[row][col] / aug[col][col];
      for (let j = col; j <= n; j++) {
        aug[row][j] -= factor * aug[col][j];
      }
    }
  }

  // Back substitution
  const x = new Array(n);
  for (let i = n - 1; i >= 0; i--) {
    x[i] = aug[i][n];
    for (let j = i + 1; j < n; j++) {
      x[i] -= aug[i][j] * x[j];
    }
    x[i] /= aug[i][i];
  }

  return x;
}
Enter fullscreen mode Exit fullscreen mode

Partial pivoting (the row-swapping step where we find maxRow) is not optional. Without it, you will encounter division by zero or severe numerical instability. Always pivot.

When systems have no solution

Consider these equations:

2x + 4y = 10
x + 2y = 7
Enter fullscreen mode Exit fullscreen mode

The first equation is twice the second, except the constants don't match (10 is not 2 times 7). These lines are parallel and never intersect. The determinant is zero, and your solver should return an indication that no solution exists.

When systems have infinite solutions 

2x + 4y = 10
x + 2y = 5
Enter fullscreen mode Exit fullscreen mode

Now the first equation is exactly twice the second. Every point on the line x + 2y = 5 is a solution. The system is "underdetermined." In practice, this means you have more unknowns than independent equations, and you need additional constraints.

Practical applications

Interpolation. Given three points, find the quadratic that passes through them. Three points give three equations with three unknowns (the coefficients a, b, c of ax^2 + bx + c).

Balancing chemical equations. Each element gives one linear equation. The unknowns are the stoichiometric coefficients.

Network flow. At each node in a flow network, the sum of inflows equals the sum of outflows. This gives one linear equation per node.

Least squares regression. The normal equations for linear regression are a system of linear equations that can be solved directly.

The numerical stability problem

Floating-point arithmetic introduces errors that accumulate through elimination steps. For a 2x2 system this is negligible. For a 100x100 system, accumulated errors can make the result meaningless if you are not careful.

Use double precision (the default in JavaScript). Always use partial pivoting. For ill-conditioned systems (where small changes in input cause large changes in output), consider using a proper numerical library with iterative refinement.

For everyday linear equation solving, whether it's a homework check, a quick engineering calculation, or testing a formula before coding it, I keep a solver at zovo.one/free-tools/linear-equation-solver. It handles systems up to any reasonable size with step-by-step solutions. Useful for verifying that your implementation produces the right answer.

I'm Michael Lip. I build free developer tools at zovo.one. 500+ tools, all private, all free.

Top comments (0)