Eigenvalue sensitivity analysis based on the transfer matrix method

For linear mechanical systems, the transfer matrix method is one of the most efficient modeling and analysis methods. However, in contrast to classical modeling strategies, the final eigenvalue problem is based on a matrix which is a highly nonlinear function of the eigenvalues. Therefore, classical strategies for sensitivity analysis of eigenvalues w.r.t. system parameters cannot be applied. The paper develops two specific strategies for this situation, a direct differentiation strategy and an adjoint variable method, where especially the latter is easy to use and applicable to arbitrarily complex chain or branched multibody systems. Like the system analysis itself, it is able to break down the sensitivity analysis of the overall system to analytically determinable derivatives of element transfer matrices and recursive formula which can be applied along the transfer path of the topology figure. Several examples of different complexity validate the proposed approach by comparing results to analytical calculations and numerical differentiation. The obtained procedure may support gradient‐based optimization and robust design by delivering exact sensitivities.


| INTRODUCTION
Linear multibody systems are composed of simple components like bodies, springs, dampers and beams, see for example, Figure 1. The transfer matrix method utilizes the property that the kinematic and kinetic properties of all these elements can be easily and exactly described by analytical transfer matrices, see Rui et al., 1 Abbas et al. 2 or Appendix A. As demonstrated by some examples in Section 6, these element transfer matrices can be combined along the topology graph of a specific multibody system to finally end up with an overall transfer equation This reduction process may be described by a n n 2 × -Boolean matrix Therefore, many eigenvalue sensitivity analysis studies are reported, for example, by Lancaster, 3 Garg,4  Further, damped mechanical vibration systems have been studied, for example, by Zmindak. 8 A major feature of all these classical sensitivity studies is that the involved matrices A or M K , are independent of the eigenvalue, whereas the coefficient matrix in Equation (2) is a highly nonlinear function of ω. Therefore, classical approaches cannot be applied to the transfer matrix method, neither for finding the eigenfrequencies ω ω p = ( ) nor for computing their sensitivities ω p ∂ /∂ k . The former problem has been solved elegantly, for example, by Bestle et al. 9 by transforming the root search problem into a minimization problem, the latter is the focus of the present paper.
Although classical sensitivity analysis, as presented by Murthy such that the last component of the re-ordered eigenvector Z P Z = c 0 0 is one. This requires also a re-ordering of the columns of U in Equation (2) resulting in Since the matrix ω U( ) is singular for eigenfrequencies ( ω U det ( ) = 0), the equations of the linear system (2) can always be re-ordered such that the last one is the dependent equation. Formally this can be achieved by pre-multiplying Equation (4) equivalent to Equation (2) but with the specific substructure Total differentiation d dp •/ k of Equation (5), respectively (6), with respect to a specific design parameter p k results in or after resorting in a linear matrix equation for the desired eigenfrequency sensitivities where the product rule and dependence ω ω p = ( ) have been applied, and derivatives are abbreviated as . An advantage of Equation (8) is that the coefficient matrix has to be computed only once whereas the right-hand side has to be adapted for the different parameter sensitivities • . k , It can be shown that this procedure works in principle, however, it has several drawbacks: (i) it requires specially normalized eigenvectors and reordering of equations according to Equation (6); (ii) for each parameter sensitivity a complete system of equations (8) has to be set up and solved; (iii) since the system transfer matrix U results from a multiplication of several element transfer matrices, the product rule will produce multiple product terms in the derivatives U ω , and U k , where identical matrix products are used several times. Thus, the procedure is computationally inefficient and not recommendable for complex models of technically relevant systems.

| ADJOINT VARIABLE METHOD
The adjoint variable approach is typically much more efficient than direct differentiation if only sensitivities of a scalar output variable (here ω) are required w.r.t. multiple input para- Scalar product with an arbitrary vector where the second term, and thus δZ, can be eliminated for a specific adjoint vector η satisfying the equation This simplifies Equation (10) to which may be solved for the variation of the eigenfrequency as By comparing this expression with the formal variation of ω ω p = ( ) for independent design variables p k , that is, the desired sensitivities can be found as The simple example in Section 6.1 demonstrates the validity of this scalar expression, especially by Equation (43), which is far simpler than solving vector Equation (8).

