-
-
Notifications
You must be signed in to change notification settings - Fork 45.9k
/
rayleigh_quotient.py
69 lines (57 loc) · 1.52 KB
/
rayleigh_quotient.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
"""
https://en.wikipedia.org/wiki/Rayleigh_quotient
"""
from typing import Any
import numpy as np
def is_hermitian(matrix: np.ndarray) -> bool:
"""
Checks if a matrix is Hermitian.
>>> import numpy as np
>>> A = np.array([
... [2, 2+1j, 4],
... [2-1j, 3, 1j],
... [4, -1j, 1]])
>>> is_hermitian(A)
True
>>> A = np.array([
... [2, 2+1j, 4+1j],
... [2-1j, 3, 1j],
... [4, -1j, 1]])
>>> is_hermitian(A)
False
"""
return np.array_equal(matrix, matrix.conjugate().T)
def rayleigh_quotient(a: np.ndarray, v: np.ndarray) -> Any:
"""
Returns the Rayleigh quotient of a Hermitian matrix A and
vector v.
>>> import numpy as np
>>> A = np.array([
... [1, 2, 4],
... [2, 3, -1],
... [4, -1, 1]
... ])
>>> v = np.array([
... [1],
... [2],
... [3]
... ])
>>> rayleigh_quotient(A, v)
array([[3.]])
"""
v_star = v.conjugate().T
v_star_dot = v_star.dot(a)
assert isinstance(v_star_dot, np.ndarray)
return (v_star_dot.dot(v)) / (v_star.dot(v))
def tests() -> None:
a = np.array([[2, 2 + 1j, 4], [2 - 1j, 3, 1j], [4, -1j, 1]])
v = np.array([[1], [2], [3]])
assert is_hermitian(a), f"{a} is not hermitian."
print(rayleigh_quotient(a, v))
a = np.array([[1, 2, 4], [2, 3, -1], [4, -1, 1]])
assert is_hermitian(a), f"{a} is not hermitian."
assert rayleigh_quotient(a, v) == float(3)
if __name__ == "__main__":
import doctest
doctest.testmod()
tests()