Highly Accurate Compact Finite Difference Schemes for Two-Point Boundary Value Problems with Robin Boundary Conditions

: In this study, a high-order compact ﬁnite difference method is used to solve boundary value problems with Robin boundary conditions. The norm is to use a ﬁrst-order ﬁnite difference scheme to approximate Neumann and Robin boundary conditions, but that compromises the accuracy of the entire scheme. As a result, new higher-order ﬁnite difference schemes for approximating Robin boundary conditions are developed in this work. Six examples for testing the applicability and performance of the method are considered. Convergence analysis is provided, and it is consistent with the numerical results. The results are compared with the exact solutions and published results from other methods. The method produces highly accurate results, which are displayed in tables and graphs.


Introduction
Boundary value problems (BVPs) are paramount in the modeling of several reallife problems with diverse applications. Boundary conditions of different types often accompany these problems. Dirichlet boundary conditions are the most common and more straightforward to implement. The Neumann boundary conditions are also not too complicated to deal with. On the other hand, the implementation of Robin boundary conditions poses significant difficulties to many researchers. As a result, BVPs with Robin boundary conditions have not received much attention. The Robin boundary condition, also referred to as the third-type boundary condition, was named after Victor Gustave Robin, who was behind its inception [1]. The Robin boundary condition is a linear combination of the solution and its derivative at the boundary point. It arises in diverse physical situations.
If a differential equation is to be solved over a domain Ω, where ∂Ω denotes the domain's boundary, the Robin boundary condition is: on the ∂Ω boundary. Our focus, in this work, is to solve second-order boundary value problems of the form y (x) = f x, y, y for a ≤ x ≤ b (2) subject to the Robin boundary conditions α a y(a) + β a y (a) = γ a , where α a , α b , β a , β b , γ a , and γ b are constants.

Quasilinearization
Before applying the CFDM, we need to linearize Equation (2) first. To this end, we use the quasilinearization (QLM) technique that was introduced by Bellman and Kalaba [31]. The QLM reduces the nonlinear BVP into a sequence of linear BVPs, which we solve iteratively until a set tolerance level is reached.
We expand the nonlinear function f (x, y, y ) in Equation (2) using Taylor's series up to the first-order terms, to obtain y = f (x, y r , y r ) + (y r+1 − y r ) ∂ f (x, y r , y r ) ∂y r + (y r+1 − y r ) ∂ f (x, y r , y r ) ∂y r , y r+1 + a 1 y r+1 + a 0 y r+1 = G(x), where a 1 = − ∂ f (x, y r , y r ) ∂y r , a 0 = − ∂ f (x, y r , y r ) ∂y r , G(x) = f (x, y r , y r ) + a 1 y r + a 0 y r .

Compact Finite Difference Schemes at Interior Nodes
In this section, we present the derivations of the sixth-order compact finite difference schemes for approximating first and second derivatives. We consider a function f (x) defined on x ∈ [a, b], where a and b are arbitrary constants. We discretize the domain into N number of nodes with a uniform step size h = (b − a)/(N − 1) and nodes x i = a + (i − 1)h, 1 ≤ i ≤ N.

Compact Finite Difference Schemes for Robin Boundary Conditions
In this section, we develop sixth-order compact finite difference schemes for approximating first and second derivatives for arbitrary functions with Robin boundary conditions. To accommodate the Robin boundary conditions and maintain a sixth-order accuracy at the boundaries, the schemes are adjusted appropriately. To that end, we use one-sided compact finite difference schemes, specifically for the points x 0 , x 1 , x N−1 , and x N . The derivatives at the interior points x 2 , x 3 , . . . , x N−2 are approximated using Equations (11) and (15).

First Derivative Approximation
We start by illustrating how the sixth-order compact schemes are formulated for first derivative approximations at the boundaries. The one-sided schemes for the first derivative at the boundary points are given as follows: At i = 1, we have The constants a i (i = 1, . . . 6) are obtained in a similar manner as in the previous section, and they are , The truncation error is At i = 2, we have 1 , (20) with the constants obtained as and At i = N − 1, we have where , and Finally, at i = N, we have Here, the constants are , and The sixth-order-accurate compact finite difference approximation for approximating the first derivative of an arbitrary function with the Robin boundary condition is obtained by combining Equation (11) with Equations (17), (20), (23), and (26), and is given by where and The first derivative f is therefore approximated by where

Second Derivative Approximation
The one-sided schemes for the second derivative at the boundary points are given as follows: At i = 1, we have The constants a i (i = 1, . . . 6) are obtained in a similar manner as in the previous section, and they are given as and At i = 2, we have where the constants are and At i = N − 1, we have and the constants are and Lastly, at i = N, we have The constants are and The sixth-order-accurate compact finite difference approximation for approximating the second derivative of an arbitrary function with the Robin boundary condition is obtained by combining Equation (15) with Equations (33), (36), (39), and (42) and is given by where The second derivative f is therefore approximated by where Therefore, to solve Equation (5), we use Equations (32) and (48) to approximate the derivatives, and so, we have where

Convergence
In this section, we discuss the convergence of the proposed method described in Section 4 to solve Equation (5).
] be vectors of the numerical solution and exact solution obtained by solving the linear system (5), respectively.
Proof. Applying the sixth-order compact schemes (29) and (45) to Equation (5), we obtain The exact solution,Ȳ, of Equation (5) is given by where (5) is given by Subtracting Equation (57) from Equation (56), we obtain where E is given by We can write B 1 and B 2 as where Substituting Equation (60) into Equation (58), we obtain Multiplying Equation (61) by h 2 , we obtain We then multiply Equation (62) by A 2 C −1 2 to obtain Solving for E in Equation (63), we obtain Recall that if · is a subordinate matrix norm and B is any matrix such that B < 1, then the matrix I + B is invertible and (66)

