Eigen values of differential operators, numerical methods

Methods for computing the eigen values and corresponding eigen functions of differential operators. Oscillations of a bounded elastic body are described by the equation

$$ \tag{1 } \frac{\partial ^ {2} \phi }{\partial t ^ {2} } = L \phi , $$

where $ L \phi $ is some differential expression. If one seeks solutions of (1) of the form

$$ \phi = T ( t) u ( x) , $$

the following equation for $ u $ is obtained:

$$ \tag{2 } L ( u) + \lambda u = 0 , $$

within a bounded domain, under certain homogeneous conditions on the boundary. The values of the parameter $ \lambda $ for which there exist non-zero solutions of (2) satisfying the homogeneous boundary conditions are called eigen values, and the corresponding solutions eigen functions. The resulting eigen value problem consists of determining the eigen values $ \lambda $ and the eigen functions corresponding to them (see [1], [2]).

The numerical solution of this problem proceeds in three steps:

1) reduction of the problem to a simpler one, for example, to an algebraic (discrete) problem (see [7], [16]);

2) specification of the precision of this discrete problem;

3) calculation of the eigen values for the discrete problem (see [3], [6] and Linear algebra, numerical methods in).

Reduction to a discrete problem.

The reduction of problem (2) to a discrete model is usually carried out using the grid method and projection methods. It is natural to require that the basic properties of the original problem, such as the self-adjointness of an operator, be preserved in its discrete analogue.

One method of reduction is the integro-interpolation method. For example, consider the problem

$$ \tag{3 } \frac{d}{dx} \left ( \frac{1}{p ( x) } \frac{du}{dx} \right ) + \lambda u = 0 ,\ 0 < x < 1 , $$

$$ \tag{4 } u ( 0) = 0 ,\ u ( 1) = 0 . $$

This problem arises, for example, in the study of transverse vibrations of a non-homogeneous string and longitudinal vibrations of a non-homogeneous rod.

On the interval $ [ 0 , 1 ] $, one introduces the uniform grid $ \overline \omega \; $ with nodes $ x _ {i} = i h $, $ i = 0 \dots N $, $ h = 1 / N $. To every node $ x _ {i} $, $ i = 1 \dots N - 1 $, is assigned an elementary domain $ S _ {i} = \{ {x } : {x _ {i} - h / 2 \leq x \leq x _ {i} + h / 2 } \} $. Integration of equation (3) over $ S _ {i} $ yields

$$ \tag{5 } \left . \begin{array}{c} \frac{w _ {i-} 1/2 - w _ {i+} 1/2 }{h} = \ \lambda \frac{1}{h} \int\limits _ {x _ {i} - h / 2 } ^ { {x } _ {i} + h / 2 } u d x , \\ w _ {i \pm 1 / 2 } = \ \left ( \frac{1}{p} \frac{du}{dx} \right ) _ {x = x _ {i} \pm h / 2 } . \\ \end{array} \right \} $$

Let

$$ u = \textrm{ const } = u _ {i} ,\ \ x _ {i} - \frac{h}{2} \leq x \leq x _ {i} + \frac{h}{2} , $$

$$ w = \textrm{ const } = w _ {i-} 1/2 ,\ x _ {i-} 1 \leq x \leq x _ {i} . $$

Then

$$ \tag{6 } w _ {i-} 1/2 = \frac{u _ {i} - u _ {i-} 1 }{h a _ {i} } ,\ \ a _ {i} = \frac{1}{h} \int\limits _ { x _ {i-} 1 } ^ { {x _ i } } p ( x) d x . $$

Substituting (6) into (5), one obtains

$$ - \frac{1}{h} \left ( \frac{v _ {i+} 1 - v _ {i} }{a _ {i+} 1 h } - \frac{v _ {i} - v _ {i-} 1 }{a _ {i} h } \right ) = \ \lambda ^ {h} v _ {i} ,\ \ i = 1 \dots N - 1 , $$

