ON THE NUMERICAL RECONSTRUCTION OF A SPRING-MASS SYSTEM FROM ITS NATURAL FREQUENCIES

Given the sequences of real numbers (λi)1 , (μi) n−1 1 , satisfying the interlacing property λi < μi < λi+1, 1 ≤ i ≤ n− 1, we present a new numerical procedure to construct a spring-mass system with eigenvalues (λi)1 , where the interlaced spectrum (μi) 1 corresponds to the modified system whose mass mr+1, 1 ≤ r ≤ n − 2, is fixed. The method, which is a modification of the fast orthogonal reduction technique, appears to be computationally less expensive than other in the literature.


Introduction
Inverse problems have important applications in physics and engineering, in areas such as circuit, control and vibration theory.The vibration analysis of an structural system is performed as an approximation to the continuous system obtained by discretizing the mechanical properties.For an undamped discrete system with n-degrees of freedom the governing equation is where the matrices M n and K n of order n are the mass matrix and the stiffness matrix, respectively.The eigenvalues λ i of the equation (1.1) are related with the natural frequencies ω i (λ i = ω 2 i ) and the eigenvectors x (i)  n = (x n,2 , . . ., x (i) n,n ) T , i = 1, 2, ...., n, represent the vibration modes.
Consider the spring-mass system consisting of n masses m i > 0 connected by n springs of stiffnesses k i > 0, as shown in figure (1).Here, the mass matrix is M n = diag{m 1 , m 2 , ..., m n } and the stiffness matrix The spring-mass system, denoted by (M n , K n ), is a real physical system if m i > 0 and k i > 0, i = 1, 2, ..., n.These properties imply that M n and K n are positive definite matrices.Moreover, as K n is also symmetric and tridiagonal, the eigenvalues are real, positive and distinct.
The structural properties of the matrices M n and K n allow us to reduce the generalized eigenvalue equation (1.1) to the standard form n = 0, where n }.Thus, the matrix J n is symmetric tridiagonal positive definite with the same eigenvalues of the system (M n , K n ).We call J n a Jacobi matrix.Clearly, the relation (1.2) implies that the diagonal elements of J n are positive, while the co-diagonal elements are negative.By the similarity transformation J = D n J n D n , where D n = diag{(−1) i , i = 1, 2, . . ., n}, the codiagonal elements become positive.Then, without loss of generality we may assume that J n has the form The natural frequencies of the system (M n , K n ) are usually computed by using the matrices M n and K n .The problem here is to identify a system from the natural frequencies.This identification problem may be solved in two ways.The first way start with a pair of matrices (the mass matrix and the stiffness matrix) and then by using an optimization method, determines the matrices M n and K n which identify the specified spectral data.The second strategy start with the spectral information to recover the matrices M n and K n , directly.
The second procedure have been studied by Gladwell [4], [5], [6] and others [10], [11], [12].Krein [8] proved that the matrices M n and K n can be uniquely reconstructed if the following information is given: The eigenvalues (λ i ) n 1 of the original system (Fig 1).The eigenvalues (µ i ) n−1 1 corresponding to the modified system whose last mass is fixed (Fig 2).An additional factor.For example, the total mass of the system m T = n i=1 m i .

The interlacing property
is a necessary and sufficient condition for the existence of a real physical system, where the set (µ i ) n−1 1 is called the interlaced spectrum.The construction of the system is performed in two parts.Firstly, the matrix J n is reconstructed by the adequate use of some standard tridiagonalization method for symmetric matrices, like Householder transformations, Givens rotations, or Lanczos method.Once J n is obtained and considering the relations in (1.2), the second part consists in to separate the masses and the stiffnesses of J n .This is possible if we consider the equation where J n was previously determined.This imply that Moreover, if the total mass of the system is known, from (1.3), we have that the last mass is given by Hence, the mass matrix M n is completely determined.Now, we may use the relations (1.2) to compute the stiffness matrix K n .In this way, we may determine the system (M n , K n ), completely.
The reconstruction of a spring-mass system may also be performed by using the interlaced spectrum corresponding to a modified system with different bounds conditions (see [4], [10]).The modified system considered in this paper consists in to fix any mass of the system (M n , K n ), other than the extreme masses.Thus, the problem we study here may be formulated as follows:

satisfying the interlacing property
It is clear that if we fix the (r + 1) th mass, the resulting modified system consists of two spring-mass systems, (M r , K r ) and (M p , K p ) , with natural frequencies (γ i ) r 1 and (ν i ) p 1 , respectively, where p = n − r − 1.The structural properties of the matrices related to the systems (M r , K r ) and (M p , K p ) allow to partition the matrix J n in the form: where the Jacobi submatrices J r and J p are given by p }.As the system (M n , K n ) can be reconstructed from the matrix J n , it is enough to reconstruct J n from the sets (λ i ) n 1 , (γ i ) r 1 , and (ν i ) p 1 to obtain (M n , K n ).Thus, the Problem 1 is reduced to: The existence and uniqueness of the solution for the Problem 2 have been considered in [1], [3], [7].Algorithms to compute the solution J n are proposed in [6] and [13].For the tridiagonalization process, they use the Lanczos algorithm and Householder transformations, respectively.The computational costs in these methods are in general of order On 3 .
In [2] and [7] is presented the fast orthogonal reduction algorithm, which reconstruct the matrix , where (µ i ) n−1 is the interlaced spectrum corresponding to the submatrix J n\1 .In section 2, we modify the fast orthogonal reduction algorithm to derive a new direct numerical procedure to solve Problem 2. This new method has only a computational cost of order On 2 .In section 3 we present some numerical examples which confirm the efficiency of the method.

