Accurate Estimations of Any Eigenpairs of N-th Order Linear Boundary Value Problems

Abstract

This paper provides a method to bound and calculate any eigenvalues and eigenfunctions of n-th order boundary value problems with sign-regular kernels subject to two-point boundary conditions. The method is based on the selection of a particular type of cone for each eigenpair to be determined, the recursive application of the operator associated to the equivalent integral problem to functions belonging to such a cone, and the calculation of the Collatz–Wielandt numbers of the resulting functions.

1. Introduction

Let J be a compact interval in $\mathbb{R}$ and let us consider the real differential operator L disconjugate on J and defined by:

$$Ly=y^{\left(n\right)}\left(x\right)+a_{n-1}\left(x\right)y^{(n-1)}\left(x\right)+\cdots +a_0\left(x\right)y\left(x\right),x\in J,$$

(1)

where $a_j\left(x\right)\in C^j\left(J\right)$ for $0\le j\le n-1$.

In this paper, we will address the eigenvalue problem:

$$Ly=\lambda r\left(x\right)y,x\in [a,b]\subseteq J,$$

(2)

subject to the boundary conditions

$$\begin{array}{c}{U}_i\left(y\right)\equiv y^{\left(k_i\right)}\left(a\right)+{\sum }_{k<k_i}{\gamma }_{ik}{y}^{\left(k\right)}\left(a\right)=0(i=1,\dots ,m),\\ \\ U_i\left(y\right)\equiv y^{\left(k_i\right)}\left(b\right)+{\sum }_{k<k_i}{\gamma }_{ik}{y}^{\left(k\right)}\left(b\right)=0(i=m+1,\dots ,n),\end{array}$$

(3)

where $1\le m\le n-1$, $n-1\ge k_1>\cdots >k_m\ge 0$, $n-1\ge k_{m+1}>\cdots >k_n\ge 0$, and $r\left(x\right)$ is a positive function piecewise continuous on $[a,b]$.

If the problem $Ly=0$ subject to the boundary conditions (3) does not have a non-trivial solution (namely, if $\lambda =0$ is not an eigenvalue of (2) and (3)), then the problem (2) and (3) is equivalent to the integral eigenvalue problem:

$$y=\lambda My,$$

(4)

where M is the operator $C[a,b]\to A{C}^{n-1}[a,b]$ defined by:

$$Mu={\int }_a^{b}G(x,t)r\left(t\right)u\left(t\right)dt,x\in [a,b],$$

(5)

and $G(x,t)$ is the Green function of the problem,

$$\begin{array}{c}LG(x,t)=0,x\in (a,t)\cup (t,b),\\ \\ U_i\left(G(x,t)\right)\equiv \frac{{\partial }^{k_i}G(a,t)}{\partial x^{k_i}}+{\sum }_{k<k_i}{\gamma }_{ik}\frac{{\partial }^{k}G(a,t)}{\partial x^k}=0(i=1,\dots ,m),\\ \\ U_i\left(G(x,t)\right)\equiv \frac{{\partial }^{k_i}G(b,t)}{\partial x^{k_i}}+{\sum }_{k<k_i}{\gamma }_{ik}\frac{{\partial }^{k}G(b,t)}{\partial x^k}=0(i=m+1,\dots ,n).\end{array}$$

(6)

Let us introduce the naming convention:

$$G\left(\begin{array}{c}{x}_1,x_2,\dots ,x_p\\ t_1,t_2,\dots ,t_p\end{array}\right)=\left|\begin{array}{ccc}G(x_1,t_1)& \dots & G(x_1,t_p)\\ \cdots & \cdots & \cdots \\ G(x_p,t_1)& \dots & G(x_p,t_p)\end{array}\right|,$$

where $|·|$ is the determinant of the matrix.

According to [1,2], $G(x,t)$ is called a sign-regular kernel of the integral Equations (4) and (5) if the following condition is satisfied for a sequence of real numbers ${ϵ}_1,{ϵ}_2,\dots$, each one taking either the value 1 or the value $-1$, for $p=1,2,\dots$:

$${ϵ}_{p}G\left(\begin{array}{c}{x}_1,x_2,\dots ,x_p\\ t_1,t_2,\dots ,t_p\end{array}\right)\ge 0\left(a<\begin{array}{c}{x}_1<x_2<\cdots <x_p\\ t_1<t_2<\cdots <t_p\end{array}<b\right).$$

(7)

A sign-regular kernel is called strongly sign-regular if, in addition, it satisfies the following two conditions:

$${ϵ}_{1}G(t,s)>0(a<t<b,a<s<b),$$

(8)

$${ϵ}_{p}G\left(\begin{array}{c}{x}_1,x_2,\dots ,x_p\\ x_1,x_2,\dots ,x_p\end{array}\right)>0(a<x_1<x_2<\cdots <x_p<b,p=1,2,\dots ).$$

(9)

If all the ${ϵ}_p$ (for $p=1,2,\dots$) are equal to $+1$, then a sign-regular kernel is called totally non-negative, whereas a strongly sign-regular kernel is called oscillatory. For the Green functions of two-point boundary problems like that of (6), sign-regularity is equivalent to strong sign-regularity (see [1]), and the condition (7) needs only to hold for $p=1,\dots ,n-1$ (see [3], Condition A).

Throughout the paper, we will assume that the Green function of (6) is sign-regular.

The interest of sign-regular Green functions resides on the Sturm–Liouville-like properties that it gives to the boundary value problem (2) and (3), namely (see [3]):

• The problem (2) and (3) has countable infinitely many different eigenvalues ${\lambda }_1,{\lambda }_2,\dots$, which are real, algebraically and geometrically simple, and can be ordered as $0<{ϵ}_0{ϵ}_1{\lambda }_1<{ϵ}_1{ϵ}_2{\lambda }_2<\cdots$, where ${ϵ}_0=1$ and ${ϵ}_1,{ϵ}_2,\dots$ are the same numbers as in (7).
• The eigenfunction ${\phi }_i\left(x\right)$ (i = 1, 2, …), corresponding to each ${\lambda }_i$, has exactly $i-1$ zeroes in $(a,b)$, all of which are simple. Moreover, the zeroes of ${\phi }_i$ and ${\phi }_{i+1}$ alternate (i = 2, 3,…). At the extremes a, b, all the eigenfunctions ${\phi }_i$ have zeroes of the order exactly imposed by the boundary conditions.
• Each non-trivial linear combination $c_r{\phi }_r+\cdots +c_l{\phi }_l$ with $r\le l$ has at least $r-1$ nodal zeroes (that is, zeroes where the function changes its sign) and at most $l-1$ zeroes in $I_1$, where $I_1$ is the interval obtained from $[a,b]$ by removing a if $k_m=0$ and b if $k_n=0$, and the zeroes which are antinodes (that is, zeroes where the function does not change its sign) are counted twice.

There is a fourth property of the eigenfunctions ${\phi }_i$, which requires the introduction of the following definition (see [4], Chapter 3, Section 5):

Definition 1.

A system of continuous functions,

$$y_1\left(x\right),y_2\left(x\right),\dots ,y_p\left(x\right),x\in I,$$

is called a Chebyshev system on the interval I if every linear combination of these functions

$$y\left(x\right)=\sum _{i=1}^p{c}_i{y}_i\left(x\right)\left(\sum _{i=1}^p{c}_i^2>0\right)$$

vanishes on the interval I at most $p-1$ times. Likewise, a sequence of (finite or infinite) functions

$$y_1\left(x\right),y_2\left(x\right),y_3\left(x\right),\dots$$

is a Markov sequence within the interval I if for every p ($p=1,2,\dots$), the functions $y_1\left(x\right),y_2\left(x\right),$ $\dots ,y_p\left(x\right)$ form a Chebyshev system on I.

According to [3], the eigenfunctions ${\phi }_1,{\phi }_2,\dots$ of the problem (2) form a Markov sequence on $[a,b]$ within $I_1$.

The problem (2)–(7) appears in the analysis of the vibrations of a loaded continuum, but the associated theory can be applied to multiple differential problems as long as the boundary conditions imply the sign-regularity of the Green function $G(x,t)$. In particular, in [5] (Appendix D), one can find multiple examples of differential problems of the type (2) in the theory of fluid dynamics, problems whose Green function satisfies (7) under some boundary conditions. Likewise, [6] contains many examples of physical and biological problems with sign regular kernels satisfying (7)–(9).

The eigenvalue problem (2) with a sign-regular kernel has been studied thoroughly in the literature. It was Kellogg [7,8,9] who first assessed symmetric totally non-negative kernels satisfying condition (7). The non-symmetric case was developed by Gantmakher and Krein in [10,11,12]. Karlin obtained new results by attacking the problem from the theory of spline interpolation and Chebyshev and Markov systems [2,13]. Other important breakthroughs were achieved by Levin and Stepanov [1], who extended the results to the sign-regular case, and Borovskikh and Pokornyi [14], who applied them to discontinuous kernels. Later, Stepanov provided necessary and sufficient conditions for the Green function of (6) to be sign-regular in [3]. The research on this topic was continued by Pokornyi and his collaborators due to its relationship with the theory of differential equations in networks [6,15]. Some more recent contributions include [16,17].

While the aforementioned papers cover multiple aspects of the theory of sign-regular kernels and the properties of the solutions of (2), none of them seems to have attempted to use them to calculate its eigenvalues and eigenfunctions, as far as the authors are aware. That will be the purpose of this paper, which will provide an iterative procedure to:

• Bound and calculate any eigenvalues ${\lambda }_i$ and
• Calculate the associated eigenfunctions ${\phi }_i$,

with as much precision as desired.

Our approach will make use of Krein–Rutman cone theory, which was also employed by many of the papers mentioned before, in the following manner:

• Defining a Banach space and a cone of functions, and picking up a function u which belongs to it. Concretely, there will be a Banach space and cone for each eigenvalue ${\lambda }_i$ to be determined.
• Calculating $M^{j}u$ iteratively, where $M^j$ is the composition of M with itself $j-1$ times, $\underset{\underset{jtimesM}{︸}}{M\circ M\circ \cdots \circ M}$.
• Calculating the so-called Collatz–Wielandt numbers of $M^{j}u$ in that cone, for different values of j. These numbers are bounds for the inverse of the eigenvalue ${\lambda }_i$, and converge to this as the iteration index j grows.
• Determining the eigenfunctions ${\phi }_i$ from $M^{j}u$.

The procedure requires a sequential calculation of the eigenpairs, that is, in order to calculate ${\lambda }_i$ and ${\phi }_i$, one has to run the process for the eigenpairs associated with eigenvalues of smaller absolute value.

For self-completeness, let us recall that, given a Banach space B, a cone $P\subset B$ is a non-empty closed set defined by the conditions:

• If $u,v\in P$, then $cu+dv\in P$ for any real numbers $c,d\ge 0$. Note that this condition implies that $0\in P$.
• If $u\in P$ and $-u\in P$, then $u=0$.

A cone P is reproducing if any $y\in B$ can be expressed as $y=u-v$ with $u,v\in P$. The existence of a cone in a Banach space B allows the definition of a partial ordering relationship in that Banach space by setting $u\le v$ if and only if $v-u\in P$. Thus, we will say that the operator M is $u_0$-positive if there exists an integer $q>0$ and a $u_0\in P$ such that for any $v\in P\setminus \{0\}$ one can find positive constants ${\delta }_1$, ${\delta }_2$ such that ${\delta }_1{u}_0\le M^{q}v\le {\delta }_2{u}_0$ (note that ${\delta }_1$ and ${\delta }_2$ will not be the same for all v). We will denote by $int\left(P\right)$ the interior of the cone P, provided that it exists. Let us note that if $M^q$ maps $P\setminus \{0\}$ into $int\left(P\right)$, then it is $u_0$-positive, with $u_0$ being any member of $int\left(P\right)$.

Following the Forster–Nagy definition [18], if $u\in P\setminus \{0\}$, the upper and lower Collatz–Wielandt numbers are defined, respectively, as:

