Numeričke metode za rešavanje nelinearnih jednačina

U ovom dokumentu ćemo dati MATLAB implementacije metoda za rešavanje nelinearnih jednačina koje su obrađivane na časovima vežbi.

Metoda polovljenja intervala:

Ovo je verovatno idejno najjednostavnija metoda za rešavanje nelinearnih jednačina i zasniva se na određivanju niza umetnutih intervala takvih da je svaki sledeći ona polovina prethodnog intervala koja sadrži tačno rešenje. Niz aproksimacija koji konvergira tačnom rešenju čine središta tako određenih intervala.

MetodaPolovljenja.png

U nastavku je data jedna verzija implementacije. Uz niz rešenja dobijenih ovom metodom, data je i prateći grafički prikaz rada metode, kao i grafik greške u logaritamskoj skali. Program staje kad vrednost funkcije u tekućoj aproksimaciji postane manja od zadate tolerancije. Greška aproksimacije se poredi sa rešenjem dobijenim korišćenjem ugrađene MATLAB funkcije fzero().

% function [x,fval] = Polovljenje(f,a,b,n,tol)

% %Implementacija metode polovljenja

% %f neka funkcija

% %a i b su krajevi intervala izolacija

% %n maksimalan broj iteracija

% %tol zadata tolerancija

%

% if (f(a)*f(b)>0)

% error('Funkcija mora biti razlicitog znaka u krajevima!');

% end

%

% xZero = fzero(f,[a,b]);

% fprintf('Resenje dobijeno MATLAB funkcijom fzero(): %.8f\n', xZero);

%

% x = (a+b)/2;

% fval = f(x);

%

% fig = figure(1);

% set(fig, 'Name', 'Prozor', 'NumberTitle', 'off');

% subplot(2,1,1)

% fplot(f,[a,b],'-r','LineWidth',1.5);

% text(xZero, -0.05, '$\xi$', 'interpreter', 'latex');

% hold on

% grid on

% title('Metoda polovljenja intervala');

% xlabel('x');

% ylabel('y');

% subplot(2,1,1)

% plot(a,0,'.b','MarkerSize',20);

% plot(b,0,'.b','MarkerSize',20);

% line([a b], [0 0],'Color','black');

% plot(x,0,'.r','MarkerSize',20);

% s = text(x, 0.05, '$x_0$', 'interpreter', 'latex');

%

% err(1) = abs(x-xZero);

% subplot(2,1,2)

% semilogy([0,n],[0,0])

% title('Greska metode');

% grid minor

% hold on

% semilogy(0,err(1),'*b');

%

% fprintf('Rezultati po iteracijama:\n');

% for i = 1:n

% if (abs(fval) < tol)

% msgbox('Postignuta je trazena tacnost');

% pause(1);

% break;

% end

% pause(1);

% s.String = '';

% if (f(a)*f(x)<0)

% subplot(2,1,1);

% plot(b,0,'.k','MarkerSize',20);

% b=x;

% plot(b,0,'.b','MarkerSize',20);

% else

% subplot(2,1,1);

% plot(a,0,'.k','MarkerSize',20);

% a=x;

% plot(a,0,'.b','MarkerSize',20);

% end

% x = (a+b)/2;

% fval = f(x);

% subplot(2,1,1)

% plot(x,0,'.r','MarkerSize',20);

% s = text(x, 0.05, ['$x_{', num2str(i),'}$'], 'interpreter', 'latex');

% err(i+1) = abs(x-xZero);

% subplot(2,1,2)

% semilogy([i-1 i], [err(i) err(i+1)],'-*b');

% fprintf('%d %f %f\n',i,x,fval);

% end

% figure(1)

% hold off

%

% fprintf('Greske kroz iteracije:\n');