| COMPUTATION OF EIGENVECTORS
Besides the derivatives U k , and U ω , , evaluation of expression (15) requires the knowledge of right eigenvector Z of system transfer matrix U according to Equation (2) and left eigenvector η according to Equation (11). Both can be obtained from a singular value decom- for example, Golub and van Loan. 11 The singular value 0 expresses the singularity of U for eigenfrequencies.
As already discussed in Section 2, eigenvectors are not unique.
Let Z and η be arbitrarily selected eigenvectors; then any other equivalent eigenvectors may be expressed as ζ Z Z = and ξ η η = with arbitrary scaling factors. With these general eigenvectors, sensitivities (15) read as Thus, all eigenvectors finally yield the same result for the sensitivities, which is why they can be selected freely and there is no need to normalize them.
The right eigenvector Z is contained in the orthogonal matrix R of SVD (16). Substitution of the decomposed U in Equation (2) yields A possible solution is the unit vector because of orthogonality property (16), the final solution is the last column of R: Analogously the left eigenvector η is the last column of the orthogonal matrix L. Substitution of the decomposed U in Equa- and orthogonality property (16) we obtain after pre-multiplication with R T the equation with the possible solution Z e =

| SENSITIVITY ANALYSIS BASED ON ELEMENT TRANSFER MATRICES
The remaining problem in utilizing the sensitivity Equation (15) is to deliver the derivatives U k , and U ω , while avoiding some of the drawbacks discussed in Section 2. A first simplification is given by the observation, that these derivatives are not required explicitly, but only as products with the eigenvector Z. A second simplification is given by the fact that reduction (3) is not necessary, but these sensitivity products may be obtained directly from the transfer matrix U all in Equation (1). To see this, let x be any of the quantities p k or ω. Then, with relations (3) we where Z B Z = all corresponds to and can be easily computed from the eigenvector (20).
The most significant simplification is, that the quantities (24) can be deduced directly from the topology graph by introducing elementspecific adjoints. Let us consider a chain consisting of, for example, three single-input-single-output elements shown in Figure 2A, which can be described by element transfer equations The symbol x may represent the eigenfrequency, any parameter of this specific element or any other design parameter, see Appendix A.
By recursive substitution for the state vectors of the chain we get The expression (24) can then be computed as This can be simplified by recursively assigning adjoints Z i x , as according to the rule There is no need to firstly compute (28) with all its magic resorting, but the assignment of adjoints Z i x , can be made directly along the transfer path according to the rule (30) starting with ≔ Z 0 In case of branched multibody systems, dummy bodies with two inputs I 1 and I 2 may be used, Abbas et al. 2 According to Figure 2B, such massless dummies are described by a transfer equation and a consistency equation where transfer matrices U U , for adjoint assignment. The use of these quantities becomes clear from the example in Section 6.3.

| APPLICATION EXAMPLES
In the following, three examples with different complexity and intensions will be provided. The first is a simple spring-mass vibrator which allows demonstrating the whole concept analytically. The second one shall validate the derivatives of the beam transfer equation, and by attaching an additional spring support, the concept for branched multibody system can be shown in the third example.

