# Chapter 10 : LU Decomposition and matrix inverse¶

## Ex10.1 Page 277¶

In :
from numpy import mat
from numpy.linalg import det
A = mat([[3,-0.1,-0.2],[0.1,7,-0.3],[0.3,-0.2,10]])
U = A#
print "A =",A
m = U[0,0]
n = U[1,0]
a = n/m#
U[1:2] = U[1:2] - U[0:1] / (m/n)#
n = U[2,0]
b = n/m

U[2:3] = U[2:3] - U[0:1] / (m/n)#
m = U[1,1]
n = U[2,1]
c = n/m#
U[2:3] = U[2:3] - U[1:2] / (m/n)#
print "\nU = \n",U
L = mat([[1,0,0],[a,1,0],[b,c,1]])
print "\nL calculated based on gauss elimination method = \n",L

A = [[  3.   -0.1  -0.2]
[  0.1   7.   -0.3]
[  0.3  -0.2  10. ]]

U =
[[  3.          -0.1         -0.2       ]
[  0.           7.00333333  -0.29333333]
[  0.           0.          10.01204188]]

L calculated based on gauss elimination method =
[[ 1.          0.          0.        ]
[ 0.03333333  1.          0.        ]
[ 0.1        -0.02712994  1.        ]]


## Ex10.2 Page 279¶

In :
from numpy import mat
A = mat([[3,-0.1,-0.2],[0.1,7,-0.3],[0.3,-0.2,10]])
B = mat([[7.85],[-19.3],[71.4]])
X = (A**-1) * B
print "X = \n",X

X =
[[ 3. ]
[-2.5]
[ 7. ]]


## Ex10.3 Page 284¶

In :
from numpy import mat,array
A = mat([[3,-0.1,-0.2],[0.1,7,-0.3],[0.3,-0.2,10]])
B = (A**-1)
L = mat([[1,0,0],[0.033333,1,0],[0.1,-0.02713,1]])
U = mat([[3,-0.1,-0.2],[0,7.0033,-0.293333],[0,0,10.012]])
for i  in range(1,4):
if i==1:
m = mat([,,])
else:
if i==2:
m = mat([,,])
else:
m = mat([,,])

d = (L**-1) * m#
x = (U**-1) * d#
B[:,i-1] = x

print "B=\n",B

B=
[[ 0.33248872  0.00494409  0.00679813]
[-0.00518174  0.14290333  0.00418348]
[-0.01007834  0.00270975  0.09988014]]


## Ex10.4 Page 291¶

In :
from numpy import mat
from __future__ import division
A = mat([[1,1/2,1/3],[1/2,1/3,1/4],[1/3,1/4,1/5]])
n = A[1,0]
A[1,:] = A[1,:]/n
n = A[2,0]
A[2,:] = A[2,:]/n
B = (A**-1)#
print "A = \n",A

m=range(1,4)
su=range(1,4)
for j in range(1,4):
a = 0#
for i in range(1,4):
m[i-1]= A[j-1,i-1]
su[j-1] = a + m[i-1]#
a = su[j-1]#

if su< su:
if su< su:
z = su
else:
z = su

else:
if su < su:
z = su#
else:
z = su#

m=range(1,4)
sm=range(1,4)
for j in  range(1,4):
a = 0#
for i in range(1,4):
m[i-1]= B[j-1,i-1]
sm[j-1]= a + abs(m[i-1])
a = sm[j-1]#

if sm< sm:
if sm< sm :
y = sm
else:
y = sm

else:
if sm< sm:
y = sm
else:
y = sm#

C = z*y#
print "\nCondition number for the matrix =\n",C

A =
[[ 1.          0.5         0.33333333]
[ 1.          0.66666667  0.5       ]
[ 1.          0.75        0.6       ]]

Condition number for the matrix =
451.2