Non-self-adjoint boundary eigenvalue problems for matrix
differential operators describe distributed non-conservative
systems with the coupled modes and appear in structural mechanics,
fluid dynamics, magnetohydrodynamics, to name a few.
\parPractical needs for optimization and rational experiment planning in modern applications
allow both the differential expression and the boundary conditions to depend
analytically on the spectral parameter and smoothly on several
physical parameters (which can be scalar or distributed). According to the ideas going back to
von Neumann and Wigner , in the multiparameter operator families, eigenvalues
with various algebraic and geometric multiplicities can be generic .
In some applications additional symmetries yield the existence of spectral meshes  in the plane ‘eigenvalue
versus parameter’ containing infinite number of nodes with the multiple eigenvalues [10, 15, 30, 38].
As it has been pointed out already by Rellich  sensitivity analysis of multiple
eigenvalues is complicated by their non-differentiability as functions of several parameters.
Singularities corresponding to the multiple eigenvalues  are related to such important effects
as destabilization paradox in near-Hamiltonian and near-reversible systems [5, 18, 25, 26, 27, 35],
geometric phase ,
reversals of the orientation of the magnetic field in MHD dynamo models , emission of sound by
rotating continua interacting with the friction pads  and other phenomena .
\parAn increasing number of multiparameter non-self-adjoint boundary eigenvalue problems and
the need for simple constructive estimates of critical parameters and eigenvalues as well as for verification of
numerical codes, require development of applicable methods, allowing one to track relatively easily
and conveniently the changes in simple and multiple eigenvalues and the corresponding eigenvectors due to variation
of the differential expression and especially due to transition from one type of boundary conditions to another one
of the original distributed problem, see e.g. [11, 15, 18, 20, 36, 37, 38].
\parA systematical study of bifurcation of eigenvalues of a non-self-adjoint linear operator due to perturbation
where is a small parameter, dates back to 1950s. Apparently, Krein was the first who derived a
formula for the splitting of a double eigenvalue with the Jordan block at the Hamiltonian resonance,
which was expressed through the generalized eigenvectors of the double eigenvalue .
In 1960 Vishik and Lyusternik and in 1965 Lidskii created a perturbation theory
for nonsymmetric matrices and non-self-adjoint differential operators allowing one to find the perturbation coefficients of
eigenvalues and eigenfunctions in an explicit form by means of the eigenelements of the unperturbed operator
[4, 6]. Classical monographs by Rellich , Kato , and Baumgärtel , mostly
focusing on the self-adjoint case, contain a detailed
treatment of eigenvalue problems linearly or quadratically dependent on the spectral parameter.
\parRecently Kirillov and Seyranian proposed a perturbation theory of multiple eigenvalues and eigenvectors
for two-point non-self-adjoint boundary eigenvalue problems with scalar differential expression and boundary conditions,
which depend analytically on the spectral parameter and smoothly on a vector of physical parameters,
and applied it to the sensitivity analysis of distributed non-conservative problems prone to dissipation-induced instabilities [16, 17, 19, 21, 20, 24, 26, 27, 38]. An extension to the case of intermediate boundary conditions with an application
to the problem of the onset of friction-induced oscillations in the moving beam was considered in .
In  this approach was applied to the study of MHD -dynamo model with idealistic boundary conditions.
\parIn the following we develop this theory further and consider boundary eigenvalue problems for linear non-self-adjoint
-th order matrix differential operators on the
interval . The coefficient matrices in the differential
expression and the matrix boundary conditions are assumed to depend
analytically on the spectral parameter and smoothly on a vector of real physical
parameters . The matrix formulation of the boundary conditions is chosen for the
convenience of its implementation in computer algebra systems
for an automatic derivation of the adjoint eigenvalue problem and perturbed eigenvalues and
eigenvectors, which is especially helpful when the order of the derivatives in the differential expression is high.
Based on the eigenelements of the unperturbed problem explicit formulae are derived describing
bifurcation of the semi-simple multiple eigenvalues (diabolical
points) as well as non-derogatory
eigenvalues (branch points, exceptional points) under small variation of the parameters
in the differential expression and in the boundary conditions.
Finally, the general technique is applied to the
investigation of the onset of oscillatory instability in rotating continua.
\parFollowing [8, 19, 22, 24, 26, 27] we consider the boundary eigenvalue problem
where . The differential expression
of the operator is
and the boundary forms are
Introducing the block matrix
and the vector
the boundary conditions can be compactly rewritten as [26, 27]
where and .
It is assumed that the matrices , , and
are analytic functions of the complex spectral parameter
and smooth functions of the real vector of physical parameters .
For some fixed vector the eigenvalue ,
to which the eigenvector corresponds, is a root of the characteristic equation obtained
after substitution of the general solution to equation into the boundary conditions (\@setrefs7) .
\parLet us introduce a scalar product 
where the asterisk denotes complex-conjugate transpose ().
Taking the scalar product of and a vector-function
and integrating it by parts yields
Lagrange formula for the case of operator matrices (cf. [8, 22, 26, 27])
with the adjoint differential expression [8, 22]
and the block matrix
where the matrices are
Extend the original matrix (cf. (\@setrefs7)) to a square matrix , which is
made non-degenerate in a neighborhood of the point and the eigenvalue
by an appropriate choice of the auxiliary matrices and
Similar transformation for the adjoint boundary conditions yields
Then, the form in (\@setrefs13) can be represented as 
so that without loss in generality we can assume [26, 27]
Differentiating the equation (\@setrefs28) we find
Hence, we obtain the formula for calculation of the matrix of the adjoint
boundary conditions and the auxiliary matrix
which exactly reproduces and extends the corresponding result of [26, 27].
\parAssume that in the neighborhood of the point
the spectrum of the boundary eigenvalue problem (\@setrefs1) is discrete. Denote
Let us consider a smooth perturbation of parameters in the form
where and is a small real number.
Then, as in the case of analytic matrix functions ,
the Taylor decomposition of the differential operator matrix and
the matrix of the boundary conditions are [19, 24, 26, 27]
where dot denotes differentiation with respect to at
and partial derivatives are evaluated at , .
Our aim is to derive explicit expressions for the leading terms in the expansions for multiple-semisimple and non-derogatory
eigenvalues and for the corresponding eigenvectors.
\parLet at the point the spectrum contain a semi-simple -fold eigenvalue
with linearly-independent eigenvectors , ,
, . Then, the perturbed eigenvalue and the eigenvector
are represented as Taylor series in
[4, 24, 28, 29, 33]
Substituting expansions (\@setrefp2) and (\@setrefp6) into (\@setrefs1) and collecting
the terms with the same powers of we derive the boundary value problems
Scalar product of (\@setrefp8) with the eigenvectors , of the adjoint boundary
With the use of the Lagrange formula (\@setrefs13), (\@setrefs28) and the boundary conditions (\@setrefp8)
the left hand side of (\@setrefp11) takes the form
Together (\@setrefp11) and (\@setrefp13) result in the equations
Assuming in the equations (\@setrefp14) the vector as a linear combination
and taking into account that
where , we arrive at the matrix eigenvalue problem (cf. [LMZ03])
The entries of the matrices and are defined by the expressions
Therefore, in the first approximation the splitting of the semi-simple eigenvalue due to variation
of parameters is ,
where the coefficients are generically distinct roots of the -th order polynomial
For the formulas (\@setrefp18) and (\@setrefp19) describe perturbation of a simple eigenvalue
The formulas (\@setrefp18),(\@setrefp19), and (\@setrefp20) generalize the corresponding results of the works  to the case
of the multiparameter non-self-adjoint boundary eigenvalue problems for operator matrices.
\parLet at the point the spectrum contain a -fold eigenvalue
with the Keldysh chain of length , consisting of the
eigenvector and the associated vectors . The vectors of the Keldysh chain solve
the following boundary value problems [8, 22]
Consider vector-functions , , , .
Let us take scalar product of the differential equation (\@setrefp21) and
the vector-function . For each
we take the scalar product of the equation (\@setrefp22) and the
vector-function . Summation of the results yields the expression
Applying the Lagrange identity (\@setrefs13), (\@setrefs28) and taking into account relation (\@setrefs29), we transform (\@setrefp23)
to the form
Equation (\@setrefp24) is satisfied in case when the vector-functions , , ,
originate the Keldysh chain of the adjoint boundary value problem, corresponding to the -fold eigenvalue [26, 27]
Taking the scalar product of equation (\@setrefp22) and the vector and employing the expressions (\@setrefs13), (\@setrefs28)
we arrive at the orthogonality conditions
Substituting into equations (\@setrefs1) the Newton-Puiseux series for the perturbed
eigenvalue and eigenvector [7, 13, 26, 27]
where , taking into account expansions (\@setrefp2) and (\@setrefp26)
and collecting terms with the same powers of , yields
boundary value problems
serving for determining the functions ,
where and , ,
are positive integers. The vector-function is a solution
of the following boundary value problem
Comparing equations (\@setrefp31) and (\@setrefp32) with the expressions (\@setrefp22)
we find the first functions in the expansions (\@setrefp26)
With the vectors (\@setrefp33) we transform the equations (\@setrefp31) and (\@setrefp32) into
Applying the expression following from the Lagrange formula
and taking into account the equations for the adjoint Keldysh chain (\@setrefp24a) and
(\@setrefp25a) yields the coefficient in (\@setrefp26). Hence, the splitting of the -fold
non-derogatory eigenvalue due to perturbation of the parameters is described by
the following expression, generalizing the results of the works [19, 24, 26, 27]
For equation (\@setrefp37) is reduced to the equation (\@setrefp20) for a simple
\parConsider a circular string of displacement
, radius , and mass per unit length that rotates with the speed and
passes at through a massless eyelet
generating a constant frictional follower force on the
string, as shown in Fig. \@setreffig1.
The circumferential tension in the string is assumed to be constant;
the stiffness of the spring supporting the eyelet is and the damping coefficient of the viscous damper is .
Introducing the non-dimensional variables and parameters
and assuming we arrive at
the non-self-adjoint boundary eigenvalue problem for a scalar differential operator 
where prime denotes differentiation with respect to .
, , , and express the speed of rotation, and damping, stiffness, and
At the nodes with
the external damper creates a circle of complex eigenvalues being a latent source of subcritical flutter instability
responsible for the emission of sound in the squealing brake and the singing wine glass 
as shown in Fig. \@setreffig3.
Non-conservative perturbation yields eigenvalues with
both the real and imaginary parts of the eigenvalue branches show a degenerate crossing,
touching at the node , Fig. \@setreffig3.
Deformation patterns of the spectral mesh and first-order approximations of the instability tongues
obtained by the perturbation theory and shown in Fig. \@setreffig2 and Fig. \@setreffig3 are in a good
qualitative and quantitative agreement with the results of numerical calculations of .
\parConsider a non-self-adjoint boundary eigenvalue problem
appearing in the theory of MHD -dynamo [14, 32, 33, 37]
with the matrices of the differential expression
and of the boundary conditions
where it is assumed that with .
For the fixed the differential expression
depends on the parameters and , while interpolates between the idealistic ()
and physically realistic () boundary conditions [14, 32, 33, 36].
\parThe matrix of the boundary conditions and auxiliary matrix
for the adjoint differential expression
follow from the formula (\@setrefs31) where the matrices , and are chosen as
In our subsequent consideration we assume that and
interpret and as perturbing parameters.
It is known  that for and the spectrum of
the unperturbed eigenvalue problem (\@setrefdm1)
forms the spectral mesh in the plane , as shown in Fig. \@setreffig4.
The eigenelements of the spectral mesh are
The branches (\@setrefdm6) intersect and
originate a double semi-simple eigenvalue with two linearly independent eigenvectors (\@setrefdm7)
at the node , where 
Taking into account that the components of the eigenfunctions of the adjoint
problems are related as and , we find from equations (\@setrefp18) and (\@setrefp19) the asymptotic formula
for the perturbed eigenvalues, originating after the splitting of the double
semi-simple eigenvalues at the nodes of the spectral mesh
When and , one of the two simple eigenvalues (\@setrefdm10) remains unshifted in all orders
of the perturbation theory with respect to the parameter : .
The sign of the first-order increment to another eigenvalue depends
on the sign of ,
which is directly determined by the Krein signature of the modes involved in the crossing .
This is in the qualitative and quantitative agreement with
the results of numerical calculations of  shown in Fig. \@setreffig4.