Local solvability and stability of an inverse spectral problem for higher-order differential operators

In this paper, we for the first time prove local solvability and stability of an inverse spectral problem for higher-order ($n>3$) differential operators with distribution coefficients. The inverse problem consists in the recovery of differential equation coefficients from $(n-1)$ spectra and the corresponding weight numbers. The proof method is constructive. It is based on the reduction of the nonlinear inverse problem to a linear equation in the Banach space of bounded infinite sequences. We prove that, under a small perturbation of the spectral data, the main equation remains uniquely solvable. Furthermore, we estimate the differences of the coefficients in the corresponding functional spaces.

We study the inverse spectral problem that consists in the recovery of the coefficients (p k ) n−2 k=0 from the eigenvalues {λ l,k } l≥1 and the weight numbers {β l,k } l≥1 of the boundary value problems L k , k = 1, n − 1, for equation (1.1) with the corresponding boundary conditions y (j) (0) = 0, j = 0, k − 1, y (s) (1) = 0, s = 0, n − k − 1. (1. 2) The main goal of this paper is to prove local solvability and stability of the inverse problem. This is the first result of such kind for arbitrary-order differential operators with distribution coefficients.

Historical background
Spectral theory of linear ordinary differential operators has a fundamental significance for mathematics and applications. For n = 2, equation (1.1) turns into the Sturm-Liouville (onedimensional Schrödinger) equation y ′′ + p 0 y = λy, which models various processes in quantum and classical mechanics, material science, astrophysics, acoustics, electronics. The third-order differential equations are applied to describing thin membrane flow of viscous liquid [1] and elastic beam vibrations [2]. The third-order spectral problems also arise in the integration of the nonlinear Boussinesq equation by the inverse scattering transform [3]. The fourth-order and the sixth-order linear differential operators appear in geophysics [4] and in vibration theory [5,6]. Therefore, the development of general mathematical methods for investigation of spectral problems for arbitrary-order linear differential operators is fundamentally important.
Inverse problems of spectral analysis consist in the reconstruction of differential operators by their spectral characteristics. Such problems have been studied fairly completely for Sturm-Liouville operators −y ′′ + q(x)y with regular (integrable) potentials q(x) (see the monographs [7][8][9][10][11] and references therein) as well as with distribution potentials of class W −1 2 (see, e.g., the papers [12][13][14][15][16][17][18] and [19] for a more extensive bibliography). The basic results for the inverse Sturm-Liouville problem were obtained by the method of Gelfand and Levitan [20], which is based on transformation operators. However, inverse problems for differential operators of higher orders n > 2 are significantly more difficult for investigation, since the Gelfand-Levitan method does not work for them. Therefore, Yurko [21][22][23] has developed the method of spectral mappings, which is based on the theory of analytic functions. The central place in this method is taken by the contour integration in the complex plane of the spectral parameter λ of some meromorphic functions (spectral mappings), which were introduced by Leibenson [24,25]. Applying the method of spectral mappings, Yurko has created the theory of inverse spectral problems for arbitrary order differential operators with regular coefficients and also with the Bessel-type singularities on a finite interval and on the half-line (see [21-23, 26, 27]). Another approach was developed by Beals et al [28,29] for inverse scattering problems on the line, which are essentially different from the spectral problems on a finite interval.
In [30], Mirzoev and Shkalikov have proposed a regularization approach for higher-order differential equations with distribution coefficients. This has motivated a number of researchers to study solutions and spectral theory for such equations (see the bibliography in the recent papers [31,32]). Inverse spectral problems for higher-order differential operators with distribution coefficients have been investigated in [32][33][34][35][36]. In particular, the uniqueness theorems were proved in [32][33][34]. The paper [35] is concerned with a reconstruction approach, based on developing the ideas of the method of spectral mappings. In [36], the necessary and sufficient conditions for the inverse problem solvability were obtained for a third-order differential equation.
In this paper, we focus on local solvability and stability of the inverse problem. These aspects for the Sturm-Liouville operators were studied in [8-10, 12, 15, 16, 18, 37-46] and many other papers. Local solvability has a fundamental significance in the inverse problem theory, especially for such problems, for which global solvability theorems are absent or contain hardto-verify conditions. Stability is important for justification of numerical methods. For the higher-order differential equation (1.1) on a finite interval, stability of the inverse problem in the uniform norm was proved by Yurko for regular coefficients p k ∈ W k 1 [0, 1] (see Theorem 2.3.2 in [23]). Furthermore, local solvability and stability in the L 2 -norm were formulated without proofs as Theorems 2.3.4 and 2.5.3 in [23]. For distribution coefficients, local solvability and stability theorems have been proved in [12] for n = 2 and in [36] for n = 3. However, to the best of the author's knowledge, there were no results in this direction for n ≥ 4.
It is worth mentioning that, for the Sturm-Liouville operators with distribution potentials of the classes W α 2 [0, 1], α > −1, Savchuk and Shkalikov [15] obtained the uniform stability of the inverse spectral problem. This result was extended by Hryniv [16] to α = −1 by another method. Recently, the uniform stability of inverse problems was also proved for some nonlocal operators (see [47][48][49]). However, the approaches of the mentioned studies do not work for higher-order differential operators. Thus, the uniform stability for them is an open problem, which is not considered in this paper.

