Реализации алгоритмов/Метод бисекции

Материал из Викиучебника — открытых книг для открытого мира
Перейти к навигации Перейти к поиску

На языке C[править]

                   
#include <stdio.h>                  // подключаем к компилятору библиотеку stdio.h;  
#include <math.h>                   // подключаем к компилятору библиотеку math.h;
#define EPS 1e-10                   // задаём точность результата 1*10^(-10)

double f(double x) {                // задаём вызываемой функции и аргументу - тип двойной точности;  
   return exp(x) - 2 - x;           // задаём описание функции f(x); 
}

int main ( ) {                      // главная часть программы; 
   double xl = 0, xr = 2, xm, xd, signfxl, signfxm; // задаём переменным тип двойной длины и начальные значения; 
   int n = 0;                       // задаём переменной тип целая и начальное значение;
   xd = xr - xl;                    // вычисляем длину отрезка;
   while ( abs(f(xl))>EPS || abs(f(xr))>EPS ) {  // пока абсолютные значения функции больше заданной точности делаем;  
      n = n + 1;                    // прибавляем 1 в счётчик числа проходов (делений на 2, итераций);
      xd = xd / 2;                  // вычисляем длину новых отрезков;  
      xm = xl + xd;                 // вычисляем значение x в середине отрезка;
      signfxl = copysign(1, f(xl)); // придаём единице знак f(xl); '''copysign is undefined!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!;'''
      signfxm = copysign(1, f(xm)); // придаём единице знак f(xm); '''copysign is undefined!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!;'''
       
      if ( signfxl != signfxm )     // узнаём, находится ли искомое приближение к корню в левой части; 
         xr = xm;                   // берём левую часть;
      else                          // иначе искомое приближение к корню находится в правой части;
         xl = xm;                   // берём правую часть;                    
   }

   printf ("Value of function: %.10lf\n", f(xm));      // выводим значение функции вблизи корня
   printf ("Left bound equal: %.10lf\n", xl );         // выводим xl
   printf ("Middle of line segment: %.10lf\n", (xl + xr) / 2);  // выводим приближение к корню
   printf ("Right bound equal: %.10lf\n", xr );        // выводим xr
   printf ("Numbers of iterations equal: %10i\n", n ); // выводим число проходов (делений на 2, итераций) n
}

В результате прогона программы на устройстве ввода-вывода должен получиться следующий вывод:

Value of function:           -0.0000000027
Left bound equal:             1.1461932193
Middle of line segment:       1.1461932203  
Right bound equal:            1.1461932212
Numbers of iterations equal:         30

На языке Matlab[править]

 
function [res, err] = bisection(fun, left, right, tol)
   if fun(left)*fun(right) > 0
      error('Значения функции на концах интервала должны быть разных знаков');
   end;
  %Деление отрезка пополам
   middle = 0.5*(left +right);
  %Итерационный цикл
   while abs(fun(middle)) > tol     
     %Нахождение нового интервала
      left = left*(fun(left)*fun(middle) < 0) + middle*(fun(left)*fun(middle) > 0);
      right = right*(fun(right)*fun(middle) < 0) + middle*(fun(right)*fun(middle) > 0);
     %Деление нового интервала пополам
      middle = 0.5*(left +right);
   end
   res = middle;
   err = abs(fun(middle));

Пример работы алгоритма для поиска корня функции y = tan(x) на интервале [1; 2] с точностью 1e-3. Результат вполне ожидаемый:

[res, err] = bisection('tan', 1, 2, 1e-3)  
res =
  
      1.5713
  
err =
  
    9.7656e-004

На языке Python[править]

import math

func_glob =  lambda x: 2 * x ** 3 - 3 * x ** 2 - 12 * x - 5
func_first = lambda x: 6 * x ** 2 - 6 * x - 12

a1, b1 = 0.0, 10.0

e = 0.001

def half_divide_method(a, b, f):
    x = (a + b) / 2
    while math.fabs(f(x)) >= e:
        x = (a + b) / 2
        a, b = (a, x) if f(a) * f(x) < 0 else (x, b)
    return (a + b) / 2


def newtons_method(a, b, f, f1):
    x0 = (a + b) / 2
    x1 = x0 - (f(x0) / f1(x0))
    while True:
        if math.fabs(x1 - x0) < e: return x1
        x0 = x1
        x1 = x0 - (f(x0) / f1(x0))

import pylab
import numpy

X = numpy.arange(-2.0, 4.0, 0.1)
pylab.plot([x for x in X], [func_glob(x) for x in X])
pylab.grid(True)
#pylab.show()

print 'root of the equation half_divide_method %s' % half_divide_method(a1, b1, func_glob),
print 'root of the equation newtons_method %s' % newtons_method(a1, b1, func_glob, func_first)

См. также[править]