Solution of a linear system


In [1]:
import sys
import re
import numpy as np
#import parser as ps

In [2]:
%ls -l
print()
%ls -l input/


total 16
-rw-r--r-- 1 root root 7419 Oct  8 16:28 Refactor.ipynb
drwxrwxrwx 2 1000 1000 4096 Jun 11 17:20 input/
-rwxrwxrwx 1 1000 1000   49 Jun  7 19:23 input.txt*

total 40
-rwxrwxrwx 1 1000 1000 51 Jun  7 19:23 eqs0*
-rwxrwxrwx 1 1000 1000 50 Jun  7 19:23 eqs1*
-rwxrwxrwx 1 1000 1000 52 Jun  7 19:23 eqs2*
-rwxrwxrwx 1 1000 1000 56 Jun  7 19:23 eqs3*
-rwxrwxrwx 1 1000 1000 56 Jun  7 19:23 eqs4*
-rwxrwxrwx 1 1000 1000 56 Jun  7 19:23 eqs5*
-rwxrwxrwx 1 1000 1000 56 Jun  7 19:23 eqs6*
-rwxrwxrwx 1 1000 1000 56 Jun  7 19:23 eqs7*
-rwxrwxrwx 1 1000 1000 51 Jun  7 19:23 eqs8*
-rwxrwxrwx 1 1000 1000 56 Jun  7 19:23 eqs9*

In [3]:
# if len(sys.argv) > 1:
#     #print sys.argv[1]
#     f = open('input/'+sys.argv[1], 'r')
# else:
#     f = open('input.txt', 'r')

f = open('input.txt', 'r')

In [4]:
# pattern = re.compile(r'(?P<a>(-?\d+/))x(?P<b>([+-]\d+/))y(?P<c>([+-]\d+/))(\=)(?P<r>[0-9]+)')

pattern = re.compile(r'(?P<a>(-?\d+))x(?P<b>([+-]\d+))y(?P<c>([+-]\d+))(\=)(?P<r>[0-9]+)')

In [5]:
_ = """
matrixA = np.array([])
matrixB = np.array([])

for line in f:
    #line = f.readline()
    #newLine = ps.expr(line).compile()
    #match = pattern.match(line) #or findall
    match = re.search(pattern, line)#.groups()
    print match, type(match)
    print pattern, type(pattern)
    
    if match:
        a = int(match.group('a'))
        b = int(match.group('b'))
        c = int(match.group('c'))
        newRow = [a, b, c]
        matrixA = np.concatenate(matrixA, newRow)    
    
        r = int(match.group('r'))
        newRow = [r]
        matrixB = np.concatenate(matrixB, newRow)
    else:
        raise Exception('Pattern not matched')
"""

In [6]:
matrixA = np.array([[2, 3, 6], [0, 4, 2], [1, 3, 5]])
matrixB = np.array([[2], [3], [1]])

print('matrixA:\n %s\n' % matrixA)
print('matrixB:\n %s' % matrixB)


matrixA:
 [[2 3 6]
 [0 4 2]
 [1 3 5]]

matrixB:
 [[2]
 [3]
 [1]]

In [7]:
%%time
detA = np.linalg.det(matrixA)

if detA == 0:
    print('Determinant is 0. Program will stop')
    quit()

matrixAinv = np.linalg.inv(matrixA)


CPU times: user 3.86 ms, sys: 33.9 ms, total: 37.7 ms
Wall time: 44.3 ms

In [8]:
print('---', 'Inverse matrix is: ')
print(matrixAinv)
checkInv = np.allclose(np.dot(matrixA, matrixAinv), np.eye(3))

print("\nIs matrixAinv ok ?", checkInv)

#x = np.linalg.solve(a, b, c)

#Check that the solution is correct:
#np.allclose(np.dot(a, x), b)


--- Inverse matrix is: 
[[ 1.4  0.3 -1.8]
 [ 0.2  0.4 -0.4]
 [-0.4 -0.3  0.8]]

Is matrixAinv ok ? True