Bounded x-axis, graph of an integrable function and line segments x=a\,\! and x=b\,\!, where a\,\! and b\,\!- integration limits (see figure).

The need to apply numerical integration can most often be caused by the absence of representation in and, therefore, the impossibility of analytically calculating the value of a certain integral over . It is also possible that the form of the antiderivative is so complex that it is faster to calculate the value of the integral numerically.

One-dimensional case

The main idea of ​​most methods of numerical integration is to replace the integrand with a simpler one, the integral of which can be easily calculated analytically. In this case, to estimate the value of the integral, formulas of the form

I \approx \sum_(i=1)^(n) w_i\, f(x_i),

where n\,\! is the number of points at which the value of the integrand is calculated. points x_i\,\! are called method nodes, numbers w_i\,\!- node weights. When the integrand is replaced by a polynomial of zero, first and second degree, the methods , and (Simpson) are obtained, respectively. Often formulas for estimating the value of the integral are called quadrature formulas.

Rectangle Method

Rectangle Method is obtained by replacing the integrand with a constant. As a constant, you can take the value of the function at any point on the segment \left\,\!. The most commonly used function values ​​are in the middle of a segment and at its ends. The corresponding modifications are called methods medium rectangles, left rectangles and right rectangles. The formula for the approximate calculation of the value of a definite integral by the method of rectangles is

I \approx f(x) (b-a),

where x=\frac(\left(a+b\right))(2), a\,\! or b\,\!, respectively.

Trapezoidal method

If we draw a straight line through the ends of the integration segment, we get trapezoidal method. From geometric considerations, it is easy to obtain

I \approx \frac(f(a)+f(b))(2) (b-a).

Parabola method

Using three points of the integration segment, we can replace the integrand with a parabola. Usually, the ends of the segment and its midpoint are used as such points. In this case, the formula is very simple

I \approx \frac(b-a)(6)\left(f(a)+4f\left(\frac(a+b)(2)\right)+f(b)\right).

Increasing Accuracy

Approximation of a function by one polynomial over the entire interval of integration, as a rule, leads to a large error in estimating the value of the integral.

To reduce the error, the integration segment is divided into parts and a numerical method is used to evaluate the integral on each of them.

As the number of partitions tends to infinity, the estimate of the integral tends to its true value for any numerical method.

The above methods allow for a simple procedure of halving the step, while at each step it is required to calculate the function values ​​only at newly added nodes. To estimate the calculation error, is used.

Gauss method

The methods described above use fixed line segment points (ends and midpoints) and are low (1, 1 and 3, respectively). If we can choose the points at which we calculate the function values f(x)\,\!, then it is possible, with the same number of calculations of the integrand, to obtain methods of a higher order of accuracy. So for two (as in the trapezoid method) calculations of the values ​​of the integrand, you can get a method no longer of the 1st, but of the 3rd order of accuracy:

I \approx \frac(b-a)(2)\left(f\left(\frac(a+b)(2) - \frac(b-a)(2\sqrt(3)) \right)+f\left( \frac(a+b)(2) + \frac(b-a)(2\sqrt(3)) \right) \right).

In general, using n\,\! points, you can get a method with an order of accuracy 2n-1\,\!. The values ​​of the nodes of the Gaussian method by n\,\! points are the roots of the Legendre polynomial of degree n\,\!.

The values ​​of the nodes of the Gaussian method and their weights are given in the reference books of special functions. The best known is the Gaussian five-point method.

Gauss-Kronrod method

The disadvantage of the Gauss method is that it does not have an easy (from a computational point of view) way to estimate the error of the obtained value of the integral. The use of Runge's rule requires the calculation of the integrand at approximately the same number of points, without giving practically any gain in accuracy, in contrast to simple methods, where the accuracy increases several times with each new partition. Kronrod proposed the following method for estimating the value of the integral

I \approx \sum_(i=1)^(n) a_i\, f(x_i) + \sum_(i=1)^(n+1) b_i\, f(y_i),

where x_i\,\!- Gauss method nodes by n\,\! points, and 3n+2\,\! parameters a_i\,\!, b_i\,\!, y_i\,\! are chosen in such a way that the order of accuracy of the method is equal to 3n+1\,\!.

Then, to estimate the error, one can use the empirical formula

\Delta = \left(200 |I - I_G|\right)^(1.5),