where $ v $ is the described grid function. The boundary conditions $ v _ {0} = 0 $, $ v _ {N} = 0 $ reduce to an algebraic eigen value problem:

$$ \tag{7 } A v = \lambda ^ {h} v , $$

where $ A $ is tri-diagonal symmetric matrix of order $ n = N - 1 $.

The variational-difference method of reduction to a discrete problem is used when the eigen value problem can be formulated as a variational one. For example, the eigen values of the problem (3), (4) are the stationary values of the functional

$$ \frac{D [ u] }{H [ u ] } ,\ \ D [ u] = \int\limits _ { 0 } ^ { 1 } \frac{1}{p(} x) \left ( \frac{du}{dx} \right ) ^ {2} \ d x ; \ H [ u] = \int\limits _ { 0 } ^ { 1 } u ^ {2} d x . $$

By changing the integrals to quadrature sums and the derivatives to difference relations, the discrete analogue of this functional has the form

$$ \frac{D ^ {h} [ v] }{H ^ {h} [ v] } ,\ \ D ^ {h} [ v] = \sum _ { i= } 1 ^ { N } \frac{1}{a _ {i} } \left ( \frac{v _ {i} - v _ {i-} 1 }{h} \right ) ^ {2} h ; $$

$$ H ^ {h} [ v] = \sum _ { i= } 0 ^ { N } v _ {i} ^ {2} h , $$

where $ a _ {i} $ is the difference analogue of the coefficient $ p( x) $, which can be calculated by formula (6). The discrete analogue of the problem (3), (4) is obtained from the necessary condition for an extremum:

$$ \frac \partial {\partial v _ {i} } ( D ^ {h} [ v] - \lambda ^ {h} H ^ {h} [ v] ) = 0 ,\ \ i = 1 \dots N - 1 . $$

Differentiation again leads to the problem (7) (see [9]).

The projection-difference method of reduction to a discrete problem consists of the following. One chooses a linearly independent system of coordinate functions $ \alpha _ {i} $, $ i = 1 \dots n $, and a linearly independent system of projection functions $ \beta _ {j} $, $ j= 1 \dots n $. One seeks approximate eigen functions in the form

$$ \widetilde{u} = \sum _ { i= } 1 ^ { n } v _ {i} \alpha _ {i} . $$

The coefficients $ v _ {i} $ and the approximate eigen values are determined by the condition

$$ \tag{8 } ( L \widetilde{u} + \widetilde \lambda \widetilde{u} , \beta _ {j} ) = 0 ,\ \ j = 1 \dots n , $$

where $ ( , ) $ is the scalar product in the Hilbert space. When the coordinate and projection systems coincide, one talks of the Bubnov–Galerkin method. If, in addition, the operator of the differential problem is self-adjoint, the method is called the Raleigh–Ritz method. In particular, for problem (3), (4), if all $ \alpha _ {i} = \beta _ {i} $ satisfy (4), condition (8) takes the form

$$ \tag{9 } \sum _ { i= } 1 ^ { n } v _ {i} \int\limits _ { 0 } ^ { 1 } \left ( \frac{1}{p(} x) \frac{d \alpha _ {i} }{dx} \frac{d \alpha _ {j} }{dx} - \widetilde \lambda \alpha _ {i} \alpha _ {j} \right ) d x = 0 , $$

$$ j = 1 \dots n. $$

In order to simplify the derivation of the algebraic problem, one chooses an almost-orthogonal system of functions $ \alpha _ {i} $.

By taking as coordinate and projection systems functions of the form

