Bisection Root Finder

Jedna z prostszych metod do znajdowania miejsc zerowych funkcji. Miejsce zerowe poszukiwane jest w zadanym obszarze (a,b). Wymagane jest (jak zresztą dla większości metod znajdowania miejsc zerowych) by $sign(f(a)) \ne sign(f(b))$. Wtedy też mamy pewność, że o ile funkcja jest ciągła (różniczkowalna) to istnieje miejsce zerowe (jedno lub więcej). Metoda ta podobnie jak inne znajduje tylko jedno miejsce zerowe. Ponieważ jest to metoda numeryczna miejsce zerowe znajdowane jest tylko z pewnym przybliżeniem - precyzją, która jest parametrem metody. Metoda ta potrafi radzić sobie z lokalnymi ekstremami w zakresie (a,b), czego nie potrafią metody oparte na podążaniu wzdłuż pochodnej.

Metoda ta to nic innego jak zmodyfikowane poszukiwanie binarne. Dzięki temu możemy z góry określić ile maksymalnie potrzebujemy kroków by dojść do zakładanej precyzji. W dowolnym momencie mamy trzy punkty $(y_1,y_2=\frac{y_1+y+2}{2},y_3)$. Znaki $y_1,y_3$ są przeciwne. W zależności od znaku $y_2$ wybieramy taki nowy zakres poszukiwań, by granice poszukiwań miały przeciwne znaki.

Jeśli stworzymy klasę która jako parametr będzie wymagała funktora (wyrażenie lambda) to otrzymamy dość uniwersalne narzędzie. Oprócz znajdowania pierwiastków wielomianów w zakresach wstępnie przybliżonych np. metodą Sturma możemy ją na przykład wykorzystać do znajdowania maksymalnego stopnia kompresji by wielkość pliku była poniżej np. 100KB.

Wszystkie metody mają różne tempo zbieżności tzn. ile iteracji jest potrzebnych by przybliżyć miejsce zerowe z zadaną precyzją. Może to przekładać się na prędkość danej metody ale nie musi (drugi czynnik to czas spędzany w pojedynczym kroku)

Przykład implementacji w C#:

public class BisectionRootFinder : RootFinder
{
    public override double FindRoot(double a_a, double a_b, 
        Func<double, double> a_func)
    {
        Debug.Assert(a_a < a_b);

        double ya = a_func(a_a);
        double yb = a_func(a_b);

        if (Math.Abs(ya) < Constants.ROOT_FINDER_PRECISION)
            return a_a;

        if (Math.Abs(yb) < Constants.ROOT_FINDER_PRECISION)
            return a_b;

        if (Math.Sign(ya) * Math.Sign(yb) > 0.0)
            return Double.PositiveInfinity;

        double x = 0;

        for (int i = 0; i < Constants.ROOT_FINDER_MAXIMUM_ITERATIONS; i++)
        {
            x = (a_a + a_b) / 2;
            double y = a_func(x);

            if (Math.Abs(x) > Constants.ROOT_FINDER_PRECISION)
            {
                if (Math.Abs(y / x) < Constants.ROOT_FINDER_PRECISION)
                    return x;
            }
            else if (Math.Abs(y) < Constants.ROOT_FINDER_PRECISION)
                return x;

            if ((Math.Sign(ya) * Math.Sign(y)) < 0)
            {
                a_b = x;
                yb = y;
            }
            else
            {
                a_a = x;
                ya = y;
            }
        }

        return x;
    }
}