[numerical computation] False-Position Method (w/ MATLAB)

False-Position Method

수치해석에서 근을 찾는 방법 Bisection Method와 비슷한 논리를 가진 방법, 가위치법(False Position Method)입니다.

False-Position Method

False-Position Method은 Bracketing Method. 이분법과 마찬가지로 새로운 근을 추정한 다음, 또 다시 새로운 구간을 정하며 근을 찾아나가는 방법입니다.

 

Bisection 과 False-Position Method의 차이점은 Bisection은 구간 끝의 두점을 계속 반으로 잘가며 근을 찾는 방법이라면, False-Position Method는 구간 끝의 두 점을 잇는 선을 그어 근을 찾는 방법입니다. 기존 Bisection Method에 비해 수렴 속도가 더 빠르다라는 장점을 갖고 있습니다.

Algorithm

폐구간 양 끝의 두점을 잇는 선을 그어 근을 찾는 방법

False-Postion를 통해 근을 계산하는 방법은 다음과 같습니다.

 

1. 탐색 구간 $[x_l, x_u]$과 임계 기준값(tolerance)를 선정합니다.
2. 선정된 구간의 양 끝점에서의 함수 값을 곱하여 음수인지 양수인지 확인합니다.
3. 가상의 직석의 x절편을 근이라 결정하고 $x_r$로 지정.
4. $f(x_r)$을 $f(x_l)$와 $f(x_r)$중 부호가 같은 쪽과 값을 비교.
5. 비교된 값의 차이가 임계값안으로 들어오면 해당되는 $x_r$을 근으로 추정.
6. 만약, 오차값이 임계 기준값 안으로 들어오지 못하면, $x_r$혹은 $x_l$을 $x_r$로 업데이트하여 4번으로 되돌아가 다시 오차값안으로 들어올때까지 iterate한다. 

Result

Bisection Method의 결과 사진(iteration)을 비교해보면 알수 있듯, 보다 빠르게 근에 수렴 하는 것을 확인 할 수 있습니다.

False-Postion Method를 활용한 근 찾기

iter 횟수 및 결과

MATLAB 코드

clc;
clear all;

x_l = 1; x_u = 2; x_r = 0;
tol_x = 0.001;
tol_y = 0.0000001;

y_l = func(x_l);
y_u = func(x_u);

if(func(x_l) * func(x_u) > 0)
	fprintf(1, "No root exists\n");
else
    iter = 0;
    fprintf(1, '\n\n');
    fprintf(1,' iter|   x_l    |    x_u   |    x_r   |    y_l    |    y_r   \n'); 
    fprintf(1,'------------------------------------------------------------\n');
    while(1)
        iter = iter + 1;
        x_r = (x_l*y_u - x_u*y_l) / (y_u - y_l);
        y_r = func(x_r);
        fprintf(1,'  %d  | %f | %f | %f | %f | %f\n',iter, x_l, x_u, x_r, y_l ,y_r);

        if(y_l * y_r < 0)
            if(abs(x_r - x_u) < tol_x || abs(y_r) < tol_y)
        	    fprintf(1,'Approximate solution c = %.7f\n', x_r);
                break;
            else
                x_u = x_r;
                y_u = y_r;
            end
        else
            if(abs(x_r - x_l) < tol_x || abs(y_r) < tol_y)
        	    fprintf(1,'Approximate solution c = %.7f\n', x_r);
                break;
            else
                x_l = x_r;
                y_l = y_r;
            end
        end
    end
end

x = linspace(-1, 3);
y = func(x);

figure;
plot([-1,3], [0,0], x,y, '-', x_r, y_r, '*', LineWidth=2);
grid on;
% % set(gca, 'XAxisLocation', 'origin');
set(gca, 'YAxisLocation', 'origin');

function y = func(x)
	y = x.^3 - x.^2 - x - 1;
end