When Double Precision Is Not Ample

80
When Double Precision Is Not Ample

Click the next button to launch an ipython notebook on Google Colab which implements the code developed in this publish:

Open In Colab

We desire to resolve the next now not new eigenvalue bid

[ mathbf{A}Psi=lambdaPsi, ]

where (mathbf{A}) is the lesp matrix, (Psi) is the eigenvector and (lambda) the eigenvalue. The lesp matrix is a tridiagonal matrix with accurate, sensitive eigenvalues, with a (5 times 5) example considered in equation (1).

[ mathbf{M}=begin{pmatrix} -5 & 2 & 0 & 0 & 0 \
frac{1}{2} & -7 & 3 & 0 & 0 \
0 & frac{1}{3} & -9 & 4 & 0 \
0 & 0 & frac{1}{4} & -11 & 5 \
0 & 0 & 0 & frac{1}{5} & -13 \
end{pmatrix} tag{1} ]

We’re going to numerically mark that the eigenvalues of a matrix (mathbf{A}) are equal to the eigenvalues of its matrix transpose, (mathbf{A}^T), without problems proved as follows

[text{det}(mathbf{A}^T – lambda mathbf{I})=text{det}((mathbf{A} – lambda mathbf{I})^T)=text{det}(mathbf{A} – lambda mathbf{I}),]

where (mathbf{I}) is the identification matrix. Whilst you homicide no longer desire to write your comprise characteristic to invent the lesp matrix then you definately are in success, there is a python library, rogues, which has one built in.

The following python code calculates the eigenvalues of (mathbf{A}) and (mathbf{A}^T) and displays them utilizing seaborn.

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
from rogues import lesp
from matplotlib import pyplot
import seaborn as sns 
from scipy.linalg import eigvals

sns.space()
palette = sns.color_palette("sparkling")

# Dimension of matrix
shadowy = 100
# salvage lesp matrix
A = lesp(shadowy)
# Transpose matrix A
AT = A.T
# Calculate eigenvalues of A
Aev = eigvals(A)
# Calculate eigenvalues of A^T
ATev = eigvals(A.T)
# Extract accurate and imaginary aspects of A
A_X = [x.real for x in Aev]
A_Y = [x.imag for x in Aev]
# Extract accurate and imaginary aspects of A^T
AT_X = [x.real for x in ATev]
AT_Y = [x.imag for x in ATev]

# Place
ax = sns.scatterplot(x=A_X, y=A_Y, color = 'gray', marker='o', heed=r'$mathbf{A}$')
ax = sns.scatterplot(x=AT_X, y=AT_Y, color = 'blue', marker='x', heed=r'$mathbf{A}^T$')
# Give axis labels
ax.space(xlabel=r'accurate', ylabel=r'imag')
# Arrangement legend
ax.legend()

pyplot.mark()

This produces the next reveal of the eigenvalues for a (100 times 100 ) matrix.

Desktop View

One thing is inappropriate, many of the eigenvalues of matrix (mathbf{A}) are no longer equal to the eigenvalues of (mathbf{A}^T), even supposing they need to be an analogous. About a of them additionally comprise advanced substances! There is now not any longer any error with the program; this discrepancy is prompted by a lack of numerical accuracy in the eigenvalue calculation due to the limitation of hardware double precision (16-digit).

The scipy eigvals characteristic calls LAPACK routine _DGEEV which first reduces the enter matrix to upper Hessenberg salvage by the employ of orthogonal similarity transformations. The QR algorithm is then at threat of additional lower the matrix to upper quasi-triangular Schur salvage, (mathbf{T}), with 1 by 1 and 2 by 2 blocks on the principle diagonal. The eigenvalues are computed from (mathbf{T}).

For most functions this could occasionally well construct very moral eigenvalues, but when a bid is ill-conditioned, or a runt alternate to the enter matrix ends in a predominant alternate to the eigenvalues; more numerically stable recommendations are required. Failing this, increased working precision is required such as quadruple (32-digit), octuple (64-digit) and even arbitrary precision. Whenever you catch yourself shedding precision in a calculation it’s always suggested to first analyse the recommendations and algorithms you are utilizing to gape if they’ll also be improved. Fully after doing here’s it dazzling to originate higher the working precision. Prolonged precision calculations retract drastically longer than double precision calculations as they are run in instrument in preference to on hardware.

Rising precision is straightforward in python with the employ of the mpmath library. Indicate that this library is incredibly monotonous for immense matrices, so is grand shunned for most functions. The following code calculates the an analogous eigenvalues as sooner than, this time at quadruple, 32-digit precision.

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
from rogues import lesp
from matplotlib import pyplot
import seaborn as sns 
from scipy.linalg import eigvals
from mpmath import *
# House precision to 32-digit
mp.dps = 32

sns.space()
palette = sns.color_palette("sparkling")

# Dimension of matrix
shadowy = 100
# Lesp matrix
A = lesp(shadowy)
# Transpose matrix A
AT = A.T
# Calculate eigenvalues of A
Aev, Eeg = mp.eig(mp.matrix(A))
# Calculate eigenvalues of A^T
ATev, ETeg = mp.eig(mp.matrix(AT))
# Extract accurate and imaginary aspects of A
A_X = [x.real for x in Aev]
A_Y = [x.imag for x in Aev]
# Extract accurate and imaginary aspects of A^T
AT_X = [x.real for x in ATev]
AT_Y = [x.imag for x in ATev]

# Place
ax = sns.scatterplot(x=A_X, y=A_Y, color = 'gray', marker='o', heed=r'$mathbf{A}$')
ax = sns.scatterplot(x=AT_X, y=AT_Y, color = 'blue', marker='x', heed=r'$mathbf{A}^T$')
# Give axis labels
ax.space(xlabel=r'accurate', ylabel=r'imag')
# Arrangement legend
ax.legend()

pyplot.mark()

This produces the next reveal of the eigenvalues for a (100 times 100 ) matrix.

Desktop View

Success! Now we comprise got proven that the eigenvalues of (mathbf{A}) are equal to the eigenvalues of (mathbf{A}^T), equation (1) and all eigenvalues are accurate. This explicit example highlights a ordinary theme within computational stories of quantum programs. Quantum simulations of atoms and molecules generally involve ill-conditioned now not new and generalized eigenvalue concerns. Numerically stable algorithms exist for such concerns, but these ceaselessly fail leaving us to brute force the calculation with increased precision to lower floating point rounding errors.

The following publish will focus on better ways to put in force increased precision in numerical calculations.

Read More

Vanic
WRITTEN BY

Vanic

“Simplicity, patience, compassion.
These three are your greatest treasures.
Simple in actions and thoughts, you return to the source of being.
Patient with both friends and enemies,
you accord with the way things are.
Compassionate toward yourself,
you reconcile all beings in the world.”
― Lao Tzu, Tao Te Ching

you're currently offline