A new method for high-order boundary value problems

*Correspondence: zhych0314@163.com 1School of Applied Science and Civil Engineering, Zhuhai Campus, Beijing Institute of Technology, Zhuhai, China Abstract This paper presents a numerical algorithm for solving high-order BVPs. We introduce the construction method of multiscale orthonormal basis inWm 2 [0, 1]. Based on the orthonormal basis, the numerical solution of the boundary value problem is obtained by finding the ε-approximate solution. In addition, the convergence order, stability, and time complexity of the method are discussed theoretically. At last, several numerical experiments show the feasibility of the proposed method.


Introduction
High-order BVPs are important mathematical models in the field of electro-magnetics, fluid mechanics, and material science. Many problems in the theory of elastic stability can be handled by BVPs [1]. It is difficult to find the analytic solutions of high-order BVPs because of the complexity of the systems, many numerical algorithms for high-order BVPs have been proposed in recent years. The multistage integration method is an important method to solve the numerical solution of high-order models by reducing the order gradually [2][3][4][5]. Ref. [6][7][8] discuss the existence of solutions to higher-order differential equations. Cao [9] solved a class of high-order fractional ordinary differential equations by the quadratic interpolation function method. The collocation method proposed by [10] and the orthonormal Bernstein polynomials method proposed by Mirzaee [11,12] can solve high-order linear complex differential equations effectively. Mirzaee et al. [13][14][15][16][17][18][19][20][21][22] proposed a variety of numerical algorithms for solving high-order integro-differential equations. Many scholars have also proposed many methods in the field of numerical solution of high-order partial differential equations [23][24][25]. Reproducing kernel space is an important Banach space which has been used in the field of numerical analysis. The reproducing kernel methods are used in the numerical solutions of high-order models, singular BVPs, and interface problems [26][27][28][29][30][31].
In this paper, we construct a set of multiscale orthonormal bases based on the idea of wavelet in the reproducing kernel space. This set of bases is orthonormal, which can improve the computational efficiency. For the numerical solution of differential equations, many literature works, such as [32,33], use the idea of ε-approximate solution. ε-approximate solution provided the stability of the algorithm, good order of convergence in the calculation method. In this article, we construct the multiscale method for the following high-order boundary value problems (BVPs): The paper is organized as follows: In Sect. 2, we construct a set of multiscale orthonormal bases in the reproducing kernel space W m 2 [0, 1]. In Sect. 3, we introduce a method to obtain the numerical solution of BVPs by finding an ε-approximate solution, and verify the existence of the ε-approximate solution. In Sect. 4, the convergence, stability, and complexity of this method are discussed. In Sect. 5, we report the numerical result obtained by the present method and compare this method with other previous methods.

Multiscale orthonormal basis
In this section, the reproducing kernel space is defined and a set of multiscale orthonormal bases is constructed. This knowledge is very useful in the following article.
In order to solve Eq. (1.1), this paper constructs a set of orthonormal bases in W m   Proof We just prove the orthogonality and completeness. First, orthonormality. Obviously, . By the definition of inner product and Eq. (2.2)-Eq. (2.4), we can obtain a = 1 k! .

ε-approximate solution of high-order BVPs
In this section, we give the ε-approximation of Eq. (1.1) and get the numerical solution of BVPs by finding the ε-approximate solution of Eq. (1.1). Put

From Eq. (3.2) and Eq. (3.3), it follows that
where M is a positive constant.
Then Eq. (1.1) is equivalent to the following equation: Zhang [32] proposed the ε-approximate theory of second-order differential equations, now we define the ε-approximate solution of Eq. (3.4) based on the idea.
is the ε-approximate solution of Eq. According to [32], the unique solution of Eq. (3.8) is the minimum point of J.

Theoretical analysis
In this section, the properties of the algorithm, such as uniform convergence, stability, and complexity, are introduced.

Convergence analysis
That is, We can obtain (c i,k ) 2 ≤ ( 1 2 ) 3i C 1 . In fact, According to Holder's inequality, Then So then where M is a constant.
where C is a constant, K(x, y) is the reproducing kernel of W m 2 .
From Theorem 4.1, u n uniformly converges to u.

Stability analysis
It is well known that if A is a reversible symmetrical matrix, then the condition number of A is where λ 1 and λ n are the maximum and minimum eigenvalues of A respectively.
Obviously, A of Eq. (3.8) is an invertible symmetric matrix. Therefore, in order to prove the stability of the algorithm, we can first prove the boundedness of the eigenvalues.
Proof According to λx = Ax, Multiply both sides of (4.2) by x i (i = 1, 2, . . . , n) and add up to get Without loss of generality, put u W m 2,0 = 1. According to the inverse operator theorem [21], To sum up From Lemma 4.1, we get That is, the presented method is stable.

Complexity analysis
Complexity analysis includes time complexity and space complexity. But ultimately, it is the time efficiency of the algorithm that matters. As long as the algorithm does not take up storage space that is unacceptable to the computer. So this part analyzes the time complexity.

Theorem 4.2 The time complexity of the algorithm is O(n 3 ).
Proof There are four steps to calculate ε-the approximate solution u n (x) of (3.1). First, the calculation of matrix A n in Eq. (3.8). The matrix A n is Set the number of multiplication required to compute Lρ k , Lρ j and B l ρ k B l ρ j as C 1 , C 2 respectively, C 1 , C 2 are constant. So each term of A n is evaluated C 1 + mC 2 times. Since A n is a symmetric matrix, we only need to consider the calculation of the main diagonal and above elements. The first row of A n is evaluated n times, the second row n -1 times, and so on, so that the total number of multiplications required in calculation of A n is n(n + 1) 2 (C 1 + mC 2 ).
Second, the calculation of vector b n in Eq. Third, solve Eq. (3.8). We solve the system by Gaussian elimination. From the mathematical knowledge, Gaussian elimination requires operations n(n + 1)(2n + 1) 6 .
Forth, calculation u n . When calculating the u n , the total number of multiplications required is n.

Numerical experiments
In this section, we give several numerical experiments to verify the effectiveness of the proposed algorithm. We denote by u n (x) the approximation to the exact solution u(x) obtained by the numerical schemes in the present work, and we measure the errors in the following sense: where n is the number of bases. C.R. represents the convergence order. All numerical experiments are computed by Mathematica 9.0.
Example 5.1 Ref. [34] mentioned that in order to get the shear deformation of sandwich beams, consider the following third-order BVP: where the physical constants are k = 5 and r = 1. The function u(x) shows the shear deformation of sandwich beams. The analytic solution of this problem is u(x) = r(k(2x -1) sinh(kx) + 2 cosh(kx) tanh( k 2 )) 2k 3 .
The numerical results are given in Table 1.

Conclusion
In summary, this study used a set of multiscale orthonormal bases to find the εapproximate solutions of higher-order BVPs. This paper not only demonstrates the convergence and stability in theory, but also demonstrates the feasibility of the method through numerical experiments. Through theoretical analysis and numerical experiments, this method can be extended to solve general linear models, such as linear integral equations, differential equations, and fractional differential equations.