Main results
In this paper, we for the first time prove local solvability and stability of an inverse problem for equation (1.1) and, in addition, obtain some sufficient conditions for global solvability. To the best of the author's knowledge, these aspects were not considered in previous studies for arbitrary-order differential equations with distribution coefficients, so our results are fundamentally novel. In order to prove the main theorems, we derive reconstruction formulas for coefficients (p k ) n−2 k=0 , which are also novel for the considered class of operators. Let us formulate the inverse problem for equation (1.1) and the corresponding local solvability and stability theorem.
Consider the problems L k given by (1.1) and (1.2). For each k ∈ {1, 2, . . . , n − 1}, the spectrum of the problem L k is a countable set of eigenvalues {λ l,k } l≥1 . They are supposed to be numbered counting with multiplicities according to the asymptotics obtained in [50]: where θ k are some constants independent of (p s ) n−2 s=0 and {κ l,k } ∈ l 2 . The weight numbers {β l,k } l≥1 k=1,n−1 will be defined in Section 3. Thus, the spectral data {λ l,k , β l,k } l≥1, k=1,n−1 are generated by the coefficients p = (p k ) n−2 k=0 of equation (1.1). In this paper, we confine ourselves to p of some class W with certain restrictions on the spectra. Definition 1.1. We say that p = (p k ) n−2 k=0 belongs to the class W if Thus, we study the following inverse problem. Problem 1.2. Given the spectral data {λ l,k , β l,k } l≥1, k=1,n−1 , find the coefficients p = (p k ) n−2 k=0 ∈ W . Problem 1.2 generalizes the classical inverse Sturm-Liouville problem studied by Marchenko [10] and Gelfand and Levitan [20] (see Example 3.3). The uniqueness for solution of Problem 1.2 follows from the results of [34,35] (see Section 3 for details). Note that, if the conditions (W-2) and (W-3) are violated, then the spectral data {λ l,k , β l,k } l≥1, k=1,n−1 do not uniquely specify the coefficients p and additional spectral characteristics are needed. In particular, for multiple eigenvalues in the case n = 2 the generalized weight numbers were defined in [51,52]. For n = 3 and λ l,1 = λ l,2 , the additional weight numbers γ l were used in [36]. For higher orders n, the situation becomes much more complicated, so in this paper we confine ourselves to the class W . Anyway, in view of the eigenvalues asymptotics (1.3), the assumptions (W-2) and (W-3) hold for all sufficiently large indices l.
Along with p, we consider an analogous vectorp = (p k ) n−2 k=0 ∈ W . We agree that, if a symbol α denotes an object related to p, then the symbolα with tilde will denote the analogous object related top. The main result of this paper is the following theorem on local solvability and stability of Problem 1.2. Theorem 1.3. Letp = (p k ) n−2 k=0 ∈ W be fixed. Then, there exists δ > 0 (which depends onp) such that, for any complex numbers {λ l,k , β l,k } l≥1, k=1,n−1 satisfying the inequality there exists a unique p = (p k ) n−2 k=0 with the spectral data {λ l,k , β l,k } l≥1, k=1,n−1 . Moreover, where the constant C depends only onp and δ. Theorem 1.3 generalizes the previous results of [12] for n = 2 and of [36] for n = 3. However, for n ≥ 4, to the best of the author's knowledge, Theorem 1.3 is the first existence result for the inverse problem solution in the case of distribution coefficients.
The proof of Theorem 1.3 is based on the constructive approach of [35]. Namely, we reduce the nonlinear inverse problem to the so-called main equation, which is a linear equation in the Banach space of bounded infinite sequences. The unique solvability of the main equation follows from the smallness of δ. Furthermore, we derive reconstruction formulas for the coefficients (p k ) n−2 k=0 in the form of infinite series. The crucial step in the proof is establishing the convergence of those series in the corresponding spaces W k−1 2 [0, 1] (including the space of generalized functions W −1 2 [0, 1]). In order to prove the convergence, we rigorously analyze the solution of the main equation and obtain the precise estimates for the Weyl solutions. Along with Theorem 1.3, we also prove Theorem 6.1 on global solvability of the inverse problem under several requirements on an auxiliary model problem.
The paper is organized as follows. In Section 2, we discuss the regularization of equation (1.1) and provide other preliminaries. In Section 3, the weight numbers are defined and the properties of the spectral data are described. In Section 4, we derive the main equation basing on the results of [35]. In Section 5, the reconstruction formulas for the coefficients (p k ) n−2 k=0 are obtained. In Section 6, we prove solvability and stability of the inverse problem. Section 7 contains concluding remarks.

