Optimization of nodes of Newton-Cotes formulas in the presence of an exponential boundary layer

The problem of numerical integration of a function of one variable with large gradients in the region of the exponential boundary layer it is studied. The problem is that the use of composite Newton-Cotes formulas on a uniform grid leads to significant errors when decreasing the small parameter ε, regardless of the number of nodes of the basic quadrature formula. In the paper it is proposed to choose nodes based on minimizing the error of the composite Newton-Cotes formula. It is proved that the minimum error is achieved on the Bakhvalov mesh, while the error of the quadrature formula becomes uniform in the small parameter ε.


Introduction
Consider the problem of numerical integration assuming that the following decomposition holds for the function u(x) : where The functions p(x) and Φ(x) do not explicitly given, α > 0, ε > 0. The coefficient α is separated from zero, the parameter ε may be close to zero, m will be set below. Throughout the paper, by C and C j we mean constants, not depending on the parameter ε and the number of grid intervals N. According to (3), the regular component p(x) has derivatives bounded to some order m, and the derivatives of the singular component Φ(x) are not uniformly bounded with respect to the parameter ε, ε ∈ (0, 1]. According to [1,2], the decomposition (2) is valid for the solution of a singularly perturbed boundary value problem: εu (x) + a 1 (x)u (x) − a 2 (x)u(x) = f (x), u(0) = A, u(1) = B,  2 where a 1 (x) ≥ α > 0, a 2 (x) ≥ 0, the functions a 1 (x), a 2 (x), f (x) are smooth enough. For small values of the parameter ε, the solution of the problem (4) has a boundary layer region of large gradients at the boundary x = 0. It corresponds to the representation (2).
According to [4,5] in the case of functions of the form (2) compound Newton-Cotes formulas on a uniform grid with a step of h have the error of the order of O(h) regardless of the number of nodes in the quadrature formulas.
In [4,5] the modification of the Newton-Cotes formulas is proposed based on the fact that the formulas become exact on the singular component Φ(x) from (2), known up to a factor. It is proved that the compound modified formula has an error of the order of O(h k−1 ) uniformly in ε, where k is the number of nodes in the grid template of the formula.
In order to achieve that the error of the Newton-Cotes formulas become uniform in the small parameter ε, in [6] it is proposed to apply these formulas on the Shishkin mesh [1]. It is proved that then the error of the compound Newton-Cotes formula with k nodes becomes of order (O (ln(N )/N ) k uniformly in the parameter ε. Note that in regular case, when the grid is uniform and the derivatives of the integrable function are bounded, the error of the Newton-Cotes formula is of the order of O 1/N k .
The monograph [3] considers the problem of constructing an optimal grid for the numerical integration of functions with known estimates of derivatives. We apply this approach for the numerical integration of functions of the form (2) with estimates of derivatives (3).

Optimization of nodes of the Newton-Cotes formula
First, we briefly describe the approach proposed in [3].
Let Ω h = {x n : x n = x n−1 + h n , n = 1, 2, . . . , N, x 0 = 0, x N = 1} be nonuniform grid of the interval [0, 1]. In [3], the grid nodes are constructed from the condition of minimizing the error of the compound quadrature formula. The trapezoid compound quadrature formula is considered for calculating the integral (1). Let The error of the trapezoid compound formula has an estimate In [3] it is assumed that the estimate |u (x)| ≤ F (x), x ∈ [0, 1] is given, where the function F (x) is continuous differentiable. Using the estimate (6) the error of the quadrature formula is minimized by choosing the nodes {x n }. A continuously differentiable strictly increasing function g(x) is constructed such that x n = g(n/N ), n = 0, 1, . . . , N, g(0) = 0, g(1) = 1. (7) In view of (7), the problem of minimizing the sum (6) reduces to minimizing the integral It remains to solve the equation (9) and to find the unknown constants from the boundary conditions (7). Note that the function F (g) in (9) can be specified up to a factor, since the right-hand side contains an unknown constant. Now we generalize the above approach [3], assuming that at each interval [x n−1 , x n ] the Newton-Cotes formula with k nodes is applied. So, we assume that to calculate the integral I n (u) the Newton-Cotes formula S k,n (u) with k equally spaced nodes on the interval [x n−1 , x n ] is applied. Then, for some constant C, the error estimate holds [3]: Then according to (10) By minimizing the sum in (12), by analogy with the case of the compound trapezoid formula, we obtain the Euler equation: Given the boundary conditions g(0) = 0, g(1) = 0, the function g(t) is found and then x n = g(n/N ), n = 0, 1, . . . , N. So, if we have the restriction (11) on the derivative u (k) (x), then to minimize the error of the compound Newton-Cotes formula based on basic formulas with k nodes used on the intervals [x n−1 , x n ], we have to construct the nodes {x n } using the solution of the equation (13). Now we proceed to the construction of the optimal grid for the numerical integration of functions of the form (2)-(3). When constructing the grid, we will assume that in (2) ε ≤ ε 0 < 1, where some constant ε 0 is separated from unity. In the case of ε ≥ ε 0 , the derivatives of the functions u(x) are uniformly bounded, and to calculate the integral (1), one can apply Newton-Cotes formulas on a uniform grid with known error estimates.
Remark 1. In the case of a uniform grid and the function u(x) of the form (2), the Newton-Cotes formula has an error of the order of O(1/N ) for small values of the parameter ε regardless of the number of nodes of the basic formula [4], [5]. According to the estimate (19), the error of the composite Newton-Cotes formula with k nodes of order O(1/N k ) is uniform in the parameter ε. This error estimate is the same in order, as in the regular case when the grid is uniform and derivatives of the functions u(x) are uniformly bounded [3].

Results of numerical experiments
Consider a function of the form (2) where Φ(x) = e −x/ε .
In tables e − m means 10 −m .
In Table 1 for various values of ε and N is given error of the trapezoid composite formula ∆ N,ε = |I(u) − S 2 (u)| in the case of a uniform grid. The calculated order of accuracy is also given As ε decreases, the accuracy of the formula decreases and the order of accuracy decreases from the second to the first.  Table 2 similarly shows the error and the calculated order of accuracy of the composite trapezoidal formula on a grid condensing according to the proposed algorithm. The second order of accuracy of the quadrature formula is confirmed, corresponding to lemma 1 for k = 2.
In Table 3 for various values of ε and N are given the error and the calculated accuracy order of the composite Simpson formula in the case of a uniform grid Ω h : where u n−1/2 = u((x n−1 + x n )/2). As ε decreases, the accuracy of the formula decreases and the order accuracy decreases from fourth to first. Table 4 shows similarly the error and the calculated order of accuracy of the Simpson formula on a grid condensing according to the proposed algorithm. The fourth order of accuracy of the quadrature formula is confirmed, which corresponds to lemma 2 for k = 4. It is known that there are three nodes in the Simpson formula because the middle node can be considered as a multiple. Note that the order of accuracy is the same as in the regular case when the integrable function does not have large gradients.

Conclusion
The problem of numerical integration of a function of one variable with large gradients is investigated. It is proposed to apply the composite Newton-Cotes formula on a grid condensing in the boundary layer. Nodes are selected based on minimizing the error of the composite quadrature formula. It is proved that the minimum error is achieved on the Bakhvalov mesh and the error of formula is uniform in parameter ε. The results of numerical experiments are consistent with the estimates obtained.