Skip to content

[Hacktoberfest] Add Rayleigh quotient iteration algorithm #3401

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
wants to merge 3 commits into from
Closed
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
57 changes: 57 additions & 0 deletions linear_algebra/src/rayleigh_quotient_iteration.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,57 @@
from typing import Tuple

import numpy as np


def rayleigh_quotient_iteration(
input_matrix: np.array,
vector: np.array,
value: float,
error_tol=1e-12,
max_iterations=100,
) -> Tuple[float, np.array]:
"""
https://en.wikipedia.org/wiki/Rayleigh_quotient_iteration

Rayleigh quotient iteration
Find an approximate eigenvector and eigenvalue of the input matrix. Which
eigenvalue, eigenvector pair is found depends on the vector and value
input. Has very rapid convergence.

>>> import numpy as np
>>> input_matrix = np.array([
... [41, 4, 20],
... [ 4, 26, 30],
... [20, 30, 50]
... ])
>>> vector = np.array([1, 1, 1])
>>> value = 30
>>> rayleigh_quotient_iteration(input_matrix, vector, value)
(33.60344292195729, array([-0.86306005, 0.45012643, 0.22915836]))
>>> value = 1
>>> rayleigh_quotient_iteration(input_matrix, vector, value)
(3.735693290153355, array([-0.23946876, -0.7641015 , 0.59900218]))
"""

N = len(input_matrix)

error = error_tol + 1
prev_eigenvalue = value
iterations = 0
while error > error_tol and iterations < max_iterations:
# Construct the iterating matrix
A = np.linalg.inv(input_matrix - np.eye(N) * value)
# Construct new eigenvector
vector = np.dot(A, vector)
# Normalize new eigenvector
vector = vector / np.linalg.norm(vector)

# Construct better approximation for eigenvalue
value = np.dot(np.conj(vector).T, np.dot(input_matrix, vector))

error = abs(value - prev_eigenvalue)
prev_eigenvalue = value

iterations += 1

return value, vector