- Email: [email protected]
An inverse eigenvalue problem for the finite element model of a vibrating rod
Accepted Manuscript An inverse eigenvalue problem for the finite element model of a vibrating rod Ying Wei, Hua Dai PII: DOI: Reference:
S0377-0427(16)00004-2 http://dx.doi.org/10.1016/j.cam.2015.12.038 CAM 10436
To appear in:
Journal of Computational and Applied Mathematics
Received date: 3 December 2013 Revised date: 10 November 2015 Please cite this article as: Y. Wei, H. Dai, An inverse eigenvalue problem for the finite element model of a vibrating rod, Journal of Computational and Applied Mathematics (2016), http://dx.doi.org/10.1016/j.cam.2015.12.038 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Manuscript Click here to view linked References
An inverse eigenvalue problem for the finite element model of a vibrating rod Ying Wei, Hua Dai∗ Department of Mathematics, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China
Abstract An inverse eigenvalue problem for the finite element model of a longitudinally vibrating rod whose one end is fixed and the other end is supported on a spring is considered. It is known that the mass and stiffness matrices are both tridiagonal for the finite element model of the rod based on linear shape functions. It is shown that the cross section areas can be determined from the spectrum of the rod. The inverse vibration problem can be recast into an inverse eigenvalue problem of a special Jacobi matrix. The necessary and sufficient conditions for the construction of a physically realizable rod with positive cross section areas are established. A numerical method is presented and an illustrative example is given. Keywords: Eigenvalue, Inverse problem, Rod, Finite element model 2000 MSC: 15A18, 65F15, 65F18
1. Introduction The process of formulating, analyzing and deriving the spectral information and, hence, inducing the dynamic behavior of a system from its priori known physical parameters, is referred to as a direct vibration problem. The inverse vibration problem, in contrast, is to validate, determine, or estimate the physical parameters of the system according to its observed or expected dynamical behavior. The inverse vibration problem is just as important as the direct vibration problem in applications. There are different kinds of inverse vibration problems depending on the type of system, the model of system and prescribed spectral data. Inverse vibration problems are studied extensively, see, for example, [1]. The discrete inverse problems in vibration may be transformed into inverse eigenvalue problems for structured matrices. See [2, 3, 4] for an overall treatment of inverse eigenvalue problems, [5, 6, 7, 8, 9, 10] for an exhaustive ∗ Corresponding
author Email addresses: [email protected] (Ying Wei), [email protected] (Hua Dai)
Preprint submitted to Elsevier
November 10, 2015
classification of, and computational procedures for the inverse eigenvalue problems. A system has a structure, a connectivity pattern, and this will be mirrored in the structures, the patterns of zero and non-zero entries in the corresponding matrix. For the problem to be more significant, it is often necessary to confine the reconstruction to certain special classes of matrices. One the other hand, lots of specially structured matrices enjoy many interesting properties. These properties not only can make the algorithms for some direct problems more efficient, for example, [11, 12], but also can help to solve inverse eigenvalue problems. The inverse eigenvalue problem for the finite element model of a vibrating rod considered in this paper is a special case of the inverse eigenvalue problems of Jacobi matrices. The free longitudinal vibrations of a thin rod of length L, cross-sectional area A(x), Young’s modulus E and density ρ are governed by the following equation d du(x) EA(x) + λρA(x)u(x) = 0, 0 ≤ x ≤ L, λ = ω 2 , (1) dx dx where ω is the natural frequency of the rod. Assume that the left end of the rod is fixed and the right end is supported on a spring having the stiffness kn+1 . In this case, the end conditions are du(x) = kn+1 u(L), (2) u(0) = 0, EA(L) dx x=L
and are said to be the fixed-elastically restrained end conditions. The fixedfree and the fixed-fixed end conditions correspond to kn+1 = 0 and kn+1 = ∞, respectively. Inverse eigenvalue problems of the continuous model of the rod are associated with inverse Sturm-Liouville problems, which have been addressed by numerous authors, see, for example, [14, 15, 16, 17, 18, 19]. Levitan [20] and Gladwell [1] expounded the results associated with the inverse Sturm-Liouville problems. Ram [21] considered an inverse mode problem for the continuous model of an axially vibrating rod, and showed that the density and axial rigidity functions are determined by two eigenvalues, their corresponding eigenfunctions and the total mass of the rod. Recently, Gao [22], Morassi [23] presented numerical methods for recovering approximately the cross-sectional area of the vibrating rods having prescribed values of the first n eigenvalues. Inverse vibration problems for the finite difference model of the vibrating rod are formulated as inverse eigenvalue problems for Jacobi matrices, which have been studied extensively, see, for example, [24, 25, 26, 27, 28]. Gladwell and Gbadeyan [29] showed that the mass and stiffness matrices of the finite difference model of the vibrating rod may be reconstructed uniquely from two sets of eigenvalues for the fixed-free and fixed-fixed boundary conditions. Gladwell [30] considered the inverse mode problem for the finite difference model of the vibrating rod, and showed that the discrete system may be constructed uniquely, apart from a scale factor, from two eigenvalues and corresponding 2
eigenvectors. However, frequently, the physical properties, Young’s modulus E and density ρ, of a homogeneous rod are constants and the cross-sectional area A(x) varies. Ram and Elishakoff [31] showed that the cross-sectional area of the finite difference model of an axially vibrating non-uniform rod can be reconstructed from one eigenpair and the total mass of the rod. Ram and Elhay [32] considered the problem of reconstructing the cross-sectional area of the finite difference model of an axially vibrating rod from one set of eigenvalues for the fixed-free boundary condition, and proposed an iterative algorithm for solving the problem. Lu et al [33, 34] recasted the problem into an inverse eigenvalue problem of a specially structured Jacobi matrix, and gave a sufficient and some necessary conditions for the inverse eigenvalue problem to have a solution, and developed a numerical method for this problem. Discretizing (1) and (2) by the finite element method, based on linear shape functions, we get Kn u = λMn u, (3) where
and
nE Kn = L
with
ρL Mn = 6n
A1 + A2 −A2
−A2
A2 + A3 .. .
..
.
..
.
..
.
2(A1 + A2 )
A2
A2
2(A2 + A3 ) .. .
..
.
An−1 + An −An
..
.
..
.
..
.
−An An + 3An+1
..
,
.
2(An−1 + An ) An
An 2An
,
kn+1 L . (4) 3nE Note that the finite element model, based on linear shape functions, of equation (1) subject to the fixed-fixed end condition is given by An+1 =
(Kn−1 − λMn−1 )u = 0, where Kn−1 and Mn−1 are the (n − 1) × (n − 1) leading principal submatrices of Kn and Mn , respectively. Ram and Gladwell [35] considered an inverse mode problem for the finite element model of the vibrating rod with the fixed-free end condition, and proposed a procedure to reconstruct a finite element model from a single eigenvalue, two eigenvectors and the total mass of the rod. Tian 3
and Dai [36] showed that the mass and stiffness matrices of the finite element model of the vibrating rod under the fixed-free and the fixed-fixed end conditions may be reconstructed from two eigenpairs and the total mass of the rod, and established the necessary and sufficient conditions for the construction of a physical realizable system with positive stiffness and masses. Lai and Ananthasuresh [37], Jim´ enez et al. [38] showed that the cross-sectional area of the finite element model of the vibrating rod can be reconstructed uniquely, apart from a scale factor, from one eigenpair. Jim´ enez et al. [39] studied an inverse eigenvalue problem of the finite element model of a vibrating rod under the fixed-free end condition, and showed that the cross-sectional areas of the rod can be reconstructed uniquely, apart from a scale factor, from half of the eigenvalues of the matrix pencil (Kn , Mn ) and half of the eigenvalues of the matrix pencil (Kn−1 , Mn−1 ). To the best of the authors’ knowledge, how to determine the cross-sectional areas of the finite element model of a vibrating rod from one set of eigenvalues for the fixed-elastically restrained end conditions has not been considered in the literature earlier. The problem can be stated as follows. Problem 1. Given the positive numbers L, E, ρ, kn+1 , and a set of distinct positive real numbers 0 < λ1 < λ2 < · · · < λn , find a set of positive values {Ai }ni=1 such that {λi }ni=1 are the eigenvalues of the matrix pencil (Kn , Mn ). The main contribution of this paper is as follows. This paper recasts Problem 1 into an equivalent inverse eigenvalue problem of a specially structured Jacobi matrix, and presents the necessary and sufficient conditions for the problem to have a solution, and provides a numerical algorithm for solving this problem. Throughout this paper, we use I to denote the identity matrix of size implied by context, ei to denote the ith column of the identity matrix, xT to denote the transpose of a vector x, det(A) and tr(A) to denote the determinant and the trace of a matrix A, respectively. 2. Inverse eigenvalue problem of the specially structured Jacobi matrix In this section, Problem 1 is formulated as an inverse eigenvalue problem of a specially structured Jacobi matrix. Let Dn = diag(A1 + A2 , A2 + A3 , · · · , An−1 + An , An + An+1 ), and
0 A2 Bn =
A2 0 .. .
A3 .. . An−1
..
. 0 An
An −2An+1 4
6n2 E − 2λρL2 . , µ = λρL2 + 6n2 E
Then
nE (Dn − Bn ) , L Substituting (5) into (3), we have Kn =
Mn =
ρL (2Dn + Bn ) . 6n
(5)
Bn u = µDn u.
(6)
Let
− 21
− 12
Jn = D n B n D n
with
and
0
b1 =
b1 0 .. .
..
.
..
.
..
.
..
.
0 bn−1
Ai+1 bi = p (Ai + Ai+1 )(Ai+1 + Ai+2 )
bn−1 a
,
1
x = Dn2 u
(i = 1, 2, · · · , n − 1),
p p 1 Dn2 = diag( A1 + A2 , · · · , An + An+1 ), a = −
(7)
2An+1 . An + An+1
Then the generalized eigenvalue problem (6) is reduced to the following standard eigenvalue problem for the Jacobi matrix Jn Jn x = µx. Let
2 2 ˜i = 6n E − 2λi ρL λ λi ρL2 + 6n2 E
(i = 1, 2, · · · , n).
(8)
˜ i }n are the eigenvalues of Jn if and only if {λi }n are It is obvious that {λ i=1 i=1 the eigenvalues of the matrix pencil (Kn , Mn ). Note that 6n2 E − 2λρL2 f (λ) = λρL2 + 6n2 E is a monotone decreasing function on the interval [0, +∞), and f (0) = 1,
lim f (λ) = −2.
λ→+∞
Then 0 < λ1 < λ2 < · · · < λn if and only if
˜1 > λ ˜2 > · · · > λ ˜n > −2. 1>λ
Let cn−1 = 1 −
2b2n−1 , 2+a
ci = 1 −
b2i ci+1 5
(i = n − 2, n − 3, · · · , 1).
From (7), we get ci =
Ai Ai + Ai+1
(i = 1, 2, · · · , n − 1).
Thus Ai > 0(i = 1, 2, · · · , n) if and only if −2 < a < 0,
and 0 < ci < 1 (i = 1, 2, · · · , n − 1).
Problem 1 can be recasted into the equivalent inverse eigenvalue problem of the specially structured Jacobi matrix Jn as follows. ˜ i }n satisfying Problem J1. Given a set of real numbers {λ i=1 ˜1 > λ ˜2 > · · · > λ ˜n > −2, 1>λ find a specially structured Jacobi matrix 0 b1 b1 0 . . . .. .. Jn = . . . ..
..
.
0 bn−1
bn−1 a
,
˜ i }n are the eigenvalues of Jn , −2 < a < 0 and such that {λ i=1 0 < ci < 1 (i = 1, 2, · · · , n − 1), where cn−1 = 1 −
2b2n−1 , 2+a
ci = 1 −
b2i ci+1
(i = n − 2, n − 3, · · · , 1).
3. The solvability of Problem 1 In this section, the necessary and sufficient conditions for the solvability of Problem J1 and Problem 1 are derived. First, we introduce some lemmas for later discussion. Lemma 1. [26] Given two sets of distinct real values ˜1 > λ ˜2 > · · · > λ ˜n, λ
λ1 > λ2 > · · · > λn ,
there exists an n × n Jacobi matrices Jn and a positive number c such that Jn ˜ i }n if and has eigenvalues {λi }ni=1 , and J˜n = Jn + cen eTn has eigenvalues {λ i=1 only if ˜1 > λ1 > λ ˜2 > · · · > λn−1 > λ ˜n > λn . λ 6
Lemma 2. [39] Let Jn be an n × n Jacobi matrix, Jn−1 be the (n − 1) × (n − 1) leading principal submatrix of Jn , φn (λ) = det(λI − Jn )
and
φn−1 (λ) = det(λI − Jn−1 ).
Then all the diagonal elements of Jn are zeros if and only if (−1)n φn (−λ) = φn (λ)
and
(−1)n−1 φn−1 (−λ) = φn−1 (λ).
Lemma 3. Given a set of distinct real values λ1 > λ2 > · · · > λn , there exists a specially structured 0 b1 b1 0 .. Jn = .
Jacobi matrix ..
.
..
.
..
.
..
.
0 bn−1
bn−1 a
whose eigenvalues are {λi }ni=1 if and only if
,
a<0
(9)
−λn > λ1 > −λn−1 > λ2 > −λn−2 > · · · > −λ1 > λn . Proof. Necessity. Let Sn = diag (−1, 1, · · · , (−1)n ). Then, Sn = SnT ,
Sn2 = I,
−SnT Jn Sn = J˜n ,
where J˜n = Jn − 2aen eTn . Thus, {−λi }ni=1 are the eigenvalues of J˜n if and only if {λi }ni=1 are the eigenvalues of Jn . By Lemma 1, we get −λn > λ1 > −λn−1 > λ2 > −λn−2 > · · · > −λ1 > λn . Sufficiency. Assume that −λn > λ1 > −λn−1 > λ2 > −λn−2 > · · · > −λ1 > λn . It follows from Lemma 1 that there exists an n × n Jacobi matrix Jn and a positive number c such that Jn has eigenvalues {λi }ni=1 , and J˜n = Jn + cen eTn has eigenvalues {−λi }ni=1 . n P λi = −tr(J˜n ). Then c = tr(J˜n ) − tr(Jn ) = −2a, and Let a = tr(Jn ) = i=1
Let
J˜n − Jn = −2aen eTn .
(10)
Jˆn = Jn − aen eTn .
(11)
7
From (11) and (10), we have Jˆn = J˜n + aen eTn .
(12)
Let Jn−1 be the (n − 1) × (n − 1) leading principal submatrix of Jn , φˆn (λ) = det(λI − Jˆn ),
φn−1 (λ) = det(λI − Jn−1 ), φn (λ) = det(λI −Jn ) =
n Y
φ˜n (λ) = det(λI − J˜n ) =
(λ−λi ),
i=1
From (11) and (12) , we have φˆn (λ) = φn (λ) + aφn−1 (λ),
(λ+λi ). (13)
i=1
and φˆn (λ) = φ˜n (λ) − aφn−1 (λ).
From (14), we obtain φn−1 (λ) = and
n Y
φ˜n (λ) − φn (λ) , 2a
(14)
(15)
φ˜n (λ) + φn (λ) . φˆn (λ) = 2
(16)
It follows from (13) that (−1)n φn (−λ) = φ˜n (λ),
(−1)n φ˜n (−λ) = φn (λ).
(17)
Substituting (17) into (15) and (16), respectively, we get and (−1)n φˆn (−λ) = φˆn (λ).
(−1)n−1 φn−1 (−λ) = φn−1 (λ),
It follows from Lemma 2 that all the diagonal elements of Jˆn are zeros. Hence, Jn is of the form described in (9). Lemma 4. Let
0
b1
b1 ∗ Jn =
0 .. .
..
.
..
.
..
.
0 bn−1
bn−1 −a 2
,
∗ Jk+1,n be the submatrix of Jn∗ which is obtained by deleting the first k columns and the first k rows of Jn∗ , ∗ φ∗k+1,n (λ) = det(λI − Jk+1,n ),
cn−1 = 1 −
2b2n−1 , 2+a
ci = 1 −
b2i ci+1 8
k = 0, 1, · · · , n − 1, ,
i = n − 2, n − 3, · · · , 1.
T ∗ n ∗ ˜ n Assume that Jn = Jn∗ + 3a 2 en en , {λi }i=1 and {λi }i=1 are the eigenvalues of Jn and Jn ,respectively, and
˜1 > λ ˜2 > · · · > λ ˜n > −2. 1>λ Then φ∗i,n (1) , i = 1, 2, · · · , n − 1; φ∗i+1,n (1) 2) 0 < ci < 1, i = 1, 2, · · · , n − 1 if and only if 1) ci =
λ∗i < 1,
i = 1, 2, · · · , n;
n Q
n n Q ˜ i ) + 3 Q (λ + λ ˜ i ); (λ − λ∗i ) = 14 (λ − λ 4 i=1 i=1 i=1 4) 0 < ci < 1, i = 1, 2, · · · , n − 1 if and only if ! n Y ˜i 1+λ 1 >− . ˜ 3 1 − λi i=1
3)
Proof. 1) Obviously, cn−1 = 1 −
2 ∗ 1+ a 2b2n−1 2 − bn−1 = φn−1,n (1) . = 2+a φ∗n,n (1) 1+ a 2
Assume that ci =
φ∗i,n (1) . φ∗i+1,n (1)
Expanding the determinant φ∗i−1,n (1) = det(I − Ji−1,n ), we have φ∗i−1,n (1) = φ∗i,n (1) − b2i−1 φ∗i+1,n (1). Thus,
φ∗i−1,n (1) φ∗i+1,n (1) b2 = 1 − b2i−1 ∗ = 1 − i−1 = ci−1 . ∗ φi,n (1) φi,n (1) ci
By mathematic induction, we have ci =
φ∗i,n (1) , φ∗i+1,n (1)
i = 1, 2, · · · , n − 1.
2) Necessity. From ci > 0, we have
i = 1, 2, · · · , n − 1,
φ∗i,n (1) > 0, φ∗i+1,n (1)
i = 1, 2, · · · , n − 1. 9
(18)
By Lemma 3, we have ˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n . −λ Then
Thus,
n−1 1 X ˜ ˜ n > −2. ˜ ˜ λi = λn + λi + λn−i > λ a= 2 i=1 i=1 n X
φ∗n,n (1) = 1 +
a > 0. 2
(19)
(20)
From (20) and (18), we have φ∗i,n (1) > 0,
i = 1, 2, · · · , n.
Therefore, I − Jn∗ is a positive matrix. Hence Sufficiency. Since
λ∗i < 1,
i = 1, 2, · · · , n.
λ∗i < 1,
i = 1, 2, · · · , n,
I − Jn∗ is a positive matrix. Therefore, φ∗i,n (1) > 0, Consequently, ci =
i = 1, 2, · · · , n.
φ∗i,n (1) > 0, φ∗i+1,n (1)
i = 1, 2, · · · , n − 1.
(21)
From (19), we get a > −2. Therefore cn−1 = 1 −
2b2n−1 < 1. 2+a
Then, for i = n − 2, n − 3, · · · , 1, ci = 1 −
b2i ci+1
< 1,
if and only if
ci+1 > 0.
(22)
Putting together (21) and (22), we have 0 < ci < 1,
i = 1, 2, · · · , n − 1.
˜ i }n are the 3) Let J˜n = Jn − 2aen eTn . From the proof of Lemma 3, {−λ i=1 ˜ eigenvalues of Jn . Let φn (λ) = det(λI − Jn ),
φ˜n (λ) = det(λI − J˜n ),
10
φ∗n (λ) = det(λI − J˜n∗ ).
Note that 1 3 ˜ 4 (λI − Jn )en + 4 (λI − Jn )en 1 (0, · · · , 0, −b , λ − a)T + 3 (0, · · · , 0, −b , λ + a)T n n 4 4 T = 0, · · · , 0, −bn , λ + a 2
=
= (λI − Jn∗ )en . Then
3 1 φn (λ) + φ˜n (λ). 4 4
φ∗n (λ) = That is
n Y
n
n
Y 1Y ˜i) + 3 ˜ i ). (λ − λ (λ + λ 4 i=1 4 i=1
(λ − λ∗i ) =
i=1
4) By Lemma 1, we have
˜ 1 > λ∗ > λ ˜2 > · · · > λ∗ > λ ˜n. 1>λ 2 n Therefore,
n Y
i=2
Then
(1 − λ∗i ) > 0.
λ∗i < 1,
if and only if 0<
n Y
i=1
i = 1, 2, · · · , n n
(1 − λ∗i ) =
which is equivalent to
n Y
i=1
n
Y 1Y ˜i) + 3 ˜i ), (1 − λ (1 + λ 4 i=1 4 i=1 ˜i 1+λ ˜ 1 − λi
!
1 >− . 3
From above discussion, we get the necessary and sufficient conditions of Problem J1 and Problem 1 as follows. Theorem 1. Problem J1 has a solution if and only if ˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n , −λ and
n Y
i=1
˜i 1+λ ˜i 1−λ 11
!
1 >− . 3
Theorem 2. Let 2 2 ˜i = 6n E − 2λi ρL , λ λi ρL2 + 6n2 E
i = 1, 2, · · · , n,
Then, Problem 1 has a solution if and only if ˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n , −λ and
n Y
i=1
˜i 1+λ ˜ 1 − λi
!
1 >− . 3
4. Numerical method and examples Let a=
n X
˜i, λ
An+1 =
i=1
and
2 2 ˜i = 6n E − 2λi ρL , λ 2 2 λi ρL + 6n E
kn+1 L , 3nE
i = 1, 2, · · · , n.
(23)
From (15), we have ˜ ˜ ˜ ˜ ˜ ˜ i ) = φn (λi ) − φn (λi ) = φn (λi ) = φn−1 (λ 2a 2a
n Q ˜i + λ ˜j ) (λ
j=1
2a
,
i = 1, 2, · · · , n.
˜ i is Thus, the last elements pi of eigenvector of Jn corresponding to λ v u Q n u s ˜ +λ ˜j ) (λ u ˜ i ) u j=1 i φn−1 (λ , i = 1, 2, · · · , n. =u pi = n u Q ˜i) φ′n (λ ˜i − λ ˜j ) t 2a (λ j=1 j6=i
There are many effective numerical methods to construct the Jacobi matrix Jn ˜ i }n are the eigenvalues of Jn and {pi }n are the last elements such that {λ i=1 i=1 ˜ i }n , such as Lanczos method, of the eigenvectors of Jn corresponding to {λ i=1 Householder transformations, and Givens transformations. For the details, see [4]. From (7), we have −An+1 (2 + a) , (24) An = a and A2i+1 Ai = 2 − Ai+1 , i = n − 1, n − 2, · · · , 1. (25) bi (Ai+1 + Ai+2 ) To summarize, the algorithm for solving Problem 1 is presented as follows. Algorithm 1. 12
Step 1. Compute a =
n P ˜i, λ
i=1
kn+1 L , and An+1 = 3nE
2 2 ˜ i = 6n E − 2λi ρL , λ λi ρL2 + 6n2 E
If
i = 1, 2, · · · , n.
˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n , −λ ! n Y ˜i 1+λ 1 >− , ˜ 3 1 − λi i=1
goto step 2, else end the algorithm. Step 2. Compute v u Q n u ˜i + λ ˜j ) (λ u u j=1 pi = u , n u Q ˜i − λ ˜j ) t 2a (λ
i = 1, 2, · · · , n.
j=1 j6=i
˜ i }n are the eigenvalues Step 3. Construct the Jacobi matrix Jn such that {λ i=1 n of Jn and {pi }i=1 are the last elements of the eigenvectors of Jn corre˜ i }n . sponding to {λ i=1 −An+1 (2 + a) and Step 4. Compute An = a Ai =
A2i+1 − Ai+1 , b2i (Ai+1 + Ai+2 )
i = n − 1, n − 2, · · · , 1.
We give a numerical example to show the performance of Algorithm 1. The computations are carried out in Matlab 7.1 in double precision arithmetic. Step 3 is completed by Householder transformations method. First, the cross-section areas {Ai }ni=1 are given. We compute the eigenvalues {λi }ni=1 of (Kn , Mn ) by the eig function in Matlab. Then, we reconstruct Ai from the eigenvalues {λi }ni=1 . Let {A˜i }ni=1 be the computed approximation of {Ai }ni=1 , A˜ − A i i ǫ1 = max |A˜i − Ai |, and ǫ2 = max . i i Ai Example 1. Given Young’s modulus of E = 2 × 1011 , the density ρ = 7800, the length L = 1, and x = i ∗ L/n,
i = 1, 2, · · · , n.
Table 1 is the result of the stiffness of the spring kn+1 = 5 × 1010 with different n.
13
Table 1. n 50 100 300 500
A(x) = x ǫ1 5.3474 × 10−12 3.7800 × 10−11 1.1989 × 10−9 7.0087 × 10−9
A(x) = x2
ǫ2 9.9803 × 10−11 8.6622 × 10−10 3.1545 × 10−8 2.1783 × 10−7
ǫ1 3.5290 × 10−12 1.1020 × 10−11 4.4111 × 10−10 2.4215 × 10−9
ǫ2 4.1077 × 10−10 1.0596 × 10−8 1.0033 × 10−6 9.6775 × 10−6
Table 1 shows that when kn+1 = 5 × 1010 , Algorithm 1 is numerical stable. Table 2 is the result of n = 100 with different kn . Table 2. kn 5 × 100 5 × 105 5 × 1010 5 × 1015 5 × 1025
A(x) = x ǫ1 2.2514 × 10−1 9.7477 × 10−7 3.7800 × 10−11 6.2021 × 10−9 6.6790 × 10−4
A(x) = x2
ǫ2 3.7001 × 10−1 1.9893 × 10−6 8.6622 × 10−10 8.1976 × 10−8 8.7655 × 10−3
ǫ1 1.9243 × 10−1 1.8172 × 10−6 1.1020 × 10−11 5.6123 × 10−9 5.7830 × 10−4
ǫ2 3.4876 × 10−1 3.8911 × 10−6 1.0596 × 10−8 2.8685 × 10−6 2.7534 × 10−1
Note that when kn+1 = 0, λi + λn+1−i = 0 and when kn+1 = ∞, λi + λn−i = 0. The accuracy of Algorithm 1 decrease when kn+1 is very large or very small. This may be because λi + λn+1−i or λi + λn−i is very small and then Step 3 is unstable. Acknowledgments The research is supported by the National Natural Science Foundation of China under grant No.11071118. The authors are grateful to the referee for his valuable comments. References [1] G.M.L. Gladwell, Inverse Problems in Vibration, 2nd ed., Kluwer Academic Publishers, Dordrecht, 2004. [2] S.Q. Zhou, H. Dai, The Algebraic Inverse Eigenvalue Problem, Henan Science and Technology Press, Zhengzhou, 1991 (in Chinese). [3] S.F. Xu, An Introduction to Inverse Algebraic Eigenvalue Problems, Peking University Press, Beijing, 1998. [4] M.T. Chu, G.H. Golub, Inverse Eigenvalue Problems, Theory, Algorithms and Applications, Oxford University Press, Oxford, 2005.
14
[5] D.L. Boley, G.H. Golub, A survey of matrix inverse eigenvalue problems, Inverse Problems 3 (1987) 595-622. [6] M.T. Chu, Inverse eigenvalue problems, SIAM Review 40 (1998) 1-39. [7] M.T. Chu, G.H. Golub, Structured inverse eigenvalue problems, Acta Numerica 12 (2002) 1-71. [8] T.J. Laffey, Inverse eigenvalue problems for matrices, Proceedings of the Royal Irish Academy 95A (1995) S81-S88. [9] C.R. Johnson, T.J. Laffey, R. Loewy, The real and the symmetric nonegative inverse eigenvalue problems are different, Proc. American Math. Soc. 124 (1996) 3647-3651. ˇ [10] T.J. Laffey, H. Smigoc, Construction of nonnegative symmetric matrices with given spectrum, Linear Algebra Appl. 421 (2007) 97-109. [11] J.S. Respondek, Numerical recipes for the high efficient inverse of the confluent Vandermonde matrices, Applied Mathematics and Computation, 218 (2011), 2044-2054. [12] J.S. Respondek, Recursive numerical recipes for the high efficient inversion of the confluent Vandermonde matrices, Applied Mathematics and Computation, 225 (2013), 718-730. [13] G. Borg, Eine Umkehrung der Sturm-Liouvilleschen Eigenwertaufgabe, Acta Math. 78 1946) 1-96. [14] M.G. Krein, Determination of the density of nonhomogeneous symmetric cord from its frequency spectrum, Proc. USSR Acad. Sci. 76 (1951) 345-348. [15] M.G. Krein, On method of efficient solution of an inverse boundary problem, Proc. USSR Acad. Sci. 95 (1954) 767-770. [16] I.M. Gel’fand, B.M. Levitan, On the determination of a differential equation from its spectral function, Izv. Akad. Nauk SSSR Ser. Mat. 15 (1951) 309360. [17] B.M. Levitan, On the determination of a SturmCLiouville equation by two spectra. Izv. Akad. Nauk SSSR Ser. Mat. 28 (1964) 63-78. [18] H. Hochstadt, The inverse Sturm-Liouville problem, Communications on Pure and Applied Mathematics 26 (1973) 715-729. [19] O.H. Hald, The inverse Sturm-Liouville problem and the Rayleigh-Ritz method, Mathematics of Computation 32 (1978) 687-705. [20] B.M. Levitan, Inverse Sturm-Liouville Problems, WSP, Zeist, 1987. [21] Y.M. Ram, An inverse mode problem for the continuous model of an axially vibrating rod, Journal of Applied Mechanics 61 (1994) 624-628. 15
[22] Q. Gao, Decent flow methods for inverse SturmCLiouville problem, Applied Mathematical Modelling 36 (2012) 4452-4465. [23] A. Morassi, Constructing rods with given natural frequencies, Mechanical Systems and Signal Processing 40 (2013) 288-300. [24] H. Hochstadt, On the construction of a Jacobi matrix from spectral data, Linear Algebra Appl. 8 (1974) 435-446. [25] O.H. Hald, Inverse eigenvalue problems for Jacobi matrices, Linear Algebra Appl. 14 (1976) 63-85. [26] C.de Boor, G.H. Golub, The numerically stable reconstruction of a Jacobi matrix from spectral data, Linear Algebra Appl. 21 (1978) 245-260. [27] W.B. Gragg, W.J. Harrod, The numerically stable reconstruction of Jacobi matrices from spectral data, Numer. Math. 44 (1984) 317-335. [28] H. Dai, The numerically stable reconstruction of Jacobi matrices from spectral data, Chinese Journal on Numerical Mathematics and Applications 12 (1990) 58-66. [29] G.M.L. Gladwell, J.A. Gbadeyan, On the inverse problem of the vibrating string or rod, Quarterly Journal of Mechanics and Applied Mathematics 38 (1985) 169-175. [30] G.M.L. Gladwell, The inverse mode problem for lumped-mass systems, Quarterly Journal of Mechanical Applied Mathematics 39 (1986) 297-307. [31] Y.M. Ram, I. Elishakoff, Reconstructing the cross-sectional area of an axially vibrating non-uniform rod from one of its mode shapes, Proc. R. Soc. Lond. A 460 (2004) 1583-1596. [32] Y. Ram, S. Elhay, Constructing the shape of a rod from elgenvalues, Comun. Numer. Methods Engng. 14 (1998) 597-608. [33] L.Z. Lu, W.W. Sun, On necessary conditions for reconstruction of a specially structured Jacobx matrix from elgenvalues, Inverse Problems 15 (1999) 977-987. [34] L.Z. Lu, M.K. Ng, On sufficient and necessary conditions for the Jacobi matrix inverse eigenvalue problem, Numer. Math. 98 (2004) 167-176. [35] Y.M. Ram, G.M.L. Gladwell, Constructing a finite element model of a vibratory rod from eigendata, Journal of Sound and Vibration 169 (1994) 229-237. [36] X. Tian, H. Dai, Inverse mode problems for the finite element model of a vibrating rod, Applied Mathematics and Computation 214 (2009) 479-486.
16
[37] E. Lai, G.K. Ananthasuresh, On the design of bars and beams for desired mode shapes, Journal of Sound and Vibration 254 (2002) 393-406. [38] R. Jim´ enez, L. Santos, N. Kuhl, J. Ega˜ na, R. Soto, An inverse eigenvalue procedure for damage detection in rods, Computers and Mathematics with Applications 47 (2004) 643-657. [39] R. Jim´ enez, L. Santos, N. Kuhl, J. Ega˜ na, The reconstruction of a specially structured Jacobi matrix with an application to damage detection in rods, Computers and Mathematics with Apphcations, 49 (2005) 1815-1823.
17
S0377-0427(16)00004-2 http://dx.doi.org/10.1016/j.cam.2015.12.038 CAM 10436
To appear in:
Journal of Computational and Applied Mathematics
Received date: 3 December 2013 Revised date: 10 November 2015 Please cite this article as: Y. Wei, H. Dai, An inverse eigenvalue problem for the finite element model of a vibrating rod, Journal of Computational and Applied Mathematics (2016), http://dx.doi.org/10.1016/j.cam.2015.12.038 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Manuscript Click here to view linked References
An inverse eigenvalue problem for the finite element model of a vibrating rod Ying Wei, Hua Dai∗ Department of Mathematics, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China
Abstract An inverse eigenvalue problem for the finite element model of a longitudinally vibrating rod whose one end is fixed and the other end is supported on a spring is considered. It is known that the mass and stiffness matrices are both tridiagonal for the finite element model of the rod based on linear shape functions. It is shown that the cross section areas can be determined from the spectrum of the rod. The inverse vibration problem can be recast into an inverse eigenvalue problem of a special Jacobi matrix. The necessary and sufficient conditions for the construction of a physically realizable rod with positive cross section areas are established. A numerical method is presented and an illustrative example is given. Keywords: Eigenvalue, Inverse problem, Rod, Finite element model 2000 MSC: 15A18, 65F15, 65F18
1. Introduction The process of formulating, analyzing and deriving the spectral information and, hence, inducing the dynamic behavior of a system from its priori known physical parameters, is referred to as a direct vibration problem. The inverse vibration problem, in contrast, is to validate, determine, or estimate the physical parameters of the system according to its observed or expected dynamical behavior. The inverse vibration problem is just as important as the direct vibration problem in applications. There are different kinds of inverse vibration problems depending on the type of system, the model of system and prescribed spectral data. Inverse vibration problems are studied extensively, see, for example, [1]. The discrete inverse problems in vibration may be transformed into inverse eigenvalue problems for structured matrices. See [2, 3, 4] for an overall treatment of inverse eigenvalue problems, [5, 6, 7, 8, 9, 10] for an exhaustive ∗ Corresponding
author Email addresses: [email protected] (Ying Wei), [email protected] (Hua Dai)
Preprint submitted to Elsevier
November 10, 2015
classification of, and computational procedures for the inverse eigenvalue problems. A system has a structure, a connectivity pattern, and this will be mirrored in the structures, the patterns of zero and non-zero entries in the corresponding matrix. For the problem to be more significant, it is often necessary to confine the reconstruction to certain special classes of matrices. One the other hand, lots of specially structured matrices enjoy many interesting properties. These properties not only can make the algorithms for some direct problems more efficient, for example, [11, 12], but also can help to solve inverse eigenvalue problems. The inverse eigenvalue problem for the finite element model of a vibrating rod considered in this paper is a special case of the inverse eigenvalue problems of Jacobi matrices. The free longitudinal vibrations of a thin rod of length L, cross-sectional area A(x), Young’s modulus E and density ρ are governed by the following equation d du(x) EA(x) + λρA(x)u(x) = 0, 0 ≤ x ≤ L, λ = ω 2 , (1) dx dx where ω is the natural frequency of the rod. Assume that the left end of the rod is fixed and the right end is supported on a spring having the stiffness kn+1 . In this case, the end conditions are du(x) = kn+1 u(L), (2) u(0) = 0, EA(L) dx x=L
and are said to be the fixed-elastically restrained end conditions. The fixedfree and the fixed-fixed end conditions correspond to kn+1 = 0 and kn+1 = ∞, respectively. Inverse eigenvalue problems of the continuous model of the rod are associated with inverse Sturm-Liouville problems, which have been addressed by numerous authors, see, for example, [14, 15, 16, 17, 18, 19]. Levitan [20] and Gladwell [1] expounded the results associated with the inverse Sturm-Liouville problems. Ram [21] considered an inverse mode problem for the continuous model of an axially vibrating rod, and showed that the density and axial rigidity functions are determined by two eigenvalues, their corresponding eigenfunctions and the total mass of the rod. Recently, Gao [22], Morassi [23] presented numerical methods for recovering approximately the cross-sectional area of the vibrating rods having prescribed values of the first n eigenvalues. Inverse vibration problems for the finite difference model of the vibrating rod are formulated as inverse eigenvalue problems for Jacobi matrices, which have been studied extensively, see, for example, [24, 25, 26, 27, 28]. Gladwell and Gbadeyan [29] showed that the mass and stiffness matrices of the finite difference model of the vibrating rod may be reconstructed uniquely from two sets of eigenvalues for the fixed-free and fixed-fixed boundary conditions. Gladwell [30] considered the inverse mode problem for the finite difference model of the vibrating rod, and showed that the discrete system may be constructed uniquely, apart from a scale factor, from two eigenvalues and corresponding 2
eigenvectors. However, frequently, the physical properties, Young’s modulus E and density ρ, of a homogeneous rod are constants and the cross-sectional area A(x) varies. Ram and Elishakoff [31] showed that the cross-sectional area of the finite difference model of an axially vibrating non-uniform rod can be reconstructed from one eigenpair and the total mass of the rod. Ram and Elhay [32] considered the problem of reconstructing the cross-sectional area of the finite difference model of an axially vibrating rod from one set of eigenvalues for the fixed-free boundary condition, and proposed an iterative algorithm for solving the problem. Lu et al [33, 34] recasted the problem into an inverse eigenvalue problem of a specially structured Jacobi matrix, and gave a sufficient and some necessary conditions for the inverse eigenvalue problem to have a solution, and developed a numerical method for this problem. Discretizing (1) and (2) by the finite element method, based on linear shape functions, we get Kn u = λMn u, (3) where
and
nE Kn = L
with
ρL Mn = 6n
A1 + A2 −A2
−A2
A2 + A3 .. .
..
.
..
.
..
.
2(A1 + A2 )
A2
A2
2(A2 + A3 ) .. .
..
.
An−1 + An −An
..
.
..
.
..
.
−An An + 3An+1
..
,
.
2(An−1 + An ) An
An 2An
,
kn+1 L . (4) 3nE Note that the finite element model, based on linear shape functions, of equation (1) subject to the fixed-fixed end condition is given by An+1 =
(Kn−1 − λMn−1 )u = 0, where Kn−1 and Mn−1 are the (n − 1) × (n − 1) leading principal submatrices of Kn and Mn , respectively. Ram and Gladwell [35] considered an inverse mode problem for the finite element model of the vibrating rod with the fixed-free end condition, and proposed a procedure to reconstruct a finite element model from a single eigenvalue, two eigenvectors and the total mass of the rod. Tian 3
and Dai [36] showed that the mass and stiffness matrices of the finite element model of the vibrating rod under the fixed-free and the fixed-fixed end conditions may be reconstructed from two eigenpairs and the total mass of the rod, and established the necessary and sufficient conditions for the construction of a physical realizable system with positive stiffness and masses. Lai and Ananthasuresh [37], Jim´ enez et al. [38] showed that the cross-sectional area of the finite element model of the vibrating rod can be reconstructed uniquely, apart from a scale factor, from one eigenpair. Jim´ enez et al. [39] studied an inverse eigenvalue problem of the finite element model of a vibrating rod under the fixed-free end condition, and showed that the cross-sectional areas of the rod can be reconstructed uniquely, apart from a scale factor, from half of the eigenvalues of the matrix pencil (Kn , Mn ) and half of the eigenvalues of the matrix pencil (Kn−1 , Mn−1 ). To the best of the authors’ knowledge, how to determine the cross-sectional areas of the finite element model of a vibrating rod from one set of eigenvalues for the fixed-elastically restrained end conditions has not been considered in the literature earlier. The problem can be stated as follows. Problem 1. Given the positive numbers L, E, ρ, kn+1 , and a set of distinct positive real numbers 0 < λ1 < λ2 < · · · < λn , find a set of positive values {Ai }ni=1 such that {λi }ni=1 are the eigenvalues of the matrix pencil (Kn , Mn ). The main contribution of this paper is as follows. This paper recasts Problem 1 into an equivalent inverse eigenvalue problem of a specially structured Jacobi matrix, and presents the necessary and sufficient conditions for the problem to have a solution, and provides a numerical algorithm for solving this problem. Throughout this paper, we use I to denote the identity matrix of size implied by context, ei to denote the ith column of the identity matrix, xT to denote the transpose of a vector x, det(A) and tr(A) to denote the determinant and the trace of a matrix A, respectively. 2. Inverse eigenvalue problem of the specially structured Jacobi matrix In this section, Problem 1 is formulated as an inverse eigenvalue problem of a specially structured Jacobi matrix. Let Dn = diag(A1 + A2 , A2 + A3 , · · · , An−1 + An , An + An+1 ), and
0 A2 Bn =
A2 0 .. .
A3 .. . An−1
..
. 0 An
An −2An+1 4
6n2 E − 2λρL2 . , µ = λρL2 + 6n2 E
Then
nE (Dn − Bn ) , L Substituting (5) into (3), we have Kn =
Mn =
ρL (2Dn + Bn ) . 6n
(5)
Bn u = µDn u.
(6)
Let
− 21
− 12
Jn = D n B n D n
with
and
0
b1 =
b1 0 .. .
..
.
..
.
..
.
..
.
0 bn−1
Ai+1 bi = p (Ai + Ai+1 )(Ai+1 + Ai+2 )
bn−1 a
,
1
x = Dn2 u
(i = 1, 2, · · · , n − 1),
p p 1 Dn2 = diag( A1 + A2 , · · · , An + An+1 ), a = −
(7)
2An+1 . An + An+1
Then the generalized eigenvalue problem (6) is reduced to the following standard eigenvalue problem for the Jacobi matrix Jn Jn x = µx. Let
2 2 ˜i = 6n E − 2λi ρL λ λi ρL2 + 6n2 E
(i = 1, 2, · · · , n).
(8)
˜ i }n are the eigenvalues of Jn if and only if {λi }n are It is obvious that {λ i=1 i=1 the eigenvalues of the matrix pencil (Kn , Mn ). Note that 6n2 E − 2λρL2 f (λ) = λρL2 + 6n2 E is a monotone decreasing function on the interval [0, +∞), and f (0) = 1,
lim f (λ) = −2.
λ→+∞
Then 0 < λ1 < λ2 < · · · < λn if and only if
˜1 > λ ˜2 > · · · > λ ˜n > −2. 1>λ
Let cn−1 = 1 −
2b2n−1 , 2+a
ci = 1 −
b2i ci+1 5
(i = n − 2, n − 3, · · · , 1).
From (7), we get ci =
Ai Ai + Ai+1
(i = 1, 2, · · · , n − 1).
Thus Ai > 0(i = 1, 2, · · · , n) if and only if −2 < a < 0,
and 0 < ci < 1 (i = 1, 2, · · · , n − 1).
Problem 1 can be recasted into the equivalent inverse eigenvalue problem of the specially structured Jacobi matrix Jn as follows. ˜ i }n satisfying Problem J1. Given a set of real numbers {λ i=1 ˜1 > λ ˜2 > · · · > λ ˜n > −2, 1>λ find a specially structured Jacobi matrix 0 b1 b1 0 . . . .. .. Jn = . . . ..
..
.
0 bn−1
bn−1 a
,
˜ i }n are the eigenvalues of Jn , −2 < a < 0 and such that {λ i=1 0 < ci < 1 (i = 1, 2, · · · , n − 1), where cn−1 = 1 −
2b2n−1 , 2+a
ci = 1 −
b2i ci+1
(i = n − 2, n − 3, · · · , 1).
3. The solvability of Problem 1 In this section, the necessary and sufficient conditions for the solvability of Problem J1 and Problem 1 are derived. First, we introduce some lemmas for later discussion. Lemma 1. [26] Given two sets of distinct real values ˜1 > λ ˜2 > · · · > λ ˜n, λ
λ1 > λ2 > · · · > λn ,
there exists an n × n Jacobi matrices Jn and a positive number c such that Jn ˜ i }n if and has eigenvalues {λi }ni=1 , and J˜n = Jn + cen eTn has eigenvalues {λ i=1 only if ˜1 > λ1 > λ ˜2 > · · · > λn−1 > λ ˜n > λn . λ 6
Lemma 2. [39] Let Jn be an n × n Jacobi matrix, Jn−1 be the (n − 1) × (n − 1) leading principal submatrix of Jn , φn (λ) = det(λI − Jn )
and
φn−1 (λ) = det(λI − Jn−1 ).
Then all the diagonal elements of Jn are zeros if and only if (−1)n φn (−λ) = φn (λ)
and
(−1)n−1 φn−1 (−λ) = φn−1 (λ).
Lemma 3. Given a set of distinct real values λ1 > λ2 > · · · > λn , there exists a specially structured 0 b1 b1 0 .. Jn = .
Jacobi matrix ..
.
..
.
..
.
..
.
0 bn−1
bn−1 a
whose eigenvalues are {λi }ni=1 if and only if
,
a<0
(9)
−λn > λ1 > −λn−1 > λ2 > −λn−2 > · · · > −λ1 > λn . Proof. Necessity. Let Sn = diag (−1, 1, · · · , (−1)n ). Then, Sn = SnT ,
Sn2 = I,
−SnT Jn Sn = J˜n ,
where J˜n = Jn − 2aen eTn . Thus, {−λi }ni=1 are the eigenvalues of J˜n if and only if {λi }ni=1 are the eigenvalues of Jn . By Lemma 1, we get −λn > λ1 > −λn−1 > λ2 > −λn−2 > · · · > −λ1 > λn . Sufficiency. Assume that −λn > λ1 > −λn−1 > λ2 > −λn−2 > · · · > −λ1 > λn . It follows from Lemma 1 that there exists an n × n Jacobi matrix Jn and a positive number c such that Jn has eigenvalues {λi }ni=1 , and J˜n = Jn + cen eTn has eigenvalues {−λi }ni=1 . n P λi = −tr(J˜n ). Then c = tr(J˜n ) − tr(Jn ) = −2a, and Let a = tr(Jn ) = i=1
Let
J˜n − Jn = −2aen eTn .
(10)
Jˆn = Jn − aen eTn .
(11)
7
From (11) and (10), we have Jˆn = J˜n + aen eTn .
(12)
Let Jn−1 be the (n − 1) × (n − 1) leading principal submatrix of Jn , φˆn (λ) = det(λI − Jˆn ),
φn−1 (λ) = det(λI − Jn−1 ), φn (λ) = det(λI −Jn ) =
n Y
φ˜n (λ) = det(λI − J˜n ) =
(λ−λi ),
i=1
From (11) and (12) , we have φˆn (λ) = φn (λ) + aφn−1 (λ),
(λ+λi ). (13)
i=1
and φˆn (λ) = φ˜n (λ) − aφn−1 (λ).
From (14), we obtain φn−1 (λ) = and
n Y
φ˜n (λ) − φn (λ) , 2a
(14)
(15)
φ˜n (λ) + φn (λ) . φˆn (λ) = 2
(16)
It follows from (13) that (−1)n φn (−λ) = φ˜n (λ),
(−1)n φ˜n (−λ) = φn (λ).
(17)
Substituting (17) into (15) and (16), respectively, we get and (−1)n φˆn (−λ) = φˆn (λ).
(−1)n−1 φn−1 (−λ) = φn−1 (λ),
It follows from Lemma 2 that all the diagonal elements of Jˆn are zeros. Hence, Jn is of the form described in (9). Lemma 4. Let
0
b1
b1 ∗ Jn =
0 .. .
..
.
..
.
..
.
0 bn−1
bn−1 −a 2
,
∗ Jk+1,n be the submatrix of Jn∗ which is obtained by deleting the first k columns and the first k rows of Jn∗ , ∗ φ∗k+1,n (λ) = det(λI − Jk+1,n ),
cn−1 = 1 −
2b2n−1 , 2+a
ci = 1 −
b2i ci+1 8
k = 0, 1, · · · , n − 1, ,
i = n − 2, n − 3, · · · , 1.
T ∗ n ∗ ˜ n Assume that Jn = Jn∗ + 3a 2 en en , {λi }i=1 and {λi }i=1 are the eigenvalues of Jn and Jn ,respectively, and
˜1 > λ ˜2 > · · · > λ ˜n > −2. 1>λ Then φ∗i,n (1) , i = 1, 2, · · · , n − 1; φ∗i+1,n (1) 2) 0 < ci < 1, i = 1, 2, · · · , n − 1 if and only if 1) ci =
λ∗i < 1,
i = 1, 2, · · · , n;
n Q
n n Q ˜ i ) + 3 Q (λ + λ ˜ i ); (λ − λ∗i ) = 14 (λ − λ 4 i=1 i=1 i=1 4) 0 < ci < 1, i = 1, 2, · · · , n − 1 if and only if ! n Y ˜i 1+λ 1 >− . ˜ 3 1 − λi i=1
3)
Proof. 1) Obviously, cn−1 = 1 −
2 ∗ 1+ a 2b2n−1 2 − bn−1 = φn−1,n (1) . = 2+a φ∗n,n (1) 1+ a 2
Assume that ci =
φ∗i,n (1) . φ∗i+1,n (1)
Expanding the determinant φ∗i−1,n (1) = det(I − Ji−1,n ), we have φ∗i−1,n (1) = φ∗i,n (1) − b2i−1 φ∗i+1,n (1). Thus,
φ∗i−1,n (1) φ∗i+1,n (1) b2 = 1 − b2i−1 ∗ = 1 − i−1 = ci−1 . ∗ φi,n (1) φi,n (1) ci
By mathematic induction, we have ci =
φ∗i,n (1) , φ∗i+1,n (1)
i = 1, 2, · · · , n − 1.
2) Necessity. From ci > 0, we have
i = 1, 2, · · · , n − 1,
φ∗i,n (1) > 0, φ∗i+1,n (1)
i = 1, 2, · · · , n − 1. 9
(18)
By Lemma 3, we have ˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n . −λ Then
Thus,
n−1 1 X ˜ ˜ n > −2. ˜ ˜ λi = λn + λi + λn−i > λ a= 2 i=1 i=1 n X
φ∗n,n (1) = 1 +
a > 0. 2
(19)
(20)
From (20) and (18), we have φ∗i,n (1) > 0,
i = 1, 2, · · · , n.
Therefore, I − Jn∗ is a positive matrix. Hence Sufficiency. Since
λ∗i < 1,
i = 1, 2, · · · , n.
λ∗i < 1,
i = 1, 2, · · · , n,
I − Jn∗ is a positive matrix. Therefore, φ∗i,n (1) > 0, Consequently, ci =
i = 1, 2, · · · , n.
φ∗i,n (1) > 0, φ∗i+1,n (1)
i = 1, 2, · · · , n − 1.
(21)
From (19), we get a > −2. Therefore cn−1 = 1 −
2b2n−1 < 1. 2+a
Then, for i = n − 2, n − 3, · · · , 1, ci = 1 −
b2i ci+1
< 1,
if and only if
ci+1 > 0.
(22)
Putting together (21) and (22), we have 0 < ci < 1,
i = 1, 2, · · · , n − 1.
˜ i }n are the 3) Let J˜n = Jn − 2aen eTn . From the proof of Lemma 3, {−λ i=1 ˜ eigenvalues of Jn . Let φn (λ) = det(λI − Jn ),
φ˜n (λ) = det(λI − J˜n ),
10
φ∗n (λ) = det(λI − J˜n∗ ).
Note that 1 3 ˜ 4 (λI − Jn )en + 4 (λI − Jn )en 1 (0, · · · , 0, −b , λ − a)T + 3 (0, · · · , 0, −b , λ + a)T n n 4 4 T = 0, · · · , 0, −bn , λ + a 2
=
= (λI − Jn∗ )en . Then
3 1 φn (λ) + φ˜n (λ). 4 4
φ∗n (λ) = That is
n Y
n
n
Y 1Y ˜i) + 3 ˜ i ). (λ − λ (λ + λ 4 i=1 4 i=1
(λ − λ∗i ) =
i=1
4) By Lemma 1, we have
˜ 1 > λ∗ > λ ˜2 > · · · > λ∗ > λ ˜n. 1>λ 2 n Therefore,
n Y
i=2
Then
(1 − λ∗i ) > 0.
λ∗i < 1,
if and only if 0<
n Y
i=1
i = 1, 2, · · · , n n
(1 − λ∗i ) =
which is equivalent to
n Y
i=1
n
Y 1Y ˜i) + 3 ˜i ), (1 − λ (1 + λ 4 i=1 4 i=1 ˜i 1+λ ˜ 1 − λi
!
1 >− . 3
From above discussion, we get the necessary and sufficient conditions of Problem J1 and Problem 1 as follows. Theorem 1. Problem J1 has a solution if and only if ˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n , −λ and
n Y
i=1
˜i 1+λ ˜i 1−λ 11
!
1 >− . 3
Theorem 2. Let 2 2 ˜i = 6n E − 2λi ρL , λ λi ρL2 + 6n2 E
i = 1, 2, · · · , n,
Then, Problem 1 has a solution if and only if ˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n , −λ and
n Y
i=1
˜i 1+λ ˜ 1 − λi
!
1 >− . 3
4. Numerical method and examples Let a=
n X
˜i, λ
An+1 =
i=1
and
2 2 ˜i = 6n E − 2λi ρL , λ 2 2 λi ρL + 6n E
kn+1 L , 3nE
i = 1, 2, · · · , n.
(23)
From (15), we have ˜ ˜ ˜ ˜ ˜ ˜ i ) = φn (λi ) − φn (λi ) = φn (λi ) = φn−1 (λ 2a 2a
n Q ˜i + λ ˜j ) (λ
j=1
2a
,
i = 1, 2, · · · , n.
˜ i is Thus, the last elements pi of eigenvector of Jn corresponding to λ v u Q n u s ˜ +λ ˜j ) (λ u ˜ i ) u j=1 i φn−1 (λ , i = 1, 2, · · · , n. =u pi = n u Q ˜i) φ′n (λ ˜i − λ ˜j ) t 2a (λ j=1 j6=i
There are many effective numerical methods to construct the Jacobi matrix Jn ˜ i }n are the eigenvalues of Jn and {pi }n are the last elements such that {λ i=1 i=1 ˜ i }n , such as Lanczos method, of the eigenvectors of Jn corresponding to {λ i=1 Householder transformations, and Givens transformations. For the details, see [4]. From (7), we have −An+1 (2 + a) , (24) An = a and A2i+1 Ai = 2 − Ai+1 , i = n − 1, n − 2, · · · , 1. (25) bi (Ai+1 + Ai+2 ) To summarize, the algorithm for solving Problem 1 is presented as follows. Algorithm 1. 12
Step 1. Compute a =
n P ˜i, λ
i=1
kn+1 L , and An+1 = 3nE
2 2 ˜ i = 6n E − 2λi ρL , λ λi ρL2 + 6n2 E
If
i = 1, 2, · · · , n.
˜n > λ ˜ 1 > −λ ˜ n−1 > λ ˜ 2 > −λ ˜ n−2 > · · · > −λ ˜1 > λ ˜n , −λ ! n Y ˜i 1+λ 1 >− , ˜ 3 1 − λi i=1
goto step 2, else end the algorithm. Step 2. Compute v u Q n u ˜i + λ ˜j ) (λ u u j=1 pi = u , n u Q ˜i − λ ˜j ) t 2a (λ
i = 1, 2, · · · , n.
j=1 j6=i
˜ i }n are the eigenvalues Step 3. Construct the Jacobi matrix Jn such that {λ i=1 n of Jn and {pi }i=1 are the last elements of the eigenvectors of Jn corre˜ i }n . sponding to {λ i=1 −An+1 (2 + a) and Step 4. Compute An = a Ai =
A2i+1 − Ai+1 , b2i (Ai+1 + Ai+2 )
i = n − 1, n − 2, · · · , 1.
We give a numerical example to show the performance of Algorithm 1. The computations are carried out in Matlab 7.1 in double precision arithmetic. Step 3 is completed by Householder transformations method. First, the cross-section areas {Ai }ni=1 are given. We compute the eigenvalues {λi }ni=1 of (Kn , Mn ) by the eig function in Matlab. Then, we reconstruct Ai from the eigenvalues {λi }ni=1 . Let {A˜i }ni=1 be the computed approximation of {Ai }ni=1 , A˜ − A i i ǫ1 = max |A˜i − Ai |, and ǫ2 = max . i i Ai Example 1. Given Young’s modulus of E = 2 × 1011 , the density ρ = 7800, the length L = 1, and x = i ∗ L/n,
i = 1, 2, · · · , n.
Table 1 is the result of the stiffness of the spring kn+1 = 5 × 1010 with different n.
13
Table 1. n 50 100 300 500
A(x) = x ǫ1 5.3474 × 10−12 3.7800 × 10−11 1.1989 × 10−9 7.0087 × 10−9
A(x) = x2
ǫ2 9.9803 × 10−11 8.6622 × 10−10 3.1545 × 10−8 2.1783 × 10−7
ǫ1 3.5290 × 10−12 1.1020 × 10−11 4.4111 × 10−10 2.4215 × 10−9
ǫ2 4.1077 × 10−10 1.0596 × 10−8 1.0033 × 10−6 9.6775 × 10−6
Table 1 shows that when kn+1 = 5 × 1010 , Algorithm 1 is numerical stable. Table 2 is the result of n = 100 with different kn . Table 2. kn 5 × 100 5 × 105 5 × 1010 5 × 1015 5 × 1025
A(x) = x ǫ1 2.2514 × 10−1 9.7477 × 10−7 3.7800 × 10−11 6.2021 × 10−9 6.6790 × 10−4
A(x) = x2
ǫ2 3.7001 × 10−1 1.9893 × 10−6 8.6622 × 10−10 8.1976 × 10−8 8.7655 × 10−3
ǫ1 1.9243 × 10−1 1.8172 × 10−6 1.1020 × 10−11 5.6123 × 10−9 5.7830 × 10−4
ǫ2 3.4876 × 10−1 3.8911 × 10−6 1.0596 × 10−8 2.8685 × 10−6 2.7534 × 10−1
Note that when kn+1 = 0, λi + λn+1−i = 0 and when kn+1 = ∞, λi + λn−i = 0. The accuracy of Algorithm 1 decrease when kn+1 is very large or very small. This may be because λi + λn+1−i or λi + λn−i is very small and then Step 3 is unstable. Acknowledgments The research is supported by the National Natural Science Foundation of China under grant No.11071118. The authors are grateful to the referee for his valuable comments. References [1] G.M.L. Gladwell, Inverse Problems in Vibration, 2nd ed., Kluwer Academic Publishers, Dordrecht, 2004. [2] S.Q. Zhou, H. Dai, The Algebraic Inverse Eigenvalue Problem, Henan Science and Technology Press, Zhengzhou, 1991 (in Chinese). [3] S.F. Xu, An Introduction to Inverse Algebraic Eigenvalue Problems, Peking University Press, Beijing, 1998. [4] M.T. Chu, G.H. Golub, Inverse Eigenvalue Problems, Theory, Algorithms and Applications, Oxford University Press, Oxford, 2005.
14
[5] D.L. Boley, G.H. Golub, A survey of matrix inverse eigenvalue problems, Inverse Problems 3 (1987) 595-622. [6] M.T. Chu, Inverse eigenvalue problems, SIAM Review 40 (1998) 1-39. [7] M.T. Chu, G.H. Golub, Structured inverse eigenvalue problems, Acta Numerica 12 (2002) 1-71. [8] T.J. Laffey, Inverse eigenvalue problems for matrices, Proceedings of the Royal Irish Academy 95A (1995) S81-S88. [9] C.R. Johnson, T.J. Laffey, R. Loewy, The real and the symmetric nonegative inverse eigenvalue problems are different, Proc. American Math. Soc. 124 (1996) 3647-3651. ˇ [10] T.J. Laffey, H. Smigoc, Construction of nonnegative symmetric matrices with given spectrum, Linear Algebra Appl. 421 (2007) 97-109. [11] J.S. Respondek, Numerical recipes for the high efficient inverse of the confluent Vandermonde matrices, Applied Mathematics and Computation, 218 (2011), 2044-2054. [12] J.S. Respondek, Recursive numerical recipes for the high efficient inversion of the confluent Vandermonde matrices, Applied Mathematics and Computation, 225 (2013), 718-730. [13] G. Borg, Eine Umkehrung der Sturm-Liouvilleschen Eigenwertaufgabe, Acta Math. 78 1946) 1-96. [14] M.G. Krein, Determination of the density of nonhomogeneous symmetric cord from its frequency spectrum, Proc. USSR Acad. Sci. 76 (1951) 345-348. [15] M.G. Krein, On method of efficient solution of an inverse boundary problem, Proc. USSR Acad. Sci. 95 (1954) 767-770. [16] I.M. Gel’fand, B.M. Levitan, On the determination of a differential equation from its spectral function, Izv. Akad. Nauk SSSR Ser. Mat. 15 (1951) 309360. [17] B.M. Levitan, On the determination of a SturmCLiouville equation by two spectra. Izv. Akad. Nauk SSSR Ser. Mat. 28 (1964) 63-78. [18] H. Hochstadt, The inverse Sturm-Liouville problem, Communications on Pure and Applied Mathematics 26 (1973) 715-729. [19] O.H. Hald, The inverse Sturm-Liouville problem and the Rayleigh-Ritz method, Mathematics of Computation 32 (1978) 687-705. [20] B.M. Levitan, Inverse Sturm-Liouville Problems, WSP, Zeist, 1987. [21] Y.M. Ram, An inverse mode problem for the continuous model of an axially vibrating rod, Journal of Applied Mechanics 61 (1994) 624-628. 15
[22] Q. Gao, Decent flow methods for inverse SturmCLiouville problem, Applied Mathematical Modelling 36 (2012) 4452-4465. [23] A. Morassi, Constructing rods with given natural frequencies, Mechanical Systems and Signal Processing 40 (2013) 288-300. [24] H. Hochstadt, On the construction of a Jacobi matrix from spectral data, Linear Algebra Appl. 8 (1974) 435-446. [25] O.H. Hald, Inverse eigenvalue problems for Jacobi matrices, Linear Algebra Appl. 14 (1976) 63-85. [26] C.de Boor, G.H. Golub, The numerically stable reconstruction of a Jacobi matrix from spectral data, Linear Algebra Appl. 21 (1978) 245-260. [27] W.B. Gragg, W.J. Harrod, The numerically stable reconstruction of Jacobi matrices from spectral data, Numer. Math. 44 (1984) 317-335. [28] H. Dai, The numerically stable reconstruction of Jacobi matrices from spectral data, Chinese Journal on Numerical Mathematics and Applications 12 (1990) 58-66. [29] G.M.L. Gladwell, J.A. Gbadeyan, On the inverse problem of the vibrating string or rod, Quarterly Journal of Mechanics and Applied Mathematics 38 (1985) 169-175. [30] G.M.L. Gladwell, The inverse mode problem for lumped-mass systems, Quarterly Journal of Mechanical Applied Mathematics 39 (1986) 297-307. [31] Y.M. Ram, I. Elishakoff, Reconstructing the cross-sectional area of an axially vibrating non-uniform rod from one of its mode shapes, Proc. R. Soc. Lond. A 460 (2004) 1583-1596. [32] Y. Ram, S. Elhay, Constructing the shape of a rod from elgenvalues, Comun. Numer. Methods Engng. 14 (1998) 597-608. [33] L.Z. Lu, W.W. Sun, On necessary conditions for reconstruction of a specially structured Jacobx matrix from elgenvalues, Inverse Problems 15 (1999) 977-987. [34] L.Z. Lu, M.K. Ng, On sufficient and necessary conditions for the Jacobi matrix inverse eigenvalue problem, Numer. Math. 98 (2004) 167-176. [35] Y.M. Ram, G.M.L. Gladwell, Constructing a finite element model of a vibratory rod from eigendata, Journal of Sound and Vibration 169 (1994) 229-237. [36] X. Tian, H. Dai, Inverse mode problems for the finite element model of a vibrating rod, Applied Mathematics and Computation 214 (2009) 479-486.
16
[37] E. Lai, G.K. Ananthasuresh, On the design of bars and beams for desired mode shapes, Journal of Sound and Vibration 254 (2002) 393-406. [38] R. Jim´ enez, L. Santos, N. Kuhl, J. Ega˜ na, R. Soto, An inverse eigenvalue procedure for damage detection in rods, Computers and Mathematics with Applications 47 (2004) 643-657. [39] R. Jim´ enez, L. Santos, N. Kuhl, J. Ega˜ na, The reconstruction of a specially structured Jacobi matrix with an application to damage detection in rods, Computers and Mathematics with Apphcations, 49 (2005) 1815-1823.
17