Numerical Examples
In this section, we consider several examples to highlight the applicability and high accuracy of the proposed algorithm in solving two-point boundary value problems subjected to Robin boundary conditions. To check the accuracy, we compute the maximum absolute error (L ∞ ) between the exact and numerical solutions, i.e., where y(x i ) andȳ(x i ) represent the CFDM solution and exact solution at the grid point x i , respectively. In addition, we compute the numerical rate of convergence, which is defined as where E 1 and E 2 are absolute errors corresponding to grids with step size h 1 and h 2 , respectively. We compare our results with published results obtained using other methods.
Exact solution : y = cos x.
We present the results of Example 1 in Table 1. We display the maximum absolute error computed using different values of grid points N. In this example, we obtain the maximum accuracy, with an error of about 10 −14 , when N = 50. When compared to the results in [9,11], our results are much more accurate. For instance, with a step size h = 0.01, the error in [9] is about 10 −10 , but our error is about 10 −14 , with a slightly larger step size. Our results are also better compared to the results of [11] obtained using Bernoulli polynomials. The agreement between the exact and numerical solutions is depicted in Figure 1, with the errors from the different values of N again shown in Figure 2.
Exact solution : y = e x . Figure 3 depicts the solution of Example 2. A good agreement between the exact and the numerical solution can be observed. In this example, N = 32 gives the maximum error of about 10 −13 , and this is shown in Table 2. Again, when compared to the results in [9], the method shows remarkable accuracy on a slightly larger step size than in [9]. The errors are also shown in Figure 4 for the different values of N. Figure 5 shows that six iterations of the quasilinearization are required to reach full convergence.
The results of Example 3 are similar to those of the previous examples. Figure 6 shows the plots of the exact and numerical solution, which are in good agreement. The maximum absolute error, computational times, and rates of convergence (ROCs) are presented in Table 3, with results from Nasir et al. [9] included. As shown in the table, the proposed method gives accurate results with fewer grid points compared to the results in [9]. Plots of the maximum absolute errors for various values of N are shown in Figure 7. Finally, Figure  8 shows the convergence plots for the various N values. It can be seen from the plot that the method reaches full convergence after six iterations.

Example 4.
Consider the nonlinear second-order boundary value problem: with boundary conditions: Exact solution : y = ln(1 + x).
For Example 4, we also compute the approximate solution using the proposed CFDM for various values of N. Likewise, the error norm, rates of convergence, and computational times are presented in Table 4. As shown in the table, for h = 0.0159, the CFDM has an error of about 10 −13 , whereas in [9], the error is 10 −11 for h = 0.01. This comparison confirms that the CFDM gives accurate results with fewer grid points. A graphical comparison of the CFDM results and the exact solution is shown in Figure 9. Plots of the absolute errors for various N values is shown in Figure 10. Lastly, Figure 11 shows the convergence plots for the various N values. The plot indicates the method reaches full convergences after six iterations.

Example 5.
Consider the nonlinear second-order boundary value problem: with boundary conditions: Exact solution : y = −2 ln cos Likewise, in Example 5, we estimate the solution of the nonlinear differential equation using the CFDM for various values of N. The error norm, rates of convergence, and execution times of the proposed method are presented in Table 5, with results from Nasir et al. [9]. The table indicates that the CFDM achieves accuracy on a slightly larger step size than in [9]. A graphical comparison of the CFDM results and the exact solution is shown in Figure 12. The errors for various values of N are shown in Figure 13. Lastly, convergence plots for the different N values are depicted in Figure 14, depicting that the method reaches full convergence after about 6 to 8 iterations, depending on the value of N.

Example 6.
Consider the nonlinear second-order boundary value problem: with boundary conditions: Lastly, a graphical comparison of the CFDM results and the exact solution is shown in Figure 15. The maximum absolute errors of Example 6 obtained using various values of N are displayed in Table 6 and also shown in Figure 16. We observe that with the same step size, we obtain better accuracy than the results in [32,34]. Lastly, convergence plots for the different N values are depicted in Figure 17, depicting that the method reaches full convergence after 8 or 9 iterations, depending on the value of N used.    In Tables 1-6, we present the computed rates of convergence for all the examples, which are consistent with the theoretical order of convergence obtained in Theorem 1.

Conclusions
In this paper, we presented a highly accurate compact finite difference method (CFDM) to solve both linear and nonlinear second-order boundary value problems subjected to Robin boundary conditions. The method utilizes sixth-order compact finite difference schemes. We successfully developed new sixth-order schemes to approximate the Robin boundary conditions, and this leads to a highly accurate method, as shown by the results. Nonlinear equations were first linearized using the quasilinearization (QLM) technique. Convergence of the CFDM was established by using the properties of the standard matrix norm. By computing the absolute error norm, and rates of convergence, the high accuracy of the CFDM was confirmed by comparing its numerical results against the numerical results of the diagonal block method presented in Nasir et al. [9] and Majid et al. [34], symmetric spline [32], Bernoulli polynomials [11], and quintic B-spline collocation [33]. We also observed that the CFDM approximate solution is in excellent agreement with the exact solution in all of the considered examples. Numerical results further confirmed that the rate of convergence of the presented CFDM is seven, consistent with the theoretical approximation.
We conclude that to solve linear and nonlinear boundary value problems subject to Robin boundary conditions, the CFDM is highly accurate and computationally efficient and a dependable method to utilize.

Conflicts of Interest:
The authors declare no conflict of interest.