Metoda bisecției#
Aceasta este una dintre cele mai utilizate metode pentru aflarea unei soluții pentru o ecuație neliniară dată. Algoritmul este printre cele mai simple, asemănându-se foarte mult cu căutarea binară.
Explicația matematică#
Fie \(f:[a,b]\subseteq\mathbb{R}\rightarrow\mathbb{R}\) o funcție continuă, definită astfel încât semnul lui \(f(a)\) să fie diferit de semnul lui \(f(b)\), adică \(f(a)\cdot f(b)\lt 0\).
Din moment ce \(f(x)\) este continuă, avem garanția că există cel puțin un punct \(x_0\in[a,b]\) astfel încât \(f(x_0)=0\), căci altfel nu putem face trecerea de la \(\text{sign}\big(f(a)\big)\) la \(\text{sign}\big(f(b)\big)\).
Așadar, condiția este necesară pentru a ne asigura că există o soluție.
Să considerăm valoarea \(c\in\mathbb{R}\) definită simplist prin media aritmetică dintre \(a\) și \(b\):
\[ \boxed{c=\frac{a+b}{2}} \]Există trei cazuri posibile în ceea ce privește studiul lui \(f\) în \(c\):
-
\(f(c)=0\) - Jackpot! Aceasta era valoarea căutată. Șansele de a nimeri însă această valoare sunt foarte mici, motiv pentru care vom fi interesați ca \(f(c)\approx 0\).
Putem interpreta această condiție în mod absolut, adică \(\Big|f(c)\Big|\lt \varepsilon\), unde \(\varepsilon\rightarrow0\), sau în mod relativ față de valorile precedente ale lui \(c\) (vom detalia ulterior);
-
Semnul lui \(f(c)\) este diferit de semnul lui \(f(a)\), ceea ce înseamnă că, datorită continuității, va exista cel puțin un punct \(\delta\in(a,c)\) astfel încât \(f(\delta)=0\).
Așadar, vom continua căutarea, restrângând-o pe intervalul \((a,c)\);
-
Semnul lui \(f(c)\) este diferit de semnul lui \(f(b)\), ceea ce înseamnă că, datorită continuității, va exista cel puțin un punct \(\delta\in(c,b)\) astfel încât \(f(\delta)=0\).
Așadar, vom continua căutarea, restrângând-o pe intervalul \((c,b)\);
Algoritmul poate fi oprit fie atunci când găsim o soluție cu eroare acceptabilă, fie atunci când se depășește un număr fixat de pași pe care ne permitem să îi executăm.
Interpretarea geometrică#
Preferăm să vizualizăm acest algoritm înainte de a continua cu partea teoretică. Așadar, să considerăm funcția:
\[ f(x)=\frac{(x-8)^3}{100}-\sin x \]Să considerăm totodată intervalul generat de \(a=0\) și \(b=16\). Să desenăm funcția noastră pe acest interval și să marcăm totodată punctul de abscisă \(c=\frac{0+16}{2}=8\):
Graficul funcției \(f(x)\)
Punctele de abscisă \(a\) (albastru), \(b\) (roșu) și \(c\) (galben) evidențiate pentru prima iterație
Deoarece \(f(c)\) și \(f(b)\) au semne diferite (adică \(f(b)\cdot f(c)\lt 0\)), vom trece în partea dreaptă a graficului, în dreapta paralelei la \(Oy\) desenată punctat. Astfel, \(a\rightarrow8\) și \(b=16\), deci \(c=\frac{8+16}{2}=12\). Să continuăm acum cu algoritmul:
Graficul funcției \(f(x)\)
Punctele de abscisă \(a\) (albastru), \(b\) (roșu) și \(c\) (galben) evidențiate pentru a doua iterație
Cum \(f(a)\) și \(f(c)\) au semne diferite (adică \(f(a)\cdot f(c)\lt 0\)), vom trece în partea stângă a graficului, adică în stânga paralelei la \(Oy\) desenată punctat. Așadar, \(a=8\) și \(b\rightarrow12\).
Vom mai desena doar o paralelă la \(Oy\) din \(c=\frac{8+12}{2}=10\):
Paralela cu \(Oy\) care trece prin \(c=10\)
Așadar, acesta este progresul nostru: \(\begin{cases} f(c_0)=f(8)\approx-0.989358\\ f(c_1)=f(12)\approx1.176573\\ f(c_2)=f(10)\approx0.624021 \end{cases}\)
Observăm că ne-am apropiat de \(p\approx9.397\), una dintre soluțiile ecuației \(f(x)=0\). Totodată, prin acest exemplu, am justificat importanța punctelor de pornire \(a\) și \(b\) - algoritmul soluționează ecuația \(f(x)=0\) cu un singur rezultat.
Numărul de iterații#
Considerând că algoritmul este lăsat să ruleze până când suntem mulțumiți de eroarea sa relativă, adică \( \left|\dfrac{c_{\text{nou}} - c_{\text{vechi}}}{c_{\text{nou}}}\right|\leq\varepsilon \), atunci algoritmul va parcurge exact \(n\in\mathbb{N}\) iterații:
\[ \boxed{n=\left\lceil\log_2\left(\frac{b-a}{\varepsilon}\right)\right\rceil} \]Intuitiv, acest fapt reiese din similitudinea algoritmului față de căutarea binară.
Ordin și rată de convergență#
Vrem să găsim o metrică mai generală a acestor metode, ceva ce poate fi folosit și pentru alți algoritmi, așadar găsim aici oportunitatea perfectă de a explica ce înțelegem prin ordinul, respectiv prin rata de convergență a unui algoritm.
Fie \((x_n)_n\in\mathbb{N}\in\mathbb{R}\) o secvență care converge către o valoare \(x^*\in\mathbb{R}\), adică \(x_n\rightarrow x^*\) atunci când \(n\rightarrow+\infty\).
În această situație, vom nota cu \(q\in[1,+\infty)\) ordinul de convergență, respectiv cu \(\mu\in[0,1]\) rata de convergență a secvenței (algoritmului), iar aceste valori pot fi determinate prin formula[11]:
\[ \boxed{\mu=\lim_{n\to+\infty}\frac{\left|x_{n+1}-x^*\right|}{\left|x_{n}-x^*\right|^q}} \]Calculul este posibil deoarece impunem \(\mu\in[0,1]\). Alternativ, pentru aflarea ordinului de convergență \(q\), se poate utiliza o aproximare practică ce adesea se aliniază cu valoarea reală[11]:
\[ \boxed{q\approx\frac{\ln\left|\dfrac{x_{k+1}-x_k}{x_k-x_{k-1}}\right|}{\ln\left|\dfrac{x_{k}-x_{k-1}}{x_{k-1}-x_{k-2}}\right|}} \,,\,\,\forall k\in\mathbb{N} \]În funcție de valoarea ordinului de convergență \(q\), spunem că:
-
Pentru \(q=1\): șirul (algoritmul) converge liniar;
-
Pentru \(q=2\): șirul (algoritmul) converge quadratic;
-
Pentru \(q=3\): șirul (algoritmul) converge cubic.
Pentru metoda bisecției#
În situația metodei bisecției, este evident că aceasta converge liniar, adică ordinul său de convergență este \(q=1\), iar eroarea continuă să se înjumătățească la fiecare pas, ceea ce implică faptul că rata de convergență este egală cu \(\mu=\frac{1}{2}\).
Dacă notăm \(c_k\) valoarea obținută pentru subintervalul \([a_k,b_k]\subseteq[a,b]\), unde \(k\in\mathbb{N}^*\), atunci devine evident că:
\[c_k=\frac{a_k+b_k}{2}\]Din formula anterioară, \(c_{k+1}=\frac{a_{k+1}+b_{k+1}}{2}\). Prin modul în care se micșorează intervalul de căutare:
-
Fie \(a_{k+1}=c_k\) și \(b_{k+1}=b_k\), deci:
\[c_{k+1}=\frac{c_k+b_{k}}{2}=\frac{\frac{a_k+b_k}{2}+b_{k}}{2}=\frac{a_k+3b_k}{4}\]Folosind această relație, obținem:
\[c_{k+1}-c_k=\frac{a_k+3b_k}{4}-\frac{a_k+b_k}{2}=\frac{b_k-a_k}{4}\] -
Fie \(b_{k+1}=c_k\) și \(a_{k+1}=a_k\). Analog, se obține:
\[c_{k+1}-c_k=\frac{3a_k+b_k}{4}-\frac{a_k+b_k}{2}=\frac{a_k-b_k}{4}\]
Ne întoarcem acum la formula lui \(q\):
\[q\approx\frac{\ln\left|\dfrac{c_{k+1}-c_k}{c_k-c_{k-1}}\right|}{\ln\left|\dfrac{c_{k}-c_{k-1}}{c_{k-1}-c_{k-2}}\right|}=\frac{\ln\left|\dfrac{a_k-b_k}{a_{k-1}-b_{k-1}}\right|}{\ln\left|\dfrac{a_{k-1}-b_{k-1}}{a_{k-2}-b_{k-2}}\right|}\]Din moment ce mereu înjumătățim intervalul pe care facem căutarea, este evident că \(\left|\frac{a_k-b_k}{a_{k-1}-b_{k-1}}\right|=\left|\frac{a_{k-1}-b_{k-1}}{a_{k-2}-b_{k-2}}\right|=2\), așadar:
\[q\approx\frac{\ln2}{\ln2}=1\]Cod ilustrativ#
Să implementăm acum metoda bisecției în Matlab sub forma funcției metoda_bisectiei
. Aceasta va primi drept parametri:
-
a
- capătul stâng de la care se pornește; -
b
- capătul drept de la care se pornește; -
f
- funcția pentru care căutăm o rădăcină; -
tol
- toleranța acceptată, echivalentul valorii minime acceptate pentru \(\frac{b-a}{2}\) până să se oprească algoritmul; -
max_iter
- numărul maxim de pași pe care îi poate executa algoritmul până se renunțe.
Codul de mai jos returnează c
, echivalentul punctului în care \(f(c)\approx0\).
function c = metoda_bisectiei(a, b, f, tol, max_iter)
if f(a) == 0
c = a;
return;
elseif f(b) == 0
c = b;
return;
end
for k = [1 : max_iter]
c = (a + b) / 2;
if c - a < tol
return
end
val_a = f(a);
val_c = f(c);
if val_c == 0
return
elseif val_a * val_c < 0
b = c;
else
a = c;
end
end
end
Generalizare pe numere complexe#
Există o metodă care generalizează algoritmul pe numere complexe, denumită metoda Lehmer–Schur, însă aceasta depășește materia descrisă în această carte. Pentru curioși, detalii suplimentare pot fi găsite aici.
Licență#
The book "Metode Numerice", written by Valentin-Ioan Vintilă, is licensed under CC BY-NC-SA 4.0