Metoda trapezelor#

De ce am folosi dreptunghiuri când trapezele se pretează mai bine pe funcții? Explorăm acum o metodă mai reușită, cunoscută drept metoda trapezelor.

Metoda trapezelor compuse#

Fie \(f:[x_0,x_n]\subset\mathbb{R}\rightarrow\mathbb{R}\) o funcție integrabilă și continuă pe întreg intervalul de definiție.

Totodată, fie mulțimea \(M=\Big\{(x_k,y_k)\in\mathbb{R}^2\,|\,y_k=f(x_k),\,k\in\overline{0,n}\Big\}\), unde \(n\in\mathbb{N}^*\) și \(x_k\lt x_{k+1}\), \(\forall k\lt n\), o mulțime de \(n+1\) puncte distincte ordonate crescător, reprezentând evaluări ale funcției \(f\).

Metoda trapezelor aproximează aria aflată sub graficul funcției \(f\) folosind \(n\) trapeze caracterizate prin înălțimea dată de perechea \((x_k,x_{k+1})\), deci \(h_k=x_{k+1}-x_k\), \(\forall k\in\overline{0,n-1}\), și prin două baze paralele cu axa \(Oy\), anume \(f(x_k)\) și \(f(x_{k+1})\).

Știind că aria unui trapez ce are bazele de lungime \(b_1\) și \(b_2\), respectiv înălțimea de lungime \(h\) are aria egală cu \(S=\frac{h(b_1+b_2)}{2}\), putem aproxima integrala folosind sume de astfel de arii:

\[ \boxed{\int_{x_0}^{x_n}f(x)\,dx\,\approx\,\frac{1}{2}\sum_{k=0}^{n-1}(x_{k+1}-x_k)(y_{k+1}+y_k)} \]

Alternativ, dacă se folosește un singur trapez pentru a calcula această arie, metoda este cunoscută drept metoda trapezului.

Metode Newton-Cotes#

Vom prezenta pe scurt conceptul de metodă Newton-Cotes. În subcapitolele ce urmează, vom extinde aceste noțiuni, mai mult drept o curiozitate, întrucât materia pe care încercăm să o predăm este totuși limitată.

O metodă de aproximare numerică a integralelor se va considera, suplimentar, o metodă tip Newton-Cotes în situația în care punctele ce alcătuiesc mulțimea \(M\) pe care se aplică algortimul sunt echidistante.

În cazul metodei trapezelor, înălțimea acestora se păstrează constantă, de exemplu \(h\in\mathbb{R}_+\). Formula se va putea rescrie astfel:

\[ \boxed{\int_{x_0}^{x_n}f(x)\,dx\,\approx\,\frac{h}{2}\sum_{k=0}^{n-1}(y_{k+1}+y_k)=\frac{h}{2}\left(y_0+2\sum_{k=1}^{n-1}y_k+y_n\right)} \]

În subcapitolele următoare, se va explica mai pe larg această noțiune.

Vizualizare#

Să considerăm funcția \(f:[0,6]\Rightarrow\mathbb{R}\), \(f(x)=\sin(x+2)+\sqrt{x}\). Vom ilustra metoda trapezelor folosind mulțimi de \(2\), \(4\), respectiv \(7\) puncte echidistante pe întreg intervalul \([0,6]\):

Un trapez (două puncte)

Trei trapeze (patru puncte)

Șase trapeze (șapte puncte)

Prin calcul matematic, se observă că \(\int_0^6\sin(x+2)+\sqrt{x}\,dx=\cos2-\cos8+4\sqrt{6}\approx9.527\). În comparație, vom arăta acum rezultatele integrării noastre:

\[ \int_0^6\sin(x+2)+\sqrt{x}\,dx\approx\begin{cases} \text{un trapez:}&13.044\\ \text{trei trapeze:}&9.104\\ \text{șase trapeze:}&9.359\\ \end{cases} \]

