A numerical method for eigensolution of locally modified systems based on the inverse power method.pdf

1、Finite Elements in Analysis and Design 45 (2009) 113-120 Contents lists available at ScienceDirect FiniteElementsinAnalysisandDesign journal homepage: Anumericalmethodforeigensolutionoflocallymodifiedsystemsbasedontheinverse powermethod Mitsuhiro Kashiwagi Department of Architecture, School of Indus

2、trial Engineering, Tokai University, 9-1-1 Toroku, Kumamoto 862-8652, Japan A R T I C L EI N F OA B S T R A C T Article history: Received 7 March 2008 Received in revised form 11 July 2008 Accepted 27 July 2008 Available online 2 September 2008 Keywords: Locally modified systems Inverse power method

3、 Eigenvalue Eigenvector Reanalysis Rapid reanalysis of eigensolutions after modification is a problem of considerable practical importance. Several methods have been developed for computing eigenvalues by using modified parts. This paper pro- poses just such a method based on inverse power method, u

4、sing the DOF of only the modified part. This approach enables exact solutions for eigenvalues to be found quickly, regardless of the magnitude of the modification, using the degree of the condensed modified part. The advantages of the proposed method are examined comparing the solutions of inverse p

5、ower method by several numerical examples. This approach will be useful for a modified system with large degrees of freedom and small modification part. 2008 Elsevier B.V. All rights reserved. 1. Introduction Vibration engineering encompasses a wide variety of fields, from mechanical vibrations to e

6、lectrical oscillations and the swaying of large structures. Dynamic analysis of eigenvalues and eigenvectors during the free vibration of such systems (kinetics) is a fundamen- tal branch of engineering, and a large number of numerical methods have been proposed for efficiently performing such analy

7、ses. Recent advances in the finite element method (FEM) and in large-scale com- puting have enabled analysts to treat many more degrees of free- dom (DOF), and increasingly versatile software is being developed to accommodate these factors. However, if any parameter is changed (e.g., shape, material

8、, initial conditions or environmental conditions) in a previously analyzed system, the entire system has to be rean- alyzed. Eigenvalue analysis is a much more arduous task than the corresponding static analysis, and the computation time for the for- mer is generally several times or even several te

9、ns of times greater than that for the latter. This can make the process quite expensive if a system undergoes multiple revisions, since analyzing the system after the revision requires just as much time, labor and money as analyzing the original system. “Reanalysis” is the term used to describe the

10、partial analysis (as opposed to a full analysis) that is performed in order to obtain data efficiently when a portion of a matrix has been altered due to a design revision or some other change. If a convenient method can be E-mail address: mkashiktmail.tokai-u.jp. 0168-874X/$-see front matter2008 El

11、sevier B.V. All rights reserved. doi:10.1016/j.finel.2008.07.009 found for performing reanalysis that utilizes the previously calcu- lated eigenvalues and eigenvectors, it will be very useful for solving eigenvalue problems involving partially altered matrices. Reanalysis methods have received inten

12、se attention, and many papers have ad- dressed the problem of changes in eigenvalues in locally modified systems19.Manyofthesepapersusedperturbationtheory1,4,7. Fox and Kapoor 1 and Rogers 7 applied first-order differential to the eigenvalues and eigenvectors of the design variables. This kind of pe

13、rturbation-based solution method is valid when there are rel- atively small adjustments to the design variables and when the ini- tial values are assumed to be unchanged. However, if the changes to the design variables are large, these approaches generally offer low accuracy even when the perturbati

14、on order is increased. Hirai et al. 2 demonstrated a way to obtain the exact eigenval- ues and eigenvectors using only the degree of the modified part, and Parazzola et al. 6,10 applied this method to theoretical solu- tions for eigenvalues and eigenvectors in damped systems. These approaches use th

15、e eigenvalues and eigenvectors of the unmodified system to determine exact solutions on the basis of a fundamental formula having the same degree as the matrix representing the modified system. They remain valid regardless of the magnitude of the modifications. If eigenvalue problems can be reanalyz

16、ed using such a condensed equation, it will mean that matrices can be sim- plified to lower degrees, which is a very effective technique for re- ducing the calculation time. This fundamental equation is nonlinear and it can be used to solve a matrix equation consisting of rational functions. It is v

17、ery time consuming to solve this kind of problem by trial and error and such a process does not indicate the order of the eigenvalues. Furthermore, there is little prospect of achieving 114M. Kashiwagi / Finite Elements in Analysis and Design 45 (2009) 113-120 a stable calculation process using the

18、NewtonRaphson method or other methods in which the differential coefficients are set to ap- propriate values for local solutions and the initial estimate must be reasonably close to the actual solution (since the process becomes unpredictable if it is not). Even if an approximate range for the solu-

19、 tion has been determined, it is difficult to find stable solutions using successive approximation methods such as reverse linear interpo- lation or polynomial approximation. Using one of these approaches may well provide some solutions, but it is well known that numer- ical calculations can provide

20、 solutions in a random order. When a solution has been found using, for example, the NewtonRaphson method, a deflation is set up and the procedure is repeated to deter- mine the next eigenvalue. Thus, conventional methods have some fundamental weaknesses for solving these nonlinear equations. Kashiw

21、agi et al. have conducted a systematic study of compressed versions of the equations published by Hirai et al. 1116. They have addressed the problem of finding all the eigenvalues in a locally mod- ified system, proposing combinations of the DurandKerner method 12 with Newtons method 16 and with rat

22、ional functions 13 to find reliable solutions for eigenvalues in low-degree systems. They have also demonstrated that the Sturm sequence is useful for certain compressed, strongly nonlinear equations 14,15, discussed how to identify regions in which Sturms eigenvalues exist, and described solutions

23、using the Sturm sequence bisection method. The Sturm sequence method is a conventional, commonly used method for cal- culating eigenvalues; generally, it is most suitable for matrices that have been transformed into tridiagonal matrices 9. The process of converting a matrix to tridiagonal form requi

24、res a considerable amount of calculation and generally represents more than half of the computation of identifying eigenvalues. The version of the Sturm se- quence method developed by Kashiwagi et al. requires identifying all of the eigenvalues and eigenvectors of the unmodified system, but does not

25、 require transforming the matrix into tridiagonal form, so it represents a unique contribution to this field. Many eigensolution problems seek just a few of the eigenvalues and their corresponding vectors from low-degree systems. Inverse power functions are an example of the iterative methods used t

26、o obtain solutions in such systems 17. The inverse power method is a basic method for obtaining solutions in eigenvalue analysis; many approaches have been based on it. For example, subspace iteration 18,19 is a very commonly used method for determining eigenval- ues. Thus, we have a nearly complete

27、 toolbox of methods for deter- mining the eigenvalues of a locally modified system, but currently no method based on inverse power functions has been proposed. This paper proposes just such a method based on inverse power method, using the DOF of only the modified part. This approach en- ables exact

28、 solutions for eigenvalues to be found quickly, regardless of the magnitude of the modification, using the degree of the con- densed modified part. When the degree of that matrix is low, the calculation time is short, especially when the inverse power func- tion is combined with a shift of the origi

29、n. All of the eigenvalues and eigenvectors of the unmodified system must be known, so for relatively small systems it is effective to determine just the first sev- eral eigensolutions (10 or less for practical systems), beginning from the smallest eigensolution. In the following sections, the theory

30、 and the algorithm of the proposed method are described and the effec- tiveness of this approach is demonstrated by numerically solving a typical eigenvalue problem. 2. Theory and algorithm for the inverse power method for locally modified systems This section describes the theory of the inverse pow

31、er method for locally modified systems and the theory of what is here termed the shifted inverse power method for locally modified systems. These theories enable the eigensolution of just the condensed version of the modified part to be exactly determined. 2.1. Inverse power method for locally modif

32、ied systems The general eigenvalue problem is as follows, assuming an n n real symmetric matrix A and a positive real symmetric matrix B: A?i=?iB?i(1) where?iis the ith eigenvalue, beginning with the smallest value and ?iis the eigenvector corresponding to eigenvalue?i. A and B take the following fo

33、rms: A = aij = a11a12.a1n a22.a2n sym. . . . . ann (2) B = bij = b11b12.b1n b22.b2n sym. . . . . bnn (3) If?is the mode matrix containing the eigenvectors for the problem expressed in Eq. (1) ?= ?1,?2, .,?n(4) then the following relationships hold: ?TA?=?= ?1o ?2 . o?n (5) ?TB?= I(6) Here, I is the

34、n n identity matrix. From Eq. (5) we obtain A1= ? i ?i?T i ?i (7) Free structural vibration and buckling are two eigenvalue problems thataretypicallyanalyzedbytheFEM.Forthefreevibrationproblem, A is the stiffness matrix K and B is the mass matrix M. In buckling, A is the stiffness matrix K and B is

35、the geometric stiffness matrix KG. Under a local modification, A is replaced by A +?A and B is replaced by B +?B. Let us assume that?A and?B are as follows: ?A = ?aij = 00 ?aii00?aij . . . 00 . . . . . . . sym.00 . . . ?ajj . . . 0 (8) M. Kashiwagi / Finite Elements in Analysis and Design 45 (2009)

36、113-120115 ?B = ?bij = 00 ?bii00?bij . . . 00 . . . . . . . sym.00 . . . ?bjj . . . 0 (9) If the zero elements are removed from?A and?B, we obtain?A and?B: ?A = ? ?aii?aij sym.?ajj ? (10) ?B = ? ?bii?bij sym.?bjj ? (11) These are condensed versions of the modified parts (mm matrices, assuming that m

37、 = 2 in Eqs. (8) and (9). Let us consider an m n Boolean matrix (shown with m = 2 in Eq. (12): E = ij ?0 100 0010 ? (12) ?A and?B are given by ?A = ET?AE(13) ?B = ET?BE(14) The general eigenvalue problem expressed in Eq. (1) can be used to express the locally modified matrix given above: A +?Ax =?B

38、+?Bx(15) Eq. (15) can be rewritten in Boolean matrix form: A + ET?AEx =?B + ET?BEx(16) Let us rewrite the above, with P as an n n matrix, Q as an n m matrix, R as an m n matrix and Imas the m m identity matrix. Then, we can state P + QR1= P1 P1QIm+ RP1Q1RP1(17) The inverse matrix of A + ET?AE is A +

39、 ET?AE1= A1 A1ETIm+?AEA1ET1?AEA1 = A1 A1ETIm+?AA11?AEA1(18) Here,A1=EA1ETand is an mm matrix composed of the elements corresponding to the modified portion of A1. By the inverse power method, if the next trial value for the eigen- value is x, then from Eqs. (16) and (18) we obtain x = A + ET?AE1B +

40、ET?BEx = I A1ETIm+?AA11?AEA1Bx + A1?Bx(19) We can also state x = ? j ?bj?j(20) ?biis?iin Eq. (25) below;?iis the approximate eigenvector obtained in the previous iteration. Considering the locality of?B (Eq. (9), we obtain A1Bx + A1?Bx = ? i ?i?T i ?i B ? j ?bj?j+ ? i ?i?T i ?i ?Bx = ? i ?bi ?i ?i+

41、? i ?T i?Bx ?i ?i = ? i ?bi+ ?T i? B x ?i ?i(21) Here, ?iand x are condensed vectors from the unmodified part ?T i = ?ii,?ij = E?i(22) x T = xi,xj = Exi(23) Then, x is given by x = I ? j ?j ? T j ?j I +?A ? k ?k ? T k ?k 1 ?AE ? i ?bi+ ? T i? B x ?i ?i = ? i ?bi+? T i? B x ?i ?i ? i ? T i I+?A?j ?j

42、? T j ?j 1 ?A?k ?bk+? T k? B x ?k ?k ?i ?i = ? i ?bi+ ? T i? B x ? T i I +?A?j ?j ? T j ?j 1 ?A?k ?bk+ ? T k? B x ?k ?k ?i ?i = ? i ?i?i(24) where ?i= ?bi+ ?T i? B x ?T i I +?A?j ?j ?T j ?j 1 ?A?k ?bk+ ?T k? B x ?k ?k ?i (25) Accordingly, we can seek the eigenvector using an expression that accounts

43、 for locality. Also, from Eq. (24), we require only the DOF of the condensed eigenvector, which can be obtained using only the condensed modified part: x = I ? j ?j ? T j ?j I +?A ? k ?k ? T k ?k 1 ?A ? i ?bi+ ? T i? B x ?i ?i = ? i ?bi+ ? T i? B x ? T i I +?A?j ?i ? T i ?j 1 ?A?k ?bk+ ? T k? B x ?k

44、 ?k ?i ?k = ? i ?i ?i(26) ?i= ?bi+ ?T i? B x ?T i I +?A?j ?j ?T j ?j 1 ?A?k ?bk+ ?T k? B x ?k ?k ?i (27) Once?iin Eq. (26) have been found and stored, the approximate eigenvector can be calculated using the DOF indicated by Eq. (24). 116M. Kashiwagi / Finite Elements in Analysis and Design 45 (2009)

45、 113-120 2.2. Shifted inverse power method for locally modified systems Let us designate the current eigenvalue?and the next eigen- value as?next. The rate of convergence of the power function is then given by |?|/|?next|. When the origin is shifted to?0, the rate of con- vergence of the inverse pow

46、er function with origin shift is given by |?0|/|?next?0|, a lower value than for the unshifted inverse power function. Thus, the convergence can potentially be greatly improved during reanalysis by selecting appropriate values of?0. In this section, we derive a shifted inverse power method for local

47、ly modified systems. The fundamental equation (Eq. (28) has the same format as a shifted inverse power method for a complete system, so it is expected to be effective for a locally modified system. This restates Eq. (1) as a general eigenvalue problem representing a shifted inverse power method for

48、a locally modified system: A +?A ?0B +?Bx = (?0)B +?Bx(28) Eq. (7) can also be rewritten as A ?0B1= ? i ?i?T i ?i?0 (29) We follow the same method as the derivation given in the previous section: x = I ? j ?j ?T j ?j?0 I + ?A ?0?B ? k ?k ?T k ?k?0 1 ?A ?0?BE ? i ?bi+ ?T i? B x ?i?0 ?i = ? i ?bi+ ?T

49、i? B x ?i?0 ?i ? i ?T i I + ?A ?0?B?j ?j ?T j ?j?0 1 ?A ?0?BE?k ?bk+ ?T k? B x ?k?0 ?k ?i?0 ?i = ? i ?bi+ ?T i? B x ?T i I + ?A ?0?B ? j ?j ?T j ?j?0 1 ?A ?0?B ? k ?bk+ ?T k? B x ?k?0 ?k ?i?0 ?i= ? i ?i?i(30) ?i= ?bi+ ?T i? B x ?T i I + ?A ?0?B ? j ?j ?T j ?j?0 1 ?A ?0?B ? k ?bk+?T k? B x ?k?0 ?k ?i

50、?0 (31) This enables us to seek an approximate eigenvector while accounting for locality. Also, from Eq. (30), the approximate eigenvector based on the condensed modified part: x = I ? j ?j ?T j ?j?0 I + ?A ?0?B ? k ?k ?T k ?k?0 1 ?A ?0?B ? i ?bi+ ?T i? B x ?i?0 ?i = ? i ?bi+ ?T i? B x ?T i I + ?A ?

51、0?B ? j ?j ?T j ?j?0 1 ?A ?0?B ? k ?bk+ ?T k? B x ?k?0 ?k ?i?0 ?i= ? i ?i ?i(32) ?i= ?bi+ ?T i? B x ?T i I + ?A ?0?B ? j ?j ?T j ?j?0 1 ?A ?0?B ? k ?bk+?T k? B x ?k?0 ?k ?i?0 (33) Once?ihave been found with Eq. (32) and saved, the approximate eigenvector can be calculated using the DOF given by Eq.

52、(30). 2.3. Initial approximate eigenvector and initial approximate eigenvalue This method of re-analysis requires previous knowledge of the eigenvectors and eigenvalues of the unmodified system. We discuss how these data can be used to generate initial approximations for the ith eigenvector xiand ei

53、genvalue?i. We start with: xi=?i?i(34) Considering the normalization condition EQ xT i B +?Bxi= 1 we find ?2 i = 1 1 + ?T i? B ?i (35) indicates that if we consider the locality of?B, we obtain xi= ? ? ? ? 1 1 + ?T i? B ?i ?i(36) The initial approximate eigenvector representing only the condensed mo

54、dified part is x i= ? ? ? ? 1 1 + ?T i? B ?i ?i(37) Also, considering the initial approximate eigenvalues?i EQ?i= xT i A +?Axi M. Kashiwagi / Finite Elements in Analysis and Design 45 (2009) 113-120117 and considering the localization of?A and?B, we find ?i= ?i+ ?T i? A ?i 1 + ?T i? B ?i (38) 2.4. E

55、igenvector normalization and eigenvalues We turn to a description of how the ith approximate eigenvector x iis normalized and how the approximate eigenvalues?iare calcu- lated. We specify x i=? x i (39) and the conditions EQ x T iB +?B x i= 1 Eqs. (24) and (30) and the localization of?B hold; we the

56、n obtain ?2= 1 ? i?2i + xT i ?B xi (40) Then, x i= ? ? ? ? 1 ? i?2i + xT i ?B xi x i (41) and the initial approximate eigenvector according to the condensed modified part is x i= ? ? ? ? 1 ? i?2i + xT i ?B xi x i (42) We specify the initial approximate eigenvalues?i EQ?i= x T iA +?A x i Also, from E

57、qs. (24) and (30) and the localization of?A, we can find ?i: ?i=?2 ? i ?2 i ?i+ x T i? A x i (43) This manipulation requires that?ibe replaced with?i. 2.5. Extraction of the already-obtained eigenvalues When finding higher-order eigenvalues, it is necessary to extract the eigensolutions that have al

58、ready been found using the inverse power function. GramSchmidt orthogonalization is used, in con- junction with Eqs. (6), (24) and (30) and the locality of?B: x i= xi i1 ? j=1 xT j B +?B xixj = xi i1 ? j=1 ? k ?jk?T kB ? l ?jl?T l + xT j?B xi xj = xi i1 ? j=1 ? k ?jk?ik+ xT j? B x i xj (44) Also, x i= x i i1 ? j=1 ? k ?jk?ik+ xT j? B x i x j (45) An equation is obtained from just the condensed modified part. If we define the?th component of?ias?il, then we can rewrite?ilas follows: ?i?=?i? i1 ? j=1 ? k ?jk?ik+ xT j? B x i ?j