The algorithm
The basic method in which our proposal is based is the fast orthogonal reduction method presented by Gragg and Harrow in [7].The idea is firstly to construct a bordered diagonal matrix A n+1 of the form where a 0 is a dummy entry, Λ n = diag{λ 1 , λ 2 , . . ., λ n }, and ) T is the first (or last) row of eigenvectors of J n .Then, we reduce A n+1 to the tridiagonal form by applying a sequence of orthogonal plane rotations in a particular order.
This orthogonal reduction process is called algorithm of Rutishauser, which reduce the matrix (2.1) to the tridiagonal form where e 1 = (1, 0, . . ., 0) T and J n is the required Jacobi matrix with eigenvalues (λ i ) n 1 .The fast orthogonal reduction method has a computational cost of O(n 2 ) operations, which is increased to O(n 3 ) if we use Lanczos or Householder transformations.A complete discussion of the fast orthogonal reduction method can be found in [2], [7].
The idea of our proposal is to construct the submatrices J r , J p , and the elements b r , b r+1 , and a r+1 , where p = n − r − 1, by using an adequate modification of the fast orthogonal reduction algorithm.Clearly, the diagonal element a r+1 can be computed by

Now let us define the bordered diagonal matrices
where u r = (u (1)  r,r , u (2)  r,r , . . ., u (r) r,r ) T and w p = (w p,1 , . . ., w p,1 ) T are respectively, the last and first row of the matrices of eigenvectors of J r and J p .The submatrices Λ r and Λ p are diagonal matrices having as their diagonal entries the eigenvalues of J r and J p , respectively.If the eigenvectors u r and w p are known, then it is possible to apply the fast orthogonal reduction algorithm to orthogonally reduce the matrices A r+1 and A p+1 to their tridiagonal forms, obtaining in this way the searched matrices J r and J p .
Consider the following notation: In general, A j\i will denote the principal submatrix of order (j − 1), obtained from A j , by deleting its i th row and column.P n (λ), Q r (λ), and S p (λ) will denote the characteristic polynomials of J n , J r , and J p , with p = n − r − 1, respectively, while Q r\r (λ) and S p\1 (λ) will denote the characteristic polynomials of J r\r , and J p\1 , respectively.
By expanding det (λI n − J n ) along its r th row we find that Setting λ = γ i and λ = ν i en (2.4) we obtain r is the orthogonal matrix of eigenvectors of J r .Then U T r (λI r − J r )U r = λI r − Λ r , and we have The left side of (2.8) is Comparing the entries in position (r, r) in both sides of (2.8) we find that Taking the limit when λ tends to γ i , we obtain ) and (2.9), we obtain The interlacing property of the eigenvalues (1.5) guarantee that the right side of (2.10) is always positive.Thus, the equations (2.5), (2.10) and (2.11) allow us to determine all the elements of the vector u r and the codiagonal element b r .
If W p = w (1)  p w (2)   p p is the matrix of eigenvectors of J p , the same procedure leads to where S p (ν i ) = r j=1,j =i (ν i − ν j ).From (2.12) and (2.7) we have Hence, from (2.5), (2.13) and (2.14), the first row of the matrix of eigenvectors of J p and the codiagonal element b r+1 are computed.
The algorithm can be summarized as follows: Algorithm 1.Given the sequences of real numbers the algorithm produces the matrices J n , J r , and J p of (1.4) 1. Compute a r+1 by (2.2).

Compute the elements of the vectors u r and w p by (2.10) and
(2.13), respectively.
4. Use the fast orthogonal reduction algorithm to reduce the matrices A r+1 and A p+1 of (2.3) to the tridiagonal form J r+1 and J p+1 , respectively, to obtain the Jacobi matrix of (1.4).
This algorithm compute the elements a r+1 , b r , and b r+1 , with approximately O(n 2 ) operations.In step 4, the algorithm produces the matrices J r and J p with about O(r 2 ) and O(p 2 ) operations, respectively.Thus, to compute the whole matrix J n we need O(n 2 ) operations.Then this procedure is computationally less expensive than the algorithms proposed in [6] and [13], which require about O(n 3 ) operations.Since it updates only a few entries of the matrix in each step (Householder transformations update all the entries in each step), the rounding errors has less influence on the process.

Numerical Results
In this Section we present two numerical examples to show the accuracy of the algorithm given in Section 2. The data used correspond to the eigenvalues of a known Jacobi matrix J n , so that the reconstructed matrix can be compared with J n .To this end, we define e a = , where a, b are the vectors whose components are diagonal and codiagonal elements of J n , and a, b are the obtained after the algorithm 3 is applied.Table 1 shows the numerical results obtained when we recover a Jacobi matrix J n of order n = 6 with diagonal elements a i = 2.0 and codiagonal elements b i = 1.0.The exact eigenvalues of J n , J r , and J p are gives by λ k = 2 − 2 cos kπ n+1 , γ k = 2 − 2 cos kπ r+1 , and ν k = 2 − 2 cos kπ n−p−1 , respectively.The same example is considered in Table 2, for n = 20.The obtained relative errors confirm a better accuracy than the examples presented in ( [13]).