Calculation of a definite integral by the method of rectangles. Numerical integration

Yekaterinburg

Calculation of a definite integral

Introduction

The task of numerical integration of functions is to calculate the approximate value of a certain integral:

based on a series of values of the integrand.( f(x) |x=x k = f(x k) = y k ).

Formulas for the numerical calculation of a single integral are called quadrature formulas, double and more multiple - cubature.

The usual technique for constructing quadrature formulas is to replace the integrand f(x) on a segment with an interpolating or approximating function g(x) of a relatively simple form, for example, a polynomial, followed by analytical integration. This leads to the presentation

Neglecting the remainder term R[f], we obtain the approximate formula

Denote by y i = f(x i) the value of the integrand at various points on . Quadrature formulas are formulas of closed type if x 0 =a, x n =b.

As an approximate function g(x), we consider the interpolation polynomial on in the form of the Lagrange polynomial:

, wherein , where is the remainder term of the Lagrange interpolation formula.

Formula (1) gives

, (2)

. (3)

In formula (2), the quantities () are called nodes, () - weights, - the error of the quadrature formula. If the weights () of the quadrature formula are calculated by formula (3), then the corresponding quadrature formula is called the quadrature formula of the interpolation type.

Summarize.

1. The weights () of the quadrature formula (2) for a given arrangement of nodes do not depend on the type of the integrand.

2. In quadrature formulas of interpolation type, the remainder term R n [f] can be represented as the value of a particular differential operator on the function f(x). For

3. For polynomials up to order n inclusive, the quadrature formula (2) is exact, i.e. . The highest degree of a polynomial for which the quadrature formula is exact is called the degree of the quadrature formula.