% disp(err');

%

% end

f = @(x)(sqrt(x+2)-2*cos(x));

a = -1.5;

b = -1;

n = 15;

tol = 1e-3;

[x,fval] = Polovljenje(f,a,b,n,tol)

Resenje dobijeno MATLAB funkcijom fzero(): -1.06666613
Rezultati po iteracijama:
1   -1.125000   0.073061
2   -1.062500   -0.005133
3   -1.093750   0.033657
4   -1.078125   0.014181
5   -1.070313   0.004503
6   -1.066406   -0.000321
Greske kroz iteracije:
    0.1833
    0.0583
    0.0042
    0.0271
    0.0115
    0.0036
    0.0003
x = -1.0664
fval = -3.2057e-04

Metoda tangenti

Njutnova metoda ili metoda tangenti ima jednostavnu geometrijsku interpretaciju: Ako nam je poznata aproksimacija rešenja , tada se sledeća aproksimacija dobija u preseku tangente na krivu u tački i -ose.

MetodaTangenti.png

Kao i za prethodnu metodu, i ovde program vraća niz aproksimacija kao i grafički prikaz rada metode i grafik greške.

% function [x,fval] = MetodaTangenti(f,df,a,b,x0,n,tol)

%

% if (f(a)*f(b)>0)

% error('Funkcija mora biti razlicitog znaka u krajevima\n');

% end

%

% xZero = fzero(f,[a,b]);

% fprintf('Resenje dobijeno MATLAB funkcijom fzero(): %.8f\n', xZero);

%

% %Grafik funkcije i metode

% fig = figure(2);

% subplot(2,1,1)

% set(fig, 'Name', 'Prozor', 'NumberTitle', 'off');

% a0 = a-(b-a)/4;

% b0 = b+(b-a)/4;

% fplot(f,[a0,b0],'-r','LineWidth',1.5);

% hold on

% text(xZero,-0.3,'$\xi$','Interpreter','latex');

% grid on

% title('Metoda tangenti');

% xlabel('x');

% ylabel('y');

% line([a0 b0], [0 0]);

% plot(a,0,'.b','MarkerSize',20);

% plot(b,0,'.b','MarkerSize',20);

%

% %Grafik greske

% err(1) = abs(x0-xZero);

% subplot(2,1,2)

% semilogy([0,n],[0,0]);

% title('Greska metode');

% hold on

% grid minor

% semilogy(0,err(1),'*b');

%

% xOld = x0;

% fval = f(x0);

% subplot(2,1,2)

% plot(x0,0,'.r','MarkerSize',20)

% s = text(x0, 0.3, '$x_0$', 'interpreter', 'latex');

% i = 1;

% pause(1);

% fprintf('Rezultati kroz iteracije:\n');

% while (i<n && fval>tol)

% subplot(2,1,1);

% s.String = '';

% if i>1

% plot(xOld,0,'.k','MarkerSize',20);

% end

% line([xOld,xOld],[0,fval],'LineStyle','--');

% pause(0.5);

% x = xOld - f(xOld)/df(xOld);

% line([xOld,x],[fval,0]);

% plot(x,0,'.r','MarkerSize',20);

% s = text(x, 0.3, ['$x_{', num2str(i),'}$'], 'interpreter', 'latex');

% fval = f(x);

% fprintf('%d %f %f\n', i, x, fval);

% err(i+1) = abs(x-xZero);

% subplot(2,1,2)

% semilogy([i-1 i], [err(i) err(i+1)],'-*b');

% xOld = x;

% i = i+1;

% pause(1);

% end

% hold off

% fprintf('Greske kroz iteracije:\n');

% disp(err');

%

% end

f = @(x)x.*exp(x)+x.^2-1;

df = @(x)(x+1).*exp(x)+2*x;

a = -2;

b = -1;

x0 = -2;

n = 10;

tol = 1e-4;

[x,fval] = MetodaTangenti(f,df,a,b,x0,n,tol);

Resenje dobijeno MATLAB funkcijom fzero(): -1.16758553
Rezultati kroz iteracije:
1   -1.339998   0.444721
2   -1.179392   0.028343
3   -1.167651   0.000155
4   -1.167586   0.000000
Greske kroz iteracije:
    0.8324
    0.1724
    0.0118
    0.0001
    0.0000

Metoda regula falsi

Metoda regula falsi se može shvatiti kao modifikacija obe prethodne metode: Sa jedne strane se sastoji u određivanju niza umetnutih intervala koji ne moraju biti polovina intervala koji im prethodi, a sa druge se može smatrati modifikacijom Njutnove metode jer predstavlja prirodno rešenje za zamenu izvoda koji se javlja u formuli te metode. Ako su nam poznate uzastopne aproksimacije rešenja , koje određuju interval izolacije, tada se naredna javlja u preseku duži određene tačkama i i -ose. Za nastavak procesa se ažurira interval izolacije.

RegulaFalsi.png

U nastavku je data varijanta implementacije ove metode. Prosleđuje se početna aproksimacija , a onda se bira novi interval izolacije između i . Pretpostavka je da su pravilno uneti ulazni podaci. Program daje i grafički prikaz metoda (izostavljen je grafik greške).

% function [x,fval] = RegulaFalsi(f,a,b,x0,n,tol)

%

% xZero = fzero(f,[a,b]);

% fprintf('Resenje dobijeno MATLAB funkcijom fzero(): %.8f\n', xZero);

%

% fig = figure(3);

% set(fig, 'Name', 'Prozor', 'NumberTitle', 'off');

% a0 = a-(b-a)/4;

% b0 = b+(b-a)/4;

% x = x0;

% fval = f(x);

%

% d = 0.1;

% fplot(f,[a0,b0],'-r','LineWidth',1.5);

% hold on

% text(xZero, -d, '$\xi$', 'interpreter', 'latex');

% grid on

% title('Metoda regula falsi');

% line([a0 b0], [0 0]);

% xlabel('x');

% ylabel('y');

% plot(a,0,'.b','MarkerSize',20);

% plot(b,0,'.b','MarkerSize',20);

% plot(x,0,'.r','MarkerSize',20);

% s = text(x, d, '$x_0$', 'interpreter', 'latex');

%

% if (fval*f(a)<0)

% leftEnd = a;

% rightEnd = x;

% fLeft = f(a);

% fRight = fval;

% else

% leftEnd = x;

% rightEnd = b;

% fLeft = fval;

% fRight = f(b);

% end

% fprintf('Rezultati po iteracijama:\n');

% fprintf('0: %f %f\n',x,fval);

% for i = 1:n

% if (abs(fval) < tol)

% break;

% end

% pause(1);

% s.String = '';

% line([leftEnd leftEnd],[0 fLeft],'LineStyle','--');

% line([rightEnd rightEnd],[0 fRight],'LineStyle','--');

% pause(0.5);

% line([leftEnd rightEnd],[fLeft,fRight],'LineStyle','--');

% if (f(leftEnd)*fval<0)

% plot(rightEnd,0,'.k','MarkerSize',20);

% rightEnd=x;

% fRight = f(rightEnd);

% plot(rightEnd,0,'.b','MarkerSize',20);

% else

% plot(leftEnd,0,'.k','MarkerSize',20);

% leftEnd=x;

% fLeft = f(leftEnd);

% plot(leftEnd,0,'.b','MarkerSize',20);

% end

% xOld = x;

% x = (leftEnd*fRight-rightEnd*fLeft)/(fRight-fLeft);

% fval = f(x);

% plot(x,0,'.r','MarkerSize',20);

% s = text(x, d, ['$x_{', num2str(i),'}$'], 'interpreter', 'latex');

% fprintf('%d: %f %f\n',i,x,fval);

% end

% figure(3)

% hold off

% end

f = @(x)(log(x+1)-2*x.^2+1);

a = -0.5;

b = 0;

x0 = 0;

n = 10;

tol = 5e-4;

[x,fval] = RegulaFalsi(f,a,b,x0,n,tol)

Resenje dobijeno MATLAB funkcijom fzero(): -0.44921528
Rezultati po iteracijama:
0:   0.000000   1.000000
1:   -0.419060   0.105670
2:   -0.447683   0.005528
3:   -0.449138   0.000278
x = -0.4491
fval = 2.7793e-04

Metoda iteracije

Ideja metode je da se iz jednačine odredi ekvivalentan oblik za funkciju takvu da je ispunjeno odakle bi se dobio iterativni niz

za neko , koji konvergira ka tačnom rešenju. Za ovu metodu postoje dve ocene greška koje ćemo koristiti:

  1. apriorna ocena greške:

2. aposteriorna ocena greške:

Prvu ocenu ćemo koristiti da odredimo broj iteracija koji nam garantuje neku zadatu tačnost, dok ćemo drugu koristiti za izlazni kriterijum:

U nastavku je data jednostavna implementacija metode uz prateći grafički prikaz. Određivanje funkcije g(x) i parametra q je netrivijalno i ostavlja se korisniku, te su i to ulazni parametri.

% function [x,fval] = Iteration(f,a,b,g,q,x0,tol)

%

% xZero = fzero(f,[a,b]);

% fprintf('Resenje dobijeno MATLAB funkcijom fzero(): %.7f\n', xZero);

% fig = figure(4);

% set(fig, 'Name', 'Prozor', 'NumberTitle', 'off');

% a0 = a-(b-a)/4;

% b0 = b+(b-a)/4;

% d = 0.05;

% fplot(g,[a0,b0],'-r','LineWidth',1.5);

% hold on

% line([a0 b0],[a0 b0]);

% L = legend('$y=g(x)$', '$y=x$');

% set(L,'Location','southeast','interpreter','latex');

% grid on

% axis equal

% title('Metoda iteracije');

% xlabel('x');

% ylabel('y');

% line([a0 b0], [0 0]);

% plot(a,0,'.b','MarkerSize',20);

% plot(b,0,'.b','MarkerSize',20);

%

% plot(x0,0,'.r','MarkerSize',20)

% s = text(x0, d, '$x_0$', 'interpreter', 'latex');

% x = g(x0);

% line([x0 x0],[0 x],'LineStyle','--');

% pause(0.5);

% line([x0 x],[x x],'LineStyle','--');

% pause(0.5);

% line([x x],[x 0],'LineStyle','--');

% s.String = '';

% plot(x0,0,'.k','MarkerSize',20)

% plot(x,0,'.r','MarkerSize',20)

% s = text(x, d, '$x_1$', 'interpreter', 'latex');

% xOld =x0;

% k = ceil(log((1-q)*tol/abs(x-x0))/log(q)); %broj iteracija koji garantuje tacnost

% fprintf('Rezultati po iteracijama:\n');

% fprintf('0: %.7f %.7f',x0,f(x0));

% fprintf('1: %.7f %.7f',x,f(x));

% err = tol*(1-q)/q;

% for i=2:k

% xOld = x;

% x = g(x);

% line([xOld xOld],[0 x],'LineStyle','--');

% pause(0.5);

% line([xOld x],[x x],'LineStyle','--');

% pause(0.5);

% line([x x],[x 0],'LineStyle','--');

% pause(0.5)

% s.String = '';

% plot(xOld,0,'.k','MarkerSize',20)

% plot(x,0,'.r','MarkerSize',20)

% s = text(x, d, ['$x_{', num2str(i),'}$'], 'interpreter', 'latex');

% fval = f(x);

% fprintf('%d: %.7f %.7f\n',i,x,fval);

% if abs(x-xOld)<err

% break;

% end

% end

% text(xZero, -d, '$\xi$','interpreter', 'latex','color','r');

% hold off

% end

f = @(x)exp(x)-2*x-2;

a = -1;

b = 0;

g = @(x)(exp(x)-2)/2;

q = 0.5;

x0 = 0;

tol = 1e-4;

x = Iteration(f,a,b,g,q,x0,tol);

Resenje dobijeno MATLAB funkcijom fzero(): -0.7680390
Rezultati po iteracijama:
0: 0.0000000    -1.0000000
1: -0.5000000    -0.3934693
2: -0.6967347    -0.1083212
3: -0.7508953    -0.0262656
4: -0.7640281    -0.0061574
5: -0.7671068    -0.0014318
6: -0.7678227    -0.0003323
7: -0.7679889    -0.0000771
8: -0.7680274    -0.0000179