|
MATH 2605: Calculus III for CS Computer Project: Jacobi Algorithm Author: Dawn Finney (gth683e) Assignment description DOWNLOAD THE ENTIRE MATLAB PROJECT I. Implementing the Jacobi Algorithm The objective of the Jacobi Algorithm is to diagonalize an N x N symmetric matrix until the diagonal holds the eigenvalues of the matrix or at least the eigenvalues approximated within some accuracy. For this project, we were asked to diagonalize a randomly generated 5 x 5 symmetric matrix until the value of the sum of the square of each off-diagonal element was less than or equal to 10-9. A. Generating the Random Matrix In order to randomly generate a 5 x 5 symmetric matrix, generate random numbers (in the case of this project number of in the set [0,9]) for the lower half of the matrix denoted below with x's: matrix = | x 0 0 0 0 |Then mirror the value onto the other side of the diagonal. For example: matrix(row,column) = 5, then matrix(column,row) = 5 as well. The values on the diagonal are in a sense "mirrored with themselves," because of how the program is ran; this does not create a problem, because the values on the diagonal are merely assigned twice to the same number both times. Ex: row=2, column=2, randNumber is some random valueB. Finding the Max Off-Diagonal Entry This portion of the algorithm is different in the two functions: jacobi() and jacobiNoMax(). In jacobi(), the program searches for the off-diagonal entry with the greatest magnitude and creates a 2 x 2 symmetric matrix using the value, a method to be discussed a little later now. However, in jacobiNoMax() the program randomly selects an off-diagonal entry and creates a 2 x 2 symmetric matrix using the value. The algorithm in jacobi() is more efficient than jacobiNoMax() because jacobi() will not select a 2 x 2 symmetric matrix that is already diagonalized, that is with 0's in its off-diagonal positions. diagonalized |3 0|jacobi() will always search with the off-diagonal entry with the greater magnitude, in other words, the greatest absolute value. In fact it is necessary to debug the jacobiNoMax() method and check if the 2 x 2 symmetric matrix is already diagonalized. If it is not done, it sometimes causes the eigenvalues to be complex values (explained later) and thus skewing the results. After the "maximum point" is determined by either methods, the program creates a 2 x 2 symmetric matrix based on the coordinates of the "maximum point." If the maximum point is defined as having the coordinate (maxRow, maxColumn) then the 2 x 2 symmetric matrix, B, should look like the following: This should hold true for N x N matrices where N >=3.C. Diagonalizing the 2 x 2 Symmetric Matrix B In order to make the formula for finding the eigenvalues a bit simpler to type and program, several variables will be set to assess the different elements within B: B = |a b|In order to diagonalize the matrix, mu must be found using the equation: det(B-(mu*I))= 0Since we know that: I = |1 0|we also know that: 0 = det | a-mu b |In this form, we can use the quadratic formula to solve for mu: mu+ = ((d+a) + sqrt((d-a)^2 + 4*b^2))/2It actually is not necessary that we solve for mu- because it can be obtained by a "flipping" of the u obtain from mu+, where the normalized U is composed of [u1 u2]/||u1||. u1 is the first column of the soon to be normalized matrix U. u1 is obtained from the first row of the u1 with the appropriate sign flipping. u = |a b|but technically we could skip the skip of actually looking for the eigenvalues if the sum of the first row is equal to the sum of the second. Ex: | 2 1 | = 3So we could assume that |1| to be u1. |1| u2 is the second column of matrix U and is obtained from u1 again with the appropriate flipping. u1 = |a|U is now the composition of u1 and u2 divided by its length. U = [u1 u2]/||u1||Since u2 is obtained from u1, their length would be the same, so it is only necessary to calculate one of them. D. Upgrading to the Givens Matrix G Now U must be upgraded into a 5 x 5 symmetric matrix G. Place U into a matrix filled with 0's. The coordinates of the U are the same as the coordinates of B and merely need to be placed back from where they came from. G(maxRow, maxCol) = U(1,2);Then we need to place 1 wherever there is a 0 on the diagonal. Example:E. Finding the More-Diagonal Matrix D and the off(D) D is equal to G^T * A * G. The off(D) is equal to the sum of the square of the each of the off-diagonal elements of D. As long as the off(D) < 10^-9, D needs to be diagonalized more, but this responsibility of checking and recalling is delegated to testJacobiSameMatrices() method. II. Actual Data Versus Theoretical Bound For the jacobi algorithm with sorting, the actual data was less than or equal to the theoretical bound, while for the algorithm without sorting, the actual data could be greater than, equal to, less than the theoretical bound at any given time except at the very end. The actual data for the algorithm with sorting was also decreasing. The algorithm without sorting did not have a clear pattern. (please see jacobiGraphWOSort.fig and jacobiGraphWSort.fig and for the actual data the graphs were based off of please see data.m)  Jacobi Algorithm WITHOUT Sorting (Same as jacobiGraphWOSort.fig)  III. Comparing the Efficiency of Sorting Versus Non-Sorting It is obvious from the number of iterations that the jacobi algorithm with sorting is more efficient than the algorithm without sorting. Each time testJacobiSameMatrices() is called it returns an average of the number of iterations for the sorting method and an average for the non-sorting based on 10 matrices. I wrote a script that called testJacobiSameMatrices() 100 times and averaged the number of iterations of the sorting method, number of iterations of the non-sorting method, also the difference between two iterations returned by each call (basically calculating how many steps non-sorting wastes). The Efficiency of Sorting Versus Non-Sorting (Based on 1000 matrices)From the data, it is safe to say that the Jacobi Algorithm is more efficient with sorting. |