Preliminaries
In this section, we explain in which sense we understand equation (1.1). For this purpose, an associated matrix and quasi-derivatives are introduced. In addition, we provide other preliminaries. We begin with some notations, which are used throughout the paper: • δ j,k is the Kronecker delta.
• C j k = k! j!(k−j)! are the binomial coefficients.
• In estimates, the same symbol C is used for various positive constants that do not depend on x, l, λ, etc.
Fix any function σ ∈ L 2 [0, 1] such that p 0 = σ ′ . Define the associated matrix F (x) = [f k,j (x)] n k,j=1 for the differential expression ℓ n (y) by the formulas All the other entries f k,j are assumed to be zero. Clearly, f k,j ∈ L 2 [0, 1]. Using the matrix function F (x), introduce the quasi-derivatives and the domain Due the the special structure of the associated matrix F (x), we have  Thus, for y ∈ D F , ℓ n (y) is a function of L 1 [0, 1] and the relation ℓ n (y) = y [n] gives the regularization of this differential expression. We call a matrix function F (x) an associated matrix of the differential expression ℓ n (y) if F (x) defines the quasi-derivatives y [k] and the domain D F so that the assertion of Lemma 2.1 holds. A function y is called a solution of equation (1.1) if y ∈ D F and ℓ n (y) = λy a.e. on (0, 1).
Following the technique of [35, Section 2], we consider along with F (x) the matrix function Using (2.1), we obtain and all the other entries f ⋆ k,j equal zero. For example, for n = 6, we have Using the matrix function F ⋆ (x), define the quasi-derivatives and the differential expression Note that, in (2.7) and (2.8), we use the quasi-derivatives defined by (2.6). Below we call a function z a solution of the differential equation if z ∈ D F ⋆ and the equality (2.9) holds a.e. on (0, 1). Throughout this paper, we always use the quasi-derivatives (2.2) for functions of D F and the quasi-derivatives (2.6) for functions of D F ⋆ . For z ∈ D F ⋆ , relations (2.5) and (2.6) imply Therefore, one can show by induction that . For z ∈ D F ⋆ and y ∈ D F , define the Lagrange bracket: (2.12) Then, the Lagrange identity holds (see [35,Section 2]): In particular, if z and y solve the equations ℓ ⋆ (z) = µz and ℓ(y) = λy, respectively, then we get Substituting (2.10) and (2.11) into (2.12), we derive the relation where all the derivatives are regular, since z, y ∈ W n−1 Remark 2.2. The associated matrix F (x) given by (2.1) regularizes the differential expression ℓ n (y) only for n ≥ 3. For n = 2, the associated matrix is constructed in a different way (see [53]): Nevertheless, the main result of this paper (Theorem 1.3) holds for n = 2 and, moreover, has been already proved in [12] for real-valued potentials. Therefore, below in the proofs, we confine ourselves to the case n ≥ 3 and use the associated matrix (2.1). For n = 2, the proofs are valid with minor modifications. Remark 2.3. For regularization of the differential expression ℓ n (y), different associated matrices can be used (see [30,32,34]). However, it has been proved in [34] that the spectral data {λ l,k , β l,k } l≥1, k=1,n−1 do not depend on the choice of the associated matrix.

