Solving a system of linear equations using matrices with NumPy

Overview:

A linear equation is of the form a₁x₁ + a₂x₂ + ... + a_nx_n + b = 0, where a₁, a₂,...a_n and b are coefficients and x₁, x₂...x_n are variables.
A linear equation can be solved by finding the value to the variables x₁, x₂ and x_n using one or more of the of several methods which include,
- Elimination
- Substitution
- Graphical method
- Using matrices
- Using determinants
When number of variable increase more than three it is easier to deal with the variables as a matrix than using the elimination and backward-substitution method.
With matrices the system of linear equations can be solved using the inverse of the matrix, LU decomposition, Gaussian elimination, Gauss-Jordan elimination, Cramer's rule which uses determinants to solve a system of linear equations and other techniques.
The function from the linalg module of the NumPy library internally uses the LPAC to solve the system of linear equations. Using LU decomposition the matrix A is factored as A = P * L * U.
- A – Matrix representing the system of equations
  
  P – Permuation
  
  L – Lower triangular matrix
  
  U – Upper triangular matrix

which is used in solving the system of equations A * X = B (Source: https://www.netlib.org/lapack/double/dgesv.f)

Example:

# Example Python program that solves a system of linear equations
# using the function solve() from the linalg module of the NumPy
# library
# System of linear equations:
# 2x + 2y = 16
# x - y = 4

import numpy

# Represent a system of linear equations using matrices
a = numpy.array([[2, 2], [1, -1]])
print("a:")
print(a)

b = numpy.ndarray(shape = (2, 1), dtype = numpy.int64)
b[0][0] = 16
b[1][0] = 4
print("b:")
print(b)

# Find the soultion for the linear equations
sol = numpy.linalg.solve(a, b)
print("x=%d"%sol[0][0])
print("y=%d"%sol[1][0])

Output:

a:
[[ 2 2]
[ 1 -1]]
b:
[[16]
[ 4]]
x=6
y=2