math4610

Tasksheet 10

Author: David Merkley

Language: Python

Links to Code:

Task 1

Description/Purpose: Power Iteration code. The largest eigenvalue was found.

Implementation/Code:

from math4610 import matrixVectorMultiply, vectorScalarMultiply, vectorDotProduct

def powerIteration(A, x0, tol=1e-10, maxiter=1e5):
    x = x0
    l0 = float('inf')
    for i in range(int(maxiter)):
        y = matrixVectorMultiply(A, x)
        v = vectorScalarMultiply(y, 1/l2(y))
        l1 = vectorDotProduct(v, matrixVectorMultiply(A, v))
        error = abs(l1 - l0)
        if (error < tol):
            return l1, v
        else:
            l0 = l1
            x = v
    print("Max Iter")
    
A = generateSquareRandomSymmetricDiagonallyDominant(100)
x = rv(100)
l1, v1 = powerIteration(A1, x1)
print("Largest Eigenvalue:", l1)

Output:

Largest Eigenvalue: 149.7580931710224

Task 2

Description/Purpose: This is the code used for inverse Power iteration. The smallest eigenvalue was found.

Implementation/Code:

from math4610 import l2, vectorScalarMultiply, vectorDotProduct, JacobiIteration

def inversePowerIteration(A, x0, tol=1e-10, maxiter=1e5):
    y = x
    x = x0
    Aiv = x
    l0 = float('inf')
    for i in range(int(maxiter)):
        y = JacobiIteration(A, x, y, tol/100, maxiter)
        v = vectorScalarMultiply(y, 1/l2(y))
        Aiv = JacobiIteration(A, v, Aiv, tol/100, maxiter)
        l1 = vectorDotProduct(v, Aiv)
        error = abs(l1 - l0)
        if error < tol:
            print(i, "iterations")
            return 1/l1, v
        else:
            l0 = l1
            x = v
    print("Max Iter")

  A = generateSquareRandomSymmetricDiagonallyDominant(100)
  x = rv(100)
  l1a, v1a = inversePowerIteration(A, x)
  print("Smallest Eigenvalue:", l1a)

Output:

Smallest Eigenvalue: 93.75590968550863

Task 3

Description/Purpose: Computed the l1 matrix norm for the square matrix.

Implementation/Code:

def l1matrix(A):
    m = len(A)
    n = len(A[0])
    sums = [0]*n
    for i in range(m):
        for j in range(n):
            sums[j] += abs(A[i][j])
    return max(sums)
    
A = generateSquareRandomSymmetricDiagonallyDominant(100)
print(l1matrix(A))

Output:

156.6866765892807

Task 4

Description/Purpose: Did the same with the infty matrix norm.

Implementation/Code:

def linftymatrix(A):
    n = len(A[0])
    m = len(A)
    sums = [0]*n
    for i in range(m):
        for j in range(n):
            sums[i] += abs(A[i][j])
    return max(sums)

A = generateSquareRandomSymmetricDiagonallyDominant(100)
print(linftymatrix(A))

Output

156.77634420476213

Task 5

Description/Purpose: Also did the same for the l2 condition norm.

Implementation/Code:

def l2ConditionNumber(A):
    b0 = rv(len(A))
    lmax = powerIteration(A, b0)[0]
    lmin = inversePowerIteration(A, b0)[0]
    return abs(lmax / lmin)

A = generateSquareRandomSymmetricDiagonallyDominant(100)
print(l2ConditionNumber(A))

Output:

1.583954211085809

Task 6

Here is the link to my Software Manual