Spectral data
In this section, we discuss the properties of the spectral characteristics for the boundary value problems L k , k = 1, n − 1. In particular, the weight numbers {β l,k } are defined as the residues of some entries of the Weyl-Yurko matrix.
For k = 1, n, denote by C k (x, λ) and Φ k (x, λ) the solutions of equation (1.1) satisfying the initial conditions C and the boundary conditions respectively. The functions {Φ k (x, λ)} n k=1 are called the Weyl solutions of equation (1.1). Let us summarize the properties of the solutions C k (x, λ) and Φ k (x, λ). For details, see [33]. The functions C k (x, λ), k = 1, n, are uniquely defined as solutions of the initial value problems, and they are entire in λ for each fixed x ∈ [0, 1] together with their quasi-derivatives C k (x, λ), j = 1, n − 1. The Weyl solutions Φ k (x, λ), k = 1, n, and their quasi-derivatives are meromorphic in λ. Furthermore, the fundamental matrices are related to each other as follows: is a unit lower-triangular matrix, whose entries under the main diagonal are meromorphic in λ, and the poles of the k-th column coincide with the zeros of ∆ k,k (λ).
On the other hand, the zeros of ∆ k,k (λ) coincide with the eigenvalues of the problem L k for equation (1.1) with the boundary conditions (1.2) for each k = 1, n − 1. Therefore, under the assumption (W-2) of Definition 1.1, all the poles of M(λ) are simple. The Laurent series at λ = λ l,k has the form where M j (λ l,k ) are the corresponding (n × n)-matrix coefficients. Define the weight matrices  The structural properties of the weight matrices N (λ l,k ) = [N s,j (λ l,k )] n s,j=1 are similar to the ones for the case of regular coefficients (see [23,35]). In view of (3.1), N s,j (λ l,k ) = 0 for s ≤ j. Moreover, under the condition (W-3), N s,j (λ l,k ) = 0 for s > j + 1. Thus, the only non-zero entries of N (λ l,k ) are N j+1,j (λ l,k ) for such j that ∆ j,j (λ l,k ) = 0. Therefore, instead of the weight matrices N (λ l,k ), it is sufficient to use the weight numbers .
Indeed, if λ l 1 ,k 1 = λ l 2 ,k 2 = · · · = λ lr,kr is a group of equal eigenvalues (of different problems L k j ), that is maximal by inclusion, we have where E i,j denotes the matrix with the unit entry at the position (i, j) and all the other entries equal zero. Hence, Proposition 3.1 implies the following corollary.

