# Solving linear equations Ax=b by Gaussian elimination
# without checking for zero-coefficients
import numpy as np
n=3
A = np.array( [[ -3., 2., -6. ],
[ 5., 7., -5. ],
[ 1., 4., -2. ] ] )
b = np.array( [ 6.,
6.,
8. ] )
x = np.array( [0]*3 ) # answer: -2, 3, 1
# Construction of upper triangular system
for k in range(n-1): # substep A-k
for i in range(k+1, n): # i-th equation
for j in range(k+1, n): # j-th term
A[i][j] = A[i][j] - (A[i][k]/A[k][k])*A[k][j] # Aij -> Aij' -> Aij'' -> ...
b[i] = b[i] - (A[i][k]/A[k][k])*b[k] # Bi -> Bi' -> Bi'' -> ...
x[n-1]= b[n-1]/A[n-1][n-1] # solution of last equation of upper triangle
# back-substitution
for k in range(n-1): # substep B-k
i=n-k-2 # i-th equation
for j in range(i+1, n): # j-th term
b[i] -= A[i][j] * x[j] # subtract j-th term in i-th equation
x[i] = b[i] / A[i][i]
# output results
print(x)