EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS BY ...
MATHEMATICS OF COMPUTATION Volume 68, Number 227, Pages 1057–1066 S 0025-5718(99)01053-4 Article electronically published on February 10, 1999
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS BY THE SHANNON-WHITTAKER SAMPLING THEOREM AMIN BOUMENIR
Abstract. We are concerned with the computation of eigenvalues of a periodic Sturm-Liouville problem using interpolation techniques in Paley-Wiener spaces. We shall approximate the Hill discriminant by sampling a few of its values and then find its zeroes which are the square roots of the eigenvalues. Computable error estimates are provided together with eigenvalue enclosures.
1. Introduction We would like to introduce a new method for computing the eigenvalues of classical regular Sturm-Liouville problems with periodic boundary conditions t ∈ [0, ω] −ϕ00 (t, µ) + Q(t)ϕ(t, µ) = µ2 ϕ(t, µ) (1.1) ϕ(0, µ) = ϕ(ω, µ) and ϕ0 (0, µ) = ϕ0 (ω, µ), where Q(t) ∈ L1 (0, ω). For the spectral theory of periodic differential equations, we shall refer to [4], [1] and [8]. Recall that in general the spectrum of (1.1) may not be simple as in the case of separated boundary conditions, and this is a major difficulty. In this work, we shall extend the method in [3], which relies on the interpolation in Paley-Wiener spaces of a certain boundary function, whose zeros are square roots of eigenvalues. In our case the boundary function turns out to be the wellknown Hill discriminant (see [4]) and only few values are needed for its approximation. The truncation error can be minimized by increasing the number of sampling points, which gives higher accuracy on the numerical side and provides eigenvalue enclosures. We shall examine a few examples, where eigenvalue enclosures and comparisons of the results with the well-known code Sleign2 are provided. In the last example our interpolation scheme is compared with an “exact” solution. One of the features of the method is the possibility of locating double eigenvalues. Indeed, by minimizing the truncation error, we can zoom in on close eigenvalues, and if they happen to be simple, then we can find two disjoint enclosures separating them. 2. Hill’s discriminant The Hill discriminant, which is at the heart of the spectral theory of the periodic Sturm-Liouville operator, is defined by ∆(µ) := ϕ1 (ω, µ) + ϕ02 (ω, µ), Received by the editor June 9, 1997 and, in revised form, October 28, 1997. 1991 Mathematics Subject Classification. Primary 34L15, 42A15. Key words and phrases. Periodic Sturm-Liouville problems, eigenvalue problem, interpolation, Shannon-Whittaker sampling theorem.
1057
c
1999 American Mathematical Society
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1058
AMIN BOUMENIR
where ϕ1 (t, µ) and ϕ2 (t, µ) are two independent solutions of (1.1) satisfying ϕ1 (0, µ) = 1 ϕ2 (0, µ) = 0 (2.1) and ϕ01 (0, µ) = 0 ϕ02 (0, µ) = 1. Recall that µ2 is an eigenvalue associated with (1.1) if and only if µ is solution of (2.2)
∆(µ) = 2
(see [4]). In order to interpolate the Hill discriminant we shall need to find a suitable representation in terms of analytic functions. From the inverse spectral theory (see [7]) ϕ1 (t, µ) and ϕ2 (t, µ) can be represented by Z t ϕ1 (t, µ) = cos tµ + (2.3) K1 (t, η) cos ηµ dη, 0 Z t sin tµ sin ηµ + dη, K2 (t, η) ϕ2 (t, µ) = µ µ 0 where K1 and K2 are the kernels of transformation operators. We recall that K1 and K2 have m + 1 locally integrable derivatives if Q has m locally integrable ∂ ∂ K2 and ∂t K1 are locally derivatives. In our case m = 0, which means that ∂t integrable and K1 and K2 are continuous functions. Thus the integrals in (2.3) are well defined. The Lyapunov function or Hill’s discriminant for (1.1) can now be described using K1 and K2 : (2.4)
∆(µ)
= =
ϕ1 (ω, µ) + ϕ02 (ω, µ) sin ωµ sin ωµ + K2 (ω, ω) 2 cos ωµ + K1 (ω, ω) µ µ Z ω ∂ ∂ sin ηµ K2 (ω, η) + K1 (ω, η) dη. + ∂t ∂t µ 0
Let us agree to denote the last integral term by Z ω ∂ ∂ sin ηµ (2.5) K2 (ω, η) + K1 (ω, η) dη, S(µ) := ∂t ∂t µ 0 and several useful constants by sin z exp (|Imz|) z ≤ c 1 + |z| , Z ω q := |Q(η)| dη 0
and
Z ξ := exp c
0
ω
η |Q(η)| dη .
We now recall that Paley-Wiener spaces are defined by Z P Wω := F (µ) entire: |F (µ)| < M eω|Imµ| and
∞
−∞
Also (see [7]) 1 K1 (ω, ω) = K2 (ω, ω) = 2
Z 0
|F (µ)| dµ < ∞ .
ω
Q(η)dη.
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
2
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS
Thus in order to recover (2.6)
∆(µ) := 2 cos(ωµ) +
sin ωµ µ
Z
1059
ω
Q(η)dη + S(µ)
0
we need to find S(µ) only. Proposition 2.1. If Q(t) ∈ L1 (0, ω), then sZ r ∞ 3 2 ω. S(µ) ∈ P Wω and E1 := |µS(µ)| dµ ≤ ξc2 q 2 2 −∞ ∂ ∂ K2 (ω, η) + ∂t K1 (ω, η) Proof. Because Q(t) is locally integrable, it follows that ∂t is also locally integrable, and applying the Riemann-Lebesgue theorem to (2.5) we deduce that S(µ) = o( µ1 ), i.e., S(µ) ∈ L2dµ (−∞, ∞). Next observe that Z ω ∂ cη exp (η |Imµ|) ∂ dη |S(µ)| ≤ ∂t K2 (ω, η) + ∂t K1 (ω, η) 1 + |ηµ| 0 Z ω ∂ cη K2 (ω, η) + ∂ K1 (ω, η) ≤ dη exp (ω |Imµ|) ∂t ∂t 1 + |ηµ| 0 ≤ C1 exp (ω |Imµ|) . Hence S(µ) ∈ P Wω . If C1 is known, then we can estimate the constant E1 , which is crucial in the computation of error bounds. To this end we shall use the following integral equations (see [10]): ϕ1 (t, µ)
(2.7)
= cos(tµ) + Tµ ϕ1 (t, µ), sin (tµ) + Tµ ϕ2 (t, µ), = µ
ϕ2 (t, µ) where
Tµ
Z
C(0, ω) −→ C(0, ω),
sin(µ(t − η)) Q(η)f (η)dη, t ∈ (0, ω), µ 0 Z ω c(ω − η) ||Tµ || ≤ |Q(η)| dη 1 + |µ| (ω − η) 0 ω ≤ cq . 1 + |µ| ω t
Tµ f (t) := (2.8)
For µ large enough, the method of successive approximation yields (see [10]) X X sin ηµ , ϕ1 (t, µ) = Tµn cos(tµ) and ϕ2 (t, µ) = Tµn µ n≥0
n≥0
and by the Gronwall lemma applied to (2.7) we obtain Z t Z t c η |Q(η)| |ϕ1 (t, µ)| ≤ exp dη ≤ exp c η |Q(η)| dη , 0 1 + |µ| η 0 |ϕ2 (t, µ)|
t exp 1 + |µ| t t . ≤ cξ 1 + |µ| t ≤ c
Z 0
t
cη |Q (η)| dη 1 + |µ| η
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1060
AMIN BOUMENIR
Applying Tµ to (2.7) yields ϕ1 (t, µ) = cos(tµ) + Tµ cos(tµ) + Tµ2 ϕ1 (t, µ),
(2.9) where
2 T ϕ1 (ω, µ) µ
≤ ≤ ≤
(2.10)
Z
2
(c qω)
2
(1 + |µ| ω) 2
(c qω)
(1 + |µ| ω)2 2
ξ (c qω)
c η |Q(η)| dη 0 1 + |µ| η Z t exp c η |Q(η)| dη t
exp
0
1 2
(1 + |µ| ω)
.
Similarly we have ϕ2 (t, µ) = ϕ02 (t, µ) =
(2.11)
where d 2 Tµ ϕ2 (ω, µ) dt
= ≤ ≤ ≤ ≤ ≤ ≤
(2.12)
sin(tµ) sin(tµ) + Tµ + Tµ2 ϕ2 (t, µ), µ µ d d sin(tµ) + Tµ2 ϕ2 (t, µ), cos(tµ) + Tµ dt µ dt
Z ω Z t sin(µ(t − η)) Q(η)ϕ2 (η, µ)dηdt cos(µ(ω − t))Q(t) µ 0 0 Z ω Z t c (t − η) |Q(η)| |ϕ2 (η, µ)| dηdt |Q(t)| 0 0 1 + |µ| (t − η) Z t Z ω η cη |Q(η)| ξc dηdt |Q(t)| 1 + |µ| η 0 0 1 + |µ| η 2 Z ω Z t η 2 ξc |Q(t)| |Q(η)| dηdt 1 + |µ| η 0 0 2 Z ω Z t ω ξc2 |Q(t)| |Q(η)| dηdt 1 + |µ| ω 0 0 2 Z ω 2 ξc2 ω |Q(t)| dt 2 1 + |µ| ω 0 1 2 2 2 1 ξc q ω 2. 2 (1 + |µ| ω)
We now use (2.9) and (2.11) to obtain ∆(µ) = =
ϕ1 (ω, µ) + ϕ02 (ω, µ) 2 cos(tµ) + Tµ cos(tµ) +
d d sin(tµ) Tµ + Tµ2 ϕ1 (t, µ) + Tµ2 ϕ2 (t, µ). dt µ dt
Observe that d sin(tµ) Tµ cos(tµ) + Tµ dt µ Z t sin(ηµ) sin(µ(t − η)) cos(ηµ) + cos(µ(t − η)) Q(η)dη = µ µ 0 Z t sin(µt) = Q(η)dη. µ 0
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS
1061
Thus a new representation of the Hill discriminant follows: (2.13) sin(µt) ∆(µ) = 2 cos(tµ) + µ
Z
t
0
Q(η)dη + Tµ2 ϕ1 (t, µ) +
d 2 T ϕ2 (t, µ). dt µ
Therefore (2.6) together with (2.13) imply (2.14)
S(µ) = Tµ2 ϕ1 (t, µ) +
d 2 T (ϕ2 (t, µ)) dt µ
and estimates (2.12) and (2.10) yield (2.15)
|S(µ)| ≤
3ξc2 q 2 ω 2 1 2 . 2 (1 + |µ| ω)
Hence it follows from (2.14) and (2.15) that S(µ) ∈ P Wω and sZ r r 3ξc2 2 2 2 3 2 2 2 q ω ω. |µS(µ)| dµ ≤ ≤ ξc q E1 = 2 3ω 3 2 It is also possible to use asymptotics of solutions (see [5]) to derive the estimates above. 3. General case We now address the question of how to apply the idea above to more general equations, such as 0 0 < x < a, − (p(x)y 0 (x, µ)) + q(x)y(x, µ) = µ2 s(x)y(x, µ), (3.1) y(0, µ) = y(a, µ), p(0)y 0 (0, µ) = p(a)y 0 (a, µ), 1 where p(x), q(x), and s(x) are real valued, p(x) , q(x), s(x) ∈ L1 (0, a) and s(x) ≥ 0. It is well known that µ2n (see [4]) is an eigenvalue of (3.1) if and only if µn is a root of
(3.2)
D(µ) := ψ1 (a, µ) + p(a)ψ20 (a, µ) = 2,
where ψ1 (x, µ) and ψ2 (x, µ) are two independent solutions of (p(x)y 0 (x, µ))0 + q(x)y(x, µ) = µ2 s(x)y(x, µ), with initial conditions ψ1 (0, µ) = 1 p(0)ψ10 (0, µ) = 0
and
0 < x < a,
ψ2 (0, µ) = 0 p(0)ψ20 (0, µ) = 1.
It is readily seen that D(µ) is an entire function of µ of order one, and the main difficulty is estimating its type. Once the type is obtained, which depends on p(x), q(x) and s(x), we can apply the above method with a few modifications. One simple way to deal with (3.1) is to transform it to the standard form (1.1) 1 and s(x) are C 2 (0, a), and s(x) ≥ δ > 0. which is possible when the coefficients p(x) Indeed, the Liouville-Green transformation allows us to recast equation (3.1) into an equivalent periodic problem (see [4]) −ϕ00 (t, µ) + Q(t)ϕ(t, µ) = µ2 ϕ(t, µ) t ∈ [0, A] (3.3) ϕ(0, µ) = ϕ(A, µ) and ϕ0 (0, µ) = Kϕ0 (A, µ),
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1062
AMIN BOUMENIR
where t(x) =
R x q s(η) 0
p(η) dη
for x ∈ [0, a] and
−1 d d p(x) {p(x)s(x)} 4 . dx dx Clearly the type of D(µ) is then A := t(a) and the boundary conditions q q s(x) s(ω)p(ω) 0 0 are obtained from y (x, µ) := ϕ (t(x), µ) p(x) , and so K := s(0)p(0) . In case s(ω)p (ω) =s(0)p (0) , then we have K = 1. • If p(x) and s(x) are not smooth, then we can still use ideas from the transformation operators for the string (see [2]) 1
3
ϕ(t(x), µ) = y(x, µ), Q(t(x)) = q(x) − p 4 (x)s− 4 (x)
−p(x)(p(x)y 0 (x, µ))0 + p(x)q(x)y(x, µ) = µ2 p(x)s(x)y(x, µ) 0 < x < a. Rx 1 Then setting t(x) := 0 p(η) dη, γ := t(a) and ϕ(t(x)) = y(x) we have (3.4)
−ϕ00 (t, µ) + Q(t)ϕ(t, µ) = µ2 ω(t)ϕ(t, µ)
0 0. N
µ∈IµN ∗
Proposition 4.1. Let µ > 0, be a simple eigenvalue of (1.1), then ∀ε > 0, ∃N0 such that ∀N ≥ N0 , ∃µN such that ε . |µ∗ − µN | < e0 inf ∆ (µ) µ∈IµN
N
e0 Proof. Given µ∗ , we can find N, such that µ∗ < Nωπ , inf ∆ N (µ) ≥ δ > 0 and µ∈IµN e ∗ ∗ e TN (µ) < ε for µ ∈ IµN . From ∆ (µ ) − ∆ (µ ) N N N ≤ TN (µ ), the mean value theorem implies ∗ e 0 (η) ≤ TN (µ∗ ), (4.3) η ∈ (µ∗ , µN ) ⊂ IµN , (µ − µN ) ∆ N e0 and since µ∗ is simple, inf ∆ N (µ) > 0, and we have µ∈IµN
|µ∗ − µN | ≤
TN (µ∗ ) ≤ e0 inf ∆ N (µ)
µ∈IµN
ε e0 inf ∆N (µ) .
µ∈IµN
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1064
AMIN BOUMENIR
As we increase N ≥ N0 , ∗
sup TN (µ) → 0,
µ∈IµN
e0 inf ∆ N (µ) → δ(N0 ), and so
µ∈IµN
|µ − µN | → 0 as N → ∞. • The second case is when we have a double eigenvalue, which is characterized by e (µ∗ ) = 0, ∆ e 0 (µ∗ ) = 0 and ∆ e 00 (µ∗ ) 6= 0. ∆ e Because of the concavity of the ∆(µ) at µ∗ , there are two cases. e – ∆N (µ) is also concave and may not have zeroes. In this case we look for solutions inside the intervals defined by (4.2). e N (µN ) has two zeroes then – If ∆ e N (µN ) + (µ∗ − µN )∆ e 0N (µN ) + 1 (µ∗ − µN )2 ∆ e 00N (η), e N (µ∗ ) = ∆ ∆ 2 e ∗ ) + (µ∗ − µN )∆ e 0N (µN ) + 1 (µ∗ − µN )2 ∆ e 00N (η), e N (µ∗ ) = ∆(µ ∆ 2 e ∗ ) = (µ∗ − µN )∆ e 0 (µN ) + 1 (µ∗ − µN )2 ∆ e 00 (η), e N (µ∗ ) − ∆(µ ∆ N N 2 h i 1 ∗ 2 e 00 ∗ 0 ∗ e (µN ) − ∆ e N (µ ) − ∆(µ e ∗ ) = 0, (µ − µN ) ∆N (η) + (µ − µN )∆ (4.4) N 2 i h e ∗ ) ≤ TN (µ∗ ). Clearly the coefficients of this quae N (µ∗ ) − ∆(µ where ∆ e 00 (η) ≤ max ∆ e 00 (η) ≤ ∆ e 00 (η), and simidratic are intervals, i.e., inf ∆ η∈Iµ
N
N
η∈Iµ
N
larly 0 ≤ TN (µ∗ ) ≤ max TN (η). Therefore we can find an upper bound η∈Iµ
e N (µN ) may have no roots for |µ∗ − µN | , using interval analysis. Since ∆ in the neighbourhood of µ∗ , the only practical way to locate double eigenvalues is the interval defined by (4.2). In practice it is more convenient to be given N and then compute the truncation error TN (µ) to find the intervals IµN defined by (4.2). 5. Examples We shall consider problems defined by (1.1), where Q(x) is integrable. Recall that N represents the number of points used in the Shannon-Whittaker sampling theorem. The positive eigenvalues are obtained by recovering the Hill discriminant from its sampled values. To this end we only need to compute the values ϕ1 (ω, µ) + ϕ02 (ω, µ) for µ = k ωπ , k = 0, 1, ..., N , which are obtained by solving the initial value problem (1.1) with Runge-Kutta 4-5, which is implemented with the symbolic manipulator Maple with precision set to 13 digits. • Example 1: Let us choose ω = 2, and 0 < x < 1/2 , 0 1 1/2 < x < 7/4 , q(x) := 2 1+x 7/4 < x < 2. 0 The square roots of the eigenvalues are listed below. N 5 13 19
µ1 0.5376360812 0.5376389426 0.5376389987
µ2 3.141592653 3.141592653 3.141592654
µ3 3.198533821 3.198655062 3.198665145
µ4 6.283185307 6.283185307 6.283185307
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
µ5 6.309309048 6.309851843 6.309892873
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS
1065
The enclosure for the first simple eigenvalue µ21 is given below. N 5 13 19
µmin 0.5063 0.5292 0.5327
µ1 0.5376360812 0.5376389426 0.5376389987
µmax 0.5700 0.5461 0.5425
• Example 2: q(x) := 2x(1 − x) and ω = 1. The enclosures are given below. N 4 7 15 20 N 4 7 15 20
µ1 min 0.576496 0.575138 0.576292 0.576496
µ2,3 min 6.29811294928 6.30116798050 6.30414770889 6.30496865578
µ1 0.576890580 0.576891292 0.576891400 0.576891408
µ2 6.30777169387 6.30769833630 6.30768224351 6.30768111508
µ1 max 0.577287 0.578648 0.577491 0.577287
µ3 6.31155085736 6.31163581017 6.31165451250 6.31165582479
µ2,3 max 6.32991780134 6.32196796225 6.31645456422 6.31520080034
For N = 20, an estimate of the truncation error is given by T N (µ) < 27.10−6 when µ ∈ (6.304968, 6.3152008) T N (µ) < 46.10−5
µ1 ∈ (0.5764, 0.5774).
1. A quick run of Sleign2, with a tolerance of 10−7 , gives the following eigenvalues Sleign2
µ=
√ λ
0.576919405 6.31198463 6.31204087
which agrees with our results. Observe that all three square roots are within our enclosures. • Example 3: In this example ω = 3 and 0 0 < x < 1, π 1 < x < 2, q(x) := 2 0 2 < x < 3. It is a simple task to compute the solutions ϕ1 (t, µ) and ϕ2 (t, µ) explicitly, e and so it follows that ∆(µ) is also obtained in a closed form. From the “exact” solution the first five eigenvalues are simple. The interpolation gives satisfactory results as can be seen from the following table. N 6 11 21 Exact
µ1 0.6629613387 0.6629773389 0.6629803405 0.6629808831
µ2 2.162705929 2.162572186 2.162548286 2.162544034
µ3 2.277446491 2.277814772 2.277880363 2.277892016
µ4 4.241533696 4.240002625 4.239796746 4.239762398
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
µ5 4.261003349 4.263173280 4.263474201 4.263524610
1066
AMIN BOUMENIR
6. Conclusion With a simple analysis, and with values of solutions of initial value problems computed at a few values of the eigenparameter, we have computed the eigenvalues of a periodic problem with a certain estimated error. The main advantage is the recovery of an approximation of the Hill discriminant whose graph shows possible double eigenvalues. By increasing the number of sampling points, the truncation error is reduced and this allows us to zoom in on very close eigenvalues. Observe that the eigenvalues in Example 2, found by the code Sleign2, are within our enclosures. One of the main features of our scheme is the estimate on the truncation error, TN (µ) and its use for the eigenvalues enclosures. It is also possible for TN (µ) to include the phase and amplitude errors made on the sampling values. Acknowledgment The author acknowledges the support of K.F.U.P.M. References [1] P. B. Bailey, W. N. Everitt and A. Zettl, Regular and singular Sturm-Liouville problem, with coupled boundary conditions, Proc. Roy. Soc. Edinburgh Sect. A, 126, 1995, 505-514. MR 97c:34037 1 D 2 , Proc. Roy. Soc. Edinburgh Sect. A. 125, [2] A. Boumenir, Transmutation operator for w(x) 1995, 85-98. MR 95m:34047 [3] A. Boumenir and B. Chanane, Eigenvalues of S-L Systems using Sampling Theory, Applicable Analysis, Vol. 62, 1996, 323-334. [4] M. S. P. Eastham, The spectral theory of periodic differential equations, Scottish Academic Press, 1973. [5] C. T. Fulton and S. A. Pruess, Eigenvalues and eigenfunctions asymptotics for regular SturmLiouville problems, J. Math. Anal. Appl, 188, 1994, 297-340. MR 97f:34071 [6] I. S. Kac and M. G. Krein, Spectral function of the string, Amer. Math. Soc. Transl. (2), 103, 1970, 19-103. [7] B. M. Levitan, Inverse Sturm-Liouville Problems, VNU Science Press, The Netherlands, 1987. MR 89b:34001 [8] M. A. Naimark, Linear Differential Operators, Part II, Ungar, 1968. MR 41:7485 [9] F. Stenger, Approximations via Whittaker’s cardinal functions, J. Approx. Theory, 17, 1976, 222-240. MR 58:1885 [10] E. C. Titchmarsh, Eigenfunctions expansion, Part II, Oxford University Press, 1961. [11] A. I. Zayed, Advances in Shannon’s Sampling Theory, CRC Press, 1993. MR 95f:94008 Department of Mathematics, College of Sciences, Sultan Qaboos University, P.O. Box 36, Alkhod 123, Muscat, Sultanate of Oman E-mail address: [email protected]
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS BY THE SHANNON-WHITTAKER SAMPLING THEOREM AMIN BOUMENIR
Abstract. We are concerned with the computation of eigenvalues of a periodic Sturm-Liouville problem using interpolation techniques in Paley-Wiener spaces. We shall approximate the Hill discriminant by sampling a few of its values and then find its zeroes which are the square roots of the eigenvalues. Computable error estimates are provided together with eigenvalue enclosures.
1. Introduction We would like to introduce a new method for computing the eigenvalues of classical regular Sturm-Liouville problems with periodic boundary conditions t ∈ [0, ω] −ϕ00 (t, µ) + Q(t)ϕ(t, µ) = µ2 ϕ(t, µ) (1.1) ϕ(0, µ) = ϕ(ω, µ) and ϕ0 (0, µ) = ϕ0 (ω, µ), where Q(t) ∈ L1 (0, ω). For the spectral theory of periodic differential equations, we shall refer to [4], [1] and [8]. Recall that in general the spectrum of (1.1) may not be simple as in the case of separated boundary conditions, and this is a major difficulty. In this work, we shall extend the method in [3], which relies on the interpolation in Paley-Wiener spaces of a certain boundary function, whose zeros are square roots of eigenvalues. In our case the boundary function turns out to be the wellknown Hill discriminant (see [4]) and only few values are needed for its approximation. The truncation error can be minimized by increasing the number of sampling points, which gives higher accuracy on the numerical side and provides eigenvalue enclosures. We shall examine a few examples, where eigenvalue enclosures and comparisons of the results with the well-known code Sleign2 are provided. In the last example our interpolation scheme is compared with an “exact” solution. One of the features of the method is the possibility of locating double eigenvalues. Indeed, by minimizing the truncation error, we can zoom in on close eigenvalues, and if they happen to be simple, then we can find two disjoint enclosures separating them. 2. Hill’s discriminant The Hill discriminant, which is at the heart of the spectral theory of the periodic Sturm-Liouville operator, is defined by ∆(µ) := ϕ1 (ω, µ) + ϕ02 (ω, µ), Received by the editor June 9, 1997 and, in revised form, October 28, 1997. 1991 Mathematics Subject Classification. Primary 34L15, 42A15. Key words and phrases. Periodic Sturm-Liouville problems, eigenvalue problem, interpolation, Shannon-Whittaker sampling theorem.
1057
c
1999 American Mathematical Society
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1058
AMIN BOUMENIR
where ϕ1 (t, µ) and ϕ2 (t, µ) are two independent solutions of (1.1) satisfying ϕ1 (0, µ) = 1 ϕ2 (0, µ) = 0 (2.1) and ϕ01 (0, µ) = 0 ϕ02 (0, µ) = 1. Recall that µ2 is an eigenvalue associated with (1.1) if and only if µ is solution of (2.2)
∆(µ) = 2
(see [4]). In order to interpolate the Hill discriminant we shall need to find a suitable representation in terms of analytic functions. From the inverse spectral theory (see [7]) ϕ1 (t, µ) and ϕ2 (t, µ) can be represented by Z t ϕ1 (t, µ) = cos tµ + (2.3) K1 (t, η) cos ηµ dη, 0 Z t sin tµ sin ηµ + dη, K2 (t, η) ϕ2 (t, µ) = µ µ 0 where K1 and K2 are the kernels of transformation operators. We recall that K1 and K2 have m + 1 locally integrable derivatives if Q has m locally integrable ∂ ∂ K2 and ∂t K1 are locally derivatives. In our case m = 0, which means that ∂t integrable and K1 and K2 are continuous functions. Thus the integrals in (2.3) are well defined. The Lyapunov function or Hill’s discriminant for (1.1) can now be described using K1 and K2 : (2.4)
∆(µ)
= =
ϕ1 (ω, µ) + ϕ02 (ω, µ) sin ωµ sin ωµ + K2 (ω, ω) 2 cos ωµ + K1 (ω, ω) µ µ Z ω ∂ ∂ sin ηµ K2 (ω, η) + K1 (ω, η) dη. + ∂t ∂t µ 0
Let us agree to denote the last integral term by Z ω ∂ ∂ sin ηµ (2.5) K2 (ω, η) + K1 (ω, η) dη, S(µ) := ∂t ∂t µ 0 and several useful constants by sin z exp (|Imz|) z ≤ c 1 + |z| , Z ω q := |Q(η)| dη 0
and
Z ξ := exp c
0
ω
η |Q(η)| dη .
We now recall that Paley-Wiener spaces are defined by Z P Wω := F (µ) entire: |F (µ)| < M eω|Imµ| and
∞
−∞
Also (see [7]) 1 K1 (ω, ω) = K2 (ω, ω) = 2
Z 0
|F (µ)| dµ < ∞ .
ω
Q(η)dη.
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
2
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS
Thus in order to recover (2.6)
∆(µ) := 2 cos(ωµ) +
sin ωµ µ
Z
1059
ω
Q(η)dη + S(µ)
0
we need to find S(µ) only. Proposition 2.1. If Q(t) ∈ L1 (0, ω), then sZ r ∞ 3 2 ω. S(µ) ∈ P Wω and E1 := |µS(µ)| dµ ≤ ξc2 q 2 2 −∞ ∂ ∂ K2 (ω, η) + ∂t K1 (ω, η) Proof. Because Q(t) is locally integrable, it follows that ∂t is also locally integrable, and applying the Riemann-Lebesgue theorem to (2.5) we deduce that S(µ) = o( µ1 ), i.e., S(µ) ∈ L2dµ (−∞, ∞). Next observe that Z ω ∂ cη exp (η |Imµ|) ∂ dη |S(µ)| ≤ ∂t K2 (ω, η) + ∂t K1 (ω, η) 1 + |ηµ| 0 Z ω ∂ cη K2 (ω, η) + ∂ K1 (ω, η) ≤ dη exp (ω |Imµ|) ∂t ∂t 1 + |ηµ| 0 ≤ C1 exp (ω |Imµ|) . Hence S(µ) ∈ P Wω . If C1 is known, then we can estimate the constant E1 , which is crucial in the computation of error bounds. To this end we shall use the following integral equations (see [10]): ϕ1 (t, µ)
(2.7)
= cos(tµ) + Tµ ϕ1 (t, µ), sin (tµ) + Tµ ϕ2 (t, µ), = µ
ϕ2 (t, µ) where
Tµ
Z
C(0, ω) −→ C(0, ω),
sin(µ(t − η)) Q(η)f (η)dη, t ∈ (0, ω), µ 0 Z ω c(ω − η) ||Tµ || ≤ |Q(η)| dη 1 + |µ| (ω − η) 0 ω ≤ cq . 1 + |µ| ω t
Tµ f (t) := (2.8)
For µ large enough, the method of successive approximation yields (see [10]) X X sin ηµ , ϕ1 (t, µ) = Tµn cos(tµ) and ϕ2 (t, µ) = Tµn µ n≥0
n≥0
and by the Gronwall lemma applied to (2.7) we obtain Z t Z t c η |Q(η)| |ϕ1 (t, µ)| ≤ exp dη ≤ exp c η |Q(η)| dη , 0 1 + |µ| η 0 |ϕ2 (t, µ)|
t exp 1 + |µ| t t . ≤ cξ 1 + |µ| t ≤ c
Z 0
t
cη |Q (η)| dη 1 + |µ| η
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1060
AMIN BOUMENIR
Applying Tµ to (2.7) yields ϕ1 (t, µ) = cos(tµ) + Tµ cos(tµ) + Tµ2 ϕ1 (t, µ),
(2.9) where
2 T ϕ1 (ω, µ) µ
≤ ≤ ≤
(2.10)
Z
2
(c qω)
2
(1 + |µ| ω) 2
(c qω)
(1 + |µ| ω)2 2
ξ (c qω)
c η |Q(η)| dη 0 1 + |µ| η Z t exp c η |Q(η)| dη t
exp
0
1 2
(1 + |µ| ω)
.
Similarly we have ϕ2 (t, µ) = ϕ02 (t, µ) =
(2.11)
where d 2 Tµ ϕ2 (ω, µ) dt
= ≤ ≤ ≤ ≤ ≤ ≤
(2.12)
sin(tµ) sin(tµ) + Tµ + Tµ2 ϕ2 (t, µ), µ µ d d sin(tµ) + Tµ2 ϕ2 (t, µ), cos(tµ) + Tµ dt µ dt
Z ω Z t sin(µ(t − η)) Q(η)ϕ2 (η, µ)dηdt cos(µ(ω − t))Q(t) µ 0 0 Z ω Z t c (t − η) |Q(η)| |ϕ2 (η, µ)| dηdt |Q(t)| 0 0 1 + |µ| (t − η) Z t Z ω η cη |Q(η)| ξc dηdt |Q(t)| 1 + |µ| η 0 0 1 + |µ| η 2 Z ω Z t η 2 ξc |Q(t)| |Q(η)| dηdt 1 + |µ| η 0 0 2 Z ω Z t ω ξc2 |Q(t)| |Q(η)| dηdt 1 + |µ| ω 0 0 2 Z ω 2 ξc2 ω |Q(t)| dt 2 1 + |µ| ω 0 1 2 2 2 1 ξc q ω 2. 2 (1 + |µ| ω)
We now use (2.9) and (2.11) to obtain ∆(µ) = =
ϕ1 (ω, µ) + ϕ02 (ω, µ) 2 cos(tµ) + Tµ cos(tµ) +
d d sin(tµ) Tµ + Tµ2 ϕ1 (t, µ) + Tµ2 ϕ2 (t, µ). dt µ dt
Observe that d sin(tµ) Tµ cos(tµ) + Tµ dt µ Z t sin(ηµ) sin(µ(t − η)) cos(ηµ) + cos(µ(t − η)) Q(η)dη = µ µ 0 Z t sin(µt) = Q(η)dη. µ 0
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS
1061
Thus a new representation of the Hill discriminant follows: (2.13) sin(µt) ∆(µ) = 2 cos(tµ) + µ
Z
t
0
Q(η)dη + Tµ2 ϕ1 (t, µ) +
d 2 T ϕ2 (t, µ). dt µ
Therefore (2.6) together with (2.13) imply (2.14)
S(µ) = Tµ2 ϕ1 (t, µ) +
d 2 T (ϕ2 (t, µ)) dt µ
and estimates (2.12) and (2.10) yield (2.15)
|S(µ)| ≤
3ξc2 q 2 ω 2 1 2 . 2 (1 + |µ| ω)
Hence it follows from (2.14) and (2.15) that S(µ) ∈ P Wω and sZ r r 3ξc2 2 2 2 3 2 2 2 q ω ω. |µS(µ)| dµ ≤ ≤ ξc q E1 = 2 3ω 3 2 It is also possible to use asymptotics of solutions (see [5]) to derive the estimates above. 3. General case We now address the question of how to apply the idea above to more general equations, such as 0 0 < x < a, − (p(x)y 0 (x, µ)) + q(x)y(x, µ) = µ2 s(x)y(x, µ), (3.1) y(0, µ) = y(a, µ), p(0)y 0 (0, µ) = p(a)y 0 (a, µ), 1 where p(x), q(x), and s(x) are real valued, p(x) , q(x), s(x) ∈ L1 (0, a) and s(x) ≥ 0. It is well known that µ2n (see [4]) is an eigenvalue of (3.1) if and only if µn is a root of
(3.2)
D(µ) := ψ1 (a, µ) + p(a)ψ20 (a, µ) = 2,
where ψ1 (x, µ) and ψ2 (x, µ) are two independent solutions of (p(x)y 0 (x, µ))0 + q(x)y(x, µ) = µ2 s(x)y(x, µ), with initial conditions ψ1 (0, µ) = 1 p(0)ψ10 (0, µ) = 0
and
0 < x < a,
ψ2 (0, µ) = 0 p(0)ψ20 (0, µ) = 1.
It is readily seen that D(µ) is an entire function of µ of order one, and the main difficulty is estimating its type. Once the type is obtained, which depends on p(x), q(x) and s(x), we can apply the above method with a few modifications. One simple way to deal with (3.1) is to transform it to the standard form (1.1) 1 and s(x) are C 2 (0, a), and s(x) ≥ δ > 0. which is possible when the coefficients p(x) Indeed, the Liouville-Green transformation allows us to recast equation (3.1) into an equivalent periodic problem (see [4]) −ϕ00 (t, µ) + Q(t)ϕ(t, µ) = µ2 ϕ(t, µ) t ∈ [0, A] (3.3) ϕ(0, µ) = ϕ(A, µ) and ϕ0 (0, µ) = Kϕ0 (A, µ),
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1062
AMIN BOUMENIR
where t(x) =
R x q s(η) 0
p(η) dη
for x ∈ [0, a] and
−1 d d p(x) {p(x)s(x)} 4 . dx dx Clearly the type of D(µ) is then A := t(a) and the boundary conditions q q s(x) s(ω)p(ω) 0 0 are obtained from y (x, µ) := ϕ (t(x), µ) p(x) , and so K := s(0)p(0) . In case s(ω)p (ω) =s(0)p (0) , then we have K = 1. • If p(x) and s(x) are not smooth, then we can still use ideas from the transformation operators for the string (see [2]) 1
3
ϕ(t(x), µ) = y(x, µ), Q(t(x)) = q(x) − p 4 (x)s− 4 (x)
−p(x)(p(x)y 0 (x, µ))0 + p(x)q(x)y(x, µ) = µ2 p(x)s(x)y(x, µ) 0 < x < a. Rx 1 Then setting t(x) := 0 p(η) dη, γ := t(a) and ϕ(t(x)) = y(x) we have (3.4)
−ϕ00 (t, µ) + Q(t)ϕ(t, µ) = µ2 ω(t)ϕ(t, µ)
0 0. N
µ∈IµN ∗
Proposition 4.1. Let µ > 0, be a simple eigenvalue of (1.1), then ∀ε > 0, ∃N0 such that ∀N ≥ N0 , ∃µN such that ε . |µ∗ − µN | < e0 inf ∆ (µ) µ∈IµN
N
e0 Proof. Given µ∗ , we can find N, such that µ∗ < Nωπ , inf ∆ N (µ) ≥ δ > 0 and µ∈IµN e ∗ ∗ e TN (µ) < ε for µ ∈ IµN . From ∆ (µ ) − ∆ (µ ) N N N ≤ TN (µ ), the mean value theorem implies ∗ e 0 (η) ≤ TN (µ∗ ), (4.3) η ∈ (µ∗ , µN ) ⊂ IµN , (µ − µN ) ∆ N e0 and since µ∗ is simple, inf ∆ N (µ) > 0, and we have µ∈IµN
|µ∗ − µN | ≤
TN (µ∗ ) ≤ e0 inf ∆ N (µ)
µ∈IµN
ε e0 inf ∆N (µ) .
µ∈IµN
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
1064
AMIN BOUMENIR
As we increase N ≥ N0 , ∗
sup TN (µ) → 0,
µ∈IµN
e0 inf ∆ N (µ) → δ(N0 ), and so
µ∈IµN
|µ − µN | → 0 as N → ∞. • The second case is when we have a double eigenvalue, which is characterized by e (µ∗ ) = 0, ∆ e 0 (µ∗ ) = 0 and ∆ e 00 (µ∗ ) 6= 0. ∆ e Because of the concavity of the ∆(µ) at µ∗ , there are two cases. e – ∆N (µ) is also concave and may not have zeroes. In this case we look for solutions inside the intervals defined by (4.2). e N (µN ) has two zeroes then – If ∆ e N (µN ) + (µ∗ − µN )∆ e 0N (µN ) + 1 (µ∗ − µN )2 ∆ e 00N (η), e N (µ∗ ) = ∆ ∆ 2 e ∗ ) + (µ∗ − µN )∆ e 0N (µN ) + 1 (µ∗ − µN )2 ∆ e 00N (η), e N (µ∗ ) = ∆(µ ∆ 2 e ∗ ) = (µ∗ − µN )∆ e 0 (µN ) + 1 (µ∗ − µN )2 ∆ e 00 (η), e N (µ∗ ) − ∆(µ ∆ N N 2 h i 1 ∗ 2 e 00 ∗ 0 ∗ e (µN ) − ∆ e N (µ ) − ∆(µ e ∗ ) = 0, (µ − µN ) ∆N (η) + (µ − µN )∆ (4.4) N 2 i h e ∗ ) ≤ TN (µ∗ ). Clearly the coefficients of this quae N (µ∗ ) − ∆(µ where ∆ e 00 (η) ≤ max ∆ e 00 (η) ≤ ∆ e 00 (η), and simidratic are intervals, i.e., inf ∆ η∈Iµ
N
N
η∈Iµ
N
larly 0 ≤ TN (µ∗ ) ≤ max TN (η). Therefore we can find an upper bound η∈Iµ
e N (µN ) may have no roots for |µ∗ − µN | , using interval analysis. Since ∆ in the neighbourhood of µ∗ , the only practical way to locate double eigenvalues is the interval defined by (4.2). In practice it is more convenient to be given N and then compute the truncation error TN (µ) to find the intervals IµN defined by (4.2). 5. Examples We shall consider problems defined by (1.1), where Q(x) is integrable. Recall that N represents the number of points used in the Shannon-Whittaker sampling theorem. The positive eigenvalues are obtained by recovering the Hill discriminant from its sampled values. To this end we only need to compute the values ϕ1 (ω, µ) + ϕ02 (ω, µ) for µ = k ωπ , k = 0, 1, ..., N , which are obtained by solving the initial value problem (1.1) with Runge-Kutta 4-5, which is implemented with the symbolic manipulator Maple with precision set to 13 digits. • Example 1: Let us choose ω = 2, and 0 < x < 1/2 , 0 1 1/2 < x < 7/4 , q(x) := 2 1+x 7/4 < x < 2. 0 The square roots of the eigenvalues are listed below. N 5 13 19
µ1 0.5376360812 0.5376389426 0.5376389987
µ2 3.141592653 3.141592653 3.141592654
µ3 3.198533821 3.198655062 3.198665145
µ4 6.283185307 6.283185307 6.283185307
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
µ5 6.309309048 6.309851843 6.309892873
EIGENVALUES OF PERIODIC STURM-LIOUVILLE PROBLEMS
1065
The enclosure for the first simple eigenvalue µ21 is given below. N 5 13 19
µmin 0.5063 0.5292 0.5327
µ1 0.5376360812 0.5376389426 0.5376389987
µmax 0.5700 0.5461 0.5425
• Example 2: q(x) := 2x(1 − x) and ω = 1. The enclosures are given below. N 4 7 15 20 N 4 7 15 20
µ1 min 0.576496 0.575138 0.576292 0.576496
µ2,3 min 6.29811294928 6.30116798050 6.30414770889 6.30496865578
µ1 0.576890580 0.576891292 0.576891400 0.576891408
µ2 6.30777169387 6.30769833630 6.30768224351 6.30768111508
µ1 max 0.577287 0.578648 0.577491 0.577287
µ3 6.31155085736 6.31163581017 6.31165451250 6.31165582479
µ2,3 max 6.32991780134 6.32196796225 6.31645456422 6.31520080034
For N = 20, an estimate of the truncation error is given by T N (µ) < 27.10−6 when µ ∈ (6.304968, 6.3152008) T N (µ) < 46.10−5
µ1 ∈ (0.5764, 0.5774).
1. A quick run of Sleign2, with a tolerance of 10−7 , gives the following eigenvalues Sleign2
µ=
√ λ
0.576919405 6.31198463 6.31204087
which agrees with our results. Observe that all three square roots are within our enclosures. • Example 3: In this example ω = 3 and 0 0 < x < 1, π 1 < x < 2, q(x) := 2 0 2 < x < 3. It is a simple task to compute the solutions ϕ1 (t, µ) and ϕ2 (t, µ) explicitly, e and so it follows that ∆(µ) is also obtained in a closed form. From the “exact” solution the first five eigenvalues are simple. The interpolation gives satisfactory results as can be seen from the following table. N 6 11 21 Exact
µ1 0.6629613387 0.6629773389 0.6629803405 0.6629808831
µ2 2.162705929 2.162572186 2.162548286 2.162544034
µ3 2.277446491 2.277814772 2.277880363 2.277892016
µ4 4.241533696 4.240002625 4.239796746 4.239762398
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
µ5 4.261003349 4.263173280 4.263474201 4.263524610
1066
AMIN BOUMENIR
6. Conclusion With a simple analysis, and with values of solutions of initial value problems computed at a few values of the eigenparameter, we have computed the eigenvalues of a periodic problem with a certain estimated error. The main advantage is the recovery of an approximation of the Hill discriminant whose graph shows possible double eigenvalues. By increasing the number of sampling points, the truncation error is reduced and this allows us to zoom in on very close eigenvalues. Observe that the eigenvalues in Example 2, found by the code Sleign2, are within our enclosures. One of the main features of our scheme is the estimate on the truncation error, TN (µ) and its use for the eigenvalues enclosures. It is also possible for TN (µ) to include the phase and amplitude errors made on the sampling values. Acknowledgment The author acknowledges the support of K.F.U.P.M. References [1] P. B. Bailey, W. N. Everitt and A. Zettl, Regular and singular Sturm-Liouville problem, with coupled boundary conditions, Proc. Roy. Soc. Edinburgh Sect. A, 126, 1995, 505-514. MR 97c:34037 1 D 2 , Proc. Roy. Soc. Edinburgh Sect. A. 125, [2] A. Boumenir, Transmutation operator for w(x) 1995, 85-98. MR 95m:34047 [3] A. Boumenir and B. Chanane, Eigenvalues of S-L Systems using Sampling Theory, Applicable Analysis, Vol. 62, 1996, 323-334. [4] M. S. P. Eastham, The spectral theory of periodic differential equations, Scottish Academic Press, 1973. [5] C. T. Fulton and S. A. Pruess, Eigenvalues and eigenfunctions asymptotics for regular SturmLiouville problems, J. Math. Anal. Appl, 188, 1994, 297-340. MR 97f:34071 [6] I. S. Kac and M. G. Krein, Spectral function of the string, Amer. Math. Soc. Transl. (2), 103, 1970, 19-103. [7] B. M. Levitan, Inverse Sturm-Liouville Problems, VNU Science Press, The Netherlands, 1987. MR 89b:34001 [8] M. A. Naimark, Linear Differential Operators, Part II, Ungar, 1968. MR 41:7485 [9] F. Stenger, Approximations via Whittaker’s cardinal functions, J. Approx. Theory, 17, 1976, 222-240. MR 58:1885 [10] E. C. Titchmarsh, Eigenfunctions expansion, Part II, Oxford University Press, 1961. [11] A. I. Zayed, Advances in Shannon’s Sampling Theory, CRC Press, 1993. MR 95f:94008 Department of Mathematics, College of Sciences, Sultan Qaboos University, P.O. Box 36, Alkhod 123, Muscat, Sultanate of Oman E-mail address: [email protected]
License or copyright restrictions may apply to redistribution; see http://www.ams.org/journal-terms-of-use