Main equation
In this section, we reduce Problem 1.2 to a linear equation in the Banach space m of bounded infinite sequences. First, we deduce an infinite system of linear equations. Second, this system is transformed to achieve the absolute convergence of the series by the method of [23]. Although we rely on the general approach of [35], the construction of the main equation is simplified because of the separation condition (W-3).
Consider the two coefficient vectors p = (p k ) n−2 k=0 andp = (p k ) n−2 k=0 of the class W . Note that the differential expressionl n (y) with the coefficientsp has the associated matrixF (x), which can be different from F (x), so the corresponding quasi-derivatives differ. The matrixF ⋆ (x) and the corresponding quasi-derivatives are defined analogously to F ⋆ (x) and (2.6), respectively.
In [35], the following infinite linear system has been obtained: Our next goal is to combine the terms in (4.3) for achieving the absolute convergence of the series. Introduce the numbers which characterize the "distance" between the spectral data {λ l,k , β l,k } l≥1, k=1,n−1 and {λ l,k ,β l,k } l≥1, k=1,n−1 of p andp, respectively. The asymptotics (1.3) and (3.3) imply that {ξ l } ∈ l 2 . In addition, define the functions w l,k (x) := l −k exp(−xl cot(kπ/n)), which characterize the growth of ϕ l,k,ε (x): the estimate |ϕ l,k,ε (x)| ≤ Cw l,k (x) holds by Lemma 7 in [35]. Pass to the new variables Analogously to ψ l,k,ε (x), defineψ l,k,ε (x). For brevity, denote v = (l, k, ε), v 0 = (l 0 , k 0 , ε 0 ), v, v 0 ∈ V . Then, relation (4.3) is transformed into where ( In view of the estimates (4.8), ψ(x),ψ(x) ∈ m and the operatorR(x) is bounded in m for each fixed x ∈ [0, 1]. Denote by I the identity operator in m. Then, the system (4.7) can be represented as a linear equation in the Banach space m: (4.9) Equation (4.9) is called the main equation of Problem 1.2. We have derived (4.9) under the assumption that {λ l,k , β l,k } and {λ l,k ,β l,k } are the spectral data of the two problems with the coefficients p = (p k ) n−2 k=0 andp = (p k ) n−2 k=0 , respectively. Anyway, the main equation (4.9) can be used for the reconstruction of p by {λ l,k , β l,k }. Indeed, one can choose an arbitraryp ∈ W , findψ(x) andR(x) by usingp, {λ l,k ,β l,k }, and {λ l,k , β l,k }, then find ψ(x) by solving the main equation. In order to find p from ψ(x), we need reconstruction formulas, which are obtained in the next section.

Reconstruction formulas
In this section, we derive formulas for recovering the coefficients (p k ) n−2 k=0 from the solution ψ(x) of the main equation (4.9). For the derivation, we use the special structure of the associated matrices F (x) and F ⋆ (x). The arguments of this section are based on formal calculations with infinite series. The convergence of those series will be rigorously studied in the next section.

Constructψ
, v ∈ V , by solving the main equation (4.9). 5. For (l, k, ε) ∈ V , determine ϕ l,k,ε (x) by (5.1) andη l,k,ε (x) by (5.4). 6. For s = n − 2, n − 3, . . . , 1, 0, find p s (x) by formula (5.9), in which t r,s (x) and T j 1 ,j 2 (x) are defined by (5.5) and (5.6), respectively. Procedure 5.1 will be used in the next section for proving Theorem 1.3. In general, there is a challenge to choose a model problemp so that the series for p s converge in the corresponding spaces. Note that steps 1-5 work for anyp ∈ W , since for them the estimate {ξ l } ∈ l 2 is sufficient. But the situation differs for step 6. In Section 6, we prove the validity of step 6 in the case {l n−2 ξ l } ∈ l 2 .

Solvability and stability
In this section, we prove the following theorem on the solvability of Problem 1.2.
Then, {λ l,k , β l,k } l≥1, k=1,n−1 are the spectral data of some (unique) p = (p k ) n−2 k=0 ∈ W . Theorem 6.1 provides sufficient conditions for global solvability of the inverse problem. Theorem 1.3 on local solvability and stability will be obtained as a corollary of Theorem 6.1. Thus, Theorem 6.1 plays an auxiliary role in this paper but also has a separate significance. The proof of Theorem 6.1 is based on Procedure 5.1. We investigate the properties of the solution ψ(x) of the main equation and prove the convergence of the series in (5.6) and (5.9) in the corresponding spaces of regular and generalized functions. This part of the proofs is the most difficult one, since the series converge in different spaces and precise estimates for the Weyl solutions are needed. Finally, we show that the numbers {λ l,k , β l,k } satisfying the conditions of Theorem 6.1 are the spectral data of the coefficients p = (p k ) n−2 k=0 reconstructed by formulas (5.9). In the end of this section, we prove Theorem 1.3.
Proceed to the proof of Theorem 6.1. Let {λ l,k , β l,k } l≥1, k=1,n−1 andp satisfy the hypotheses (S-1)-(S-5). We emphasize that {λ l,k , β l,k } are not necessarily the spectral data corresponding to some p. We have to prove this.
Proof. The assertion of the lemma immediately follows from the definition (5.4), the estimates (6.7) for the weight numbers, and Corollary 6.3.
Proceed to the investigation of the convergence for the series t r,s (x) and T j 1 ,j 2 (x) in the reconstruction formulas (5.9). We rely on Proposition 6.6, which (due to our notations) readily follows from Lemma 8 in [35] and its proof. Proposition 6.6 ( [35]). The following statements hold.
The constants a j 1 ,j 2 ,l,k are explicitly found in the proof of Lemma 8 in [35]. However, we do not provide them here in order not to introduce many additional notations. Moreover, explicit formulas for a j 1 ,j 2 ,l,k are not needed in the proofs.
Below, similarly to the series in Proposition 6.6, we consider the series T j 1 ,j 2 (x) with the brackets: Moreover, we agree that we understand the summation of several series T j 1 ,j 2 (x) (in particular, in (5.5) and in (5.9)) in the sense Using Lemmas 6.4, 6.5 and Proposition 6.6, we get the lemma on the convergence of T j 1 ,j 2 (x).
Lemma 6.7. The following statements hold.
Since p s (x) =p s (x) + S 1 (x) + S 2 (x) + S 3 (x), we arrive at the assertion of the lemma.
Thus, we have obtained the vector p = (p k ) n−2 k=0 by the reconstruction formulas. It remains to prove the following lemma. Lemma 6.9. The spectral data of p coincide with {λ l,k , β l,k } l≥1, k=1,n−1 .
Proof. The proof is based on the approximation approach which has been considered in [36] for n = 3 in detail. The proof for higher orders n is similar, so we omit the technical details and just outline the main idea.
Then, by using {λ N l,k , β N l,k } instead of {λ l,k , β l,k }, one can construct the main equation (6.13) analogously to (4.9). It can be shown that equation (6.13) is uniquely solvable for sufficiently large values of N. Therefore, one can find the functions ϕ N l,k,ε (x), (l, k, ε) ∈ V , by using the solution ψ N (x) similarly to (5.1). Then, construct the functions {Φ N k (x, λ)} n k=1 and p N = (p N k ) n−2 k=0 analogously to (5.2) and (5.9), respectively. The advantage of the "truncated" data {λ N l,k , β N l,k } over {λ l,k , β l,k } is that the series in (5.2) and (5.9) are finite, so one can show by direct calculations that {Φ N k (x, λ)} n k=1 are the Weyl solutions of equation (1.1) with the coefficients p N and deduce that {λ N l,k , β N l,k } are the spectral data of p N . At this stage, the assumptions (S-1)-(S-3) of Theorem 6.1 are crucial. If these assumptions do not hold, then one needs additional data to recover p N . Next, using (5.9), we show that Furthermore, it can be shown that the spectral data depend continuously on the coefficients p = (p k ) n−2 k=0 . This concludes the proof. Let us summarize the arguments of this section in the proof of Theorem 6.1.
Proof of Theorem 1.3. Let us show that, if data {λ l,k , β l,k } l≥1, k=1,n andp satisfy the conditions of Theorem 1.3 for sufficiently small δ > 0, then they also satisfy the hypothesis of Theorem 6.1.