[{"@context":"http:\/\/schema.org\/","@type":"BlogPosting","@id":"https:\/\/wiki.edu.vn\/en\/wiki2\/jacobi-eigenvalue-algorithm-wikipedia\/#BlogPosting","mainEntityOfPage":"https:\/\/wiki.edu.vn\/en\/wiki2\/jacobi-eigenvalue-algorithm-wikipedia\/","headline":"Jacobi eigenvalue algorithm – Wikipedia","name":"Jacobi eigenvalue algorithm – Wikipedia","description":"before-content-x4 In numerical linear algebra, the Jacobi eigenvalue algorithm is an iterative method for the calculation of the eigenvalues and","datePublished":"2016-03-27","dateModified":"2016-03-27","author":{"@type":"Person","@id":"https:\/\/wiki.edu.vn\/en\/wiki2\/author\/lordneo\/#Person","name":"lordneo","url":"https:\/\/wiki.edu.vn\/en\/wiki2\/author\/lordneo\/","image":{"@type":"ImageObject","@id":"https:\/\/secure.gravatar.com\/avatar\/c9645c498c9701c88b89b8537773dd7c?s=96&d=mm&r=g","url":"https:\/\/secure.gravatar.com\/avatar\/c9645c498c9701c88b89b8537773dd7c?s=96&d=mm&r=g","height":96,"width":96}},"publisher":{"@type":"Organization","name":"Enzyklop\u00e4die","logo":{"@type":"ImageObject","@id":"https:\/\/wiki.edu.vn\/wiki4\/wp-content\/uploads\/2023\/08\/download.jpg","url":"https:\/\/wiki.edu.vn\/wiki4\/wp-content\/uploads\/2023\/08\/download.jpg","width":600,"height":60}},"image":{"@type":"ImageObject","@id":"https:\/\/wikimedia.org\/api\/rest_v1\/media\/math\/render\/svg\/4611d85173cd3b508e67077d4a1252c9c05abca2","url":"https:\/\/wikimedia.org\/api\/rest_v1\/media\/math\/render\/svg\/4611d85173cd3b508e67077d4a1252c9c05abca2","height":"","width":""},"url":"https:\/\/wiki.edu.vn\/en\/wiki2\/jacobi-eigenvalue-algorithm-wikipedia\/","wordCount":18151,"articleBody":"     (adsbygoogle = window.adsbygoogle || []).push({});before-content-x4In numerical linear algebra, the Jacobi eigenvalue algorithm is an iterative method for the calculation of the eigenvalues and eigenvectors of a real symmetric matrix (a process known as diagonalization). It is named after Carl Gustav Jacob Jacobi, who first proposed the method in 1846,[1] but only became widely used in the 1950s with the advent of computers.[2]       (adsbygoogle = window.adsbygoogle || []).push({});after-content-x4Table of ContentsDescription[edit]Convergence[edit]Algorithm[edit]Notes[edit]Example[edit]Applications for real symmetric matrices[edit]Generalizations[edit]References[edit]Further reading[edit]External links[edit]Description[edit]Let        (adsbygoogle = window.adsbygoogle || []).push({});after-content-x4S{displaystyle S} be a symmetric matrix, and G=G(i,j,\u03b8){displaystyle G=G(i,j,theta )} be a Givens rotation matrix. Then:S\u2032=GSG\u22a4{displaystyle S’=GSG^{top },}is symmetric and similar to        (adsbygoogle = window.adsbygoogle || []).push({});after-content-x4S{displaystyle S}.Furthermore, S\u2032{displaystyle S^{prime }} has entries:Sii\u2032=c2Sii\u22122scSij+s2SjjSjj\u2032=s2Sii+2scSij+c2SjjSij\u2032=Sji\u2032=(c2\u2212s2)Sij+sc(Sii\u2212Sjj)Sik\u2032=Ski\u2032=cSik\u2212sSjkk\u2260i,jSjk\u2032=Skj\u2032=sSik+cSjkk\u2260i,jSkl\u2032=Sklk,l\u2260i,j{displaystyle {begin{aligned}S’_{ii}&=c^{2},S_{ii}-2,sc,S_{ij}+s^{2},S_{jj}\\S’_{jj}&=s^{2},S_{ii}+2sc,S_{ij}+c^{2},S_{jj}\\S’_{ij}&=S’_{ji}=(c^{2}-s^{2}),S_{ij}+sc,(S_{ii}-S_{jj})\\S’_{ik}&=S’_{ki}=c,S_{ik}-s,S_{jk}&kneq i,j\\S’_{jk}&=S’_{kj}=s,S_{ik}+c,S_{jk}&kneq i,j\\S’_{kl}&=S_{kl}&k,lneq i,jend{aligned}}}where s=sin\u2061(\u03b8){displaystyle s=sin(theta )} and c=cos\u2061(\u03b8){displaystyle c=cos(theta )}.Since G{displaystyle G} is orthogonal, S{displaystyle S} and S\u2032{displaystyle S^{prime }} have the same Frobenius norm ||\u22c5||F{displaystyle ||cdot ||_{F}} (the square-root sum of squares of all components), however we can choose \u03b8{displaystyle theta } such that Sij\u2032=0{displaystyle S_{ij}^{prime }=0}, in which case S\u2032{displaystyle S^{prime }} has a larger sum of squares on the diagonal:Sij\u2032=cos\u2061(2\u03b8)Sij+12sin\u2061(2\u03b8)(Sii\u2212Sjj){displaystyle S’_{ij}=cos(2theta )S_{ij}+{tfrac {1}{2}}sin(2theta )(S_{ii}-S_{jj})}Set this equal to 0, and rearrange:tan\u2061(2\u03b8)=2SijSjj\u2212Sii{displaystyle tan(2theta )={frac {2S_{ij}}{S_{jj}-S_{ii}}}}if Sjj=Sii{displaystyle S_{jj}=S_{ii}}\u03b8=\u03c04{displaystyle theta ={frac {pi }{4}}}In order to optimize this effect, Sij should be the off-diagonal element with the largest absolute value, called the pivot.The Jacobi eigenvalue method repeatedly performs rotations until the matrix becomes almost diagonal. Then the elements in the diagonal are approximations of the (real) eigenvalues of S.Convergence[edit]If  p=Skl{displaystyle p=S_{kl}}  is a pivot element, then by definition  |Sij|\u2264|p|{displaystyle |S_{ij}|leq |p|}  for  1\u2264i,j\u2264n,i\u2260j{displaystyle 1leq i,jleq n,ineq j} . Let \u0393(S)2{displaystyle Gamma (S)^{2}} denote the sum of squares of all off-diagonal entries of S{displaystyle S}. Since S{displaystyle S} has exactly  2N:=n(n\u22121){displaystyle 2N:=n(n-1)}  off-diagonal elements, we have  p2\u2264\u0393(S)2\u22642Np2{displaystyle p^{2}leq Gamma (S)^{2}leq 2Np^{2}}  or 2p2\u2265\u0393(S)2\/N{displaystyle 2p^{2}geq Gamma (S)^{2}\/N} . Now \u0393(SJ)2=\u0393(S)2\u22122p2{displaystyle Gamma (S^{J})^{2}=Gamma (S)^{2}-2p^{2}}. This implies\t\u0393(SJ)2\u2264(1\u22121\/N)\u0393(S)2{displaystyle Gamma (S^{J})^{2}leq (1-1\/N)Gamma (S)^{2}}   or  \u0393(SJ)\u2264(1\u22121\/N)1\/2\u0393(S){displaystyle Gamma (S^{J})leq (1-1\/N)^{1\/2}Gamma (S)},i.e. the sequence of Jacobi rotations converges at least linearly by a factor  (1\u22121\/N)1\/2{displaystyle (1-1\/N)^{1\/2}}  to a diagonal matrix.A number of N{displaystyle N} Jacobi rotations is called a sweep; let S\u03c3{displaystyle S^{sigma }} denote the result. The previous estimate yields\u0393(S\u03c3)\u2264(1\u22121N)N\/2\u0393(S){displaystyle Gamma (S^{sigma })leq left(1-{frac {1}{N}}right)^{N\/2}Gamma (S)},i.e. the sequence of sweeps converges at least linearly with a factor  \u2248 e1\/2{displaystyle e^{1\/2}} .However the following result of Sch\u00f6nhage[3] yields locally quadratic convergence. To this end let S have m distinct eigenvalues  \u03bb1,...,\u03bbm{displaystyle lambda _{1},…,lambda _{m}}  with multiplicities  \u03bd1,...,\u03bdm{displaystyle nu _{1},…,nu _{m}}  and let  d > 0 be the smallest distance of two different eigenvalues. Let us call a number ofNS:=n(n\u22121)2\u2212\u2211\u03bc=1m12\u03bd\u03bc(\u03bd\u03bc\u22121)\u2264N{displaystyle N_{S}:={frac {n(n-1)}{2}}-sum _{mu =1}^{m}{frac {1}{2}}nu _{mu }(nu _{mu }-1)leq N}Jacobi rotations a Sch\u00f6nhage-sweep. If  Ss{displaystyle S^{s}}  denotes the result then\u0393(Ss)\u2264n2\u22121(\u03b32d\u22122\u03b3),\u03b3:=\u0393(S){displaystyle Gamma (S^{s})leq {sqrt {{frac {n}{2}}-1}}left({frac {gamma ^{2}}{d-2gamma }}right),quad gamma :=Gamma (S)}  .Thus convergence becomes quadratic as soon as\u0393(S),mn\u22121{displaystyle m_{1},,dots ,,,m_{n-1}}  with the property that mi{displaystyle m_{i}}  is the index of the largest element in row i, (i = 1, …, n\u00a0\u2212\u00a01) of the current S. Then the indices of the pivot (k, l)  must be one of the pairs (i,mi){displaystyle (i,m_{i})}. Also the updating of the index array can be done in O(n) average-case complexity: First, the maximum entry in the updated rows k and l can be found in O(n) steps. In the other rows i, only the entries in columns k and l change. Looping over these rows, if mi{displaystyle m_{i}} is neither k nor l, it suffices to compare the old maximum at mi{displaystyle m_{i}} to the new entries and update mi{displaystyle m_{i}} if necessary. If mi{displaystyle m_{i}} should be equal to k or l and the corresponding entry decreased during the update, the maximum over row i has to be found from scratch in O(n) complexity. However, this will happen on average only once per rotation. Thus, each rotation has  O(n) and one sweep O(n3) average-case complexity, which is equivalent to one matrix multiplication. Additionally the  mi{displaystyle m_{i}}  must be initialized before the process starts, which can be done in n2 steps.Typically the Jacobi method converges within numerical precision after a small number of sweeps. Note that multiple eigenvalues reduce the number of iterations since NS \u2502Sm\u00a0indm\u2502 then m\u00a0:= k endif    endfor    k\u00a0:= m; l\u00a0:= indm; p\u00a0:= Skl   \u00a0! calculate c = cos \u03c6, s = sin \u03c6    y\u00a0:= (el\u2212ek)\/2; d\u00a0:= \u2502y\u2502+\u221a(p2+y2)    r\u00a0:= \u221a(p2+d2); c\u00a0:= d\/r; s\u00a0:= p\/r; t\u00a0:= p2\/d    if y emthen            m\u00a0:= l        endif    endfor    if k \u2260 m then        swap em,ek        swap Em,Ekendifendfor4. The algorithm is written using matrix notation (1 based arrays instead of 0 based).5. When implementing the algorithm, the part specified using matrix notation must be performed simultaneously.6. This implementation does not correctly account for the case in which one dimension is an independent subspace.  For example, if given a diagonal matrix, the above implementation will never terminate, as none of the eigenvalues will change.  Hence, in real implementations, extra logic must be added to account for this case.Example[edit]LetS=(4\u22123060\u221235\u221230300\u221267542060\u22126751620\u22121050\u221235420\u22121050700){displaystyle S={begin{pmatrix}4&-30&60&-35\\-30&300&-675&420\\60&-675&1620&-1050\\-35&420&-1050&700end{pmatrix}}}Then jacobi produces the following eigenvalues and eigenvectors after 3 sweeps (19 iterations)\u00a0:e1=2585.25381092892231{displaystyle e_{1}=2585.25381092892231}E1=(0.0291933231647860588\u22120.3287120557631889970.791411145833126331\u22120.514552749997152907){displaystyle E_{1}={begin{pmatrix}0.0291933231647860588\\-0.328712055763188997\\0.791411145833126331\\-0.514552749997152907end{pmatrix}}}e2=37.1014913651276582{displaystyle e_{2}=37.1014913651276582}E2=(\u22120.1791862905354548260.741917790628453435\u22120.100228136947192199\u22120.638282528193614892){displaystyle E_{2}={begin{pmatrix}-0.179186290535454826\\0.741917790628453435\\-0.100228136947192199\\-0.638282528193614892end{pmatrix}}}e3=1.4780548447781369{displaystyle e_{3}=1.4780548447781369}E3=(\u22120.5820756994972376500.3705021850670930580.5095786345017996260.514048272222164294){displaystyle E_{3}={begin{pmatrix}-0.582075699497237650\\0.370502185067093058\\0.509578634501799626\\0.514048272222164294end{pmatrix}}}e4=0.1666428611718905{displaystyle e_{4}=0.1666428611718905}E4=(0.7926082911637635850.4519231209015997940.3224163985818249920.252161169688241933){displaystyle E_{4}={begin{pmatrix}0.792608291163763585\\0.451923120901599794\\0.322416398581824992\\0.252161169688241933end{pmatrix}}}Applications for real symmetric matrices[edit]When the eigenvalues (and eigenvectors) of a symmetric matrix are known, the followingvalues are easily calculated.Singular valuesThe singular values of a (square) matrix A are the square roots of the (non-negative) eigenvalues of ATA{displaystyle A^{T}A}. In case of a symmetric matrix S we have of STS=S2{displaystyle S^{T}S=S^{2}}, hence the singular values of S are the absolute values of the eigenvalues of S2-norm and spectral radiusThe 2-norm of a matrix A is the norm based on the Euclidean vectornorm, i.e. the largest value \u2016Ax\u20162{displaystyle |Ax|_{2}} when x runs through all vectors with \u2016x\u20162=1{displaystyle |x|_{2}=1}. It is the largest singular value of A. In case of a symmetric matrix it is the largest absolute value of its eigenvectors and thus equal to its spectral radius.Condition numberThe condition number of a nonsingular matrix A is defined as  cond(A)=\u2016A\u20162\u2016A\u22121\u20162{displaystyle {mbox{cond}}(A)=|A|_{2}|A^{-1}|_{2}}. In case of a symmetric matrix it is the absolute value of the quotient of the largest and smallest eigenvalue. Matrices with large condition numbers can cause numerically unstable results: small perturbation can result in large errors. Hilbert matrices are the most famous ill-conditioned matrices. For example, the fourth-order Hilbert matrix has a condition of 15514, while for order 8 it is 2.7\u00a0\u00d7\u00a0108.RankA matrix A has rank r if it has r columns that are linearly independent while the remaining columns are linearly dependent on these. Equivalently, r is the dimension of the range of\u00a0A. Furthermore it is the number of nonzero singular values.In case of a symmetric matrix r is the number of nonzero eigenvalues. Unfortunately because of rounding errors numerical approximations of zero eigenvalues may not be zero (it may also happen that a numerical approximation is zero while the true value is not). Thus one can only calculate the numerical rank by making a decision which of the eigenvalues are close enough to zero.Pseudo-inverseThe pseudo inverse of a matrix A is the unique matrix X=A+{displaystyle X=A^{+}} for which AX and XA are symmetric and for which AXA = A,  XAX = X  holds. If A is nonsingular, then ‘A+=A\u22121{displaystyle A^{+}=A^{-1}}.When procedure jacobi (S, e, E) is called, then the relation S=ETDiag(e)E{displaystyle S=E^{T}{mbox{Diag}}(e)E} holds where Diag(e)  denotes the diagonal matrix with vector e on the diagonal. Let e+{displaystyle e^{+}} denote the vector where  ei{displaystyle e_{i}} is replaced by 1\/ei{displaystyle 1\/e_{i}} if  ei\u22640{displaystyle e_{i}leq 0} and by 0  if  ei{displaystyle e_{i}} is (numerically close to) zero. Since matrix E is orthogonal, it follows that the pseudo-inverse of S is given by   S+=ETDiag(e+)E{displaystyle S^{+}=E^{T}{mbox{Diag}}(e^{+})E}.Least squares solutionIf matrix A  does not have full rank, there may not be a solution of the linear system  Ax = b. However one can look for a vector x for which  \u2016Ax\u2212b\u20162{displaystyle |Ax-b|_{2}} is minimal. The solution is x=A+b{displaystyle x=A^{+}b}. In case of a symmetric matrix S as before, one has x=S+b=ETDiag(e+)Eb{displaystyle x=S^{+}b=E^{T}{mbox{Diag}}(e^{+})Eb}.Matrix exponentialFrom S=ETDiag(e)E{displaystyle S=E^{T}{mbox{Diag}}(e)E} one finds exp\u2061S=ETDiag(exp\u2061e)E{displaystyle exp S=E^{T}{mbox{Diag}}(exp e)E} where  exp\u00a0e  is the vector where ei{displaystyle e_{i}} is replaced by exp\u2061ei{displaystyle exp e_{i}}. In the same way, f(S)  can be calculated in an obvious way for any (analytic) function f.Linear differential equationsThe differential equation  x’\u00a0 =\u00a0Ax, x(0) = a has the solution  x(t) = exp(t\u00a0A)\u00a0a. For a symmetric matrix S, it follows that  x(t)=ETDiag(exp\u2061te)Ea{displaystyle x(t)=E^{T}{mbox{Diag}}(exp te)Ea}. If  a=\u2211i=1naiEi{displaystyle a=sum _{i=1}^{n}a_{i}E_{i}} is the expansion of a  by the eigenvectors of S, then x(t)=\u2211i=1naiexp\u2061(tei)Ei{displaystyle x(t)=sum _{i=1}^{n}a_{i}exp(te_{i})E_{i}}.Let  Ws{displaystyle W^{s}} be the vector space spanned by the eigenvectors of S which correspond to a negative eigenvalue and Wu{displaystyle W^{u}} analogously for the positive eigenvalues. If  a\u2208Ws{displaystyle ain W^{s}} then limt\u00a0\u221ex(t)=0{displaystyle {mbox{lim}}_{t infty }x(t)=0} i.e. the equilibrium point 0 is attractive to x(t). If a\u2208Wu{displaystyle ain W^{u}} then limt\u00a0\u221ex(t)=\u221e{displaystyle {mbox{lim}}_{t infty }x(t)=infty }, i.e. 0 is repulsive to  x(t). Ws{displaystyle W^{s}} and Wu{displaystyle W^{u}} are called stable and unstable manifolds for S. If a has components in both manifolds, then one component is attracted and one component is repelled. Hence x(t) approaches Wu{displaystyle W^{u}} as t\u2192\u221e{displaystyle tto infty }.Generalizations[edit]The Jacobi Method has been generalized to complex Hermitian matrices, general nonsymmetric real and complex matrices as well as block matrices.Since singular values of a real matrix are the square roots of the eigenvalues of the symmetric matrix S=ATA{displaystyle S=A^{T}A} it can also be used for the calculation of these values. For this case, the method is modified in such a way that S must not be explicitly calculated which reduces the danger of round-off errors. Note that JSJT=JATAJT=JATJTJAJT=BTB{displaystyle JSJ^{T}=JA^{T}AJ^{T}=JA^{T}J^{T}JAJ^{T}=B^{T}B}  with  B:=JAJT{displaystyle B,:=JAJ^{T}} .The Jacobi Method is also well suited for parallelism.References[edit]Further reading[edit]Press, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007), “Section 11.1. Jacobi Transformations of a Symmetric Matrix”, Numerical Recipes: The Art of Scientific Computing (3rd\u00a0ed.), New York: Cambridge University Press, ISBN\u00a0978-0-521-88068-8Rutishauser, H. (1966). “Handbook Series Linear Algebra: The Jacobi method for real symmetric matrices”. Numerische Mathematik. 9 (1): 1\u201310. doi:10.1007\/BF02165223. MR\u00a01553948.Sameh, A.H. (1971). “On Jacobi and Jacobi-like algorithms for a parallel computer”. Mathematics of Computation. 25 (115): 579\u2013590. doi:10.1090\/s0025-5718-1971-0297131-6. JSTOR\u00a02005221. MR\u00a00297131.Shroff, Gautam M. (1991). “A parallel algorithm for the eigenvalues and eigenvectors of a general complex matrix”. Numerische Mathematik. 58 (1): 779\u2013805. CiteSeerX\u00a010.1.1.134.3566. doi:10.1007\/BF01385654. MR\u00a01098865.Veseli\u0107, K. (1979). “On a class of Jacobi-like procedures for diagonalising arbitrary real matrices”. Numerische Mathematik. 33 (2): 157\u2013172. doi:10.1007\/BF01399551. MR\u00a00549446.Veseli\u0107, K.; Wenzel, H. J. (1979). “A quadratically convergent Jacobi-like method for real matrices with complex eigenvalues”. Numerische Mathematik. 33 (4): 425\u2013435. doi:10.1007\/BF01399324. MR\u00a00553351.Yousef Saad: “Revisiting the (block) Jacobi subspace rotation method for the symmetric eigenvalue problem”, Numerical Algorithms, vol.92 (2023), pp.917-944. https:\/\/doi.org\/10.1007\/s11075-022-01377-w .External links[edit]     (adsbygoogle = window.adsbygoogle || []).push({});after-content-x4"},{"@context":"http:\/\/schema.org\/","@type":"BreadcrumbList","itemListElement":[{"@type":"ListItem","position":1,"item":{"@id":"https:\/\/wiki.edu.vn\/en\/wiki2\/#breadcrumbitem","name":"Enzyklop\u00e4die"}},{"@type":"ListItem","position":2,"item":{"@id":"https:\/\/wiki.edu.vn\/en\/wiki2\/jacobi-eigenvalue-algorithm-wikipedia\/#breadcrumbitem","name":"Jacobi eigenvalue algorithm – Wikipedia"}}]}]