$$\overline{r}(M,u)=inf\{w\in \mathbb{R}:Mu\le wu\},\underline{r}(M,u)=sup\{w\in \mathbb{R}:wu\le Mu\}.$$

(10)

They are called upper and lower Collatz–Wielandt numbers as they extend the estimates for the spectral radius of a non-negative matrix given by L. Collatz [19] and H. Wielandt [20]. We will write them too as $\overline{r}(M,u,P)$ and $\underline{r}(M,u,P)$ when we want to make an explicit reference to the concrete cone P in which they are calculated.

The properties of $\underline{r}(M,u)$ and $\overline{r}(M,u)$ and their relationship with the spectral radius of the operator M have been studied by several authors, starting with Marek [21,22], Forster and Nagy [18], who corrected some previous mistakes from Marek, and Marek again [23]. The concept has been extended to multiple types of operators, Banach spaces and cones. The references [24,25,26,27], include a good account of recent results.

The use of $\underline{r}(M,M^{j}u)$ and $\overline{r}(M,M^{j}u)$ to bound and estimate the principal eigenvalue of a boundary value problem of the type (3) and (4) seems to date from Webb [28], who applied it to define conditions for the existence of solutions to non-linear boundary value problems. Chang [25] proved it for other 1-homogeneous non-linear differential problems like p-Laplace systems, calling it the power method, so it is possible that it was known and used before. Later, the authors used it in [29,30,31] to determine the solvability of boundary value problems and in [32] to bound and estimate the principal eigenvalue of boundary value problems including higher derivatives, for which some results on the sign of the derivatives of the Green function were needed. However, to the knowledge of the authors, it has never been applied to determine the value of other eigenpairs apart from the principal one.

This paper took inspiration from an algorithm used to calculate eigenvalues of an oscillatory matrix described in [4] (Appendix 1), which is also based on the use of the Collatz–Wielandt numbers in different cones.

The organization of the paper is as follows. In Section 2, the main results will be presented. In particular, Section 2.1 will introduce the Banach space and the cones to be used for each ${\lambda }_i$, and will prove the convergence of the Collatz–Wielandt numbers of $M^{j}u$ in such cones. The subsection Section 2.2 will yield a procedure to find a function u belonging to each cone, and will show how to simplify the calculation of the Collatz–Wielandt numbers. Section 3 will give an example of how to apply the previous theory to calculate several eigenpairs of a boundary value problem. Finally, Section 4 will draw some conclusions.

2. Main Results

2.1. The Procedure to Calculate ${\lambda }_p$ and ${\phi }_p$

2.1.1. Some Preliminaries

The fact that $a_j\left(x\right)\in C^j[a,b]$ allows definition of the operator $L^{\ast }$ adjoint to L, namely [33] (Chapter 11, Section 1):

$$L^{\ast }y={(-1)}^n{y}^{\left(n\right)}\left(x\right)+{(-1)}^{n-1}{\left(a_{n-1}\left(x\right)y\left(x\right)\right)}^{(n-1)}+\cdots +a_0\left(x\right)y\left(x\right),x\in [a,b].$$

(11)

Accordingly, let us consider the eigenvalue problem adjoint to (2)

$$L^{\ast }y=\tilde{\lambda }r\left(x\right)y,x\in [a,b];U_i^{\ast }\left(y\right)=0,i=1,\dots ,n,$$

(12)

where $U_i^{\ast }$ are the boundary conditions adjoint to $U_i$ (see [33] (Chapter 11, Theorem 3.1) for a definition of adjoint boundary conditions). From the properties of the adjoint eigenvalue problems, it is well known [33] (Chapter 12, Theorem 5.2) that the eigenvalues ${\tilde{\lambda }}_i$ of (12) are the complex conjugates of those of (2) and (3) (that is $\widetilde{{\lambda }_i}={\lambda }_i^{\ast }$), and that its eigenfunctions ${\psi }_i$ form a biorthogonal system with ${\phi }_i$, namely that:

$$〈{\psi }_i,{\phi }_j〉={\int }_a^b{\psi }_i\left(x\right){\phi }_j\left(x\right)r\left(x\right)dx={\delta }_{ij}.$$

(13)

Given that the eigenvalues ${\lambda }_i$ are real, obviously ${\tilde{\lambda }}_i$ are also real and ${\tilde{\lambda }}_i={\lambda }_i$.

Next, given a set of p functions $y_i\left(x\right)\in C[a,b]$ for $i=1,\dots ,p$, let us introduce the notation:

$$\Delta \left(\begin{array}{cccc}{y}_1& y_2& \dots & y_p\\ x_1& x_2& \dots & x_p\end{array}\right)=\left|\begin{array}{ccc}{y}_1\left(x_1\right)& \dots & y_1\left(x_p\right)\\ \cdots & \cdots & \cdots \\ y_p\left(x_1\right)& \dots & y_p\left(x_p\right)\end{array}\right|,(x_1,\dots x_p)\in {[a,b]}^p\subset {\mathbb{R}}^p.$$

(14)

The determinant (14) has very interesting symmetric properties. In particular, its value in ${[a,b]}^p$ is determined by its value in the simplex $\Omega =\{(x_1,\dots ,x_p):a\le x_1<\cdots <x_p\le b\}$, as can be easily shown using the properties of determinants. We will use this simplex frequently in the rest of the paper, together with the related simplex ${\Omega }^{\ast }=\{(x_1,\dots ,x_p)\in \Omega :x_1,\cdots ,x_p\in I_1\}$.

Let $u\in C[a,b]$. We will denote by ${\Delta }_p(u;x_1,\dots ,x_p)$ (or simply ${\Delta }_p\left(u\right)$) the function:

$$\Delta \left(\begin{array}{ccccc}{\phi }_1& {\phi }_2& \dots & {\phi }_{p-1}& u\\ x_1& x_2& \dots & x_{p-1}& x_p\end{array}\right),(x_1,\dots x_p)\in {[a,b]}^p\subset {\mathbb{R}}^p,$$

(15)

where ${\phi }_i$ are the eigenfunctions of (2). As before, the value of ${\Delta }_p\left(u\right)$ is given by its value in the simplex $\Omega$.

Now we are in a position to define the Banach spaces and cones needed by our method. Thus, for each index p, we will define the Banach space $B_p$ as:

$$B_p=\{y\in C[a,b]:〈{\psi }_i,y〉=0,i=1,\dots ,p-1\},$$

(16)

subject to the sup norm $\parallel y\parallel =sup\{|f\left(x\right)|,x\in [a,b]\}$. Since $C[a,b]$ is complete with regards to the sup norm and the functional $〈{\psi }_i,y〉$ is linear and bounded for each ${\psi }_i$, and therefore continuous, it is straightforward to show that $B_p$ is also complete and therefore a proper Banach space.

Likewise, the cone $P_p$ will be defined by:

$$P_p=\{y\in B_p:{\Delta }_p\left(y\right)\ge 0,(x_1,\dots x_p)\in \overline{\Omega }\},$$

(17)

where $\overline{\Omega }$ is the closure of $\Omega$.

Similarly, the Banach space ${\mathcal{B}}_p$ will be defined as:

$$\begin{array}{c}{\mathcal{B}}_p=\{y\in A{C}^{max(n-1,p)}[a,b]:〈{\psi }_i,y〉=0,i=1,\dots ,p-1;\\ \\ U_i\left(y\right)\equiv y^{\left(k_i\right)}\left(a\right)+{\sum }_{k<k_i}{\gamma }_{ik}{y}^{\left(k\right)}\left(a\right)=0,i=1,\dots ,m;\\ \\ U_i\left(y\right)\equiv y^{\left(k_i\right)}\left(b\right)+{\sum }_{k<k_i}{\gamma }_{ik}{y}^{\left(k\right)}\left(b\right)=0,i=m+1,\dots ,n\},\end{array}$$

(18)

where $A{C}^m[a,b]$ is the space of functions whose m-th derivative is absolutely continuous in $[a,b]$. This space will be endowed with the norm $\parallel y\parallel =max\{sup\{|f^{\left(i\right)}\left(x\right)|,x\in [a,b]\},i=0,\dots ,max(n-1,p)\}$, for $p=0,1,\dots$. As before, it is straightforward to show that ${\mathcal{B}}_p$ is a Banach space.

Finally, linked to ${\mathcal{B}}_p$, the cone ${\mathcal{P}}_p$ will be given by:

$${\mathcal{P}}_p=\{y\in {\mathcal{B}}_p:{\Delta }_p\left(y\right)\ge 0,(x_1,\dots x_p)\in \overline{\Omega }\}.$$

(19)

Lemma 1.

The cones $P_p$ and ${\mathcal{P}}_p$ are actually cones.

Proof. 

From property 3 of the eigenfunctions of (2), (13) and the definition of ${\Delta }_p$ in (15), it is clear that either ${\phi }_p$ or $-{\phi }_p$ belong to both cones, so they are not empty. From the properties of the determinants one also has:

$${\Delta }_p(cy+dz)=c{\Delta }_p\left(y\right)+d{\Delta }_p\left(z\right),$$

so if $c,$ $d>0$ and ${\Delta }_p\left(y\right),{\Delta }_p\left(z\right)\ge 0$ for $\mathbf{x}\in \overline{\Omega }$, then obviously ${\Delta }_p(cy+dz)\ge 0$ for $\mathbf{x}\in \overline{\Omega }$. Last, but not least, if ${\Delta }_p\left(y\right)\ge 0$ and ${\Delta }_p(-y)\ge 0$ for $\mathbf{x}\in \overline{\Omega }$, then $\Delta \left(y\right)\equiv 0$ in $\overline{\Omega }$. This is only possible if y is a linear combination of ${\phi }_i$ for $i=1,\dots ,p-1$. As the definition of $B_p$ and ${\mathcal{B}}_p$ requires y to be orthogonal to the adjoint eigenfunctions ${\psi }_i$, $i=1,\dots ,p-1$, (13) leaves $y\equiv 0$ in $\overline{\Omega }$ as the only alternative. This completes the proof. □

2.1.2. The Operator $M_p$ and Its Properties in the Cones

Let us introduce the operator $M_p$, defined by:

$$M_p={ϵ}_{p-1}{ϵ}_{p}M.$$

(20)

The operator $M_p$ has some interesting properties in the cone $P_p$, such as, for instance, its positive character.

Theorem 1.

The operator $M_p$ maps $P_p$ into itself.

Proof. 

From (5), (6) and (20), it is clear that $M_p$ maps $C[a,b]$ into $A{C}^{n-1}[a,b]\subset C[a,b]$, and incidentally the resulting function satisfies the boundary conditions (3). Moreover,

$$〈{\psi }_i,M_{p}u〉={ϵ}_{p-1}{ϵ}_p{\int }_a^b{\psi }_i\left(x\right)r\left(x\right){\int }_a^{b}G(x,t)u\left(t\right)r\left(t\right)dtdx=$$

$${ϵ}_{p-1}{ϵ}_p{\int }_a^{b}u\left(t\right)r\left(t\right){\int }_a^{b}G(x,t){\psi }_i\left(x\right)r\left(x\right)dxdt.$$

The Green function of the adjoint problem (12) is exactly $G^{\ast }(x,t)=G(t,x)$ (see for instance [33], Chapter 11, Theorem 4.2), which yields

$$〈{\psi }_i,M_{p}u〉={ϵ}_{p-1}{ϵ}_p{\int }_a^{b}u\left(t\right)r\left(t\right){\int }_a^b{G}^{\ast }(t,x){\psi }_i\left(x\right)r\left(x\right)dxdt={ϵ}_{p-1}{ϵ}_p{\int }_a^{b}u\left(t\right)r\left(t\right)\frac{{\psi }_i\left(t\right)}{{\lambda }_i}dt=0,$$

for any $u\in B_p$. That implies that $〈{\psi }_i,M_{p}u〉=0$ for $i=1,\dots ,p-1$, and therefore $M_p$ maps $B_p$ into itself.

