[Scientific / technical calculation by Python] Solving (generalized) eigenvalue problem using numpy / scipy, using library

Introduction

It is very common to face the eigenvalue problem of linear algebra in many science and engineering fields. Here is a summary of how to solve the ** (generalized) eigenvalue problem ** using Python's numpy and scipy. The methods are summarized in the addendum.

The generalized eigenvalue problem is a question of finding a vector $ \ mathbf {x} $ and a (complex) numerical value $ \ lambda $ that satisfy the following equations, using $ A $ and $ B $ as matrices. ** **

$ A \ \mathbf{x} = \lambda \ B \ \mathbf{x} \ {\tag 1}$

If $ B $ is the identity matrix $ E $, the problem is $ A \ \mathbf{x} = \lambda \ \mathbf{x} \ {\tag 2}$ Will be. This problem is called the ** standard eigenvalue problem **.

** There are power method and QR method as numerical calculation methods for eigenvalue problems [1], but this paper does not touch on them and only describes the basic usage of the library [2]. ** **

Postscript: 2017 9/12

● I posted an article about numerical solution by power method. [Scientific / technical calculation by Python] Numerical solution of eigenvalue problem of matrix by power method, linear algebra


Contents

First, consider the above equation (2), which is a standard eigenvalue problem. Use the ** numpy.linalg.eig ** method.

** (1) For run columns ** ** (2) For complex matrix **

Next, consider the above equation (1) type. Use ** scipy.linalg.eig **.

** (3) Solving the general eigenvalue problem **

If ** $ B $ is used as the identity matrix, the problems (1) and (2) can be solved. ** **

☆ ** The eigenvectors obtained to normalize the eigenvectors in the execution column calculation codes (1) and (2) are divided by their norms. ** ** In other words

for i in range(len(eig_vec)):
    eig_vec [i] = eig_vec[i]/np.linalg.norm(eig_vec[i])

And said.

In the case of a complex matrix (code (2)), the normalized solution is output from the beginning.


Code (1): Eigenvalue / eigenvector of execution column

import numpy as np
"""
Calculation of eigenvalues and eigenvectors of execution columns
"""
#A_matrix=np.matrix([[2,4],[0,4]]) #Generate matrix A
A_matrix=np.array([[2,4],[0,4]]) #Generate matrix A


eig_val, eig_vec =np.linalg.eig(A_matrix) #Eig each eigenvalue and eigenvector_val, eig_Store in vec

for i in range(len(eig_vec)):
    eig_vec [i] = eig_vec[i]/np.linalg.norm(eig_vec[i])
    
#Output of calculation result
print (eig_val) 
print(eig_vec)

Result (1)

[2. 4.] # 2 eigenvalues (2 and 4)

[[ 0.74535599 0.66666667] [ 0. 1. ]] The corresponding eigenvector. It is standardized to 1.


Code (2): Eigenvalue / eigenvector of complex matrix

It can be calculated in exactly the same way as an execution column. Note that the complex number $ a + bi $ is expressed as a + bj in python by changing i to j.


import numpy as np
"""
Calculation of eigenvalues and eigenvectors of complex matrix
"""
#A_matrix=np.matrix([[1+4j,2j],[2,4j]]) #Complex numbers are in matrix elements
A_matrix=np.array([[1+4j,2j],[2,4j]]) #Complex numbers are in matrix elements
eig_val, eig_vec =np.linalg.eig(A_matrix) #Calculation of eigenvalues and eigenvectors
print (eig_val)
print(eig_vec)

Result (2)

eigenvalue [ 1.95907589+5.37073062j -0.95907589+2.62926938j]

** Eigenvector (complex vector with norm normalized to 1) ** [[ 0.76703668+0.j -0.36782319-0.52570033j], [ 0.52570033-0.36782319j 0.76703668+0.j ]]


Code (3): Generalized eigenvalue problem

import numpy as np
import scipy.linalg 


A = np.array([[-32,16,0],[16,-32,16],[0,16,-32]])
B = np.array([[-2, 16, 13], [16, -5, 16], [0, 1, -32]]) #Matrix B
#B=np.identity(3)


eig_val,eig_vec =  scipy.linalg.eig(A,B)  #Solve the general solidification value problem

for i in range(len(eig_vec)): #Normalization
    eig_vec [i] = eig_vec[i]/np.linalg.norm(eig_vec[i])