Consider special cases of formulas (2) and (3): the method of rectangles, trapezoids, parabolas (Simpson's method). The names of these methods are due to the geometric interpretation of the corresponding formulas.

Rectangle method

The definite integral of the function of the function f(x): is numerically equal to the area of the curvilinear trapezoid bounded by the curves y=0, x=a, x=b, y=f(x) (Figure 1).

Rice. 1 Area under the curve y=f(x) To calculate this area, the entire integration interval is divided into n equal subintervals of length h=(b-a)/n. The area under the integrand is approximately replaced by the sum of the areas of the rectangles, as shown in Figure (2).

Rice. 2 The area under the curve y=f(x) is approximated by the sum of the areas of the rectangles
The sum of the areas of all rectangles is calculated by the formula

The method represented by formula (4) is called the left box method, and the method represented by formula (5) is called the right box method:

The error in calculating the integral is determined by the value of the integration step h. The smaller the integration step, the more accurately the integral sum S approximates the value of the integral I. Based on this, an algorithm is built to calculate the integral with a given accuracy. It is considered that the integral sum S represents the value of the integral I with an accuracy of eps, if the difference in absolute value between the integral sums and calculated with the step h and h/2, respectively, does not exceed eps.

To find a definite integral using the method of middle rectangles, the area bounded by lines a and b is divided into n rectangles with the same bases h, the heights of the rectangles will be the points of intersection of the function f(x) with the midpoints of the rectangles (h/2). The integral will be numerically equal to the sum of the areas of n rectangles (Figure 3).

Rice. 3 The area under the curve y=f(x) is approximated by the sum of the areas of the rectangles

n is the number of partitions of the segment .

Trapezoidal method

To find a definite integral using the trapezoid method, the area of a curvilinear trapezoid is also divided into n rectangular trapezoids with heights h and bases y 1, y 2, y 3,..y n, where n is the number of the rectangular trapezoid. The integral will be numerically equal to the sum of the areas of rectangular trapezoids (Figure 4).

Rice. 4 The area under the curve y=f(x) is approximated by the sum of the areas of rectangular trapezoids.

n is the number of partitions

(6)

The error of the trapezoid formula is estimated by the number

The error of the trapezoid formula decreases faster with growth than the error of the rectangle formula. Therefore, the trapezoid formula allows you to get more accuracy than the rectangle method.

Simpson formula

If for each pair of segments we construct a polynomial of the second degree, then integrate it on the segment and use the additivity property of the integral, then we obtain the Simpson formula.

In Simpson's method for calculating the definite integral, the entire integration interval is divided into subintervals of equal length h=(b-a)/n. The number of partition segments is an even number. Then, on each pair of adjacent subintervals, the subintegral function f(x) is replaced by a Lagrange polynomial of the second degree (Figure 5).

Rice. 5 The function y=f(x) on the segment is replaced by a polynomial of the 2nd order

Consider the integrand on the interval . Let us replace this integrand with a second-degree Lagrange interpolation polynomial coinciding with y= at the points :

We integrate on the segment .:

We introduce a change of variables:

Given the replacement formulas,

After integrating, we get the Simpson formula:

The value obtained for the integral coincides with the area of a curvilinear trapezoid bounded by the axis , straight lines , and a parabola passing through the points. On the segment, Simpson's formula will look like:

In the parabola formula, the value of the function f (x) at odd split points x 1, x 3, ..., x 2 n -1 has a coefficient of 4, at even points x 2, x 4, ..., x 2 n -2 - coefficient 2 and at two boundary points x 0 \u003d a, x n \u003d b - coefficient 1.

The geometric meaning of Simpson's formula: the area of a curvilinear trapezoid under the graph of the function f(x) on a segment is approximately replaced by the sum of the areas of the figures lying under the parabolas.

If the function f(x) has a continuous derivative of the fourth order, then the absolute value of the error of the Simpson formula is no more than

where M is the largest value on the segment. Since n 4 grows faster than n 2 , the error of Simpson's formula decreases with increasing n much faster than the error of the trapezoid formula.

We calculate the integral

This integral is easy to calculate:

Let's take n equal to 10, h=0.1, calculate the values of the integrand at the partition points, as well as half-integer points .

According to the formula of middle rectangles, we get I straight = 0.785606 (the error is 0.027%), according to the trapezoid formula I trap = 0.784981 (the error is about 0.054. When using the method of right and left rectangles, the error is more than 3%.

To compare the accuracy of the approximate formulas, we calculate once again the integral

but now by the Simpson formula for n=4. We divide the segment into four equal parts with points x 0 \u003d 0, x 1 \u003d 1/4, x 2 \u003d 1/2, x 3 \u003d 3/4, x 4 \u003d 1 and calculate approximately the values of the function f (x) \u003d 1 / ( 1+x) at these points: y 0 =1.0000, y 1 =0.8000, y 2 =0.6667, y 3 =0.5714, y 4 =0.5000.

By Simpson's formula, we get

Let us estimate the error of the result obtained. For the integrand f(x)=1/(1+x) we have: f (4) (x)=24/(1+x) 5 , whence it follows that on the segment . Therefore, we can take M=24, and the result error does not exceed 24/(2880× 4 4)=0.0004. Comparing the approximate value with the exact one, we conclude that the absolute error of the result obtained by the Simpson formula is less than 0.00011. This is in accordance with the error estimate given above and, in addition, indicates that the Simpson formula is much more accurate than the trapezoid formula. Therefore, the Simpson formula for the approximate calculation of definite integrals is used more often than the trapezoid formula.

Comparison of methods for accuracy

Let's compare the methods in terms of accuracy, for this we will calculate the integral of the functions y=x, y=x+2, y=x 2 , at n=10 and n=60, a=0, b=10. The exact value of the integrals is respectively: 50, 70, 333.(3)

Table 1

Table 1 shows that the most accurate is the integral found by the Simpson formula, when calculating the linear functions y=x, y=x+2, accuracy is also achieved by the methods of middle rectangles and the trapezoid method, the method of right rectangles is less accurate. Table 1 shows that with an increase in the number of partitions n (an increase in the number of integrations), the accuracy of the approximate calculation of the integrals increases

Assignment for laboratory work

1) Write programs for calculating a definite integral using methods: middle, right rectangles, trapezoid and Simpson's method. Perform integration of the following functions:

on a segment with a step , ,

3. Perform a variant of an individual task (table 2)

Table 2 Individual task options

	Function f(x)	Segment of integration

2) Conduct a comparative analysis of the methods.

Calculation of a definite integral: Guidelines for laboratory work in the discipline "Computational mathematics" / comp. I.A. Selivanova. Yekaterinburg: GOU VPO USTU-UPI, 2006. 14 p.

The guidelines are intended for students of all forms of education of the specialty 230101 - "Computers, complexes, systems and networks" and bachelors of the direction 230100 - "Computer science and computer technology". Compiled by Selivanova Irina Anatolyevna

Graphic image:

Let us calculate the approximate value of the integral. To assess the accuracy, we use the calculation by the method of left and right rectangles.

Calculate the step when splitting into 10 parts:

The split points of the segment are defined as.

We calculate the approximate value of the integral using the formulas of the left rectangles:

0.1(0.6288+0.6042+0.5828+0.5642+0.5479+0.5338+0.5214+0.5105+0.5008+0.4924)0.5486

We calculate the approximate value of the integral using the formulas of the right rectangles:

0.1(0.6042+0.5828+0.5642+0.5479+0.5338+0.5214+0.5105+0.5008+0.4924+0.4848)0.5342

Solution of a boundary value problem for an ordinary differential equation by the sweep method.

For an approximate solution of an ordinary differential equation, the sweep method can be used.

Consider a linear d.p.

y""+p(x)y"+q(x)y=f(x) (1)

with two-point linear boundary conditions

Let's introduce the notation:

The sweep method consists of a "forward move", in which the coefficients are determined:

After performing the "forward move", they proceed to the execution of the "reverse move", which consists in determining the values of the desired function according to the formulas:

Using the sweep method, compose a solution to the boundary value problem for an ordinary differential equation with accuracy; Step h=0.05

2; A=1; =0; B=1.2;

The Dirichlet problem for the Laplace equation by the grid method

Find a continuous function u(x, y) that satisfies the Laplace equation inside a rectangular region

and taking on the boundary of the region given values, i.e.

where f l , f 2 , f 3 , f 4 are given functions.

Introducing the notation, we approximate the partial derivatives and at each internal grid node by the second-order central difference derivatives

and replace the Laplace equation with a finite difference equation

The error of replacing a differential equation with a difference one is .

Equations (1) together with the values at the boundary nodes form a system of linear algebraic equations for the approximate values of the function u(x, y) at the grid nodes. This system has the simplest form when:

When obtaining grid equations (2), the scheme of nodes shown in Fig. 1 was used. 1. The set of nodes used to approximate the equation at a point is called a template.

Picture 1

The numerical solution of the Dirichlet problem for the Laplace equation in a rectangle consists in finding the approximate values of the desired function u(x, y) at the internal nodes of the grid. To determine the quantities, it is required to solve the system of linear algebraic equations (2).

In this paper, it is solved by the Gauss--Seidel method, which consists in constructing a sequence of iterations of the form

(the superscript s denotes the iteration number). For , the sequence converges to the exact solution of system (2). As a condition for the termination of the iterative process, one can take

Thus, the error of the approximate solution obtained by the grid method consists of two errors: the error of approximating the differential equation by difference; error resulting from the approximate solution of the system of difference equations (2).

It is known that the difference scheme described here has the property of stability and convergence. The stability of the scheme means that small changes in the initial data lead to small changes in the solution of the difference problem. Only such schemes make sense to apply in real calculations. The convergence of the scheme means that when the grid step tends to zero (), the solution of the difference problem tends in a certain sense to the solution of the original problem. Thus, by choosing a sufficiently small step h, one can solve the original problem arbitrarily exactly.

Using the grid method, compose an approximate solution of the Dirichlet problem for the Laplace equation in the square ABCD with vertices A(0;0) B(0;1) C(1;1) D(1;0); step h=0.02. When solving the problem, use the iterative Libman averaging process until an answer is obtained with an accuracy of 0.01.

1) Calculate the values of the function on the sides:

1. On the AB side: according to the formula. u(0;0)=0 u(0;0.2)=9.6 u(0;0.4)=16.8 u(0;0.6)=19.2 u(0;0.8)=14.4 u(0;1)=0
2. BC side=0
3. On the side CD=0

4. On the AD side: by the formula u(0;0)=0 u(0.2;0)=29.376 u(0.4;0)=47.542 u(0.6;0)=47.567 u(0.8;0)=29.44 u(1;0)=0
2) To determine the values of the function at the internal points of the region using the grid method, we replace the given Laplace equation at each point with a finite-difference equation according to the formula

Using this formula, we will make an equation for each interior point. As a result, we obtain a system of equations.

The solution of this system is performed by the Liebman-type iterative method. For each value, we compose a sequence that we build up to convergence in hundredths. Let us write down the relations with the help of which we will find the elements of all sequences:

For calculations using these formulas, it is necessary to determine the initial values that can be found in any way.

3) To obtain the initial approximate solution of the problem, we assume that the function u(x,y) is uniformly distributed along the horizontals of the region.

First, consider a horizontal line with boundary points (0;0.2) and (1;0.2).

Let us denote the desired values of the function at internal points through.

Since the segment is divided into 5 parts, the measurement step of the function

Then we get:

Similarly, we find the values of the function at the interior points of other horizontals. For a horizontal, with boundary points (0;0.4) and (1;0.4) we have

For a horizontal with boundary points (0;0.6) and (1;0.6) we have

Finally, we find the values for the horizontal with the boundary points (0;0.8) and (1;0.8).

We will present all the obtained values in the following table, which is called the null pattern:

Estimation of the remainder term of the formula:

Service assignment. The service is intended for online calculation of a definite integral using the formula of rectangles.

Instruction. Enter the integrand f(x) , click Solve. The resulting solution is saved in a Word file. A solution template is also created in Excel. Below is a video instruction.

Function entry rules

Examples
≡ x^2/(x+2)
cos 2 (2x+π) ≡ (cos(2*x+pi))^2
≡ x+(x-1)^(2/3) This is the simplest quadrature formula for calculating the integral, which uses one value of the function

(8.5.1)
where ; h=x 1 -x 0 .
Formula (8.5.1) is the central formula for rectangles. Let's calculate the remainder. Let us expand the function y=f(x) at the point ε 0 into a Taylor series:
(8.5.2)
where ; . We integrate (8.5.2):

(8.5.3)

In the second term, the integrand is odd, and the limits of integration are symmetric with respect to the point ε 0 . Therefore, the second integral is equal to zero. Thus, from (8.5.3) it follows .
Since the second factor of the integrand does not change sign, then by the mean value theorem we obtain , where . After integration, we get . (8.5.4)
Comparing with the remainder of the trapezoid formula, we see that the error of the rectangle formula is two times less than the error of the trapezoid formula. This result is true if in the formula of rectangles we take the value of the function at the midpoint.
We obtain the formula of rectangles and the remainder term for the interval . Let the grid x i =a+ih, i=0,1,...,n, . Consider the grid ε i =ε 0 +ih, i=1,2,..,n, ε 0 =a-h/2. Then . (8.5.5)
Residual term .
Geometrically, the formula of rectangles can be represented by the following figure:

If the function f (x) is given in a table, then either the left-hand formula of rectangles is used (for a uniform grid)

or the right-hand formula of rectangles

.
The error of these formulas is estimated through the first derivative. For the interval, the error is

; .
After integration, we get .

Example. Calculate the integral for n=5:
a) according to the trapezoid formula;
b) according to the formula of rectangles;
c) according to the Simpson formula;
d) according to the Gauss formula;
e) according to the Chebyshev formula.
Calculate the error.
Solution. For 5 integration nodes, the grid step will be 0.125.
When solving, we will use the table of function values. Here f(x)=1/x.

	x		f(x)
x0	0.5	y0	2
x1	0.625	y1	1.6
x2	0.750	y2	1.33
x3	0.875	y3	1.14
x4	1.0	y4	1

a) trapezoid formula:
I=h/2×;
I=(0.125/2)×= 0.696;
R= [-(b-a)/12]×h×y¢¢(x);
f¢¢(x)=2/(x 3).
The maximum value of the second derivative of the function on the interval is 16: max (f¢¢(x)), xн=2/(0.5 3)=16, therefore
R=[-(1-0.5)/12]×0.125×16=- 0.0833;
b) formula of rectangles:
for the left-hand formula I=h×(y0+y1+y2+y3);
I=0.125×(2+1.6+1.33+1.14)= 0.759;
R=[(b-a)/6]×h 2×y¢¢(x);
R=[(1-0.5)/6]×0.125 2×16= 0.02;
c) Simpson's formula:
I=(2h/6)×(y0+y4+4×(y1+y3)+2×y2);
I=(2×0.125)/6×(2+1+4×(1.6+1.14)+2×1.33)= 0.693;
R=[-(b-a)/180]×h 4×y (4) (x);
f(4)(x)=24/(x5)=768;
R=[-(1-0.5)/180]×(0.125) 4×768 = - 5.2 e-4;
d) Gauss formula:
I=(b-a)/2×;
x i =(b+a)/2+t i (b-a)/2
(A i , t i - table values).

				t (n=5)		A (n=5)
x1	0.9765	y1	1.02	t1	0.90617985	A 1	0.23692688
x2	0.8846	y2	1.13	t2	0.53846931	A2	0.47862868
x3	0.75	y3	1.33	t3	0	A 3	0.56888889
x4	0.61	y4	1.625	t4	-0.53846931	A4	0.47862868
x5	0.52	y5	1.91	t5	-0.90617985	A5	0.23692688

I=(1-0.5)/2×(0.2416+0.5408+0.7566+0.7777+0.4525)= 0.6923;
e) Chebyshev formula:
I=[(b-a)/n] ×S f(x i), i=1..n,
x i =(b+a)/2+[ t i (b-a)]/2 - necessary reduction of the integration interval to the interval [-1;1].
For n=5

t1	0.832498
t2	0.374541
t3	0
t4	-0.374541
t5	-0.832498

