Example 2. Using the parabola formula with an accuracy of 0.0001, calculate $\int_{0}^{1} \frac{d x}{1+x^{2}}$
Solution. We will use formula (13.10). The inequality $\left|R_{n}\right| \leq \varepsilon$ will be satisfied when

$$
\frac{(b-a)^{5} \cdot M}{180 \cdot(2 n)^{4}} \leq \varepsilon, \text { i.e., when } 2 n \geq \sqrt[4]{\frac{(b-a)^{5} \cdot M}{180 \varepsilon}}
$$

We find the value $M=\max _{a \leq x \leq b}\left|f^{\mathrm{IV}}(x)\right|$. Since

$$
\begin{gathered}
f(x)=\frac{1}{1+x^{2}}, f^{\prime}(x)=-\frac{2 x}{\left(1+x^{2}\right)^{2}}, f^{\prime \prime}(x)=-2 \cdot \frac{1-3 x^{2}}{\left(1+x^{2}\right)^{3}} \\
f^{\prime \prime \prime}(x)=\frac{24\left(x-x^{3}\right)}{\left(1+x^{2}\right)^{4}}, f^{\mathrm{IV}}(x)=\frac{24\left(5 x^{4}-10 x^{2}+1\right)}{\left(1+x^{2}\right)^{5}}
\end{gathered}
$$

then $M=\max _{0 \leq x \leq 1}\left|f^{\mathrm{IV}}(x)\right|=24$ (The maximum modulus of the fourth derivative is found using $\left.f^{\mathrm{v}}(x)=-240 x \frac{3 x^{4}-10 x^{2}+3}{\left(1+x^{2}\right)^{6}}\right)$. Since $a=0$, $b=1, M=24, \varepsilon=0.0001$, then

$$
2 n \geq \sqrt[4]{\frac{(1-0)^{5} \cdot 24}{180 \cdot 0.0001}} \approx 6.07
$$

Since $2 n$ is an integer and even, we can take $2 n=8$. Noting that $h=\frac{b-a}{2 n}=\frac{1-0}{8}=0.125, y_{0}=1, y_{8}=0.5$, we construct a table of values $y_{k}=\frac{1}{1+x_{k}^{2}}$ (Table 13.3).

By formula (13.8) we get

$$
\begin{gathered}
\int_{0}^{1} \frac{d x}{1+x^{2}} \approx \frac{h}{3}\left[y_{0}+4\left(y_{1}+y_{3}+y_{5}+y_{7}\right)+2\left(y_{2}+y_{4}+y_{6}\right)+y_{8}\right]= \\
=\frac{1}{24}(1+4 \cdot 3.14672+2 \cdot 2.38118+0.5)=0.785395
\end{gathered}
$$

Remark 1. The value of the integral obtained differs from the value found using the Newton-Leibniz formula only in the sixth decimal place:

$$
\int_{0}^{1} \frac{d x}{1+x^{2}}=\left.\operatorname{arctg} x\right|_{0} ^{1}=\frac{\pi}{4} \approx 0.785398
$$

Remark 2. In the error estimate, the value of the maximum modulus of the fourth derivative of the integrand function is involved. Finding this value can sometimes be difficult. A very convenient way to estimate the accuracy of the result is a control calculation using the same formula but with a doubled step (with the number of intervals halved).

| $\mathrm{k}$ | $x_{k}$ | $1+x_{k}^{2}$ | $y_{k}$ (k-odd) | $y_{k}$ (k-even) |
| :---: | :---: | :---: | :---: | :---: |
| 1 | 0.125 | 1.01563 | 0.98461 |  |
| 2 | 0.250 | 1.06250 |  | 0.94118 |
| 3 | 0.375 | 1.14063 | 0.87670 |  |
| 4 | 0.500 | 1.25000 |  | 0.80000 |
| 5 | 0.625 | 1.39063 | 0.71910 |  |
| 6 | 0.750 | 1.56250 |  | 0.64000 |
| 7 | 0.875 | 1.76563 | 0.56637 |  |
| 8 | 1 | 2 |  |  |
| $\Sigma$ |  |  | 3.14678 | 2.38118 |

Table 13.3

## Problems

Using the parabola formula, with $2 n=10$, compute the definite integrals for the given values of the parameter $p$:

1. $\int_{0}^{1} \frac{d x}{1+p x^{2}}$; 1) $p=1$; 2) $p=2$; 3) $p=3$; 4) $p=4$; 5) $p=5$.
2. $\int_{0}^{1} \frac{d x}{x^{2}+p^{2}}$; 1) $p=1$; 2) $p=2$; 3) $p=3$; 4) $p=4$; 5) $p=5$.

## Answers

1. 2) 0.785398 ; 2) 0.667853 ; 3) 0.588826 ; 4) 0.536078 ; 5) 0.488999 . 2. 1) 0.785398 ; 2) 0.231824 ; 3) 0.107250 ; 4) 0.061245 ; 5) 0.039479 .

## Chapter 14. Approximate Solution of Differential Equations

## § 14.1. Euler's Method

Suppose we need to solve the Cauchy problem: find the solution of the differential equation

$$
y^{\prime}=f(x, y)
$$

satisfying the initial condition $y\left(x_{0}\right)=y_{0}$.

In the numerical solution of equation (14.1), the problem is posed as follows: in the points $x_{0}, x_{1}, x_{2}, \ldots, x_{n}$, find approximations $y_{k}(k=0,1,2, \ldots, n)$ for the values of the exact solution $y\left(x_{k}\right)$. The difference $\Delta x_{k}=x_{k+1}-x_{k}$ is often taken to be constant $h$ and is called the step of the grid, then

$$
x_{k}=x_{0}+k h(k=0,1,2, \ldots, n)
$$

Euler's method for solving the specified Cauchy problem is based on the direct replacement of the derivative with a difference ratio using the approximate formula $\frac{\Delta y}{\Delta x}=f(x, y)$, where $\Delta y=y(x+h)-y(x)$, $\Delta x=(x+h)-x=h$.

The approximate values $y_{k}$ at the points $x_{k}=x_{0}+k h$ are computed using the formula

$$
y_{k+1}=y_{k}+h f\left(x_{k}, y_{k}\right)(k=0,1,2, \ldots, n)
$$