| Spring-mass vibrator
The system in Figure 1A combining spring and lumped mass is considered only in the horizontal y-direction. From its topology graph we can immediately see and thus, by using element transfer Equations (75) and (77), find the overall system transfer Equation (1) as all all M S y for the free end, which may be expressed by Equation (3) as This leads to the reduced state transfer Equation (2 In this very simple case, the eigenfrequency ω can be found to parameters m and k yields its sensitivities which will serve as reference values later: To validate sensitivity Equation (15), we first may compute the eigenvectors. By substituting the eigenfrequency (39) in Equation (38), the (right) eigenvector may be obtained: Similarly, the left eigenvector results from Equation (11) as With these eigenvectors and the derivatives of system transfer matrix (38) the required terms for feeding Equation (15) may be computed: This results in sensitivities which is identical with the analytical results (39) obtained by direct differentiation, and thus validates Equation (15).
This calculation based on the overall system transfer matrix U is for demonstration only, because for more complex technical problems U cannot be calculated symbolically. Therefore, the recommended concept is based on eigenvector (40), derivatives (82) and (83) of element transfer matrices and adjoints drawn from the topology graph in Figure 1A according to rule (30): Since U Z Z ≡

| Simply supported beam
The next example shown in Figure 1B  (46) Exactly the same values are obtained by evaluating the equations  from a SVD of the reduced system transfer matrix with abbreviations (79). The eigenvector is extracted from Z.

| Simply supported beam with elastic support
The third example in Figure 1C comes already closer to practical problems, since the branching enforces a typical substructure of the overall transfer equation (1). By recursive substitution of element transfer relations (25) and (31) from right to left in the topology graph of Figure 1C, we find for the transfer equation (51) has to be formulated for the dummy element. Both together may be summarized as According to Equation (24) we need for sensitivity analysis the This resorting allows to recursively introduce adjoints Z k x , according to rules (30) and (33) for each intermediate state vector Z k along the transfer path in Figure 1C as where we have to set ≔ Z 0 for the formally missing terms in Z x 2, and Z x 4, . With these adjoints, the lower part of Equation (54) yields which is consistent with the rule (34) for dummy elements.  Figure 1C and the eigenvector T to deliver quantities (56) and (57) for being used as substitutes (54) in sensitivity formula (15):  To check the result, we may use forward differences based on eigenfrequency evaluations at design point and after disturbing the interesting design parameter by some small parameter variation ε. The result in Table 1 for ε = 10 −4 , however, is rather disappointing for sensitivities w.r.t. EI and k questioning the procedure. Therefore, additionally central differences are evaluated for  ε = 5 10 −5 resulting in sensitivities now consistent with those calculated with the adjoint method, and thus validating the procedure (58).
This reveals one of the major drawbacks of numerical differentiation often used in combination with gradient-based optimization algorithms. If the parameter perturbation ε is too small, numerical errors in ω-evaluation will be amplified too much; if ε is too large, approximation errors of numerical differences come into play, Bestle. 12 In Figure 3 Thus, the proposed method has not only the advantage of being exact, but is also several orders of magnitude faster than numerical difference approximations.

| EXTENSION TO DAMPED SYSTEMS
The concept in Section 3 is also applicable to damped systems where the overall transfer matrix is complex and the real eigenfrequency has to be substituted by the complex eigenvalue   λ σ iω i = ± , = −1 . For showing this, the spring in Figure 1A is substituted by an element with spring and damper in parallel, which is governed by the force law q q k y y d y y = = ( − ) + (̇−̇) The spring-damper transfer matrix and their derivatives are For damped vibrations, also the transfer relations (75) and (82) of the lumped mass need to be adapted to the changed modal transformation resulting in Analogously to Equation (36) we find the overall system transfer which validates the applicability of the concept also for damped systems.

| CONCLUSIONS
The paper develops two different strategies for computing eigenvalue sensitivities of multibody systems modeled with the transfer matrix method, where only the adjoint method is efficient enough to be applied to complex mechanical systems. Simple rules are found for assigning adjoints to each state vector along the transfer path of the topology graph, which are then combined with the left eigenvector of the reduced system transfer matrix in a simple, explicit formula for the desired eigenvalue sensitivities. The rules can be deduced element-wise and require only derivatives of the analytically given element transfer matrices. The procedure is easy to use, exact and computationally efficient by using vector relations instead of matrix multiplications.
Although demonstration examples in the paper are restricted to selected multibody system elements only, the library of element transfer matrices and their derivatives can be easily extended to bodies and spatial systems.

ACKNOWLEDGMENTS
The paper is dedicated to late Laith Abbas who was professor at the Nanjing University of Science and Technology (NUST) and challenged the author by many interesting discussions regarding the transfer matrix method. This study of eigenvalue sensitivities was initialized by Prof. Xiaoting Rui during a research stay at his Institute of Launch Dynamics supported by NUST in 2019, which is greatly appreciated.

CONFLICT OF INTEREST
The author declares that there are no conflict of interest.

DATA AVAILABILITY STATEMENT
Data are available from the corresponding author upon reasonable request.