fork download
copy 
  1. import numpy as np
  2.  
  3. def GENP(A, b):
  4.     '''
  5.     Gaussian elimination with no pivoting.
  6.     % input: A is an n x n nonsingular matrix
  7.     %        b is an n x 1 vector
  8.     % output: x is the solution of Ax=b.
  9.     % post-condition: A and b have been modified. 
  10.     '''
  11.     n =  len(A)
  12.     if b.size != n:
  13.         raise ValueError("Invalid argument: incompatible sizes between A & b.", b.size, n)
  14.     for pivot_row in xrange(n-1):
  15.         for row in xrange(pivot_row+1, n):
  16.             multiplier = A[row][pivot_row]/A[pivot_row][pivot_row]
  17.             #the only one in this column since the rest are zero
  18.             A[row][pivot_row] = multiplier
  19.             for col in xrange(pivot_row + 1, n):
  20.                 A[row][col] = A[row][col] - multiplier*A[pivot_row][col]
  21.             #Equation solution column
  22.             b[row] = b[row] - multiplier*b[pivot_row]
  23.     print A
  24.     print b
  25.     x = np.zeros(n)
  26.     k = n-1
  27.     x[k] = b[k]/A[k,k]
  28.     while k >= 0:
  29.         x[k] = (b[k] - np.dot(A[k,k+1:],x[k+1:]))/A[k,k]
  30.         k = k-1
  31.     return x
  32.  
  33. def GEPP(A, b):
  34.     '''
  35.     Gaussian elimination with partial pivoting.
  36.     % input: A is an n x n nonsingular matrix
  37.     %        b is an n x 1 vector
  38.     % output: x is the solution of Ax=b.
  39.     % post-condition: A and b have been modified. 
  40.     '''
  41.     n =  len(A)
  42.     if b.size != n:
  43.         raise ValueError("Invalid argument: incompatible sizes between A & b.", b.size, n)
  44.     # k represents the current pivot row. Since GE traverses the matrix in the upper 
  45.     # right triangle, we also use k for indicating the k-th diagonal column index.
  46.     for k in xrange(n-1):
  47.         #Choose largest pivot element below (and including) k
  48.         maxindex = abs(A[k:,k]).argmax() + k
  49.         if A[maxindex, k] == 0:
  50.             raise ValueError("Matrix is singular.")
  51.         #Swap rows
  52.         if maxindex != k:
  53.             A[[k,maxindex]] = A[[maxindex, k]]
  54.             b[[k,maxindex]] = b[[maxindex, k]]
  55.         for row in xrange(k+1, n):
  56.             multiplier = A[row][k]/A[k][k]
  57.             #the only one in this column since the rest are zero
  58.             A[row][k] = multiplier
  59.             for col in xrange(k + 1, n):
  60.                 A[row][col] = A[row][col] - multiplier*A[k][col]
  61.             #Equation solution column
  62.             b[row] = b[row] - multiplier*b[k]
  63.     print A
  64.     print b
  65.     x = np.zeros(n)
  66.     k = n-1
  67.     x[k] = b[k]/A[k,k]
  68.     while k >= 0:
  69.         x[k] = (b[k] - np.dot(A[k,k+1:],x[k+1:]))/A[k,k]
  70.         k = k-1
  71.     return x
  72.  
  73. if __name__ == "__main__":
  74.     A = np.array([[2., 4., -2., -2.], [1., 2., 4., -3.], [-3., -3., 8., -2.], [-1., 1., 6., -3.]])
  75.     b = np.array([-4., 5., 7., 7.])
  76.     print GENP(np.copy(A), np.copy(b))
  77.     print GEPP(A,b)# your code goes here
Success #stdin #stdout #stderr 0.08s 23792KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
[[ 2.   4.  -2.  -2. ]
 [ 0.5  0.   5.  -2. ]
 [-1.5  inf -inf  inf]
 [-0.5  inf  nan  nan]]
[ -4.   7. -inf  nan]
[nan nan nan nan]
[[-3.00000000e+00 -3.00000000e+00  8.00000000e+00 -2.00000000e+00]
 [-6.66666667e-01  2.00000000e+00  3.33333333e+00 -3.33333333e+00]
 [-3.33333333e-01  5.00000000e-01  5.00000000e+00 -2.00000000e+00]
 [ 3.33333333e-01  1.00000000e+00  8.88178420e-17  1.00000000e+00]]
[7.         0.66666667 7.         4.        ]
[1. 2. 3. 4.]
stderr
copy 
prog:16: RuntimeWarning: divide by zero encountered in double_scalars
prog:16: RuntimeWarning: invalid value encountered in double_scalars
https://ideone.com/bC7YwB
language:
Python (cpython 2.7.16)
created:
1 year ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!