Pierwiastki funkcji

Metoda siecznych ( ang. secant method ) jest ulepszeniem metody regula falsi. Warunki początkowe są zatem podobne:

Na osi OX wybieramy pewien przedział [a,b] o następujących własnościach:

Funkcja jest określona w przedziale [a,b]

Funkcja jest ciągła w przedziale [a,b]

Na krańcach przedziału [a,b] funkcja ma różne znaki

Dodatkowo w przedziale [a,b] pierwsza pochodna funkcji jest różna od zera. Warunek ten gwarantuje nam, iż w tym przedziale funkcja nie posiada minimum/maksimum lokalnego, a co za tym idzie jej sieczna nie będzie równoległa do osi OX, ponieważ algorytm bazuje na znajdowaniu punktów przecięcia siecznych z osią OX.

W przedziale [a,b] spełniającym warunki początkowe metody wybieramy dwa punkty początkowe x₁ i x₂, które leżą w tym przedziale. Będą to dwa pierwsze przybliżenia pierwiastka:

Z punktów wykresu o współrzędnych (x₁,f(x₁)) i (x₂,f(x₂)) prowadzimy sieczną i wyznaczamy punkt x₃, który jest punktem przecięcia tej siecznej z osią OX:

Do wyznaczenia punktu x₃ stosujemy ten sam wzór, co w metodzie regula falsi:

Sprawdzamy, czy w wyznaczonym punkcie wartość funkcji jest dostatecznie bliska zeru. Jeśli tak, to kończymy z wynikiem x₃. Jeśli nie, to za nowe przybliżenia przyjmujemy dwa ostatnie x₂ i x₃ i kontynuujemy obliczenia wg tej samej metody:

Wzory możemy uogólnić:

Metoda siecznych jest zwykle szybko zbieżna do pierwiastka funkcji. Jednak po wybraniu złych punktów początkowych może się zdarzyć, iż nie będzie ona zbieżna. Dlatego należy zastosować licznik kolejnych przybliżeń, a po przekroczeniu zadanej ich liczby algorytm powinien zatrzymać się z błędem.

Istnieją dwa warunki poprawnego zakończenia obliczeń:

1. Gdy wartość funkcji w wyznaczonym punkcie x_n jest dostatecznie bliska zeru.
2. Gdy odległość dwóch ostatnio wyliczonych przybliżeń pierwiastka jest mniejsza od założonej dokładności wyznaczania pierwiastka.

Algorytm metody siecznych

Dane wejściowe:

εx	–	dokładność przybliżenia pierwiastka
εy	–	dokładność przybliżenia do zera
x1,x2	–	początkowe przybliżenia pierwiastka
f	–	funkcja, której pierwiastek będzie obliczany
n	–	maksymalna liczba aproksymacji

Dane wyjściowe

xn

–

przybliżony pierwiastek funkcji

Zmienne pomocnicze

f1,f2	–	wartości funkcji w punktach x1 i x2
fn	–	wartość funkcji w punkcie xn

Lista kroków

K01:	f1 ← f ( x1 ) f2 ← f ( x2 )	Obliczamy wartości funkcji w podanych punktach startowych
K02:	Jeśli n = 0, to zakończ z błędem inaczej n ← n - 1	Sprawdzamy licznik obiegów
K03:		Wyznaczamy punkt przecięcia siecznej z osią OX
K04:	fn  ← f ( xn )	Obliczamy wartość funkcji w wyznaczonym punkcie przecięcia
K05:	Jeśli ( | fn | < εy ) lub (| xn - x2 | < εx), to zakończ	Sprawdzamy warunek zakończenia
K06:	x1 ← x2, f1 ← f2 x2 ←xn, f2 ← fn	Uaktualniamy dwa ostatnie przybliżenia
K07:	Idź do kroku K02	Kontynuujemy obliczenia w następnym obiegu pętli.

Poniższy program jest przykładową realizacją algorytmu siecznych. Program oblicza pierwiastek funkcji:

Wykres tej funkcji w przedziale [-3,1] jest następujący:

Z wykresu wynika, iż pierwiastek znajduje się w przedziale [-1,1 , -1]. Jako przybliżenia początkowe wybieramy krańce tego przedziału.

// Pierwiastek funkcji - metoda siecznych
// (C)2019 mgr Jerzy Wałaszek
// Metody numeryczne
//------------------------------------------

#include <iostream>
#include <cmath>
#include <iomanip>

using namespace std;

// Tutaj definiujemy funkcję, której pierwiastek jest wyliczany
//-------------------------------------------------------------

double f(double x)
{
  return sin(x*x-x+1/3.0)+0.5*x;
}

// Tutaj definiujemy parametry początkowe

double epsx = 1e-14; // Dokładność wyznaczania pierwiastka
double epsy = 1e-14; // Dokładność wyznaczania zera
double x1   = -1.1;  // Punkty startowe
double x2   = -1.0;
int n       = 64;    // Maksymalna liczba obiegów

// Program główny
//---------------

int main()
{
    // Zmienne

    double f1,f2,fn,xn;
    bool result = false; // Zmienna informująca o poprawnym zakończeniu metody

    setlocale(LC_ALL,"");

    cout << setprecision(15) << fixed;
    cout << "Obliczanie przybliżonego pierwiastka funkcji metodą siecznych" << endl
         << "-------------------------------------------------------------" << endl << endl;

    // Obliczamy i zapamiętujemy wartości funkcji na krańcach przedziału [a,b]

    f1 = f(x1);
    f2 = f(x2);

    while(--n)
    {
        // Obliczamy punkt przecięcia siecznej z osią OX

        xn = (f1 * x2 - f2 * x1) / (f1 - f2);

        // Obliczamy wartość funkcji w punkcie przecięcia

        fn = f(xn);

        // sprawdzamy warunki zakończenia

        if((fabs(fn) < epsy) || (fabs(xn - x2) < epsx))
        {
            result = true;
            break;
        }

        // Uaktualniamy przybliżenia

        x1 = x2; f1 = f2;
        x2 = xn; f2 = fn;
    }

    if(!result) cout << "Przekroczono liczbę obiegów !!!" << endl << endl;

    cout << "Pierwiastek        xn = " << setw(18) << xn << endl
         << "Wartość funkcji f(xn) = " << setw(18) << fn << endl
         << "Dokładność dla x epsx = " << setw(18) << epsx << endl
         << "Dokładność dla y epsy = " << setw(18) << epsy << endl
         << "Liczba obiegów        = " << 64 - n;

    cout << endl << endl;

    return 0;
}

Porównaj otrzymany wynik z dwiema poprzednimi metodami: