김경종
31 KiB
elements needed to discretize the penetration depth, we use
\Delta x = \frac {4}{N} \sqrt {a t _ {\min}} \tag {9.127}
Typically, for two-node elements, we use N = 6 to 10, and this element size would be used throughout the mesh.
Next, the time step \Delta t for the integration scheme must be chosen. Assume that we use two-node elements and a lumped heat capacity matrix. If the Euler forward method is employed, the stability limit dictates the time step to be (see Exercise 9.27),
\Delta t \leq \frac {(\Delta x) ^ {2}}{2 a} \tag {9.128}
Whereas using the trapezoidal rule or Euler backward method we can employ
\Delta t = \frac {(\Delta x) ^ {2}}{a} \tag {9.129}
or an even larger time step.
When using the implicit methods, \alpha = 1 or \frac{1}{2} , it is, however, frequently effective to use higher-order (parabolic) elements and the consistent heat capacity matrix because significantly better accuracy might be achieved. In this case \Delta x is the distance between adjacent nodes.
These considerations are also directly useful in two- and three-dimensional analysis. In these cases we use the Euler forward method with the low-order elements (four-node quadrilateral elements in two-dimensional analysis and eight-node brick elements in three-dimensional analysis) and a lumped heat capacity matrix, and a mesh as close as possible to a uniform mesh. The time step is given by (9.128), where \Delta x is now the smallest distance between any nodes. On the other hand, with the trapezoidal rule or Euler backward method, we usually use the parabolic elements (nine-node elements in two-dimensional analysis and 27-node elements in three-dimensional analysis), the consistent heat capacity matrix, and a time step of magnitude (9.129), where \Delta x is again the smallest distance between any nodes.
Finally, let us demonstrate the application of the \alpha -method of time integration in two examples.
EXAMPLE 9.15: Develop the equations to be solved in a nonlinear transient heat transfer analysis using the Euler backward method and the full Newton-Raphson method.
The governing equations in a nonlinear heat transfer analysis using an implicit time integration procedure are (see Section 7.2.3)
{ } ^ { t + \Delta t } \mathbf { C } ^ { ( i ) } { } ^ { t + \Delta t } \dot { \boldsymbol { \theta } } ^ { ( i ) } + { } ^ { t + \Delta t } \mathbf { K } ^ { ( i - 1 ) } \Delta \boldsymbol { \theta } ^ { ( i ) } = { } ^ { t + \Delta t } \tilde { \mathbf { Q } } ^ { ( i - 1 ) } \tag {a}
where t^{+\Delta t}\tilde{\mathbf{Q}}^{(t-1)} is a vector of nodal point heat flows corresponding to time t + \Delta t and iteration (i - 1) . Using the Euler backward method, we have
{ } ^ { t + \Delta t } \dot { \boldsymbol { \theta } } ^ { ( i ) } = \frac { { } ^ { t + \Delta t } \boldsymbol { \theta } ^ { ( i - 1 ) } + \Delta \boldsymbol { \theta } ^ { ( i ) } - { } ^ { t } \boldsymbol { \theta } } { \Delta t }
and hence (a) reduces to
\left(^ {t + \Delta t} \mathbf {K} ^ {(i - 1)} + \frac {1}{\Delta t} ^ {t + \Delta t} \mathbf {C} ^ {(i)}\right) \Delta \boldsymbol {\theta} ^ {(i)} = ^ {t + \Delta t} \tilde {\mathbf {Q}} ^ {(i - 1)} - ^ {t + \Delta t} \mathbf {C} ^ {(i)} {} ^ {t + \Delta t} \dot {\boldsymbol {\theta}} ^ {(i - 1)} \tag {b}
where
{ } ^ { t + \Delta t } \dot { \boldsymbol { \theta } } ^ { ( i - 1 ) } = \frac { { } ^ { t + \Delta t } \boldsymbol { \theta } ^ { ( i - 1 ) } - { } ^ { t } \boldsymbol { \theta } } { \Delta t }
For solution the relation (b) is further linearized (corresponding to a full Newton-Raphson iteration) using
\left(t + \Delta t \mathbf {K} ^ {(i - 1)} + \frac {1}{\Delta t} t + \Delta t \mathbf {C} ^ {(i - 1)}\right) \Delta \boldsymbol {\theta} ^ {(i)} = t + \Delta t \tilde {\mathbf {Q}} ^ {(i - 1)} - t + \Delta t \mathbf {C} ^ {(i - 1)} t + \Delta t \dot {\boldsymbol {\theta}} ^ {(i - 1)}
EXAMPLE 9.16: Referring to Example 9.15, develop the equations to be solved in incompressible fluid flow without heat transfer using the Euler backward method.
The governing finite element equations in fluid flow are given in (7.74) and (7.75), which we can restate as
\left[ \begin{array}{l l} \mathbf {M} & \mathbf {0} \\ \mathbf {0} & \mathbf {0} \end{array} \right] \left[ \begin{array}{l} \dot {\mathbf {V}} \\ \dot {\mathbf {P}} \end{array} \right] + \left[ \begin{array}{l l} \mathbf {K} & \mathbf {K} _ {p} \\ \mathbf {K} _ {p} ^ {T} & \mathbf {0} \end{array} \right] \left[ \begin{array}{l} \mathbf {V} \\ \mathbf {P} \end{array} \right] = \left[ \begin{array}{l} \mathbf {R} \\ \mathbf {0} \end{array} \right] \tag {a}
where the vector V denotes all velocity degrees of freedom and P lists all pressure degrees of freedom.
We note that the fully incompressible flow condition results in the zero diagonal elements corresponding to the pressure degrees of freedom. Hence, an implicit time integration is necessary. In the Euler backward integration we use (now showing the superscripts denoting time and iteration number),
\begin{array}{l} \left(\left[ \begin{array}{c c} \mathbf {K} & \mathbf {K} _ {p} \\ \mathbf {K} _ {p} ^ {T} & \mathbf {0} \end{array} \right] \right| _ {\text { evaluated at } t + \Delta t \mathbf {V} (i - 1), t + \Delta t \mathbf {P} (i - 1)} + \frac {1}{\Delta t} \left[ \begin{array}{c c} \mathbf {M} & \mathbf {0} \\ \mathbf {0} & \mathbf {0} \end{array} \right] \left[ \begin{array}{c} \Delta \mathbf {V} ^ {(i)} \\ \Delta \mathbf {P} ^ {(i)} \end{array} \right] \\ = \left[\begin{array}{c}t + \Delta t \mathbf {R}\\\mathbf {0}\end{array}\right] - \left( \right.\left[\begin{array}{l l}\mathbf {M}&\mathbf {0}\\\mathbf {0}&\mathbf {0}\end{array}\right]\left[\begin{array}{l}t + \Delta t \dot {\mathbf {V}} (i - 1)\\t + \Delta t \dot {\mathbf {P}} (i - 1)\end{array}\right] + \left[\begin{array}{l l}\mathbf {K}&\mathbf {K} _ {p}\\\mathbf {K} _ {p} ^ {T}&\mathbf {0}\end{array}\right]\left. \right| _ {\text {evaluated at}} \left[\begin{array}{l}t + \Delta t \mathbf {V} (i - 1)\\t + \Delta t \mathbf {P} (i - 1)\end{array}\right]\left. \right) \tag {b} \\ \end{array}
The equations in (b) correspond simply to a solution by successive substitution because the coefficient matrices are not tangent matrices. A Newton-Raphson type of iteration would have to be developed as described in Sections 6.3.1 and 8.4.1.
In practice, the simple iteration in (b) frequently works quite well, particularly if the time step \Delta t is sufficiently small. However, the right-hand-side vector must then be very efficiently calculated by minimizing the multiplications actually required.
In the preceding example solutions, we considered nonlinear systems with an implicit time integration scheme. Hence, as pointed out in Section 9.5.2, it is important to iterate for the solution. In heat transfer analysis we can also employ the explicit Euler forward method, in which case no iteration is performed, but the time step \Delta t must be smaller than 2/\lambda_{n} , where \lambda_{n} is the largest eigenvalue of the problem K\phi = \lambda C\phi , with K and C changing in nonlinear analysis (see Exercise 9.25). This time step can be estimated from (9.128) or the matrices of the individual finite elements (see Example 9.13 and Exercise 9.27). On the other hand, in incompressible fluid flow solutions, an unconditionally stable implicit time integration method must be used for the pressure equations, whereas implicit or explicit integration can be performed for the velocity and temperature equations (see Section 9.2.5 for an example of coupling integration schemes, and Exercise 9.28).
9.6.2 Exercises
9.25. The governing heat transfer equations of a general linear finite element system are
\mathbf {C} \dot {\boldsymbol {\theta}} + \mathbf {K} \boldsymbol {\theta} = \mathbf {Q}
\boldsymbol {\theta} | _ {\text { time } 0} = ^ {0} \boldsymbol {\theta} \tag {a}
(i) Consider the unforced system with \mathbf{Q} = \mathbf{0} , assume \theta = \phi e^{-\lambda t} , and develop the eigenproblem
\mathbf {K} \boldsymbol {\phi} = \lambda \mathbf {C} \boldsymbol {\phi} \tag {b}
(ii) Use the eigensolutions of (b) and show how to solve the equations in (a) by mode superposition.
(iii) Assume now the specific case
\mathbf {C} = \left[ \begin{array}{l l} \frac {1}{2} & 0 \\ 0 & 1 \end{array} \right]; \quad \mathbf {K} = \left[ \begin{array}{l l} 2 & - 1 \\ - 1 & 1 \end{array} \right]
\mathbf {Q} = \left[ \begin{array}{c} 1 \\ 0 \end{array} \right]
and obtain the solution using the mode superposition proposed in part (ii).
9.26. Assume that for the general system in Exercise 9.25 a mode superposition solution is performed using only the modes corresponding to the smallest p eigenvalues (i.e., using only the equations corresponding to \lambda_1, \ldots, \lambda_p ). Show how the error in the solution of the governing finite element equations could be evaluated and also develop a correction scheme similar to the static correction used in the dynamic analysis of structural systems [see (9.52)].
9.27. Assume that the Euler forward method is used in the solution of the governing heat transfer equation C\dot{\theta} + K\theta = Q . Show how the critical time step can be evaluated from the individual element matrices. (Hint: See Example 9.13.) Apply your result to evaluate the critical time step for the finite element model of the one-dimensional heat transfer problem shown and compare your result with the value given in (9.128).
text_image
θL x L ∞ θR ∞
Uniform slab of material constants k, \rho c ; initial temperature \theta_{i} ; suddenly \theta_{L} and \theta_{R} are applied
text_image
θL x θR
Model of n one-dimensional 2-node elements, unit cross-sectional area
9.28. Consider the governing finite element equations of incompressible transient fluid flow (7.74) to (7.76). Propose a time integration scheme in which the pressure equations are integrated implicitly and the velocity and temperature equations are integrated explicitly (see Section 9.2.5).
9.29. Use a finite element program to solve for the transient response of the mathematical model shown in the figure in Exercise 9.27. Choose a reasonable finite element discretization and select a time integration scheme and an appropriate time step \Delta t . Show that your results are reasonably accurate. In this analysis, k = 0.10, \rho c = 0.01 , \theta_{i} = 70 , \theta_{R} = 70 , \theta_{L} = 400 , L = 10.
9.30. Proceed as in Exercise 9.29 but for the two-dimensional mathematical model of a corner where the surface temperature \theta^{s} is suddenly applied at time 0^{+} .
text_image
Region discretized, unit thickness θ^S = 50 k = 35.0 ρc = 100.0 0θ = 0 0.75 y x 0.75 θ^S = 50
9.31. Use a computer program to solve for the transient response of a fluid between two rotating cylinders. Use the same geometries and material properties as in Exercise 7.28. Assume that the cylinders are initially at rest and that the full rotational velocities are developed in the time interval 0 to 1.
Compare your results with the analytically calculated steady-state solution (see, for example, F. M. White [A]).
9.32. Use a computer program to solve for the transient response of a fluid between two plates; see the figure for the data to be used. The bottom plate is at rest, and the top plate starts from rest and with a linear increase reaches a steady-state velocity V. The fluid is subjected to a pressure gradient.
Compare your results with the analytically calculated steady-state solution (see, for example, H. Schlichting [A]).
text_image
V dp/dy = -1 V = 60 h = 3 ρ = 1 μ = 0.01
Preliminaries to the Solution of Eigenproblems
10.1 INTRODUCTION
In various sections of the preceding chapters we encountered eigenproblems and the statement of their solutions. We did not at that time discuss how to obtain the required eigenvalues and eigenvectors. It is the purpose of this and the next chapter to describe the actual solution procedures used to solve the eigenproblems of interest. Before presenting the algorithms, we discuss in this chapter some important basic considerations for the solution of eigenproblems.
First, let us briefly summarize the eigenproblems that we want to solve. The simplest problem encountered is the standard eigenproblem,
\mathbf {K} \boldsymbol {\phi} = \lambda \boldsymbol {\phi} \tag {10.1}
where K is the stiffness matrix of a single finite element or of an element assemblage. We recall that K has order n, and for an element assemblage the half-bandwidth m_{K} (i.e., the total bandwidth is 2m_{K} + 1 ), and that K is positive semidefinite or positive definite. There are n eigenvalues and corresponding eigenvectors satisfying (10.1). The ith eigenpair is denoted as (\lambda, \phi_{i}) , where the eigenvalues are ordered according to their magnitudes:
0 \leq \lambda_ {1} \leq \lambda_ {2} \cdot \cdot \cdot \leq \lambda_ {n - 1} \leq \lambda_ {n} \tag {10.2}
The solution for p eigenpairs can be written
\mathbf {K} \boldsymbol {\Phi} = \boldsymbol {\Phi} \boldsymbol {\Lambda} \tag {10.3}
where \Phi is an n \times p matrix with its columns equal to the p eigenvectors and \Lambda is a p \times p diagonal matrix listing the corresponding eigenvalues. As an example, (10.3) may represent the solution to the lowest p eigenvalues and corresponding eigenvectors of \mathbf{K} , in which case \Phi = [\phi_1, \ldots, \phi_p] and \Lambda = \text{diag}(\lambda_i), i = 1, \ldots, p . We recall that if \mathbf{K} is positive definite,
\lambda_{i} > 0, i = 1, \ldots, n , and if \mathbf{K} is positive semidefinite, \lambda_{i} \geq 0, i = 1, \ldots, n , where the number of zero eigenvalues is equal to the number of rigid body modes in the system.
The solution of the eigenvalue problem in (10.1) is, for example, sought in the evaluation of an element stiffness matrix or in the calculation of the condition number of a structure stiffness matrix. We discussed in Section 4.3.2 that the representation of the element stiffness matrix in its canonical form (i.e., in the eigenvector basis) is used to evaluate the effectiveness of the element. In this case all eigenvalues and vectors of K must be calculated. On the other hand, to evaluate the condition number of a stiffness matrix, only the smallest and largest eigenvalues are required (see Section 8.2.6).
Before proceeding to the generalized eigenproblems, we should mention that other standard eigenproblems may also need to be solved. For example, we may require the eigenvalues of the mass matrix M, in which case M replaces K in (10.1). Similarly, we may want to solve for the eigenvalues of a conductivity or heat capacity matrix in heat flow analysis (see Section 7.2).
A very frequently considered eigenproblem is the one to be solved in vibration mode superposition analysis (see Section 9.3). In this case we consider the generalized eigenproblem,
\mathbf {K} \phi = \lambda \mathbf {M} \phi \tag {10.4}
where K and M are, respectively, the stiffness matrix and mass matrix of the finite element assemblage. The eigenvalues \lambda_{i} and eigenvectors \phi_{i} are the free vibration frequencies (radians/second) squared, \omega_{i}^{2} , and corresponding mode shape vectors, respectively. The properties of K are as discussed above. The mass matrix may be banded, in which case its half-bandwidth m_{M} is equal to m_{K} , or M may be diagonal with m_{ii} \geq 0 ; i.e., some diagonal elements may possibly be zero. A banded mass matrix, obtained in a consistent mass analysis, is always positive definite, whereas a lumped mass matrix is positive definite only if all diagonal elements are larger than zero. In general, a diagonal mass matrix is positive semidefinite.
In analogy to (10.3), the solution for p eigenvalues and corresponding eigenvectors of (10.4) can be written
\mathbf {K} \boldsymbol {\Phi} = \mathbf {M} \boldsymbol {\Phi} \boldsymbol {\Lambda} \tag {10.5}
where the columns in \Phi are the eigenvectors and \Lambda is a diagonal matrix listing the corresponding eigenvalues.
Of course, the generalized eigenproblem in (10.4) reduces to the standard eigenproblem in (10.1) if M is an identity matrix. In other words, the eigenvalues and eigenvectors in (10.3) can also be thought of as frequencies squared and vibration mode shapes of the system when unit mass is specified at each degree of freedom. Corresponding to the possible eigenvalues in the solution of (10.1), the generalized eigenproblem in (10.4) has eigenvalues \lambda_{i} \geq 0, i = 1, \ldots, n , where the number of zero eigenvalues is again equal to the number of rigid body modes in the system.
Two additional generalized eigenproblems should be mentioned briefly. A second problem is solved in linearized buckling analysis, in which case we consider (see Section 6.8.2)
^ \prime \mathbf {K} \phi = \lambda^ {\prime - \Delta^ {\prime}} \mathbf {K} \phi \tag {10.6}
where t^{-\Delta t}\mathbf{K} and \mathbf{K} are the stiffness matrices corresponding to times (i.e., load levels) t - \Delta t and t , respectively.
A third generalized eigenproblem is encountered in heat transfer analysis, where we consider the equation
\mathbf {K} \boldsymbol {\phi} = \lambda \mathbf {C} \boldsymbol {\phi} \tag {10.7}
where K is the heat conductivity matrix and C is the heat capacity matrix. The eigenvalues and eigenvectors are the thermal eigenvalues and mode shapes, respectively. The solution of (10.7) is required in heat transfer analysis using mode superposition (see Exercise 9.25). The matrices K and C in (10.7) are positive definite or positive semidefinite, so that the eigenvalues of (10.7) are \lambda_{i} \geq 0 , i = 1, \ldots, n .
In this and the next chapter we discuss the solution of the eigenproblems K\phi = \lambda\phi and K\phi = \lambda M\phi in (10.1) and (10.4). These eigenproblems are encountered frequently in practice. However, it should be realized that all algorithms to be presented are also applicable to the solution of other eigenproblems, provided they are of the same form and the matrices satisfy the appropriate conditions of positive definiteness, semidefiniteness, and so on. For example, to solve the problem in (10.7), the mass matrix M simply needs to be replaced by the heat capacity matrix C, and the matrix K is the heat conductivity matrix.
Considering the actual computer solution of the required eigenproblems, we recall that in the introduction to equation solution procedures in static analysis (see Section 8.1), we observed the importance of using effective calculation procedures. This is even more so in eigensystem calculations because the solution of eigenvalues and corresponding eigenvectors requires, in general, much more computer effort than the solution of static equilibrium equations. A particularly important consideration is that the solution algorithms must be stable, which is more difficult to achieve in eigensolutions.
A variety of eigensystem solution methods have been developed and are reported in the literature (see, for example, J. H. Wilkinson [A]). Most of the techniques have been devised for rather general matrices. However, in finite element analysis we are concerned with the solution of the specific eigenproblems summarized above, in which each of the matrices has specific properties such as being banded, positive definite, and so on. The eigensystem solution algorithms should take advantage of these properties in order to make a more economical solution possible.
The objective in this chapter is to lay the foundation for a thorough understanding of effective eigensolution methods. This is accomplished by first discussing the properties of the matrices, eigenvalues, and eigenvectors of the problems of interest and then presenting some approximate solution techniques. The actual solution methods recommended for use are presented in Chapter 11.
10.2 FUNDAMENTAL FACTS USED IN THE SOLUTION OF EIGENSYSTEMS
Before the working of any eigensystem solution procedure can be properly studied, it is necessary first to thoroughly understand the different properties of the matrices and eigenvalues and eigenvectors considered. In particular, we will find that all solution methods are, in essence, based on these fundamental properties. We therefore want to summarize in this section the important properties of the matrices and their eigensystems, although some of
the material has already been presented in other sections of the book. As pointed out in Section 10.1, we consider the eigenproblem K\phi = \lambda M\phi , which reduces to K\phi = \lambda\phi when M = I, but the observations made are equally applicable to other eigenproblems of interest.
10.2.1 Properties of the Eigenvectors
It was stated that the solution of the generalized eigenproblem \mathbf{K}\boldsymbol{\phi} = \lambda \mathbf{M}\boldsymbol{\phi} yields n eigenvalues \lambda_1, \ldots, \lambda_n , ordered as shown in (10.2), and corresponding eigenvectors \boldsymbol{\phi}_1, \ldots, \boldsymbol{\phi}_n . Each eigenpair (\lambda_i, \boldsymbol{\phi}_i) satisfies (10.4); i.e.; we have
\mathbf {K} \boldsymbol {\phi} _ {i} = \lambda_ {i} \mathbf {M} \boldsymbol {\phi} _ {i}; \quad i = 1, \dots , n \tag {10.8}
The significance of (10.8) should be well understood. The equation says that if we establish a vector \lambda_i\mathbf{M}\phi_i and use it as a load vector \mathbf{R} in the equation \mathbf{KU} = \mathbf{R} , then \mathbf{U} = \boldsymbol{\phi}_i . This thought may immediately suggest the use of static solution algorithms for the calculation of an eigenvector. We will see later that the \mathbf{LDL}^T decomposition algorithm is indeed an important part of eigensolution procedures.
The equation in (10.8) also shows that an eigenvector is defined only within a multiple of itself; i.e., we also have
\mathbf {K} (\alpha \boldsymbol {\phi} _ {i}) = \lambda_ {i} \mathbf {M} (\alpha \boldsymbol {\phi} _ {i}) \tag {10.9}
where \alpha is a nonzero constant. Therefore, with \phi_{i} being an eigenvector, \alpha\phi_{i} is also an eigenvector, and we say that an eigenvector is defined only by its direction in the n-dimensional space considered. However, in our discussion we refer to the eigenvectors \phi_{i} as satisfying (10.8) and also the relation \phi_{i}^{T}M\phi_{i}=1 , which fixes the lengths of the eigenvectors, i.e., the absolute magnitude of the elements in each eigenvector. However, we may note that the eigenvectors are still defined only within a multiplier of -1.
An important relation which the eigenvectors satisfy is that of M-orthonormality; i.e., we have
\boldsymbol {\phi} _ {i} ^ {T} \mathbf {M} \boldsymbol {\phi} _ {j} = \delta_ {i j} \tag {10.10}
where \delta_{ij} is the Kronecker delta. This relation follows from the orthonormality of the eigenvectors of standard eigenproblems (see Section 2.5) and is discussed further in Section 10.2.5. Premultiplying (10.8) by \phi_{j} transposed and using the condition in (10.10), we obtain
\boldsymbol {\phi} _ {i} ^ {T} \mathbf {K} \boldsymbol {\phi} _ {j} = \lambda_ {i} \delta_ {i j} \tag {10.11}
meaning that the eigenvectors are also K-orthogonal. When using the relations in (10.10) and (10.11) it should be kept in mind that the M- and K-orthogonality follow from (10.8) and that (10.8) is the basic equation to be satisfied. In other words, if we believe that we have an eigenvector and eigenvalue, then as a check we should substitute them into (10.8) (see Example 10.3).
So far we have made no mention of multiple eigenvalues and corresponding eigenvectors. It is important to realize that in this case the eigenvectors are not unique but that we can always choose a set of M-orthonormal eigenvectors which span the subspace that corresponds to a multiple eigenvalue (see Section 2.5). In other words, assume that \lambda_{i} has multiplicity m (i.e., \lambda_{i} = \lambda_{i+1} = \cdots = \lambda_{i+m-1} ); then we can choose m eigenvectors
\phi_{i},\ldots,\phi_{i+m-1} that span the m-dimensional subspace corresponding to the eigenvalues of magnitude \lambda_{i} and which satisfy the orthogonality relation in (10.10) and (10.11). However, the eigenvectors are not unique; instead, the eigenspace corresponding to \lambda_{i} is unique. We demonstrate the results by means of some examples.
EXAMPLE 10.1: The stiffness matrix and mass matrix of a two degree of freedom system are
\mathbf {K} = \left[ \begin{array}{c c} 5 & - 2 \\ - 2 & 2 \end{array} \right]; \quad \mathbf {M} = \left[ \begin{array}{c c} \frac {5}{4} & 0 \\ 0 & \frac {1}{3} \end{array} \right]
It is believed that the two eigenpairs of the problem \mathbf{K}\phi = \lambda \phi are
(d _ {1}, \mathbf {v} _ {1}) = (1, \left[ \begin{array}{l} \frac {1}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} \end{array} \right]); \quad (d _ {2}, \mathbf {v} _ {2}) = (6, \left[ \begin{array}{l} \frac {2}{\sqrt {5}} \\ - \frac {1}{\sqrt {5}} \end{array} \right]) \tag {a}
and the two eigenpairs of the problem \mathbf{K}\phi = \lambda \mathbf{M}\phi are
\left(g _ {1}, \mathbf {w} _ {1}\right) = (2, \left[ \begin{array}{l} \frac {4}{5} \\ 1 \end{array} \right]); \quad \left(g _ {2}, \mathbf {w} _ {2}\right) = (1 2, \left[ \begin{array}{l} \frac {2}{5} \\ - 2 \end{array} \right]) \tag {b}
Verify that we indeed have in (a) and (b) the eigensolutions of the problems \mathbf{K}\boldsymbol{\phi} = \lambda \boldsymbol{\phi} and \mathbf{K}\boldsymbol{\phi} = \lambda \mathbf{M}\boldsymbol{\phi} , respectively.
Consider first the problem K\phi = \lambda\phi . The values given in (a) are indeed the eigensolution if they satisfy the relation in (10.8) with M = I and, to fix the lengths of the vectors, the orthonormality relations in (10.10) with M = I. Substituting into (10.3), which expresses the relation in (10.8) for all eigenpairs, we have
\left[ \begin{array}{c c} 5 & - 2 \\ - 2 & 2 \end{array} \right] \left[ \begin{array}{c c} \frac {1}{\sqrt {5}} & \frac {2}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} & - \frac {1}{\sqrt {5}} \end{array} \right] = \left[ \begin{array}{c c} \frac {1}{\sqrt {5}} & \frac {2}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} & - \frac {1}{\sqrt {5}} \end{array} \right] \left[ \begin{array}{l l} 1 & 0 \\ 0 & 6 \end{array} \right]
\text { or } \quad \left[ \begin{array}{c c} \frac {1}{\sqrt {5}} & \frac {1 2}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} & - \frac {6}{\sqrt {5}} \end{array} \right] = \left[ \begin{array}{c c} \frac {1}{\sqrt {5}} & \frac {1 2}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} & - \frac {6}{\sqrt {5}} \end{array} \right]
Evaluating (10.10), we obtain
\mathbf {v} _ {1} ^ {T} \mathbf {v} _ {1} = \left[ \begin{array}{l l} \frac {1}{\sqrt {5}} & \frac {2}{\sqrt {5}} \end{array} \right] \left[ \begin{array}{l} \frac {1}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} \end{array} \right] = 1
\mathbf {v} _ {2} ^ {T} \mathbf {v} _ {1} = \left[ \frac {2}{\sqrt {5}} - \frac {1}{\sqrt {5}} \right] \left[ \begin{array}{l} \frac {1}{\sqrt {5}} \\ \frac {2}{\sqrt {5}} \end{array} \right] = 0
\mathbf {v} _ {2} ^ {T} \mathbf {v} _ {2} = \left[ \frac {2}{\sqrt {5}} - \frac {1}{\sqrt {5}} \right] \left[ \begin{array}{c} \frac {2}{\sqrt {5}} \\ - \frac {1}{\sqrt {5}} \end{array} \right] = 1
\mathbf {v} _ {1} ^ {T} \mathbf {v} _ {2} = \mathbf {v} _ {2} ^ {T} \mathbf {v} _ {1} = 0
Therefore, the relations in (10.3) and (10.10) are satisfied, and we have \lambda_1 = d_1 , \lambda_2 = d_2 , \phi_1 = \mathbf{v}_1 , and \phi_2 = \mathbf{v}_2 .
To check whether we have in (b) the eigensolution to K\phi = \lambda M\phi , we proceed in an analogous way. Substituting into (10.5), we have
\left[ \begin{array}{c c} 5 & - 2 \\ - 2 & 2 \end{array} \right] \left[ \begin{array}{c c} \frac {4}{5} & \frac {2}{5} \\ 1 & - 2 \end{array} \right] = \left[ \begin{array}{c c} \frac {5}{4} & 0 \\ 0 & \frac {1}{5} \end{array} \right] \left[ \begin{array}{c c} \frac {4}{5} & \frac {2}{5} \\ 1 & - 2 \end{array} \right] \left[ \begin{array}{c c} 2 & 0 \\ 0 & 1 2 \end{array} \right]
or \left[ \begin{array}{cc}2 & 6\\ \frac{2}{5} & -\frac{24}{5} \end{array} \right] = \left[ \begin{array}{cc}2 & 6\\ \frac{2}{5} & -\frac{24}{5} \end{array} \right]
and evaluating (10.10), we obtain
\mathbf {w} _ {1} ^ {T} \mathbf {M} \mathbf {w} _ {1} = \left[ \begin{array}{l l} 4 & 1 \\ \frac {4}{5} & 1 \end{array} \right] \left[ \begin{array}{l l} \frac {5}{4} & 0 \\ 0 & \frac {1}{5} \end{array} \right] \left[ \begin{array}{l} \frac {4}{5} \\ 1 \end{array} \right] = 1
\mathbf {w} _ {2} ^ {T} \mathbf {M} \mathbf {w} _ {1} = \left[ \begin{array}{l l} \frac {2}{5} & - 2 \end{array} \right] \left[ \begin{array}{l l} \frac {5}{4} & 0 \\ 0 & \frac {1}{5} \end{array} \right] \left[ \begin{array}{l} \frac {4}{5} \\ 1 \end{array} \right] = 0
\mathbf {w} _ {2} ^ {T} \mathbf {M} \mathbf {w} _ {2} = \left[ \begin{array}{l l} 2 & - 2 \end{array} \right] \left[ \begin{array}{l l} \frac {5}{4} & 0 \\ 0 & \frac {1}{5} \end{array} \right] \left[ \begin{array}{l} \frac {2}{5} \\ - 2 \end{array} \right] = 1
\mathbf {w} _ {1} ^ {T} \mathbf {M} \mathbf {w} _ {2} = \mathbf {w} _ {2} ^ {T} \mathbf {M} \mathbf {w} _ {1} = 0
Hence, the relations in (10.5) and (10.10) are satisfied and we have
\lambda_ {1} = g _ {1}; \quad \lambda_ {2} = g _ {2}; \quad \phi_ {1} = \mathbf {w} _ {1}; \quad \phi_ {2} = \mathbf {w} _ {2}
EXAMPLE 10.2: Consider the eigenproblem
\mathbf {K} \boldsymbol {\phi} = \lambda \boldsymbol {\phi} \quad \text { with } \mathbf {K} = \left[ \begin{array}{c c c} 2 & & \\ & 2 & \\ & & 3 \end{array} \right]
and show that the eigenvectors corresponding to the multiple eigenvalue are not unique.
The eigenvalues of \mathbf{K} are \lambda_1 = 2 , \lambda_2 = 2 , and \lambda_3 = 3 , and a set of eigenvectors is
\boldsymbol {\phi} _ {1} = \left[ \begin{array}{l} 1 \\ 0 \\ 0 \end{array} \right]; \quad \boldsymbol {\phi} _ {2} = \left[ \begin{array}{l} 0 \\ 1 \\ 0 \end{array} \right]; \quad \boldsymbol {\phi} _ {3} = \left[ \begin{array}{l} 0 \\ 0 \\ 1 \end{array} \right] \tag {a}
where \phi_3 is unique. These values could be checked as in Example 10.1. However, any linear combinations of \phi_1 and \phi_2 given in (a) that satisfy the orthonormality conditions in (10.10) with