Let's find x values and function values at these points:

x1	0,958	f(x1)	1,043
x2	0,844	f(x2)	1,185
x3	0,75	f(x3)	1,333
x4	0,656	f(x4)	1,524
x5	0,542	f(x5)	1,845

The sum of the function values is 6.927.
I=(1-0.5)/5×6.927=0.6927.

In general left rectangle formula on the segment as follows (21) :

In this formula x 0 =a, x n =b, since any integral in general looks like: (see the formula 18 ).

h can be calculated using the formula 19 .

y 0 ,y 1 ,...,y n-1 x 0 , x 1 ,...,x n-1 (x i =x i-1 +h).

Formula of right rectangles.

In general right rectangle formula on the segment as follows (22) :

In this formula x 0 =a, x n =b(see formula for left rectangles).

h can be calculated using the same formula as in the formula for the left rectangles.

y 1 ,y 2 ,...,y n are the values of the corresponding function f(x) at the points x 1 , x 2 ,...,x n (x i =x i-1 +h).

Medium Rectangle Formula.

In general middle rectangle formula on the segment as follows (23) :

Where x i =x i-1 +h.

In this formula, as in the previous ones, h is required to multiply the sum of the values of the function f (x), but not just by substituting the corresponding values x 0 ,x 1 ,...,x n-1 into the function f(x), and adding to each of these values h/2(x 0 +h/2, x 1 +h/2,..., x n-1 +h/2) and then only substituting them into the given function.

h can be calculated using the same formula as in the formula for left rectangles." [ 6 ]

In practice, these methods are implemented as follows:

Mathcad ;

excel .

Mathcad ;

excel .

In order to calculate the integral using the formula of average rectangles in Excel, you must perform the following steps:

Continue working in the same document as when calculating the integral using the formulas of the left and right rectangles.

Enter the text xi+h/2 in cell E6, and f(xi+h/2) in cell F6.

Enter the formula =B7+$B$4/2 in cell E7, copy this formula by dragging to the range of cells E8:E16

Enter the formula =ROOT(E7^4-E7^3+8) in cell F7, copy this formula by pulling to the range of cells F8:F16

Enter the formula =SUM(F7:F16) in cell F18.

Enter the formula =B4*F18 in cell F19.

Enter the text of averages in cell F20.

As a result, we get the following:

Answer: the value of the given integral is 13.40797.

Based on the obtained results, we can conclude that the formula for the middle rectangles is the most accurate than the formulas for the right and left rectangles.

1. Monte Carlo method

"The main idea of the Monte Carlo method is to repeat random tests many times. A characteristic feature of the Monte Carlo method is the use of random numbers (numerical values of some random variable). Such numbers can be obtained using random number generators. For example, the Turbo Pascal programming language has standard function random, whose values are random numbers uniformly distributed on the interval . This means that if you divide the specified segment into a certain number of equal intervals and calculate the value of the random function a large number of times, then approximately the same number of random numbers will fall into each interval. In the basin programming language, a similar sensor is the rnd function. In spreadsheet MS Excel, the function RAND returns a uniformly distributed random number greater than or equal to 0 and less than 1 (changes when recalculated)" [ 7 ].

In order to calculate it, you need to use the formula () :

Where (i=1, 2, …, n) are random numbers lying in the interval .

To obtain such numbers based on a sequence of random numbers x i uniformly distributed in the interval , it is enough to perform the transformation x i =a+(b-a)x i .

In practice, this method is implemented as follows:

In order to calculate the integral by the Monte Carlo method in Excel, you must perform the following steps:

In cell B1, enter the text n=.

In cell B2, enter the text a=.

In cell B3, enter the text b=.

Enter the number 10 in cell C1.

Enter the number 0 in cell C2.

In cell C3, enter the number 3.2.

In cell A5, enter I, in B5 - xi, in C5 - f (xi).

Cells A6:A15 fill with numbers 1,2,3, ..., 10 - since n=10.

Enter the formula =RAND()*3.2 in cell B6 (numbers are generated in the range from 0 to 3.2), copy this formula by pulling into the range of cells B7:B15.

Enter the formula =ROOT(B6^4-B6^3+8) into cell C6, copy this formula by dragging it into the range of cells C7:C15.

Enter the text "sum" in cell B16, "(b-a)/n" in B17, and "I=" in B18.

Enter the formula =SUM(C6:C15) in cell C16.

Enter the formula =(C3-C2)/C1 in cell C17.

Enter the formula =C16*C17 in cell C18.

As a result, we get:

Answer: the value of the given integral is 13.12416.