Next, we will prove that ${\Delta }_p\left(M_{p}u\right)\ge 0$ in $\overline{\Omega }$ when ${\Delta }_p\left(u\right)\ge 0$ in the same set. We will show first that, in fact, ${\Delta }_p\left(M_{p}u\right)\ne 0$ for $\mathbf{x}\in {\Omega }^{\ast }$ in that case. Thus, let us assume that, on the contrary, ${\Delta }_p\left(M_{p}u\right)=0$ at ${\mathbf{x}}^{\prime }\in {\Omega }^{\ast }$, that is, there exist $(x_1^{\prime },\dots ,x_p^{\prime })\in {\Omega }^{\ast }$ such that ${\Delta }_p(M_{p}u;x_1^{\prime },\dots ,x_p^{\prime })=0$. From (14), it follows that there is a linear combination of ${\phi }_1,\dots ,{\phi }_{p-1},M_{p}u$:

$$y\left(x\right)=c_1{\phi }_1\left(x\right)+\cdots +c_{p-1}{\phi }_{p-1}\left(x\right)+c_p{M}_{p}u\left(x\right),$$

(21)

which vanishes at least at $x=x_1^{\prime },\dots ,x_p^{\prime }\in I_1$.

In order for the Green function (6) to be sign-regular, it is necessary that the equation $Lz=0$ is disconjugate on $[a,b]$, that is, no solution of such an equation can have n zeroes in $[a,b]$ (see, for instance, [3], p. 1690). In that case, a result from Polya [34] allows factoring L in first order differential operators as follows:

$$L_{0}y={\rho }_{0}y,$$

$$L_{i}y={\rho }_i{\left(L_{i-1}y\right)}^{\prime },i=1,\dots ,n,$$

(22)

$$andLy=L_{n}y,$$

where ${\rho }_i>0$, ${\rho }_0{\rho }_1\cdots {\rho }_n=1$ and ${\rho }_i\in C^{n-i}[a,b]$.

Let us assume that the number of zeroes of y of (21) is finite in $(a,b)$. Following [3] (Section 5), let $S(d_0,\dots ,d_n)$ denote the number of sign changes in the sequence $d_0,\dots ,d_n$ of non-zero real numbers and let $\sigma \left(f\right)$ be the number of sign changes of $f\left(x\right)$ in $(a,b)$. For the function y, we define the one-sided limit:

$$S(y,x^{-})=\underset{\xi \to x^{-}}{lim}S(L_{0}y\left(\xi \right),L_{1}y\left(\xi \right),\dots ,L_{n}y\left(\xi \right)).$$

(23)

Using (22) and Rolle’s theorem, Stepanov proved [3] (Lemma 5.1) that:

$$\sigma \left(Ly\right)\ge \sigma \left(y\right)+S(y,b^{-})-S(y,a^{-}).$$

(24)

However, if one reviews the proof of such a lemma, one can easily conclude that in fact:

$$\sigma \left(Ly\right)\ge z\left(y\right)+S(y,b^{-})-S(y,a^{-}),$$

(25)

where $z\left(y\right)$ is the number of zeroes of y in $(a,b)$. Stepanov used (24) in several lemmata of the same paper [3] (Lemmata 6.2, 6.3 and 6.5) to prove the sufficiency of his [3] (Theorem 1.3) for the Green function (6) to be sign regular on $[a,b]$, theorem who Stepanov had previously ([3], p. 1713) shown also to be necessary for the sign regularity. Such lemmata essentially proved

$$\sigma \left(Ly\right)\ge \sigma \left(y\right),$$

and one can repeat their same arguments, using inequality (25) instead and mutatis mutandi, to show that

$$\sigma \left(Ly\right)\ge z\left(y\right)\ge p.$$

(26)

Since [3] (Theorem 1.3) is also a necessary condition, any sign regular Green function must fulfil (26) for y having p zeroes in $(a,b)$. Recalling the form of y in (21), then the function:

$$Ly\left(x\right)=r\left(x\right)\left(c_1{\lambda }_1{\phi }_1\left(x\right)+\cdots +c_{p-1}{\lambda }_{p-1}{\phi }_{p-1}\left(x\right)+c_p{ϵ}_p{ϵ}_{p-1}u\left(x\right)\right),$$

(27)

must change its sign at least p times in $(a,b)$, exactly at the same points as

$$v\left(x\right)=c_1{\lambda }_1{\phi }_1\left(x\right)+\cdots +c_{p-1}{\lambda }_{p-1}{\phi }_{p-1}\left(x\right)+c_p{ϵ}_p{ϵ}_{p-1}u\left(x\right),$$

(28)

given that $r\left(x\right)$ is piecewise continuous and positive. Let $(x_1^{″},\dots ,x_p^{″})$ be such points and let us build the matrix:

$$A=\left(\begin{array}{cccc}{\phi }_1\left(x_1^{″}\right)& \dots & {\phi }_1\left(x_{p-1}^{″}\right)& {\phi }_1\left(x_p^{″}\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_{p-1}\left(x_1^{″}\right)& \dots & {\phi }_{p-1}\left(x_{p-1}^{″}\right)& {\phi }_{p-1}\left(x_p^{″}\right)\\ u\left(x_1^{″}\right)& \dots & u\left(x_{p-1}^{″}\right)& u\left(x_p^{″}\right)\end{array}\right),$$

(29)

whose determinant $\left|A\right|$ is obviously zero. Given that ${\phi }_1,\dots ,{\phi }_{p-1}$ is a Chebyshev system on $I_1$, the matrix

$$\left(\begin{array}{ccc}{\phi }_1\left(x_1^{″}\right)& \dots & {\phi }_1\left(x_{p-1}^{″}\right)\\ \cdots & \cdots & \cdots \\ {\phi }_{p-1}\left(x_1^{″}\right)& \dots & {\phi }_{p-1}\left(x_{p-1}^{″}\right)\end{array}\right),$$

(30)

has a range of $p-1$, and therefore A must also have a range of $p-1$ (the difference between both matrices is one row and one column). That means that the null space N of A, namely the subspace of column vectors $C\in {\mathbb{R}}^p$ such that $CA=0$, has dimension 1, and the vector composed by the coefficients of ${\phi }_1\left(x\right),\dots ,{\phi }_{p-1}\left(x\right),u\left(x\right)$ in (27) belongs to N.

Expanding the determinant:

$$\left|\begin{array}{cccc}{\phi }_1\left(x_1^{″}\right)& \dots & {\phi }_1\left(x_{p-1}^{″}\right)& {\phi }_1\left(x\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_{p-1}\left(x_1^{″}\right)& \dots & {\phi }_{p-1}\left(x_{p-1}^{″}\right)& {\phi }_{p-1}\left(x\right)\\ u\left(x_1^{″}\right)& \dots & u\left(x_{p-1}^{″}\right)& u\left(x\right)\end{array}\right|$$

(31)

along its last column, one obtains a linear combination of ${\phi }_1\left(x\right),\dots ,{\phi }_{p-1}\left(x\right),u\left(x\right)$ that vanishes at $x_1^{″},\dots ,x_p^{″}\in I_1$ and which equals ${\Delta }_p(u;x_1^{″},\dots ,x_{p-1}^{″},x)$ for any $x\in (a,b)$. Therefore, the column vector composed by the coefficients of that linear combination (namely, the cofactors of the last column of (31)) must also belong to the subspace N and be a multiple of the coefficients of ${\phi }_1\left(x\right),\dots ,{\phi }_{p-1}\left(x\right),u\left(x\right)$ in (27), that is, $v\left(x\right)$ of (28) must be a multiple of ${\Delta }_p(u;x_1^{″},\dots ,x_{p-1}^{″},x)$. As $v\left(x\right)$ changed its sign at $x=x_p^{″}$, ${\Delta }_p(u;x_1^{″},\dots ,x_{p-1}^{″},x)$ must also change its sign at $x=x_p^{″}$. A similar conclusion can be obtained if y in (21) has infinitely many zeroes in $(a,b)$. Thus, we have proven that if ${\Delta }_p\left(M_{p}u\right)$ vanishes at ${\mathbf{x}}^{\prime }\in {\Omega }^{\ast }$, then ${\Delta }_p\left(u\right)$ must change sign at least at a point ${\mathbf{x}}^{″}\in {\Omega }^{\ast }$, so u cannot belong to $P_p$. In summary, if $u\in P_p$, ${\Delta }_p\left(M_{p}u\right)\ne 0$ in ${\Omega }^{\ast }$.

It remains to be proven that ${\Delta }_p\left(M_{p}u\right)$ and ${\Delta }_p\left(u\right)$ have the same sign, namely that ${\Delta }_p\left(M_{p}u\right)>0$ in ${\Omega }^{\ast }$ when ${\Delta }_p\left(u\right)>0$ in the same set. This follows from the expression (see, for instance, [4], Chapter 4, Section 3, Equation (62)),

$$\begin{array}{c}\Delta \left(\begin{array}{ccccc}M{\phi }_1& M{\phi }_2& \dots & M{\phi }_{p-1}& Mu\\ x_1& x_2& \dots & x_{p-1}& x_p\end{array}\right)\\ \\ ={\int }_a^b\cdots {\int }_a^{t_2}G\left(\begin{array}{cccc}{x}_1,& x_2,& \dots ,& x_p\\ t_1,& t_2,& \dots ,& t_p\end{array}\right)\Delta \left(\begin{array}{ccccc}{\phi }_1& {\phi }_2& \dots & {\phi }_{p-1}& u\\ t_1& t_2& \dots & t_{p-1}& t_p\end{array}\right)r\left(t_1\right)\cdots r\left(t_p\right)d{t}_1\cdots d{t}_p.\end{array}$$

(32)

The left hand side of (32) is exactly $\frac{{ϵ}_p{\Delta }_p\left(M_{p}u\right)}{|{\lambda }_1\cdots {\lambda }_{p-1}|}$. From the sign-regularity of $G(x,t)$, (7) and (32), one has that the signs of ${\Delta }_p\left(M_{p}u\right)$ and ${\Delta }_p\left(u\right)$ coincide in ${\Omega }^{\ast }$. This completes the proof. □

Although the Theorem 1 can be used to obtain some information about the nature of the eigenvalues ${\lambda }_p$, it does not provide any indication about the relationship between the Collatz–Wielandt numbers and ${\lambda }_p$. A first step towards that direction can be made if we can find a solid cone $K_p$ contained in $P_p$ ($P_p$ is not solid, as per its definition), which is mapped by $M_p$ into itself, as the next theorem will show. For that, we need to introduce the notion of weak irreducibility (see [35], Definition 7.5):

Definition 2.

Let P be a solid cone. We say that $M_p$ is weakly irreducible, if the boundary $\partial P$ of P contains no eigenvectors of $M_p$ pertaining to nonnegative eigenvalues.

Theorem 2.

Let $K_p$ be a solid cone such that $K_p\subset P_p$. If $M_p$ maps $K_p$ into itself, then ${\phi }_p\in int\left(K_p\right)$ and for any $u\in K_p\setminus \{0\}$ one has

$$\underline{r}(M_p,M_p^{j}u,K_p)\le \frac{1}{|{\lambda }_p|},j=1,\dots ,$$

(33)

$$\underset{j\to \infty }{lim}\underline{r}(M_p,M_p^{j}u,K_p)=\underset{j\to \infty }{lim}\overline{r}(M_p,M_p^{j}u,K_p)=\frac{1}{|{\lambda }_p|},$$

(34)

and

$$\underset{j\to \infty }{lim}{\left|{\lambda }_p\right|}^j{M}_p^{j}u=f\left(u\right){\phi }_p,$$

(35)

where $f\left(u\right)$ is a non-zero linear functional dependent on u and ${\phi }_p$.

Proof. 

Let us first prove that $M_p$ is weakly irreducible in $K_p$.

From the property 2 of the sign regular problems (see the Introduction) any non-trivial linear combination of ${\phi }_1,\dots ,{\phi }_{p-1},{\phi }_i$ with $i>p$, where the coefficients of ${\phi }_1,\dots ,{\phi }_{p-1}$ are zero, must have $i-1$ zeroes in $I_1$. That implies that ${\Delta }_p\left({\phi }_i\right)$ must vanish in at least $i-p$ points in ${\Omega }^{\ast }$. Using an argument similar to that used in the Theorem 1, one has that ${\Delta }_p\left({\phi }_i\right)$ must change its sign in ${\Omega }^{\ast }$, so ${\phi }_i\notin P_p$ (ergo ${\phi }_i\notin K_p$) for $i>p$, and therefore $M_p$ is weakly irreducible in $K_p$ according to the Definition 2.

Next, in the Banach space $B_p$ it is clear that $r\left(M_p\right)=\frac{1}{|{\lambda }_p|}>\frac{1}{|{\lambda }_i|}$ for $i>p$. We can apply [35] (Theorem 7.7) to conclude that $\frac{1}{|{\lambda }_p|}$ is a simple eigenvalue of $M_p$ with an eigenvector ${\phi }_p\in int\left(K_p\right)$. From here and [25] (Lemma 1.16), one gets (33). Likewise, following the same reasoning as in [32] (Theorems 6 and 7), one can get to (35) and, noting that for some $j_0>0$, $M_p^{j}u\in int\left(K_p\right)$ for all $j\ge j_0$, also to (34). □

2.1.3. The Cone ${\mathcal{P}}_p$

The previous theorem does not offer any hints for finding the solid cone $K_p$, nor does it indicate any relationship between $\overline{r}(M_p,M_p^{j}u,K_p)$ and ${\lambda }_p$ beyond the fact that the upper Collatz–Wielandt number converges to $\frac{1}{|{\lambda }_p|}$. To determine such a relationship, the solid cone $K_p$ must be such that $M_p$ maps it (excluding the zero element) to its interior. As it turns out, under certain conditions, the cone ${\mathcal{P}}_p$ defined in (19) is solid and satisfies that property with regards to $M_p$.

To establish that, let us start by identifying the interior of ${\mathcal{P}}_p$. Although one could be tempted to think that $int\left({\mathcal{P}}_p\right)$ is merely composed by the functions $y\in {\mathcal{B}}_p$ such that ${\Delta }_p\left(y\right)>0$ in ${\Omega }^{\ast }$, in the end this is only a necessary condition as one must pay attention to the value of ${\Delta }_p\left(y\right)$ in the vicinity of the closure of ${\Omega }^{\ast }$ where in fact ${\Delta }_p\left(y\right)$ vanishes, namely when $x_1$ is close to a if $a\notin I_1$, when $x_p$ is close to b if $b\notin I_1$, and when several values $x_i$ converge simultaneously to the same point $x^{\ast }$.

The next Lemma will give the value of ${\Delta }_p\left(y\right)$ in the latter case.

Lemma 2.

Let us suppose that ${\phi }_i,u\in A{C}^l[a,b]$. If $x_{i+1},\dots ,x_{i+l}\in (x_i-\delta ,x_i+\delta )$ for $\delta >0$, then

$$\begin{array}{c}{\Delta }_p(u;x_1,\dots ,x_p)\\ ={\left.\frac{1}{{\prod }_{j=1}^{l}j!}\frac{{\partial }^l}{\partial x_{i+l}^l}\frac{{\partial }^{l-1}}{\partial x_{i+l-1}^{l-1}}\cdots \frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+1}}\right|}_{x_i=\cdots =x_{i+l}}{\prod }_{i\le k<j\le i+l}(x_j-x_k)+o\left({\delta }^{\frac{l(l+1)}{2}}\right).\end{array}$$

