Eigenvalue Computation of Regular 4th Order Sturm-Liouville Problems

In this paper we present and test a numerical method for computing eigenvalues of 4th order Sturm-Liouville (SL) differential operators on finite intervals with regular boundary conditions. This method is a 4th order shooting method based on Magnus expansions (MG4) which use MG4 shooting as the integrator. This method is similar to the SLEUTH (Sturm-Liouville Eigenvalues Using Theta Matrices) method of Greenberg and Marletta which uses the 2nd order Pruess method (also known as the MG2 shooting method) for the integrator. This method often achieves near machine precision accuracies, and some comparisons of its performance against the well-known SLEUTH software package are presented.


Introduction
In this paper 1 we consider the self-adjoint differential operators which arise from the 4th order differential equation where 1 2 , A A and 1 2 , B B are any choice of real, 2 2 × matrices satisfying the T T 1 1 2 2

. A A A A
I + = (5) and The above boundary conditions can be shown to be equivalent to the general forms of boundary conditions used by Everitt [2] (in his PhD dissertation on 4th order Sturm-Liouville problems), Fulton [3] and many others.
The domain of the maximal operator 1 L associated with the Equation (1) where the bilinear concomitant is defined as : .
x y x z x y x z x y x z x y z s x y x z x y x z x y x z x Using this definition, and the boundary conditions (2) and (3), it can be shown that the operators 1

Selection of Test Problems
To investigate the performance of the method, we make the following selection of test problems. These problems are the square of a 2nd order SL problem.
1) The square of the 2nd order Bessel equation and ( ) 2) The square of the 2nd order Modified Harmonic Oscillator equation and ( ) ( ) ( ) 3) The square of the 2nd order equation 4) The square of the 2nd order Coffey-Evan equation sec sec sec tan . 16 2 Problem 5, the Legendre squared equation, arises from changes of variables to the non-LNF form discussed in [31].

The MG4 Shooting Method Associated with the 4th Order Sturm-Liouville Equation
In this section we describe an implementation of the MG4 shooting technique for the 4th order SL Equation (1) on regular intervals with ( ) ( ) , s x q x continuous. The Equation (1) can be converted to the 1st order system (Atkinson [32], pp. 323-324) Remark 2.1 Currently the most reliable software package for eigenvalues and eigenfunctions of the 4th order Sturm-Liouville equation with regular endpoints is ACM Algorithm 775: SUBROUTINE SLEUTH, produced by L. Greenberg and M. Marletta in 1997 [14], which is available from NETLIB at ORNL. The SLEUTH (Sturm-Liouville Eigenvalues Using Theta Matrices) code employs an MG2 approximation for the solution ( ) , (see [33], p. 283), on each mesh interval of a Hamiltonian system similar to the system (1) (the order of the derivatives in the above ( ) , Y x λ vector being slightly different). The SLEUTH code is based on using formulas of Greenberg [15] [16] [34] [35] for the number of eigenvalues less than λ , so the eigenvalue algorithm is quite different than the method we are proposing here for eigenvalue computation for (1) using (1).
For the IVP, trix is also basic to Magnus methods for (1). We introduce the following lemma and the definitions of the Lie-group and Lie-algebra (see [36] The Magnus methods originate (see [37]) with the observation that an analyt- can be written as, (see [33], and where the square brackets denote the matrix commutator and are defined as: The MG4 method is a well known 4th order method obtained by truncation of the above Magnus series, together with evaluation of the A matrix in (1) at two gaussian points 1 A and 2 A : For the Hamiltonian system (1), we put ( ) where the square bracket [ ] 2 1 , A A denotes the matrix commutator and is de- Case 2: Case 3: is a solution of which remains in the Lie Group, We consider the SL problem for Equation (1) with the following choices of Dirichlet boundary conditions at the left and right endpoints (compare (2) and (3)).
y a y a s a y a y a l y A A y a y a y a l y where and The left boundary conditions are implemented by fixing initial conditions for two solutions The two-dimensional subspace, By fixing (57) the two-dimensional space ℑ is fixed by the 4 2 × matrix (57). The constants in the (57) matrix were chosen to ensure that the boundary conditions at x a = (43) was satisfied, so we know that that is, that both In this case the eigenfunction, which is unique up to a constant multiple, is For multiplicity two, we need to require and in this case

Description of the MG4 Algorithm
We obtained the 4 4 × transfer matrix by doing the following steps: 1) Calculate the eigenvalues and the eigenvectors of A  . 2) Diagonalize where P denotes the matrix of eigenvectors of A  , and D denotes the diagonal matrix of the four eigenvalues of A  (12). , , : det , 0. , ,

3) Put
The computation is performed using an initial uniform mesh, applying bisection method with initial upper and lower bounds for a given eigenvalue n λ , and then doubling the number of mesh points by bisecting the mesh to generate a Richardson h 4 -extrapolation table over successively bisected meshes. Then the extrapolated eigenvalue is selected when the eigenvalue extrapolation error satisfies a tolerance test.

Description of h 4 -Richardson's Extrapolation
Putting where we have taken : 0 we have from (13)

Computing Large Eigenvalues by Using the Invariant-Imbedding Variables
In a manner similar to Greenberg and Marletta in their SLEUTH code (see [14],

Numerical Results and Discussion
In this section we give some numerical outputs for each of the 5 test problems in Squaring the self-adjoint operator corresponding to (1)-(5) gives the 4th order self-adjoint operator corresponding to the 4th order problems in Tables 1-5, respectively. Accordingly, the eigenvalues of the problems in Tables 1-5

Conclusion
In this paper we have presented the MG4 algorithm of Iserles et al. [33] [40], for eigenvalue computation of regular 4th order Sturm-Liouville problems. Applying the change to "Invariant-Imbedding" variables in a manner similar to Greenberg and Marletta in their SLEUTH code [14] provides good stabilization for the MG4 algorithm, and this resulted in its capability for achieving high accuracy, very often on the order of machine precision. The MG4 algorithm appears to be competitive to the SLEUTH code.