Ch-07 (523180), страница 4
Текст из файла (страница 4)
According to Equation 8, the orientations of the x-, y-, and zaxes can be determined by using the rows of [Q]. For example, a unit vector I' alongx-axis is equal to 0.9623i + 0.3277j–0.3860k where i, j, and k are unit vectors alongx-, y-, and z-axes, respectively. That means, the angles between x-axis and x-, y-,and z-axes by be obtained easily as:θx′x = cos−1 (.8623) = 30.42 o , θx′y = cos−1 (.3277) = 70.87oandθx′z = cos−1 ( −.3860) = 112.7oSimilarly, we can find yx = cos–1(.4921) = 60.52°, yy = cos–1(-.7218) = 136.2°,and yz = cos –1 (.4867) = 60.88°, and zx = cos –1 (.1192) = 83.15°, zy =cos–1(.6096) = 52.44°, and zz = cos–1(.7837) = 38.40°.To verify the above assertions, MATLAB is again applied to obtain:© 2001 by CRC Press LLCThe arc-cosine function acos of MATLAB is employed above and pi ( =3.14159) also has been used for converting the results in radians into degrees.7.5 PROGRAM EIGENVIT — ITERATIVE SOLUTIONOF THE EIGENVALUE AND EIGENVECTORThe program EigenvIt is designed to iteratively solve for the largest eigenvalue inmagnitude max and its associated eigenvector {V} of a given square matrix [A].Eigenvalue and eigenvector problems are discussed in the programs CharacEq,EigenODE, and EigenVec.
Since the eigenvector {V} satisfies the matrix equation:[A]{V} = λ{V}(1)which indicates that if we make a successful guess of {V} then when it is multipliedby the matrix [A] the product should be a scaled version of {V} and the scalingvector is the eigenvalue. Of course, it is not easy to guess correctly what this vector{V} is. But, we may devise a successive guessing scheme and hope for eventualconvergence toward the needed solution.
In order to make the procedure betterorganized, let us use normalized vectors, that is, to require all guessed eigenvectorsto have a length equal to unity.The iterative scheme may be written as, for k = 0,1,2,…[A]{V( k )} = λ( k ) {V( k +1)}(2)where k is the iteration counter and {V(0)} is the initial guess of the eigenvector for[A].
The iteration is to be terminated when the differences in every components(denoted in lower case of V) of {V(k + 1)} and {V(k)} are sufficiently small. Or,mathematicallyv(i k +1 − v(i k ) < ε(3)for i = 1,2,…,N and N being the order of [A]. in Equation 3 is a predeterminedtolerance of accuracy. As can be mathematically proven.7 this iterative process leadsto the largest eigenvalue in magnitude, max, and its associated eigenvectors. If it isthe smallest eigenvalue in magnitude, min, and its associated eigenvector {V} thatare necessary to be found, the iterative procedure can also be applied but instead of© 2001 by CRC Press LLC[A] its inverse [A]–1 should be utilized. The iteration should use the equation, fork = 0,1,2,…[A]−1{V( k )} = α( k ) {V( k +1)}(4)When max is found, the required smallest eigenvalue in magnitude is to becomputed as:λ min = 1 α max(5)The program EigenvIt has been developed following the concept explainedabove.
Both QuickBASIC and FORTRAN versions are made available and listedbelow along with sample applications.QUICKBASIC VERSION© 2001 by CRC Press LLCSample ApplicationFORTRAN VERSION© 2001 by CRC Press LLCSample ApplicationMATLAB APPLICATIONA EigenvIt.m file can be created and added to MATLAB m files for iterativesolution of the eigenvalue largest in magnitude and its associated normalized eigenvector of a given square matrix. It may be written as follows:© 2001 by CRC Press LLCThe arguments of EigenvIt are explained in the comment statements which inMATLAB start with a character %.
As illustrations of how this function can beapplied, two examples are given below.Notice that the first attempt allows 10 iterations but the answer has been obtainedafter 8 trials whereas the second attempt allowing only two trials fails to convergeand V is printed as a blank vector.In fact, the iteration can be carried out without a M file. To resolve the aboveproblem, MATLAB commands can be repeatedly entered by interactive operationsas follows:© 2001 by CRC Press LLC© 2001 by CRC Press LLCNotice that the command format compact makes the printout lines to be closelyspaced.
The iteration converges after 7 trials. This interactive method of continuousiteration for finding the largest eigenvalue and its associated eigenvector is easy tofollow, but the repeated entering of the statement “>> Ntry = Ntry + 1…V =V/Lambda” in the interactive execution is a cumbersome task. To circumvent thissituation, one may enter:>> A = [2,0,3;0,5,0;3,0,2]; V = [1;0;0]; format compact>> for Ntry = 1:100, Ntry, V = A*V; Lambda = sqrt(V’*V), V = V/Lambda, pause, endThe pause command enables each iterated result to be viewed. To continue thetrials, simply press any key.
The total number of trials is arbitrarily limited at 100;the actual need is 7 trials as indicated by the above printed results. Viewer canterminate the iteration by pressing the <Ctrl> and <Break> keys simultaneously aftersatisfactorily seeing the 7th, converged results of Lambda = 5 and V = [0.7071; 0;0.7071] being displayed on screen.To iteratively find the smallest eigenvector and its associated, normalized eigenvector of [A] according to Equations 4 and 5, we apply the iterative method to [A]–1by entering>> Ainv = inv(A); V = [1;0;0];>> for Ntry = 1:100, Ntry, V = Ainv*V; Lambda = 1/sqrt(V’*V), V = V*Lambda, pause, end© 2001 by CRC Press LLCThe resulting display is:For saving space, the last four iterations are listed in the column on the right.The components of the last two iterated normalized eigenvectors are numericallyequal but differ in sign, it suggests that the eigenvalue is actually equal to –1.0000instead of 1.0000.MATHEMATICA APPLICATIONSThe While command of Mathematica can be effectively applied here for iteration of eigenvalue and eigenvector.
It has two arguments, the first is a testingexpression when the condition is true the statement(s) specified in the second argument should then be implemented. When the condition is false, the While statement© 2001 by CRC Press LLCshould be terminated. To obtain the largest eigenvalue in magnitude, we proceeddirectly with a given matrix as follows:Input[1]: = A = {{2,0,3},{0,5,0},{3,0,2}}; V = {1,0,0}; I = 0;Input[2]: = While[i<8, I = I + 1; VN = A.V; Lambda = Sqrt[VN.VN];Print["Iteration # ",I," ","Lambda = ",N[Lambda]," {V} =",V = N[VN/Lambda]]]Output[2]: =IterationIterationIterationIterationIterationIterationIterationIteration########12345678Lambda =Lambda =Lambda =Lambda =Lambda =Lambda =Lambda =Lambda =3.60555 {V} = {0.5547, 0, 0.83205}4.90682 {V} = {0.734803, 0, 0.67828}4.99616 {V} = {0.701427, 0, 0.712741}4.99985 {V} = {0.708237, 0, 0.709575}4.99999 {V} = {0.70688, 0, 0.707333}5.
{V} = {0.707152, 0, 0.707062}5. {V} = {0.707098, 0, 0.707016}5. {V} = {0.707109, 0, 0.707105}Notice that the testing condition is whether the running iteration counter I isless than 8 or not. This setup enables up to eight iterations to be conducted. Thefunction N in Input[2] requests the value of the variable inside the brackets to begiven in numerical form. For example, when the value of Lambda is displayed assqrt[4], it will be displayed as 2.00000 if the input is N[Lambda]. VA·VB is the dotproduct of VA and VB. In Input[2] some sample printouts of the character stringsspecified inside a pair of parentheses are also demonstrated.For iterating the smallest eigenvalue in magnitude, we work on the inverse ofthis given matrix as follows:Input[3]: = Ainv = Inverse[A]; V = {1,0,0}; I = 0;Input[4]: = While[i<8, I = I + 1; VN = Ainv.V; Lambda = Sqrt[VN.VN];Print["Iteration # ",I," ","Lambda = ",N[Lambda]," {V} = ",V = N[VN/Lambda]]]Output[4]: =IterationIterationIterationIterationIterationIterationIterationIteration© 2001 by CRC Press LLC########12345678Lambda =Lambda =Lambda =Lambda =Lambda =Lambda =Lambda =Lambda =0.72111 {V} = {–0.5547, 0, 0.83205}0.981365 {V} = {0.734803, 0, –0.67828}0.999233 {V} = {–0.701427, 0, 0.712741}0.999969 {V} = {0.708237, –0, 0.709575}0.999999 {V} = {–0.70688, 0, 0.707333}1.
{V} = {0.707152, 0, –0.707062}1. {V} = {–0.707098, 0, 0.707016}1. {V} = {0.707109, 0, –0.707105}Inverse is a Mathematica function which inverts a specified matrix. The aboveprintout of eight iterations shows that {V} continues to change its sign. This is anindication that the eigenvalue carries a minus sign.7.6 PROBLEMSPROGRAMS EIGENODE.STBANDEIGENODE.VIB1. Apply the program EigenODE.Stb for the cases N = 6, 7, and 8 to obtainthe coefficient matrix [C] in the standard eigenvalue problem of ([C]λ[I]){Y} = {0}.2. Apply the program CharacEq to obtain the characteristic equations oforder 6, 7, and 8 for the matrices [C] derived in Problem 1.3.
Apply the program Bairstow to find the roots for the characteristic equationsderived in Problem 2. If necessary, change this program to allowinteractive input of the u and v values for the guessing quadratic factor2 + u + v. This enhancement will be helpful if the iteration fails toconverge.4. The program EigenODE.Vib has been arranged for solving the generalproblem of having N masses, m1–mN, in series connected by N + 1 springswith stiffnesses k1–kN + 1. Apply it for the case when the three massesshown in Figure 1 are connected by four springs with the fourth springattached to the ground. Use m1 = m2 = m3 = 1 and k1 = k2 = k3 = k4 = 10.5.