$$ \alpha _ {i} ( x) - \left \{ \begin{array}{cl} 0 , & 0 \leq x \leq x _ {i-} 1 , \\ \frac{x - x _ {i-} 1 }{h} , & x _ {i-} 1 \leq x \leq x _ {i} , \\ \frac{x _ {i+} 1 - x }{h} , & x _ {i} \leq x \leq x _ {i+} 1 , \\ 0 , & x _ {i+} 1 \leq x \leq 1 , \\ \end{array} \right .$$

where $ x _ {i} = i h $, it follows from (9) that

$$ \frac{v _ {i} - v _ {i-} 1 }{a _ {i} h } - \frac{v _ {i+} 1 - v _ {i} }{a _ {i+} 1 h } + $$

$$ - \lambda ^ {h} [ ( \rho v ) _ {i-} 1 + ( \rho ^ {*} v ) _ {i} + ( \rho v ) _ {i+} 1 ] = 0 , $$

where

$$ a _ {i} = \frac{1}{h} \int\limits _ { x _ {i-} 1 } ^ { {x _ i } } \frac{dx}{p(} x) , $$

$$ \rho ^ {i} = \int\limits _ {x _ {i} } ^ { {x _ i+} 1 } \alpha _ {i} \alpha _ {i+} 1 d x ,\ \rho _ {i} ^ {*} = \int\limits _ {x _ {i-} 1 } ^ { {x _ i+} 1 } \alpha _ {i} ^ {2} d x . $$

Thus, together with the boundary conditions, one obtains a generalized eigen value problem:

$$ A v = \lambda ^ {h} D v ,\ \ h = \frac{1}{N} . $$

Here $ A $ and $ D $ are tri-diagonal symmetric matrices of order $ n= N - 1 $.

These methods also yield discrete models of other equations. For example, for a rod:

$$ \frac{d ^ {2} }{d x ^ {2} } \left ( k \frac{d ^ {2} u }{d x ^ {2} } \right ) = \ \lambda r u ; $$

for a membrane:

$$ \frac \partial {\partial x } \left ( k _ {1} \frac{\partial u }{\partial x } \right ) + \frac \partial {\partial y } \left ( k _ {2} \frac{\partial u }{\partial y } \right ) + \lambda r u = 0 ; $$

and for a plate:

$$ \frac{\partial ^ {2} }{\partial x ^ {2} } \left ( k _ {11} \frac{\partial ^ {2} u }{\partial x ^ {2} } + k _ {12} \frac{\partial ^ {2} u }{\partial y ^ {2} } \right ) + 2 \frac{\partial ^ {2} }{\partial x \partial y } \left ( k _ {33} \frac{\partial ^ {2} u }{\partial x \partial y } \right ) + $$

$$ + \frac{\partial ^ {2} }{\partial y ^ {2} } \left ( k _ {21} \frac{\partial ^ {2} u }{\partial x ^ {2} } + k _ {22} \frac{\partial ^ {2} u }{\partial y ^ {2} } \right ) = \lambda r u . $$

The eigen vectors corresponding to $ \lambda _ {k} ^ {h} $ satisfy the homogeneous system of algebraic equations:

$$ ( A - \lambda _ {k} ^ {h} E ) v _ {k} = 0 . $$

The problem of finding all eigen values and eigen vectors of a matrix $ A $ is called the complete problem of eigen values, and the problem of finding some eigen values of $ A $ is called the partial problem of eigen values.

In the case of algebraic systems corresponding to a given problem, the latter problem arises most frequently. The application of traditional solution methods requires a considerable amount of calculation, in view of the poor separation of the eigen values of $ A $. In this case it is more effective to use modified gradient methods with spectrally equivalent operators (see [17]- [19]) and multi-grid methods (see [20]).

References

[1]	S. Gould, "Variational methods for eigenvalue problems" , Univ. Toronto Press (1957)
[2]	L. Collatz, "Eigenwertaufgaben mit technischen Anwendungen" , Akad. Verlagsgesell. Geest u. Portig K.-D. (1949)
[3]	D.K. Faddeev, V.N. Faddeeva, "Computational methods of linear algebra" , Freeman (1963) (Translated from Russian)
[4]	V.V. Voevodin, "Numerical methods of algebra" , Moscow (1966) (In Russian)
[5]	J.H. Wilkinson, "The algebraic eigenvalue problem" , Oxford Univ. Press (1969)
[6]	B.N. Parlett, "The symmetric eigenvalue problem" , Prentice-Hall (1980)
[7]	V.K. Saul'ev, "On estimating the error in finding eigenfunctions by the finite difference method" Vychisl. Mat. , 1 (1957) pp. 87–115 (In Russian)
[8]	A.N. Tikhonov, A.A. Samarskii, "The Sturm–Liouville difference problem" USSR Comput. Math. Math. Phys. , 1 : 4 (1961) pp. 939–961 Zh. Vychisl. Mat. i Mat. Fiz. , 1 : 5 (1961) pp. 784–805
[9]	A.A. Samarskii, "Theorie der Differenzverfahren" , Akad. Verlagsgesell. Geest u. Portig K.-D. (1984) (Translated from Russian)
[10]	Shou Hao, "The Sturm–Liouville difference problem for a fourth order equation with discontinuous coefficients" USSR Comput. Math. Math. Phys. , 3 : 6 (1963) pp. 1383–1406 Zh. Vychisl. Mat. i Mat. Fiz. , 3 : 6 (1963) pp. 1014–1031
[11]	V.G. Prikazchikov, "High-accuracy homogeneous difference schemes for the Sturm–Liouville problem" USSR Comput. Math. Math. Phys. , 9 : 2 pp. 76–106 Zh. Vychisl. Mat. Mat. Fiz. , 9 : 2 (1969) pp. 315–336
[12]	V.G. Prikazchikov, "A difference eigenvalue problem for a second-order operator with mixed boundary conditions" USSR Comput. Math. Math. Phys. , 22 : 3 (1982) pp. 165–172 Zh. Vychisl. Mat. Mat. Fiz. , 22 : 3 (1982) pp. 655–662
[13]	V.G. Prikazchikov, A.N. Khimich, "The eigenvalue difference problem for the fourth order elliptic operator with mixed boundary conditions" USSR. Comput. Math. Math. Phys. , 25 : 5 (1985) pp. 137–144 Zh. Vychisl. Mat. Mat. Fiz. , 25 (1985)
[14]	B.N. Bublik, "The numerical solution of problems in the dynamics of plates and hulls" , Kiev (1969) (In Russian)
[15]	G.J. Fix, "An analyse of the finite element method" , Prentice-Hall (1973)
[16]	P.G. Ciarlet, "The finite element method for elliptic problems" , North-Holland (1978)
[17]	B.A. Samokish, "The method of steepest descent in the eigen elements problem for semi-bounded operators" Izv. Vyssh. Uchebn. Zaved. Mat. , 5 : 6 (1958) pp. 105–114 (In Russian)
[18]	E.G. D'yakonov, M.Yu. Orekhov, "Minimization of the computational labor in determining the first eigen values of differential operators" Math. Notes , 27 (1980) pp. 382–391 Mat. Zametki , 27 : 5 (1980) pp. 795–812
[19]	V.G. Prikazchikov, Differentsial'nye Uravn. , 22 (1986) pp. 1268–1271
[20]	W. Hackbusch, "On the computation of approximate eigen values and eigen functions of elliptic operators of means of a multi-grid method" SIAM. J. Num. Anal. , 16 : 2 (1979) pp. 201–215

References

[a1]	P.B. Bailey, M.K. Gordon, L.F. Shampine, "Automatic solution of the Sturm–Liouville problem" ACM Trans. Math. Software , 4 (1978) pp. 193–208
[a2]	H.O. Kreiss, "Difference approximation for boundary and eigen values problems for ordinary differential equations" Math. Comp. , 26 (1972) pp. 605–624
[a3]	J.S. Bramley, S.C.R. Dennis, "The calculation of eigenvalues for the stationary perturbation of Poiseuille flow using initial value methods" Math. Anal. Appl. , 101 (1984) pp. 30–38