where I_G\,\!- the value of the integral, estimated by the Gauss method according to n\,\! points. Libraries [

numerical integration formula programming

Introduction

1. Numerical integration methods

2. Quadrature formulas

3. Automatic selection of the integration step

Conclusion

Bibliographic list

Introduction

The purpose of the abstract is to study and comparative analysis of methods for numerical integration of functions; implementation of these methods in the form of machine programs in a high-level language and practical solution of numerical integration problems on a computer.

When solving engineering problems, it often becomes necessary to calculate the values ​​of a certain integral of the form

. (1)

If the function is continuous on the segment [ a , b] and its antiderivative can be determined through a known function, then the calculation of such an integral is carried out according to the Newton-Leibniz formula:

.

In engineering problems, it is rarely possible to obtain the value of the integral in an analytical form. In addition, the function f (x) can be given, for example, by a table of experimental data. Therefore, in practice, to calculate a definite integral, special methods are used, which are based on the interpolation apparatus.

The idea behind these methods is as follows. Instead of calculating the integral using formula (1), the values ​​of the function are first calculated f (x i) = y i at some nodes x i Î[ a , b]. Then the interpolation polynomial is chosen P (x) passing through the obtained points ( x i , y i), which is used in calculating the approximate value of the integral (1):

.

When implementing this approach, the numerical integration formulas take the following general form:

, (2) - interpolation nodes, Ai are some coefficients, R– residual term characterizing the error of the formula. Note that formulas of the form (2) are called quadrature formulas.

The geometric meaning of numerical integration is to calculate the area of ​​a curvilinear trapezoid bounded by the graph of the function f (X), an abscissa axis and two straight lines x = a and x = b. An approximate calculation of the area leads to the rejection of the residual term in the quadrature formulas R characterizing the error of the method, which is additionally superimposed by the computational error.

1. Numerical integration methods

In applied research, it often becomes necessary to calculate the value of a definite integral

As is known from the course of mathematics, the analytical calculation of the integral can not be carried out in all cases. And even in the case when it is possible to find the analytical form of this integral, the calculation procedure gives an approximate result, so the problem of the approximate value of this integral arises.

The essence of the approximate calculation consists in two operations: 1. in choosing a finite number instead of n; 2. in point selection

in the corresponding segment.

Depending on choice

we get various formulas for calculating the integral: Formulas for left and right rectangles (5), (6) (5) (6)

Trapezium formula:

Simpson formula

b, a - ends of the considered segment.

To compare the results of calculation by the above numerical integration formulas, we calculate the following integral in 3 ways, dividing the segment into 6 equal segments: h=

According to the formula of the left rectangles:

According to the trapezoid formula:

According to Simpson's formula:

And the result obtained analytically is equal to

=1

Therefore, we can conclude that the numerical method of integration according to the Simpson formula is more accurate, but is used in the general case when dividing the segment being separated into an even number of intervals.

2. Quadrature formulas

Rectangle formulas are the simplest quadrature formulas. Let us split the integration interval [ a, b] on the P equal parts length

. Note that the value h is called the integration step. At split points X 0 = a ,X 1 = a + h , ..., x n = b note the ordinates y 0 ,y 1 ,…,y n crooked f (x), i.e. compute i = f (x i), x i = a+ ih = x i -1 + h (i =). On each segment of length h construct a rectangle with sides h and y i, where i =, i.e. by the values ​​of the ordinates calculated at the left ends of the segments. Then the area of ​​the curvilinear trapezoid, which determines the value of the integral (1), can be approximately represented as the sum of the areas of the rectangles (Fig. 1). From here we get the formula of rectangles:
. (3)

If, when calculating the integral sum, we take the values ​​of the function f (x) not at the left, but at the right ends of segments of length h, which is shown in fig. 1 with a dotted line, then we get the second version of the rectangle formula:

. (4)

The third variant of the rectangles formula can be obtained by using the values ​​of the function f (x) calculated at the midpoint of each segment of length h(Fig. 2):

. (5)

Formulas (3), (4) and (4) are called the formulas of the left, right and central rectangles, respectively.

Simpson formula. We divide the integration interval into 2 n equal parts length

. On each segment [ x i , x i+2] the integrand f (X) is replaced by a parabola passing through the points ( x i , y i), (x i +1 , y i +1), (x i +2 , y i+2). Then the approximate value of the integral is determined by the Simpson formula: . (7)

When calculating on a computer, the following formula is more convenient:

Simpson's method is one of the most widely known and used methods of numerical integration, it gives exact values ​​of the integral when integrating polynomials up to the third order inclusive.

Newton's formula. The approximate value of the integral according to Newton's formula is calculated as follows:

where the number of segments of the partition is a multiple of three, i.e., is 3 n. When developing computer programs, it is more convenient to use the equivalent formula:

Newton's method gives exact values ​​of the integral when integrating polynomials up to the fourth order inclusive.

3. Automatic selection of the integration step

As a result of the calculation by formulas (3) - (8), an approximate value of the integral is obtained, which may differ from the exact one by a certain amount, called the integration error. The error is determined by the remainder formula R, different for each of the integration methods. If it is required to calculate the value of the integral with an error not exceeding e, then it is necessary to choose such an integration step h to satisfy the inequality R (h) £e. In practice, automatic value selection is used h, which ensures the achievement of the specified error. First calculate the value of the integral I (n), dividing the integration interval into P sections, then the number of sections is doubled and the integral is calculated I (2n). The calculation process is continued until the condition becomes true.

Numerical Integration

The main questions discussed at the lecture:

2. Newton-Cotes quadrature formulas

3. Formulas of rectangles

4. Trapezoidal formula

5. Simpson formula

6. Quadrature formulas of Gauss

7. Monte Carlo Method

1. Statement of the problem of numerical integration

It is required to calculate a definite integral of the form , and the function can be given both in the form of a formula and in the form of a table.

Newton-Cotes quadrature formulas

,
where - Cotes coefficients.
These formulas give different representations for a different number n of partition segments on the same integration segment.

Rectangle formulas

Let it be required to calculate the integral .
If the integration segment is large enough, then you need to break it into smaller segments of equal length, where n is the number of segments, and replacing the curvilinear trapezoid with a rectangle on each of the segments, calculate the areas of these rectangles. Then the resulting areas must be added up, and this amount will be taken as the approximate value of the desired integral.
As for the construction of rectangles, they can be built in different ways: you can draw a perpendicular to the intersection with the curve f (x) from the right end of each segment (Fig. 1), you can - from the left end (Fig. 2)

Rice. one

Rice. 2

Depending on this, the formulas for calculating are somewhat different and are called the formulas of rectangles with right or left ordinates, respectively:

(formula of "right" rectangles)

(formula of "left" rectangles)
There is also a formula for "middle" rectangles: , for which the construction of rectangles is carried out through the midpoints of each of the segments of the partition:

· Trapezoidal formula

· Simpson formula

Replacing on each segment of the partition a part of the curve y = f(x) onto a parabolic curve, calculating the areas of the resulting figures and summing them, we obtain the Simpson formula:

·

· Quadrature formulas of Gauss

Traditionally, when obtaining quadrature Gaussian formulas in the original integral, a change of variable is performed, translating the integral over the segment into the integral over the segment [-1; one]:

.
Then .
We will use linear interpolation of the integrand.
If instead of the segment [-1; 1] to take moving nodes t1, t2 as interpolation nodes, then you need to choose these values ​​so that the area of ​​the trapezoid bounded from above by the straight line passing through the points A1 (t1, φ(t1)) and A2 (t2, φ(t2)) was equal to the integral of any polynomial of some highest degree.
Assuming that this is a polynomial of the third degree, we calculate t1, t2, which turn out to be equal to and , differing only in the numbering of values.
Further, breaking the integration segment into n parts, applying the idea described above to each of them, we can obtain the Gauss formula:

Numerical Integration

Numerical Integration(historical name: (numeric) quadrature ) - calculation of the value of a definite integral (usually approximate). Numerical integration is understood as a set of numerical methods for finding the value of a certain integral.

Numerical integration is applied when:

In these two cases, it is impossible to calculate the integral using the Newton-Leibniz formula. It is also possible that the form of the antiderivative is so complex that it is faster to calculate the value of the integral numerically.

One-dimensional case

The main idea of ​​most methods of numerical integration is to replace the integrand with a simpler one, the integral of which can be easily calculated analytically. In this case, to estimate the value of the integral, formulas of the form

where is the number of points at which the value of the integrand is calculated. The points are called the nodes of the method, the numbers are the weights of the nodes. When the integrand is replaced by a polynomial of zero, first and second degree, the methods of rectangles, trapezoids and parabolas (Simpson) are obtained, respectively. Often formulas for estimating the value of the integral are called quadrature formulas.

A special case is the method for constructing integral quadrature formulas for uniform grids, known as Cotes formulas. The method is named after Roger Coates. The main idea of ​​the method is to replace the integrand with some kind of interpolation polynomial. After taking the integral, we can write

where the numbers are called Cotes coefficients and are calculated as integrals of the corresponding polynomials in the original interpolation polynomial for the integrand with the value of the function at the node ( is the grid step; is the number of grid nodes, and the node index is ). The term is the error of the method, which can be found in different ways. For odd, the error can be found by integrating the error of the interpolation polynomial of the integrand.

Special cases of Cotes formulas are: rectangle formulas (n=0), trapezoid formulas (n=1), Simpson formula (n=2), Newton formula (n=3), etc.

Rectangle Method

Let it be required to determine the value of the integral of the function on the interval . This segment is divided by points into equal segments of length Denote by the value of the function at points Next, we make up the sums Each of the sums is an integral sum for on and therefore approximately expresses the integral

If the given function is positive and increasing, then this formula expresses the area of ​​a stepped figure made up of "incoming" rectangles, also called the formula of left rectangles, and the formula

expresses the area of ​​a stepped figure consisting of "outgoing" rectangles, also called the formula of right rectangles. The shorter the length of the segments into which the segment is divided, the more accurate the value calculated by this formula of the desired integral.

Obviously, it is worth counting on greater accuracy if we take the point in the middle of the gap as a reference point for finding the height. As a result, we get the formula for the middle rectangles:

Given the a priori greater accuracy of the last formula with the same volume and nature of calculations, it is called the formula of rectangles

Trapezoidal method

If the function on each of the partial segments is approximated by a straight line passing through the final values, then we obtain the trapezoid method.

The area of ​​the trapezoid on each segment:

Approximation error on each segment:

where

The full formula for trapezoids in the case of dividing the entire integration interval into segments of the same length:

where

Trapezoidal formula error:

where

Parabola method (Simpson's method)

Using three points of the integration segment, we can replace the integrand with a parabola. Usually, the ends of the segment and its midpoint are used as such points. In this case, the formula is very simple

.

If we divide the integration interval into equal parts, then we have

Increasing Accuracy

Approximation of a function by one polynomial over the entire interval of integration, as a rule, leads to a large error in estimating the value of the integral.

To reduce the error, the integration segment is divided into parts and a numerical method is used to evaluate the integral on each of them.

As the number of partitions tends to infinity, the estimate of the integral tends to its true value for analytic functions for any numerical method.

The above methods allow for a simple procedure of halving the step, while at each step it is required to calculate the function values ​​only at newly added nodes. To estimate the calculation error, the Runge rule is used.

Gauss method

The methods described above use fixed segment points (ends and middle) and have a low order of accuracy (1 - right and left rectangle methods, 2 - middle rectangle and trapezoid methods, 3 - parabola (Simpson) method). If we can choose the points at which we calculate the values ​​of the function , then we can obtain methods of a higher order of accuracy with the same number of calculations of the integrand. So for two (as in the trapezoid method) calculations of the values ​​of the integrand, you can get a method no longer of the 2nd, but of the 3rd order of accuracy:

.

In general, using points, you can get a method with an order of accuracy. The node values ​​of the Gauss method by points are the roots of the Legendre polynomial of degree .

The values ​​of the nodes of the Gaussian method and their weights are given in the reference books of special functions. The best known is the Gaussian five-point method.

Gauss-Kronrod method

The disadvantage of the Gauss method is that it does not have an easy (from a computational point of view) way to estimate the error of the obtained value of the integral. The use of Runge's rule requires the calculation of the integrand at about the same number of points, without giving practically any gain in accuracy, in contrast to simple methods, where the accuracy increases several times with each new partition. Kronrod proposed the following method for estimating the value of the integral

,

where are the nodes of the Gauss method by points, and the parameters , , are chosen in such a way that the order of accuracy of the method is equal to .

Then, to estimate the error, you can use the empirical formula:

,

where is the approximate value of the integral obtained by the Gauss method over points. The gsl and SLATEC libraries for calculating definite integrals contain routines using the Gauss-Kronrod method for 15, 21, 31, 41, 51 and 61 points. The library uses the Gauss-Kronrod method for 15 points.

Chebyshev method

Integration under infinite limits

To integrate over infinite limits, you need to introduce a non-uniform grid, the steps of which increase as you go to infinity, or you can make such a change of variables in the integral, after which the limits will be finite. One can proceed in a similar way if the function is singular at the ends of the integration interval

Monte Carlo Methods

Figure 3 Numerical integration of a function by the Monte Carlo method

To determine the area under the function graph, you can use the following stochastic algorithm:

For a small number of dimensions of an integrable function, the performance of Monte Carlo integration is much lower than the performance of deterministic methods. However, in some cases, when the function is implicitly specified, but it is necessary to determine the area specified in the form of complex inequalities, the stochastic method may be more preferable.

Runge-Kutta Methods

spline method

Multivariate case

In small dimensions, one can also apply quadrature formulas based on interpolation polynomials. However, in higher dimensions, these methods become unacceptable due to the rapid increase in the number of grid points and/or the complex boundary of the region. In this case, the Monte Carlo method is applied. Random points are generated in our area and the function values ​​in them are averaged. You can also use a mixed approach - divide the area into several parts, in each of which (or only in those where the integral cannot be calculated due to a complex boundary) apply the Monte Carlo method.

Literature

1. Kahaner D., Moler K., Nash S. Numerical methods and software (translated from English). M.: Mir, 2001, 575 p.