(36)

Proof. 

Noting that $1+2+\cdots +l=\frac{l(l+1)}{2}$, Taylor’s formula for multivariate functions allows expressing ${\Delta }_p(u;x_1,\dots ,x_p)$ when $x_{i+1},\dots ,x_{i+l}$ are in a neighborhood of $x_i$ of radius $\delta$, as:

$$\begin{array}{c}{\Delta }_p(u;x_1,\dots ,x_p)={\left.{\Delta }_p\left(u\right)\right|}_{x_i=\cdots =x_{i+l}}\\ \\ +{\left.\sum _{j=1}^l\frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+j}}\right|}_{x_i=\cdots =x_{i+l}}(x_{i+j}-x_i)+{\left.\frac{1}{2}\sum _{j=1}^l\sum _{k=1}^l\frac{\partial }{\partial x_{i+j}}\frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+k}}\right|}_{x_i=\cdots =x_{i+l}}(x_{i+j}-x_i)(x_{i+k}-x_i)+\cdots +\\ \\ +\frac{1}{\frac{l(l+1)}{2}!}\underset{\frac{l(l+1)}{2}summationsymbols}{\underbrace{{\displaystyle \sum _{j=1}^l}{\displaystyle \sum _{k=1}^l}\cdots {\displaystyle \sum _{q=1}^l}}}{\left.\frac{\partial }{\partial x_{i+j}}\cdots \frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+q}}\right|}_{x_i=\cdots =x_{i+l}}(x_{i+j}-x_i)\cdots (x_{i+q}-x_i)+o\left({\delta }^{\frac{l(l+1)}{2}}\right).\end{array}$$

(37)

From the properties of the determinants, (14) and (15), it is clear that all the terms in (37) where the order of the partial derivatives of two different $x_j$ coincide, vanish, which yields:

$${\Delta }_p(u;x_1,\dots ,x_p)$$

$$={\left.\sum _{(j_1,\dots ,j_l)\in K}\frac{1}{{\prod }_{j=1}^{l}j!}\frac{{\partial }^{j_l}}{\partial x_{i+l}^{j_l}}\frac{{\partial }^{j_{l-1}}}{\partial x_{i+l-1}^{j_{l-1}}}\cdots \frac{{\partial }^{j_1}{\Delta }_p\left(u\right)}{\partial x_{i+1}^{j_1}}\right|}_{x_i=\cdots =x_{i+l}}{(x_{i+l}-x_i)}^{j_l}\cdots {(x_{i+1}-x_i)}^{j_1}+o\left({\delta }^{\frac{l(l+1)}{2}}\right),$$

(38)

where K is the set of all permutations of the indexes $(1,2,\dots ,l)$.

Let us denote by $s(j_1,j_2,\dots ,j_l)$ the signature of the l-tuple $(j_1,j_2,\dots ,j_l)$ (the signature of a tuple is defined to be $+1$ whenever the reordering $(1,2,\dots ,l)$ can be achieved by successively interchanging two entries of the tuple an even number of times, and $-1$ whenever it can be achieved by an odd number of such interchanges). As the different partial derivatives appearing in (38) are continuous in $x_{i+1},\dots ,x_{i+l}$ (${\phi }_i,u\in C^l[a,b]$ by hypothesis) and are calculated at the same point $x_{i+j}=x_i$, we can exchange their order just by taking into account the impact of such a change in the determinant (14) (it is an exchange of rows), which leads us to:

$${\left.\frac{{\partial }^{j_l}}{\partial x_{i+l}^{j_l}}\frac{{\partial }^{j_{l-1}}}{\partial x_{i+l-1}^{j_{l-1}}}\cdots \frac{{\partial }^{j_1}{\Delta }_p\left(u\right)}{\partial x_{i+1}^{j_1}}\right|}_{x_i=\cdots =x_{i+l}}={\left.{(-1)}^{s(j_1,j_2,\dots ,j_l)}\frac{{\partial }^l}{\partial x_{i+l}^l}\frac{{\partial }^{l-1}}{\partial x_{i+l-1}^{l-1}}\cdots \frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+1}}\right|}_{x_i=\cdots =x_{i+l}}.$$

(39)

Equations (38) and (39) give:

$${\Delta }_p(u;x_1,\dots ,x_p)$$

$$={\left.\frac{1}{{\prod }_{j=1}^{l}j!}\frac{{\partial }^l}{\partial x_{i+l}^l}\frac{{\partial }^{l-1}}{\partial x_{i+l-1}^{l-1}}\cdots \frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+1}}\right|}_{x_i=\cdots =x_{i+l}}\sum _{(j_1,\dots ,j_l)\in K}{(-1)}^{s(j_1,j_2,\dots ,j_l)}{(x_{i+1}-x_i)}^{j_1}\cdots {(x_{i+l}-x_i)}^{j_l}+o\left({\delta }^{\frac{l(l+1)}{2}}\right).$$

(40)

The expression within the sum in (40) has exactly the form of the Vandermonde determinant

$$\left|\begin{array}{ccc}{x}_{i+1}-x_i& \cdots & x_{i+l}-x_i\\ \cdots & \cdots & \cdots \\ {(x_{i+1}-x_i)}^l& \cdots & {(x_{i+l}-x_i)}^l\end{array}\right|,$$

whose value, as it is well known, is ${\prod }_{i\le k<j\le i+l}(x_j-x_k)$. From here and (40), one gets (36). □

Given that $a\le x_1<\dots <x_p\le b$ for $(x_1,\dots ,x_p)\in \Omega$, the consequence of the Lemma 2 is that, if $\alpha$ is the lowest derivative which the boundary conditions on a do not specify to vanish, $\beta$ is the lowest derivative which the boundary conditions on b do not specify to vanish, and ${\phi }_i,u\in A{C}^p[a,b]$, the interior of ${\mathcal{P}}_p$ is defined by:

$$\begin{array}{c}int\left({\mathcal{P}}_p\right)=\left\{y\in {\mathcal{B}}_p:{\Delta }_p\left(y\right)>0,(x_1,\dots x_p)\in {\Omega }^{\ast },\right.\\ \\ {\left.\frac{{\partial }^{\alpha }{\Delta }_p\left(u\right)}{\partial x_1^{\alpha }}\right|}_{x_1=a}>0,(a,x_2,\dots ,x_p)\in \Omega ,ifa\notin I_1,\\ \\ {\left.{(-1)}^{\beta }\frac{{\partial }^{\beta }{\Delta }_p\left(u\right)}{\partial x_p^{\beta }}\right|}_{x_p=b}>0,(x_1,x_2,\dots ,b)\in \Omega ,ifb\notin I_1,\\ \\ \left.{\left.\frac{{\partial }^l}{\partial x_{i+l}^l}\frac{{\partial }^{l-1}}{\partial x_{i+l-1}^{l-1}}\cdots \frac{\partial {\Delta }_p\left(u\right)}{\partial x_{i+1}}\right|}_{x_i=\cdots =x_{i+l}=x^{\ast }}>0,i=1,\dots ,p,l=1,\dots ,p-i\right\}.\end{array}$$

(41)

Remark 1.

By the definition of $M_p$ and $G(x,t)$, it is clear that ${\phi }_i,M_{p}u\in A{C}^{n-1}[a,b]$. However, if $p>n-1$, one cannot grant that ${\phi }_i,M_{p}u\in A{C}^p[a,b]$, or even the mere existence of $int\left({\mathcal{P}}_p\right)$, without imposing extra conditions on u, r and the coefficients $a_j$ of L in (1). The next theorem will display some sufficient conditions for that.

Theorem 3.

Let us suppose that either $p<n$ or $p\ge n$, $r\left(x\right)$, $a_j\left(x\right)\in A{C}^{p-n}[a,b]$ for $j=0,\dots ,n-1$. Let q be the lowest integer greater than 1 such that $q·n>p$. Then ${\mathcal{P}}_p$ is solid and $M_p^q$ maps $P_p\setminus \{0\}$ into $int\left({\mathcal{P}}_p\right)$. In addition, if $u\in P_p\setminus \{0\}$, then,

$$\underline{r}(M_p,M_p^{j}u)\le \frac{1}{|{\lambda }_p|}\le \overline{r}(M_p,M_p^{j}u),j=q,\dots ,$$

(42)

$$\underset{j\to \infty }{lim}\underline{r}(M_p,M_p^{j}u)=\underset{j\to \infty }{lim}\overline{r}(M_p,M_p^{j}u)=\frac{1}{|{\lambda }_p|},$$

(43)

and

$$\underset{j\to \infty }{lim}{\left|{\lambda }_p\right|}^j{M}_p^{j}u=f\left(u\right){\phi }_p,$$

(44)

where $f\left(u\right)$ is a non-zero linear functional dependent on u and ${\phi }_p$.