Se observă că, odată cu creșterea numărului de trapeze utilizate, crește și acuratețea rezultatului nostru. Menționăm că, în acest exemplu, pentru a obține o eroare mai mică de \(0.1\), sunt necesare \(9\) trapeze, iar pentru o eroare mai mică de \(0.01\), se folosesc \(44\).

Estimarea erorii#

Să ne focalizăm pe metoda de tip Newton-Cotes. Dacă am particulariza formula erorii discutată anterior, am obține următoarea relație:

\[ \boxed{\varepsilon(f,h)\leq\frac{x_n-x_0}{12} h^2\cdot\max_{\xi\in[x_0,x_n]}\Big\{\big|f''(\xi)\big|\Big\}} \]

Demonstrație

Prin definiție, \(h=\dfrac{x_n-x_0}{n}\), deci \(x_k=x_0+(k-1) h\).

Fie \(g_k(t)=\dfrac{1}{2}t\Big[f(x_k)+f(x_k+t)\Big]-\displaystyle \int_{x_k}^{x_{k+1}}f(x)\,dx\) definită astfel încât \(\big|g_k(t)\big|\) să reprezinte eroarea metodei Newton-Cotes pe intervalul \([x_k,x_k+h]\). Atunci:

\[ \frac{dg_k}{dt}=\frac{1}{2}\Big[f(x_k) + f(x_k+t)\Big] + \frac{1}{2}t\cdot f'(x_k+t)-f(x_k+t) \]

Prin derivare repetată, se obține:

\[ \frac{d^2g_k}{dt^2}=\frac{1}{2}t\cdot f''(x_k+t) \]

\[ -f''(\xi)\leq f''(x_k+t)\leq f''(\xi)\Leftrightarrow -\frac{f''(\xi)\cdot t}{2} \leq g_k''(t)\leq \frac{f''(\xi)\cdot t}{2} \]

Deoarece \(g_k(0)=g'_k(0)=0\), se obține că:

\[ \int_0^tg_k''(x)\,dx=g_k'(t)\,\,;\,\,\int_0^tg_k'(x)\,dx=g_k(t) \]

Integrând de două ori în inegalitatea de mai sus, se observă că:

\[ -\frac{f''(\xi)\cdot t^2}{4} \leq g_k'(t)\leq \frac{f''(\xi)\cdot t^2}{4}\Leftrightarrow-\frac{f''(\xi)\cdot t^3}{12} \leq g_k(t)\leq \frac{f''(\xi)\cdot t^3}{12} \]

Însumând toate aceste erori și considerând \(t=h\), se obține că:

\[ \varepsilon(f,h)=\sum_{k=0}^{n-1}g_k(h)\leq n\cdot\frac{h^3}{12}\cdot \big|f''(\xi)\big|=\frac{x_n-x_0}{12}h^2\cdot\max\Big\{\big|f''(x)\big|\Big\} \]

S-a ajuns așadar la formula ce trebuia dovedită.

Cod ilustrativ#

Funcția trapez va conține o implementare în Matlab pentru metoda trapezelor compuse, în urma explicațiilor prezentate mai sus. În acest context, folosim:

a - capătul stâng al integralei de calculat;
b - capătul drept al integralei de calculat;
n - numărul de puncte ales - va seta pasul \(h=\frac{b-a}{n}\);
f - funcția ce trebuie integrată.

În această situație, funcția va returna I, valoarea integralei căutate, adică \(I\approx\int_a^bf(x)\,dx\).

function I = trapez(a, b, n, f)
    h = (b - a) / n; 
    s = 0;
    for i = [1 : n-1]
        s += f(a + i * h);
    end
    I = h * (f(a) + f(b) + 2 * s) / 2;
end

(Bonus) Scrieri alternative#

Metoda trapezelor compuse poate fi rescrisă matriceal sub diverse forme ce permit o rezolvare mai rapidă în unele situații a problemei.

Fie \(f:[x_0,x_n]\subset\mathbb{R}\rightarrow\mathbb{R}\) o funcție integrabilă și continuă pe întreg intervalul de definiție.

Formula inițială#

Pornind de la formula inițială, observăm scrierea sub formă de produs scalar între două componente vectoriale:

\[ \int_{x_0}^{x_n}f(x)\,dx\,\approx\frac{1}{2}\begin{bmatrix} x_1-x_0 & x_2-x_1 & x_3-x_2 & \dots & x_n-x_{n-1} \end{bmatrix}\begin{bmatrix} y_1+y_0\\y_2+y_1\\y_3+y_2\\\vdots\\y_n+y_{n-1} \end{bmatrix} \]

Rescrierea formulei#

Formula inițială se poate transforma astfel:

\[\begin{split} \int_{x_0}^{x_n}f(x)\,dx & \approx \frac{1}{2}\sum_{i=0}^{n-1}\left(x_{i+1}y_{i+1}+x_{i+1}y_i - x_iy_{i+1} - x_iy_i\right) \\ & = \frac{1}{2}\left[\sum_{i=1}^{n}x_{i}y_{i} + \sum_{i=0}^{n-1}x_{i+1}y_{i} - \sum_{i=0}^{n-1}x_{i}y_{i+1} - \sum_{i=0}^{n-1}x_{i}y_{i}\right] \\ & = \frac{1}{2}\left[x_ny_n - x_0y_0 + \sum_{i=0}^{n-1}x_{i+1}y_{i} - \sum_{i=0}^{n-1}x_{i}y_{i+1}\right]=S \end{split}\]

În funcție de \(x\)#

Putem desfășura și grupa termenii acestei ecuații în funcție de \(x\):

\[\begin{split} S & = \frac{1}{2}\Big(x_ny_n - x_0y_0 + x_1y_0 + x_2y_1 + \dots + x_ny_{n-1} - x_0y_1 - x_1y_2 - \dots - x_{n-1}y_n\Big) \\ &= \frac{1}{2}\Big[x_0(-y_0-y_1)+x_1(y_0-x_2)+x_2(y_1-y_3)+\dots+x_{n-1}(y_{n-2}-y_n)+x_n(y_{n-1}+y_n)\Big]\\ &= \frac{1}{2}\begin{bmatrix} x_0&x_1&x_2&\dots&x_{n-1}&x_n \end{bmatrix} \begin{bmatrix} -y_0-y_1\\ y_0-y_2\\ y_1-y_3\\ \vdots\\ y_{n-2}-y_n\\ y_{n-1}+y_{n} \end{bmatrix} \end{split}\]

În funcție de \(y\)#

Alternativ, putem desfășura și grupa termenii acestei ecuații în funcție de \(y\):

\[\begin{split} S & = \frac{1}{2}\Big(x_ny_n - x_0y_0 + x_1y_0 + x_2y_1 + \dots + x_ny_{n-1} - x_0y_1 - x_1y_2 - \dots - x_{n-1}y_n\Big) \\ &= \frac{1}{2}\Big[y_0(x_1-x_0)+y_1(x_2-x_0)+y_2(x_3-x_1)+\dots+y_{n-1}(x_n-x_{n-2})+y_n(x_n-x_{n-1})\Big]\\ &= \frac{1}{2}\begin{bmatrix} y_0&y_1&y_2&\dots&y_{n-1}&y_n \end{bmatrix} \begin{bmatrix} x_1-x_0\\ x_2-x_0\\ x_3-x_1\\ \vdots\\ x_n-x_{n-2}\\ x_n-x_{n-1} \end{bmatrix} \end{split}\]

Concluzii#

Metoda trapezelor compuse este ușor de înteles, rapid de implementat și oferă garanția că soluția va converge către rezultatul corect odată cu creșterea numărului de puncte de la care se pornește. Dezavantajul acestei metode este eroarea, în continuare foarte mare.

Subcapitolele următoare vor propune soluții mai elegante pentru calcularea integralelor definite.

The book "Metode Numerice", written by Valentin-Ioan Vintilă, is licensed under CC BY-NC-SA 4.0