print (eig_val)
print(eig_vec)

Eigenvalues: ** Eigenvalues are generally complex numbers, so they are displayed as complex numbers using j as a pure imaginary number ** [-0.81267990+0.j 2.48526372+0.j 1.00000000+0.j]

array([-0.81267990+0.j,  2.48526372+0.j,  1.00000000+0.j])

Eigenvector:

[[ -7.43016156e-01  -5.06563678e-01   4.37401683e-01]
 [ -6.38461105e-01   7.69654089e-01   2.24401124e-16]
 [ -2.11406258e-01  -2.50015624e-01  -9.44880724e-01]]

Addendum: numpy and scipy linalg methods

1.eig: Calculate eigenvalues and eigenvectors of general matrices
2.eigh: Hermitian / real symmetric matrix eigenvalue / eigenvector calculation
3.eigvals: Calculation of eigenvalues of general matrix(Do not calculate eigenvectors)
4.eigvalsh: Calculation of eigenvalues of Hermitian / real symmetric matrix(Do not calculate eigenvectors)

** In particular, the calculation of eigh for symmetric matrices is better than the eig used for general matrices. Very fast. ** **


References

There seems to be various texts about numerical linear algebra. ..

[1] F. Chatlan, ["Eigenvalues of the matrix"](https://www.amazon.co.jp/%E8%A1%8C%E5%88%97%E3%81%AE%E5%9B% BA% E6% 9C% 89% E5% 80% A4-% E6% 9C% 80% E6% 96% B0% E3% 81% AE% E8% A7% A3% E6% B3% 95% E3% 81% A8 % E5% BF% 9C% E7% 94% A8-F-% E3% 82% B7% E3% 83% A3% E3% 83% 88% E3% 83% A9% E3% 83% B3 / dp / 443171037X / ref = sr_1_1? s = books & ie = UTF8 & qid = 1504668098 & sr = 1-1), Springer Fairlark Tokyo, 2003., Masatake Mori et al. ["Iwanami Lecture Applied Mathematics <12 [Method 2] Linear Calculation"](https://www.amazon.co.jp/%E5%B2%A9%E6%B3%A2%E8%AC% 9B% E5% BA% A7-% E5% BF% 9C% E7% 94% A8% E6% 95% B0% E5% AD% A6% E3% 80% 8812% E3% 80% 89% E3% 80% 94 % E6% 96% B9% E6% B3% 951% E3% 80% 95-% E6% 95% B0% E5% 80% A4% E8% A8% 88% E7% AE% 97% E3% 81% AE% E5% 9F% BA% E7% A4% 8E% EF% BC% 8F% E3% 80% 94% E6% 96% B9% E6% B3% 952% E3% 80% 95-% E7% B7% 9A% E5 % BD% A2% E8% A8% 88% E7% AE% 97-% E6% AD% A3% E6% AD% A6 / dp / 4000108026), Iwanami Shoten, 1997.

[2] scipy website: https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.eig.html

https://docs.scipy.org/doc/scipy-0.18.1/reference/tutorial/linalg.html Note that ** scipy.linalg is always linked to BLAS / LAPACK (numpy is an option), so it is stated that it is basically faster than numpy.linalg. ** ** The part is reprinted below.

"scipy.linalg vs numpy.linal scipy.linalg contains all the functions in numpy.linalg. plus some other more advanced ones not contained in numpy.linalg

Another advantage of using scipy.linalg over numpy.linalg is that it is always compiled with BLAS/LAPACK support, while for numpy this is optional. Therefore, the scipy version might be faster depending on how numpy was installed.

Therefore, unless you don’t want to add scipy as a dependency to your numpy program, use scipy.linalg instead of numpy.linalg"

Recommended Posts

[Scientific / technical calculation by Python] Solving (generalized) eigenvalue problem using numpy / scipy, using library
[Scientific / technical calculation by Python] Solving simultaneous linear equations, numerical calculation, numpy
[Scientific / technical calculation by Python] Inverse matrix calculation, numpy
[Scientific / technical calculation by Python] 3rd order spline interpolation, scipy
[Scientific / technical calculation by Python] Basic operation of arrays, numpy
[Scientific / technical calculation by Python] Monte Carlo integration, numerical calculation, numpy
[Scientific / technical calculation by Python] Numerical integration, trapezoidal / Simpson law, numerical calculation, scipy
[Scientific / technical calculation by Python] Sum calculation, numerical calculation
[Scientific / technical calculation by Python] 1D-3D discrete fast Fourier transform, scipy
[Scientific / technical calculation by Python] 2D random walk (drunk walk problem), numerical calculation
[Scientific / technical calculation by Python] List of usage of (special) functions used in physics by using scipy
[Scientific / technical calculation by Python] Numerical solution of eigenvalue problem of matrix by power method, numerical linear algebra
[Scientific / technical calculation by Python] Fitting by nonlinear function, equation of state, scipy
[Scientific / technical calculation by Python] Calculation of matrix product by @ operator, python3.5 or later, numpy
[Scientific / technical calculation by Python] Solving ordinary differential equations, mathematical formulas, sympy
[Scientific / technical calculation by Python] Histogram, visualization, matplotlib
[Scientific / technical calculation by Python] Lagrange interpolation, numerical calculation
[Scientific / technical calculation by Python] Solving second-order ordinary differential equations by Numerov method, numerical calculation
[Scientific / technical calculation by Python] Solving one-dimensional Newton equation by the 4th-order Runge-Kutta method
[Scientific / technical calculation by Python] Logarithmic graph, visualization, matplotlib
[Scientific / technical calculation by Python] Polar graph, visualization, matplotlib
[Scientific / technical calculation by Python] Solving the boundary value problem of ordinary differential equations in matrix format, numerical calculation
[Scientific / technical calculation by Python] Numerical solution of second-order ordinary differential equations, initial value problem, numerical calculation
[Scientific / technical calculation by Python] Numerical solution of one-dimensional harmonic oscillator problem by velocity Verlet method
[Scientific / technical calculation by Python] Vector field visualization example, electrostatic field, matplotlib
Scientific / technical calculation by Python] Drawing and visualization of 3D isosurface and its cross section using mayavi
[Scientific / technical calculation by Python] Solving one-dimensional Schrodinger equation in stationary state by shooting method (1), well-type potential, quantum mechanics
[Scientific / technical calculation by Python] Drawing animation of parabolic motion with locus, matplotlib
[Scientific / technical calculation by Python] Numerical calculation to find the value of derivative (differential)
[Scientific / technical calculation by Python] Analytical solution to find the solution of equation sympy
[Scientific / technical calculation by Python] Plot, visualize, matplotlib 2D data with error bars
[Scientific / technical calculation by Python] Drawing of 3D curved surface, surface, wireframe, visualization, matplotlib
Speeding up numerical calculation using NumPy / SciPy: Application 2
Speeding up numerical calculation using NumPy / SciPy: Application 1
[Scientific / technical calculation by Python] Solving the Schrodinger equation in the steady state in the 3D isotropic harmonic oscillator potential by the matrix method, boundary value problem, quantum mechanics
[Scientific / technical calculation by Python] Solving one-dimensional Schrodinger equation in stationary state by shooting method (2), harmonic oscillator potential, quantum mechanics
[Scientific / technical calculation by Python] Plot, visualization, matplotlib of 2D data read from file
[Scientific / technical calculation by Python] Drawing, visualization, matplotlib of 2D (color) contour lines, etc.
Why what? Deep Learning Scientific Calculation Library Numpy Edition
[Science / technical calculation by Python] Numerical solution of first-order ordinary differential equations, initial value problem, numerical calculation
[Scientific / technical calculation by Python] List of matrices that appear in Hinpan in numerical linear algebra
[Scientific and technical calculation by Python] Drawing of fractal figures [Sierpinski triangle, Bernsley fern, fractal tree]
[Scientific / technical calculation by Python] Wave "beat" and group velocity, wave superposition, visualization, high school physics
python numpy array calculation
Age calculation using python
Calculation problem using mod
Speeding up numerical calculation using NumPy / SciPy: Picking up fallen ears
[Numerical calculation method, python] Solving ordinary differential equations by Eular method
Install python3 and scientific calculation library on Ubuntu (virtualenv + pip)
[Scientific / technical calculation by Python] Monte Carlo simulation of thermodynamics of 2D Ising spin system by Metropolis method
[Scientific / technical calculation by Python] Numerical solution of two-dimensional Laplace-Poisson equation by Jacobi method for electrostatic potential, elliptic partial differential equation, boundary value problem