Proof. 

From (1) and (2), one has:

$${\left(M_p^{q}u\right)}^{\left(n\right)}\left(x\right)={ϵ}_{p-1}{ϵ}_{p}r\left(x\right)M_p^{q-1}u\left(x\right)-\sum _{i=0}^{n-1}{a}_i\left(x\right){\left(M_p^{q}u\right)}^{\left(i\right)}\left(x\right),x\in (a,b).$$

(45)

Therefore, in order for ${\left(M_p^{q}u\right)}^{\left(n\right)}\left(x\right)$ to belong to $A{C}^{p-n}[a,b]$, it suffices that $M_p^{q-1}u\left(x\right)$, $r\left(x\right)$, $a_j\left(x\right)\in A{C}^{p-n}[a,b]$ for $j=0,\dots ,n-1$, which is granted by the hypotheses and the fact that $q·n>p$ and $M_{p}u\in A{C}^{n-1}[a,b]$. With this, following the same steps as in the Theorem 1 it is straightforward to show that $M_p^q$ maps $B_p$ into ${\mathcal{B}}_p$, and that ${\Delta }_p\left(M_p^{q}u\right)>0$ for $\mathbf{x}\in {\Omega }^{\ast }$ provided that $u\in P_p\setminus \{0\}$, which covers the conditions of the first line of (41).

Let us focus now on the condition of the second line of (41), related to the derivative of ${\Delta }_p\left(M_p^{q}u\right)$ at a when $a\notin I_1$. From (32), one has:

$$\begin{array}{c}\frac{{ϵ}_p}{|{\lambda }_1\cdots {\lambda }_{p-1}|}\frac{{\partial }^{\alpha }{\Delta }_p(M_p^{q}u;a,x_2,\dots ,x_p)}{\partial x_1^{\alpha }}=\frac{{\partial }^{\alpha }\Delta }{\partial x_1^{\alpha }}\left(\begin{array}{ccccc}M{\phi }_1& M{\phi }_2& \dots & M{\phi }_{p-1}& M\left(M_p^{q-1}\right)u\\ a& x_2& \dots & x_{p-1}& x_p\end{array}\right)\\ \\ ={\int }_a^b\cdots {\int }_a^{t_2}\left|\begin{array}{ccc}\frac{{\partial }^{\alpha }G(a,t_1)}{\partial x_1^{\alpha }}& \dots & \frac{{\partial }^{\alpha }G(a,t_p)}{\partial x_1^{\alpha }}\\ G(x_2,t_1)& \dots & G(x_2,t_p)\\ \cdots & \cdots & \cdots \\ G(x_p,t_1)& \dots & G(x_p,t_p)\end{array}\right|\Delta \left(\begin{array}{ccccc}{\phi }_1& {\phi }_2& \dots & {\phi }_{p-1}& M_p^{q-1}u\\ t_1& t_2& \dots & t_{p-1}& t_p\end{array}\right)r\left(t_1\right)\cdots r\left(t_p\right)d{t}_1\cdots d{t}_p\\ \\ ={\int }_a^b\cdots {\int }_a^{t_2}\left|\begin{array}{ccc}\frac{{\partial }^{\alpha }G(a,t_1)}{\partial x_1^{\alpha }}& \dots & \frac{{\partial }^{\alpha }G(a,t_p)}{\partial x_1^{\alpha }}\\ G(x_2,t_1)& \dots & G(x_2,t_p)\\ \cdots & \cdots & \cdots \\ G(x_p,t_1)& \dots & G(x_p,t_p)\end{array}\right|\Delta (M_p^{q-1}u;t_1,\dots ,t_p)r\left(t_1\right)\cdots r\left(t_p\right)d{t}_1\cdots d{t}_p.\end{array}$$

(46)

As $q>1$, ${\Delta }_p\left(M_p^{q-1}u\right)>0$ in ${\Omega }^{\ast }$ according to the Theorem 1. For that reason, the key to grant the positivity of the $\alpha$-th partial derivative at a lies on the value of the determinant of the matrix:

$$\left(\begin{array}{ccc}\frac{{\partial }^{\alpha }G(a,t_1)}{\partial x_1^{\alpha }}& \dots & \frac{{\partial }^{\alpha }G(a,t_p)}{\partial x_1^{\alpha }}\\ G(x_2,t_1)& \dots & G(x_2,t_p)\\ \cdots & \cdots & \cdots \\ G(x_p,t_1)& \dots & G(x_p,t_p)\end{array}\right).$$

(47)

Let us denote by $K_p(t_1;x_2,\dots ,x_p)$ the matrix:

$$K_p(t_1;x_2,\dots ,x_p)=\left(\begin{array}{cccc}\frac{{\partial }^{\alpha }G(a,t_1)}{\partial x_1^{\alpha }}& \frac{{\partial }^{\alpha }G(a,x_2)}{\partial x_1^{\alpha }}& \dots & \frac{{\partial }^{\alpha }G(a,x_p)}{\partial x_1^{\alpha }}\\ G(x_2,t_1)& G(x_2,x_2)& \dots & G(x_2,x_p)\\ \cdots & \cdots & \cdots & \cdots \\ G(x_p,t_1)& G(x_p,x_p)& \dots & G(x_p,x_p)\end{array}\right),$$

(48)

whose determinant we will write as $|K_p|$. Using Taylor’s formula, when $x_1$ is in the neighborhood of a one has:

$$G\left(\begin{array}{c}{x}_1,x_2,\dots ,x_p\\ t_1,x_2,\dots ,x_p\end{array}\right)=\left|\begin{array}{cccc}\frac{{\partial }^{\alpha }G(a,t_1)}{\partial x_1^{\alpha }}& \frac{{\partial }^{\alpha }G(a,x_2)}{\partial x_1^{\alpha }}& \dots & \frac{{\partial }^{\alpha }G(a,x_p)}{\partial x_1^{\alpha }}\\ G(x_2,t_1)& G(x_2,x_2)& \dots & G(x_2,x_p)\\ \cdots & \cdots & \cdots & \cdots \\ G(x_p,t_1)& G(x_p,x_2)& \dots & G(x_p,x_p)\end{array}\right|\frac{{(x_1-a)}^{\alpha }}{\alpha !}+o\left((x_1-a){}^{\alpha }\right),$$

so that the matrix $K_p(t_1;x_2,\dots ,x_p)$ must also be sign regular with ${ϵ}_p\left|K_p(t_1;x_2,\dots ,x_p)\right|\ge 0$ for $a<t_1<x_2<\cdots <x_p<b$ as per (7). We will prove in fact that

$${ϵ}_p\left|K_p\right|>0,a<t_1<x_2<\cdots x_p<b.$$

(49)

We will proceed by induction, following the ideas of [1]. Thus, from [1] (Equation (12.12)), we know that:

$${ϵ}_1\frac{{\partial }^{\alpha }G(a,t)}{\partial x^{\alpha }}>0,a<t<b.$$

(50)

Let us assume that ${ϵ}_{p-1}\left|K_{p-1}(t_1^{\prime };x_2^{\prime },\dots ,x_{p-1}^{\prime })\right|>0$ for $a<t_1^{\prime }<x_2^{\prime }<\cdots <x_{p-1}^{\prime }<b$ and ${ϵ}_p\left|K_p(t_1^{\prime };x_2^{\prime },\dots ,x_p^{\prime })\right|=0$ for $a<t_1^{\prime }<x_2^{\prime }<\cdots <x_{p-1}^{\prime }<x_p^{\prime }<b$. If we introduce an additional pair $(x^{\ast },x^{\ast })$ such that $a<t_1^{\prime }<x^{\ast }<x_2^{\prime }$, by the same argument as before on the sign regularity of $G(x,t)$ one must have that the matrix $K_{p+1}(t_1^{\prime };x^{\ast },x_2^{\prime },\dots ,x_p^{\prime })$ must be sign regular too with ${ϵ}_{p+1}\left|K_{p+1}(t_1^{\prime };x^{\ast },x_2^{\prime },\dots ,x_p^{\prime })\right|\ge 0$.

$K_{p+1}(t_1^{\prime };x^{\ast },x_2^{\prime },\dots ,x_p^{\prime })$ is therefore a $(p+1)\times (p+1)$ sign regular matrix whose first row is composed of terms of the type $\frac{{\partial }^{\alpha }G(a,t_1^{\prime })}{\partial x^{\alpha }}$ and $\frac{{\partial }^{\alpha }G(a,x_j^{\prime })}{\partial x^{\alpha }}$ while the rest of rows form the matrix:

$$\left(\begin{array}{ccccc}G(x^{\ast },t_1^{\prime })& G(x^{\ast },x^{\ast })& G(x^{\ast },x_2^{\prime })& \dots & G(x^{\ast },x_p^{\prime })\\ G(x_2^{\prime },t_1^{\prime })& G(x_2^{\prime },x^{\ast })& G(x_2^{\prime },x_2^{\prime })& \dots & G(x_2^{\prime },x_p^{\prime })\\ \cdots & \cdots & \cdots & \cdots & \cdots \\ G(x_p^{\prime },t_1^{\prime })& G(x_p^{\prime },x^{\ast })& G(x_p^{\prime },x_2^{\prime })& \dots & G(x_p^{\prime },x_p^{\prime })\end{array}\right),$$

(51)

which is $p\times (p+1)$ and sign regular, and whose last p columns are linearly independent (their determinant does not vanish as per (9)). Accordingly its range is p, and the range of $K_{p+1}(t_1^{\prime };x^{\ast },x_2^{\prime },\dots ,x_p^{\prime })$ must be at least p too.

In the same way, one can find that the range of $K_p(t_1^{\prime };x_2^{\prime },\dots ,x_p^{\prime })$, whose determinant vanishes, is $p-1$, with its last $p-1$ rows linearly independent. Since the range of $K_{p-1}(t_1^{\prime };x_2^{\prime },\dots ,x_{p-1}^{\prime })$ is $p-1$ by the induction hypothesis, one can apply [1] (Lemma 2) to conclude that the range of $K_{p+1}(t_1^{\prime };x^{\ast },x_2^{\prime },\dots ,x_p^{\prime })$ equals $p-1$, contradicting the previous assertion. Therefore, $|K_p(t_1^{\prime };x_2^{\prime },\dots ,x_p^{\prime })|$ cannot be zero for $a<t_1^{\prime }<x_2^{\prime }<\cdots <x_p^{\prime }<b$ and, due to its sign regularity, ${ϵ}_p\left|K_p(t_1^{\prime };x_2^{\prime },\dots ,x_p^{\prime })\right|>0$. By continuity, the matrix (47) must have a determinant of sign ${ϵ}_p$ for $t_i\in (x_i-\delta ,x_i+\delta )$, $i=2,\dots ,p$. From here and (46) one gets ${\left.\frac{{\partial }^{\alpha }{\Delta }_p\left(u\right)}{\partial x_1^{\alpha }}\right|}_{x_1=a}>0$ for $(a,x_2,\dots ,x_p)\in \Omega$.

The condition of the third line of (41) with respect to the derivatives at b, if $b\notin I_1$, can be proven in the same way.

As for the last condition of (41), let us assume that:

$${\left.\frac{{\partial }^l}{\partial x_{i+l}^l}\frac{{\partial }^{l-1}}{\partial x_{i+l-1}^{l-1}}\cdots \frac{\partial {\Delta }_p\left(M_p^{q}u\right)}{\partial x_{i+1}}\right|}_{x_i=\cdots =x_{i+l}=x^{\ast }}=0,$$

(52)

for an ${\mathbf{x}}^{\ast }\in {\Omega }^{\ast }$. From here, (14) and (15) one has that there exists a linear combination of ${\phi }_1,\dots ,{\phi }_{p-1},M_p^{q}u$,

$$w\left(x\right)=d_1{\phi }_1+\cdots +d_{p-1}{\phi }_{p-1}+d_p{M}_p^{q}u,$$

with p zeroes on $I_1$, counting their multiplicities (there must be a zero multiple of order $l+1$ at $x^{\ast }$). Using a similar argument to that of the Theorem 1, one has that the function

$$Lw\left(x\right)=r\left(x\right)\left(d_1{\lambda }_1{\phi }_1+\cdots +d_{p-1}{\lambda }_{p-1}{\phi }_{p-1}+d_p{ϵ}_{p-1}{ϵ}_p{M}_p^{q-1}u\right)$$

(53)

must change its sign at least p times in $(a,b)$, exactly at the same points as:

$$y\left(x\right)=d_1{\lambda }_1{\phi }_1+\cdots +d_{p-1}{\lambda }_{p-1}{\phi }_{p-1}+d_p{ϵ}_{p-1}{ϵ}_p{M}_p^{q-1}u.$$

(54)

Let these points be $x_1^{″},\dots ,x_p^{″}$. This means that the function ${\Delta }_p(M_p^{q-1}u;x_1^{″},\dots ,x_{p-1}^{″},x)$ must change its sign at $x=x_p^{″}$ and therefore u cannot belong to $P_p$. This completes the proof that $M_p^q(P_p\setminus \{0\})\subset int\left({\mathcal{P}}_p\right)$, that is, that $M_p$ is $u_0$-positive in ${\mathcal{P}}_p$. Equations (42)–(44) follow now from [32] (Theorems 6 and 7). □

Remark 2.

The Theorem 3 shows that the Collatz–Wielandt numbers $\underline{r}(M_p,M_p^{j}u)$ and $\overline{r}(M_p,M_p^{j}u)$ are lower and upper bounds of the inverse of the eigenvalue $|{\lambda }_p|$, which converge to it as the iteration index j grows. Therefore, they allow determining ${\lambda }_p$ with as much accuracy as desired, as the error in the approximation is bounded by the difference $\overline{r}(M_p,M_p^{j}u)-\underline{r}(M_p,M_p^{j}u)$.

To clarify how to calculate them, let us recall that, from (10), (14), (15) and (19), for the cone ${\mathcal{P}}_p$, they can be expressed as:

$$\overline{r}(M_p,M_p^{j}u)=inf\{w\in \mathbb{R}:{\Delta }_p\left(M_p^{j+1}u\right)\le w{\Delta }_p\left(M_p^{j}u\right),x\in \Omega \},j\ge q,$$

(55)

and

$$\underline{r}(M_p,M_p^{j}u)=sup\{w\in \mathbb{R}:w{\Delta }_p\left(M_p^{j}u\right)\le {\Delta }_p\left(M_p^{j+1}u\right),x\in \Omega \},j\ge q,$$

(56)

provided that ${\Delta }_p\left(u\right)\ge 0$ in Ω. Therefore, their calculation requires a comparison of two functions in the simplex $\Omega \subset {[a,b]}^p$.

We will close this subsection by showing that, in practice, the function u does not need to be orthogonal to the adjoint eigenfunctions ${\psi }_i$ for $i=1,\dots ,p-1$ for Equations (42)–(43) to be valid.

Theorem 4.

Under the conditions of the Theorem 3, let us suppose that $v\in C[a,b]$ and it is not identically zero. If ${\Delta }_p\left(v\right)\ge 0$ in Ω, then its Collatz–Wielandt numbers satisfy (42)–(43).

Proof. 

If we decompose v as:

$$v=\sum _{i=1}^{p-1}〈{\psi }_i,v〉{\phi }_i+u,$$

(57)

it follows that $u\in B_p$ and

$${\Delta }_p\left(u\right)=-\sum _{i=1}^{p-1}〈{\psi }_i,v〉{\Delta }_p\left({\phi }_i\right)+{\Delta }_p\left(v\right)={\Delta }_p\left(v\right)\ge 0,$$

for $\mathbf{x}\in \Omega$, that is, $u\in P_p$. Likewise,

$$M_p^{j}v={\left({ϵ}_p{ϵ}_{p-1}\right)}^j\sum _{i=1}^{p-1}\frac{〈{\psi }_i,v〉}{{\lambda }_i^j}{\phi }_i+M_p^{j}u,$$

(58)

so ${\Delta }_p\left(M_p^{j}v\right)={\Delta }_p\left(M_p^{j}u\right)$ and therefore $\underline{r}(M_p,M_p^{j}u)=\underline{r}(M_p,M_p^{j}v)$ and $\overline{r}(M_p,M_p^{j}v)=\overline{r}(M_p,M_p^{j}u)$, for $j=1,\dots$. □

Remark 3.

The property that the function u of the previous theorem fails to hold due to the lack of orthogonality with ${\psi }_i$, $i=1,\dots ,p-1$, is precisely (44), since the term in (58) associated to $\frac{1}{{\lambda }_1}$ gets bigger in absolute value than the rest of terms as the iteration index j grows.

2.1.4. The Calculation of the Adjoint Eigenfunctions ${\psi }_i$

The application of the method described in the Remark 2 for different values of p requires knowledge of the eigenfunctions ${\phi }_i$, $i=1,\dots ,p-1$. Although the Theorem 4 stated that one can start the iteration $M_p^{j}v$ with a function v not orthogonal to the adjoint eigenfunctions ${\psi }_i$, $i=1,\dots ,p-1$, such an orthogonality is necessary in order to use (44) to determine ${\phi }_p$, and employ the latter in the calculation of ${\lambda }_i$, ${\phi }_i$ for $i>p$. This implies that knowledge on ${\psi }_i$ must also be obtained as p increases.

The process to obtain ${\psi }_p$ is very similar to that followed for ${\phi }_p$. To start with, the sign regularity of $G(x,t)$ ensures the sign regularity of $G^{\ast }(x,t)$, where $G^{\ast }(x,t)$ is the Green function of the adjoint problem,

$$L^{\ast }y=0,x\in [a,b];U_i^{\ast }\left(y\right)=0,i=1,\dots ,n.$$

(59)

This is due to the fact that $G^{\ast }(x,t)=G(t,x)$ (see [33], Chapter 11, Theorem 4.2), so that, if $G(x,t)$ satisfies (7)–(9), it is immediately shown that these conditions hold for $G(t,x)$ too.

Next, one has to define the Banach space,

$$B_p^{\ast }=\{y\in C[a,b]:〈{\phi }_i,y〉=0,i=1,\dots ,p-1\},$$

(60)

subject to the sup norm $\parallel y\parallel =sup\{|f\left(x\right)|,x\in [a,b]\}$, the cone $P_p^{\ast }$

$$P_p^{\ast }=\{y\in B_p^{\ast }:{\Delta }_p^{\ast }\left(y\right)\ge 0,(x_1,\dots x_p)\in \overline{\Omega }\},$$

(61)

where ${\Delta }_p^{\ast }\left(y\right)$ is defined by:

$${\Delta }_p^{\ast }\left(y\right)=\Delta \left(\begin{array}{ccccc}{\psi }_1& {\psi }_2& \dots & {\psi }_{p-1}& y_p\\ x_1& x_2& \dots & x_{p-1}& x_p\end{array}\right),(x_1,\dots x_p)\in {[a,b]}^p\subset {\mathbb{R}}^p,$$

(62)

the Banach space ${\mathcal{B}}_p^{\ast }$

$${\mathcal{B}}_p^{\ast }=\{y\in A{C}^{max(n-1,p)}[a,b]:〈{\phi }_i,y〉=0,i=1,\dots ,p-1;U_i^{\ast }\left(y\right)=0,i=1,\dots ,n\},$$

(63)

and the cone ${\mathcal{P}}_p^{\ast }$

$${\mathcal{P}}_p^{\ast }=\{y\in {\mathcal{B}}_p^{\ast }:{\Delta }_p^{\ast }\left(y\right)\ge 0,(x_1,\dots x_p)\in \overline{\Omega }\}.$$

(64)

The operator $M_p^{\ast }$ is defined by:

$$M_p^{\ast }w\left(x\right)={ϵ}_p{ϵ}_{p-1}{\int }_a^{b}G(t,x)w\left(t\right)r\left(t\right)dt,x\in [a,b],$$

(65)

given that the Green function of the adjoint problem (12) is exactly $G(t,x)$.

However, the conditions on $a_j\left(x\right)$ required for ${\mathcal{P}}_p^{\ast }$ to be solid, and for $M_p^{\ast }$ to map $P_p^{\ast }\setminus \{0\}$ into $int\left({\mathcal{P}}_p^{\ast }\right)$ are stronger, as the next theorem will show:

Theorem 5.

Let us suppose that either $p<n$ and $a_j\left(x\right)\in A{C}^j[a,b]$ for $j=0,\dots ,n-1$, or $p\ge n$, $r\left(x\right)\in A{C}^{p-n}[a,b]$ and $a_j\left(x\right)\in A{C}^{p+j-n}[a,b]$ for $j=0,\dots ,n-1$. Let q be the lowest integer greater than 1 such that $q·n>p$. Then ${\mathcal{P}}_p^{\ast }$ is solid and ${\left(M_p^{\ast }\right)}^q$ maps $P_p^{\ast }\setminus \{0\}$ into $int\left({\mathcal{P}}_p^{\ast }\right)$. In addition, if $w\in P_p^{\ast }\setminus \{0\}$, then,

$$\underline{r}(M_p^{\ast },{\left(M_p^{\ast }\right)}^{j}w)\le \frac{1}{|{\lambda }_p|}\le \overline{r}(M_p^{\ast },{\left(M_p^{\ast }\right)}^{j}w),j=q,\dots ,$$

(66)

$$\underset{j\to \infty }{lim}\underline{r}(M_p^{\ast },{\left(M_p^{\ast }\right)}^{j}w)=\underset{j\to \infty }{lim}\overline{r}(M_p^{\ast },{\left(M_p^{\ast }\right)}^{j}w)=\frac{1}{|{\lambda }_p|},$$

(67)

and

$$\underset{j\to \infty }{lim}{\left|{\lambda }_p\right|}^j{\left(M_p^{\ast }\right)}^{j}w=g\left(w\right){\psi }_p,$$

(68)

where $g\left(w\right)$ is a non-zero linear functional dependent on w and ${\psi }_p$.

Proof. 

The proof is essentially the same as that of the Theorem 3, changing the references to ${\phi }_i$, ${\Delta }_p$ and $M_p$, to ${\psi }_i$, ${\Delta }_p^{\ast }$ and $M_p\ast$, respectively. □

Remark 4.

The Theorem 5 implies that, in practice, $r\left(x\right)\in A{C}^{p-n}[a,b]$ and $a_j\left(x\right)\in A{C}^{p+j-n}[a,b]$ for $j=0,\dots ,n-1$, in order to be able to apply the procedure to determine ${\lambda }_p$, ${\phi }_p$ and ${\psi }_p$ for values of p higher than n.

2.2. Some Practical Considerations for the Application of the Procedure

2.2.1. The Selection of the Starting Function u

A key aspect in the procedure described before is the selection of a function u belonging to the cone $P_p$. Whereas the continuity in $[a,b]$ is a condition easy to meet, and for the orthogonality with the adjoint eigenfunctions ${\psi }_i$ one can follow a similar approach to that of the Theorem 4, the satisfaction of the condition ${\Delta }_p\left(u\right)\ge 0$ in $\Omega$ is not so straightforward.

A possible solution to this problem goes by interpolating u by means of linear splines. Thus, let us assume a partition $\{{\dot{x}}_l\}$ of $[a,b]$, with t points and a mesh size $h=max\{{\dot{x}}_{l+1}-{\dot{x}}_l,l=1,\dots ,t-1\}$. The linear spline $\widehat{u}$ is defined by:

$$\widehat{u}\left(x\right)=u\left({\dot{x}}_l\right)\frac{{\dot{x}}_{l+1}-x}{{\dot{x}}_{l+1}-{\dot{x}}_l}+u\left({\dot{x}}_{l+1}\right)\frac{x-{\dot{x}}_l}{{\dot{x}}_{l+1}-{\dot{x}}_l},x\in ({\dot{x}}_l,{\dot{x}}_{l+1}).$$

(69)

The linear spline $\widehat{u}\left(x\right)$ defines a function continuous on $[a,b]$, whose interpolation error $e\left(x\right)=|u\left(x\right)-\widehat{u}\left(x\right)|$ in each subinterval, if $u\in C^2[{\dot{x}}_l,{\dot{x}}_{l+1}]$, is given by:

$$e\left(x\right)=\frac{({\dot{x}}_{l+1}-x)(x-{\dot{x}}_l)}{2}{u}^{\prime \prime }\left(\xi \right),x,\xi \in ({\dot{x}}_l,{\dot{x}}_{l+1}),$$

(70)

with $\xi$ depending on x, and can be bounded by:

$$\parallel e\parallel =\frac{{({\dot{x}}_{l+1}-{\dot{x}}_l)}^2}{8}max\left\{\left|u^{″}\left(\xi \right)\right|,\xi \in ({\dot{x}}_l,{\dot{x}}_{l+1})\right\}.$$

(71)

From (71), it follows that if the size h of the mesh is small the interpolation error will also be small.

The advantage of the use of the splines is that it allows the reduction of the determination of ${\Delta }_p\left(u\right)$ to calculations over the t points of the mesh, that is, to vectors composed by the values of ${\phi }_i$, $i=1,\dots ,p-1$, and u at the points ${\dot{x}}_l$. In [4] (Chapter 5, Section 3), one can find a couple of Lemmata, Lemma 2 and Lemma 3, which allow constructing a vector $\{v\left({\dot{x}}_l\right)\}\in {\mathbb{R}}^t$, which forms a Markov system of vectors with the vectors $\{{\phi }_i\left({\dot{x}}_l\right)\}\in {\mathbb{R}}^t$, $i=1,\dots ,p-1$. To do so, it suffices to select (more or less randomly) $p-1$ values for the points $v\left({\dot{x}}_1\right),\dots v\left({\dot{x}}_{p-1}\right)$ and pick up a value $v\left({\dot{x}}_p\right)$, such that the following inequality holds:

$$\left|\begin{array}{cccc}{\phi }_1\left({\dot{x}}_1\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_1\right)& v\left({\dot{x}}_1\right)\\ {\phi }_1\left({\dot{x}}_2\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_2\right)& v\left({\dot{x}}_2\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_1\left({\dot{x}}_p\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_p\right)& v\left({\dot{x}}_p\right)\end{array}\right|>0,$$

(72)

which is always possible if ${\dot{x}}_1\in I_1$, given that:

$$\left|\begin{array}{ccc}{\phi }_1\left({\dot{x}}_1\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_1\right)\\ {\phi }_1\left({\dot{x}}_2\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_2\right)\\ \cdots & \cdots & \cdots \\ {\phi }_1\left({\dot{x}}_{p-1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{p-1}\right)\end{array}\right|>0,$$

(73)

as ${\phi }_1,\dots ,{\phi }_{p-1}$ form a Chebyshev system on $I_1$. If ${\dot{x}}_1\notin I_1$, then the determinant (72) will vanish regardless of the value of ${\dot{x}}_i$, $i=2,\dots ,p$, so we must start the process by calculating ${\dot{x}}_{p+1}$ such that:

$$\left|\begin{array}{cccc}{\phi }_1\left({\dot{x}}_2\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_2\right)& v\left({\dot{x}}_2\right)\\ {\phi }_1\left({\dot{x}}_3\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_3\right)& v\left({\dot{x}}_3\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_1\left({\dot{x}}_{p+1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{p+1}\right)& v\left({\dot{x}}_{p+1}\right)\end{array}\right|>0.$$

The values of the following coefficients $v\left({\dot{x}}_j\right)$ can be determined using a similar inequality:

$$\left|\begin{array}{cccc}{\phi }_1\left({\dot{x}}_i\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_i\right)& v\left({\dot{x}}_i\right)\\ {\phi }_1\left({\dot{x}}_{i+1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{i+1}\right)& v\left({\dot{x}}_{i+1}\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_1\left({\dot{x}}_{i+p-1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{i+p-1}\right)& v\left({\dot{x}}_{i+p-1}\right)\end{array}\right|>0.$$

(74)

Once all entries $v\left(\dot{x_l}\right)$ have been selected, one has to make the spline defined by $\{v\left({\dot{x}}_l\right)\}$ orthogonal to the adjoint eigenfunction functions ${\psi }_i$, also constructed as splines. The way to do it is by calculating the inner vector product:

$$〈{\psi }_i,v〉=\sum _{l=1}^{t-1}{\int }_{x_l}^{x_{l+1}}\left(v\left({\dot{x}}_l\right)\frac{{\dot{x}}_{l+1}-x}{{\dot{x}}_{l+1}-{\dot{x}}_l}+v\left({\dot{x}}_{l+1}\right)\frac{x-{\dot{x}}_l}{{\dot{x}}_{l+1}-{\dot{x}}_l}\right)\left({\psi }_i\left({\dot{x}}_l\right)\frac{{\dot{x}}_{l+1}-x}{{\dot{x}}_{l+1}-{\dot{x}}_l}+{\psi }_i\left({\dot{x}}_{l+1}\right)\frac{x-{\dot{x}}_l}{{\dot{x}}_{l+1}-{\dot{x}}_l}\right)r\left(x\right)dx,$$

for $i<p$, and putting

$$u\left(\dot{x_l}\right)=v\left(\dot{x_l}\right)-\sum _{i=1}^{p-1}〈{\psi }_i,v〉{\phi }_i\left({\dot{x}}_l\right),l=1,\dots ,t.$$

(75)

From (69) and (75), one can build the linear spline $\widehat{u}\left(x\right)$ to be used as the starting point for the calculations of $M_p^j\widehat{u}$ and ${\Delta }_p\left(M_p^j\widehat{u}\right)$.

A similar procedure can be used to find a starting function w for the calculations of ${\left(M_p^{\ast }\right)}^j\widehat{w}$ and ${\Delta }_p^{\ast }\left({\left(M_p^{\ast }\right)}^j\widehat{w}\right)$ in the adjoint problem, taking into consideration $\{{\psi }_i\left({\dot{x}}_l\right)\}$ and the orthogonality of w and ${\phi }_i$, $i=1,\dots ,p-1$.

2.2.2. How to Simplify the Calculation of the Collatz-Wielandt Numbers

The Lemmata of [4] that simplified the determination of the starting function $\widehat{u}$ are based on a Theorem by Fekete ([4], Theorem 8), which relates the minors of matrices whose columns are long vectors (in our case, $\{{\phi }_i\left({\dot{x}}_l\right)\}$, $i=1,\dots ,p-1$, and $\{M_p^{j}u\left({\dot{x}}_l\right)\}$) with the minors made by consecutive entries of these vectors (see (74)). Given that the calculation of the Collatz–Wielandt numbers requires basically finding $\omega$ such that ${\Delta }_p(M_p^{j+1}u-\omega M_p^{j}u;\dot{\mathbf{x}})$ is greater or lower than zero for all $\dot{\mathbf{x}}\in \Omega$, we can also exploit this property to reduce the number of combination of points ${\dot{\mathbf{x}}}_j$ in the simplex $\Omega$, where that determinant has to be calculated.

Thus, $\underline{r}(M_p,M_p^{j}u)$ will be given by the supremum of $\omega$ such that:

$$\left|\begin{array}{cccc}{\phi }_1\left({\dot{x}}_i\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_i\right)& M_p^{j+1}u\left({\dot{x}}_i\right)-\omega M_p^{j}u\left({\dot{x}}_i\right)\\ {\phi }_1\left({\dot{x}}_{i+1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{i+1}\right)& M_p^{j+1}u\left({\dot{x}}_{i+1}\right)-\omega M_p^j\left({\dot{x}}_{i+1}\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_1\left({\dot{x}}_{i+p-1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{i+p-1}\right)& M_p^{j+1}u\left({\dot{x}}_{i+p-1}\right)-\omega M_p^j\left({\dot{x}}_{i+p-1}\right)\end{array}\right|>0,i=1,\dots ,t-p+1,$$

(76)

whereas $\overline{r}(M_p,M_p^{j}u)$ will be given by the infimum of $\omega$ such that:

$$\left|\begin{array}{cccc}{\phi }_1\left({\dot{x}}_i\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_i\right)& M_p^{j+1}u\left({\dot{x}}_i\right)-\omega M_p^{j}u\left({\dot{x}}_i\right)\\ {\phi }_1\left({\dot{x}}_{i+1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{i+1}\right)& M_p^{j+1}u\left({\dot{x}}_{i+1}\right)-\omega M_p^j\left({\dot{x}}_{i+1}\right)\\ \cdots & \cdots & \cdots & \cdots \\ {\phi }_1\left({\dot{x}}_{i+p-1}\right)& \cdots & {\phi }_{p-1}\left({\dot{x}}_{i+p-1}\right)& M_p^{j+1}u\left({\dot{x}}_{i+p-1}\right)-\omega M_p^j\left({\dot{x}}_{i+p-1}\right)\end{array}\right|<0,i=1,\dots ,t-p+1.$$

(77)

3. An Example

Let us consider the problem:

$$y^{\left(4\right)}-\lambda \frac{1}{x}y=0,x\in [1,2],y\left(1\right)=y^{\prime \prime }\left(1\right)=y^{\left(3\right)}\left(1\right)=y\left(2\right)=0,$$

(78)

which matches the structure (1)–(3) since $r\left(x\right)=\frac{1}{x}>0$. The Green function of the associated problem:

$$\frac{{\partial }^{4}G(x,t)}{\partial x^4}=0,x\in [1,2],G(1,t)=\frac{{\partial }^{2}G(1,t)}{\partial x^2}=\frac{{\partial }^{3}G(1,t)}{\partial x^3}=G(2,t)=0,$$

(79)

is defined by:

$$G(x,t)=\left\{\begin{array}{cc}-\frac{{(2-t)}^3(x-1)}{6},\hfill & x\in [1,t),\hfill \\ -\frac{{(2-t)}^3(x-1)}{6}+\frac{{(x-t)}^3}{6},\hfill & x\in (t,2].\hfill \end{array}\right.$$

(80)

In order to apply the procedure described in this paper, one must verify that $G(x,t)$ is in fact a sign-regular kernel. In [3], one can find several theorems (Theorems 1.1–1.3) that provide algorithmically effective conditions for such an identification. However, in this case it is easier to use the following theorem of Kalafati–Gantmakher–Krein–Karlin (see for instance [13], Theorem 4 or [1], Theorem 8):

Theorem 6.

Suppose that the operator L on $[a,b]$ has the form (22) and the boundary conditions are:

$$\begin{array}{c}{\sum }_{k=1}^n{\alpha }_{ik}{L}_{k-1}y\left(a\right)=0,i=1,\dots ,m,\\ \\ {\sum }_{k=1}^n{\beta }_{jk}{L}_{k-1}y\left(b\right)=0,j=1,\dots ,n-m.\end{array}$$

(81)

Suppose also that all the non-zero minors of order m of the matrix

$${A=\parallel (-1)}^k{\alpha }_{ik}\parallel ,i=1,\dots ,m,k=1,\dots ,n,$$

have the same sign and the same holds for the minors of order $(n-m)$ of the matrix

$$B=\parallel {\beta }_{jk}\parallel ,j=1,\dots ,n-m,k=1,\dots ,n.$$

Then ${(-1)}^{n-m}G(x,t)$ is an oscillatory kernel provided that the boundary value problem is non-singular.

The problem (78) has the form (22) with ${\rho }_0={\rho }_1={\rho }_2={\rho }_3=1$, so that $L_{i}y\left(x\right)=y^{\left(i\right)}\left(x\right)$. The resulting matrix A is

$$A=\left(\begin{array}{cccc}-1& 0& 0& 0\\ 0& 0& -1& 0\\ 0& 0& 0& 1\end{array}\right),$$

(82)

whereas B is just the matrix

$$B=\left(\begin{array}{cccc}1& 0& 0& 0\end{array}\right).$$

(83)

The only non-zero minor of order $m=3$ of A is

$$\left|\begin{array}{ccc}-1& 0& 0\\ 0& -1& 0\\ 0& 0& 1\end{array}\right|=1,$$

(84)

whereas the only non-zero minor of order $n-m=1$ of B equals 1. In addition, the homogeneous boundary conditions are poised in Elias’ sense [36] (that is, the number of boundary conditions set on derivatives of an order lower than t is at least t, for $t=1,\dots ,n$). This implies that $\lambda =0$ is not an eigenvalue and the problem is not singular [36] (Lemma 10.3). Accordingly, one can apply the Kalafati–Gantmakher–Krein–Karlin theorem and conclude that $-G(x,t)$ of (79) is an oscillatory kernel, that is, ${ϵ}_i={(-1)}^i$ for all i.

Given that the coefficients of (78) are infinitely continuously differentiable, one can apply the method described in previous sections to determine all eigenfunctions and eigenvalues. As an example, we will calculate ${\lambda }_1$ and ${\lambda }_2$, as well as the corresponding eigenfunctions ${\phi }_1$ and ${\phi }_2$, and the adjoint eigenfunctions ${\psi }_1$ and ${\psi }_2$.

The operator $M_p$ can be calculated as:

$$\begin{array}{c}{M}_{p}u=\frac{x-1}{6}{\int }_1^x\frac{{(2-t)}^3}{t}u\left(t\right)dt-\frac{x^3}{6}{\int }_1^x\frac{u\left(t\right)}{t}dt+\frac{x^2}{2}{\int }_1^{x}u\left(t\right)dt\\ \\ -\frac{x}{2}{\int }_1^{x}tu\left(t\right)dt+\frac{1}{6}{\int }_1^x{t}^{2}u\left(t\right)dt+\frac{x-1}{6}{\int }_x^2\frac{{(2-t)}^3}{t}u\left(t\right)dt,x\in [1,2].\end{array}$$

(85)

The Green function $G^{\ast }(x,t)$ of the problem adjoint to (78) is linked to $G(x,t)$ of (80) by $G^{\ast }(x,t)=G(t,x)$. Therefore,

$$G^{\ast }(x,t)=\left\{\begin{array}{cc}-\frac{{(2-x)}^3(t-1)}{6},\hfill & t\in [1,x),\hfill \\ -\frac{{(2-x)}^3(t-1)}{6}-\frac{{(x-t)}^3}{6},\hfill & t\in (x,2]\hfill \end{array}\right.$$

(86)

and

$$\begin{array}{c}{M}_p^{\ast }u=\frac{{(2-x)}^3}{6}{\int }_1^2\frac{t-1}{t}u\left(t\right)dt+\frac{x^3}{6}{\int }_x^2\frac{u\left(t\right)}{t}dt-\frac{x^2}{2}{\int }_x^{2}u\left(t\right)dt\\ \\ +\frac{x}{2}{\int }_x^{2}tu\left(t\right)dt-\frac{1}{6}{\int }_x^2{t}^{2}u\left(t\right)dt,x\in [1,2].\end{array}$$

(87)

Since $p\le 2$ and $n=4$, in all cases it will suffice to use $q=1$ ($q·n=1\times 4=4>3\ge p$).

The execution of the procedure, using a partition of $[1,2]$ with size $h={10}^{-6}$, gives the Collatz–Wielandt numbers associated to $p=1$ and $p=2$, which are displayed in

Table 1. Calculation of Collatz–Wielandt numbers for the first and second eigenvalue.

j	$\underline{\mathit{r}}({\mathit{M}}_1,{\mathit{M}}_1^{\mathit{j}}\mathit{u})$	$\overline{\mathit{r}}({\mathit{M}}_1,{\mathit{M}}_1^{\mathit{j}}\mathit{u})$	$\underline{\mathit{r}}({\mathit{M}}_2,{\mathit{M}}_2^{\mathit{j}}\mathit{u})$	$\overline{\mathit{r}}({\mathit{M}}_2,{\mathit{M}}_2^{\mathit{j}}\mathit{u})$
1	0.0059445	0.0074422	0.0000013	0.0006649
2	0.0061165	0.0062012	0.0001568	0.000287
3	0.006121	0.0061244	0.0001863	0.0002122
4	0.0061212	0.0061213	0.0001902	0.0001961
5	0.0061212	0.0061212	0.0001907	0.0001917
6	0.0061212	0.0061212	0.0001912	0.0001912

. According to them, one deduces that ${\lambda }_1=-163.36711$ and ${\lambda }_2=-5229.8041$.

The resulting eigenfunctions ${\phi }_1$ and ${\phi }_2$, as well as the adjoint eigenfunctions ${\psi }_1$ and ${\psi }_2$, are shown in the Figure 1. They have been normalized to sup norm.

It is worth remarking that two phenomena observed during the numerical calculations:

• The calculation of ${\Delta }_p$ is very sensitive to rounding errors when $x_i,\dots ,x_{i+p}$ are close to the extremes a or b, if there are homogeneous boundary conditions set at these. The reason for that is that the values of ${\Delta }_p$ are zero or almost zero there. In these points of the partition, it makes sense to replace the calculation of ${\Delta }_p$ by the calculation of the equivalent determinant composed by the lowest derivatives of ${\phi }_i$ and $M_p^{j}u$, which do not vanish at the extreme.
• If u is not exactly orthogonal to ${\psi }_i$, $i=1,\dots ,p-1$, beyond a certain iteration it can happen that in the decomposition of $M^{j}u$ as a sum of terms of the form $\frac{1}{{\lambda }_i^j}〈u,{\psi }_i〉{\phi }_i$, the terms associated with those $i<p$ for which $〈u,{\psi }_i〉\ne 0$ start to get a size similar to that of the term $\frac{1}{{\lambda }_p^j}〈u,{\psi }_p〉{\phi }_p$, as anticipated by Remark 3. Further iterations will make $M_p^{j}u$ diverge from ${\phi }_p$ and get closer to the eigenfunctions ${\phi }_i$, $i<p$, for which $〈u,{\psi }_j〉{\phi }_j\ne 0$. The precise orthogonality is therefore key for the accuracy of the method.

4. Discussion

The results presented in this paper allow finding the n smallest eigenvalues (and their associated eigenfunctions) of boundary value problems with sign-regular Green functions, as well as the following ones provided that certain conditions on the functions $a_j\left(x\right)$ of L and $r\left(x\right)$ (namely, the absolutely continuity of their derivatives) are met.

The procedure is sequential in the sense that it requires running it for the first $p-1$ eigenvalues in order to use it to calculate the p-th one.

For each p, it can be summarized in the following algorithm, which assumes the knowledge of the $p-1$ previous eigenfunctions ${\phi }_i$ and the $p-1$ previous adjoint eigenfunctions ${\psi }_i$:

• Calculate $q>1$ such that $q·n>p$;
• Select $u\in P_p\setminus \{0\}$ using the process described in the Section 2.2.1, so that u is orthogonal to ${\psi }_i$, $i=1,\dots ,p-1$;
• Calculate $M_p^{j}u$ for $j\ge q$;
• Calculate the Collatz–Wielandt numbers using (76) and (77) and the considerations described in the Section 2.2.2. These will be bounds for the inverse of the absolute value of the eigenvalue ${\lambda }_p$, whose sign is determined by ${ϵ}_{p-1}$ and ${ϵ}_p$ of (7), the error in the calculation being given by the difference $\overline{r}(M_p,M_p^{j}u)-\underline{r}(M_p,M_p^{j}u)$. Due to the convergence of the Collatz-Wielandt numbers as j increases, the eigenvalue ${\lambda }_p$ can be estimated with as much accuracy as desired;
• The quotient $\frac{M_p^{j}u}{\overline{r}{(M_p,M_p^{j}u)}^j}$ will converge to ${\phi }_p$ as the iteration index j grows;
• Select $w\in P_p^{\ast }\setminus \{0\}$ using the process described in the Section 2.2.1, so that w is orthogonal to ${\phi }_i$, $i=1,\dots ,p-1$;
• Calculate ${\left(M_p^{\ast }\right)}^{j}w$ for $j\ge q$. The quotient $\frac{{\left(M_p^{\ast }\right)}^{j}w}{\overline{r}{(M_p,M_p^{j}u)}^j}$ will converge to ${\psi }_p$ as the iteration index j grows;
• Once obtained, this ${\psi }_p$ will have to be normalized by dividing it by $〈{\phi }_p,{\psi }_p〉$ so as to satisfy (13).

The method, however, has also some limitations, mainly:

• For $p>n-1$ it requires some absolutely continuity conditions on $r\left(x\right)$ and $a_j\left(x\right)$ in order to be applicable, in addition to those required for the existence of the adjoint problem;
• The method depends on the accuracy of the calculation of the $p-1$ previous eigenfunctions, given that the determinant ${\Delta }_p\left(u\right)$ depends on them. For greater values of p one can expect more accumulated errors in ${\phi }_i$, $i<p$, and potentially bigger errors in ${\lambda }_p$ and ${\phi }_p$;
• As the size of the determinant ${\Delta }_p\left(M_p^{j}u\right)$ grows with p, the computations become more complicated as p increases. The use of optimized algorithms for the calculation of the determinants is key to reduce this problem;
• In a practical scenario, the calculation of the Collatz–Wielandt numbers needs to be performed on a mesh of the simplex ${\Omega }^{\ast }$, as described in the Section 2.2.2. This raises some questions about the validity of these numbers in other points of the simplex. The problem can be addressed by avoiding the use of the supremum and infimum in the calculation of the Collatz–Wielandt numbers, so that the difference between $M_p^{j+1}u$ and $\omega M_p^{j}u$ is not zero at any points of the mesh, but a proper analysis on the effect of the interpolation error needs to be performed.

In any case, we the authors believe it can be a practical alternative for the calculation of eigenpairs, especially for lower values of p, and also a source for later work, since an aspect not stressed in this paper is that this approach also allows the determination of the existence of the eigenvalues ${\lambda }_p$, the Markov character of the sequence ${\phi }_i$ and, with the right conditions on $a_j\left(x\right)$ and $r\left(x\right)$, the algebraical and geometrical simplicity of each eigenvalue and their different absolute value than the others. These properties are widely known from the previous literature (the reason our focus has been more practical, on the determination of ${\lambda }_p$), but in this sense it is worth highlighting that the approach used here differs from that used in the classical papers [1,4,14], to give a few examples. These based their analysis on expressions such as (32), whose iterative application p times leads to strictly positive kernels, and applied it to the cone of positive functions in $C{[a,b]}^p$, making use of Krein–Rutman cone theory and other classical results of Schur. While this approach has the advantage of not imposing extra conditions on $a_j\left(x\right)$ and $r\left(x\right)$, it does not lend itself easily to work with cone interiors, which are key for the calculation of the Collatz–Wielandt numbers (or rather, to their relationship with the eigenvalue ${\lambda }_p$), as Theorems 2 and 3 show.

To complete this paper, let us mention several areas of interest for future research:

• Explore ways of extending the procedure to the case $r\left(x\right)=0$ in a set of points of $[a,b]$. The effect of this is that the set $I_1$, on which the eigenfunctions ${\phi }_i$ form a Chebyshev system, does not contain these points where $r\left(x\right)$ vanishes, complicating the extensions of some of the results presented here;
• Analyze the effect of the interpolation error committed in the calculation of each eigenvalue ${\lambda }_p$ by performing the calculation of Collatz–Wielandt numbers only in the points of the mesh $\{{\dot{x}}_l\}$;
• Simplify or categorize the conditions defined by Stepanov for the sign-regularity of the Green function [3] so that their validation does not always require the calculation of Wronskians of the solutions of $Ly=0$ under certain boundary conditions. This would allow an easier identification of sign-regular problems, where the procedures of this manuscript can be applied;
• Last but not least, we have made use of the cone ${\mathcal{P}}_p$ as it allows the fixing of conditions for such a cone to be solid and for $M_p$ to map $P_p$ into its interior. However, this does not exclude the existence of other solid cones $K_p$ on which to apply Theorem 2. It would be very interesting to find some examples of these, in order to relax the hypothesis that ${\mathcal{P}}_p$ demands on r and $a_j$.