Nelinearne jednačine

19. zadatak

Neka je funkcija zadata eksplicitno komandnim M-fajlom $funkcija.m$

$funkcija.m$

​x
f = @(x) cos(x) - x;
​
% Test: Njutn(0,0,1e-5)
% Tačna vrednost: fzero(f,0)
% Za poziv Njutn(0,0,1e-500) biće dostignut maksimalan broj iteracija

• Napisati M-fajl sa funkcijom $x=Njutn(x0,y,tol)$ koja za unete argumente $x0,y$ i $tol$ vraća rešenje jednačine $f(x)=y$ (gde je $x0$ početna vrednost iterativnog procesa) izračunato Njutnovom metodom sa tačnošću $tol$ . Kriterijum zaustavljanja je $|x_n-x_{n-1}|<tol$. Broj iteracija ograničiti na 100. U slučaju da je dostignut maksimalan broj iteracija štampati odgovrajuću poruku. (Pretpostavka je da su ispunjeni svi uslovi za primenu metode.)

xxxxxxxxxx
function x = Njutn(x0, y, tol)
  
  % Ulazni argumenti su:
  %     x0 - početna vrednost iterativnog niza
  %     y - desna strana jednačine f(x)=y
  %     tol - tačnost
  % Izlazni argument:
  %     x - vrednost argumenta za koju je f(x)=y
​
% Pomoćni fajl:
%     funkcija.m - sadrzi anonimnu funkciju f
​
funkcija;
​
% Njutnova metoda rešava jednačinu oblika F(x)=0
f = @(x) f(x) - y;
​
F = sym(f); % konvertujemo funkciju u simboličku
F1 = diff(F); % izvod simboličke funkcije
f1 = matlabFunction(F1); % konvertujemo simboličku funkciju F1 u anonimnu funkciju
​
iter = 0;
while iter < 100
    x1 = x0 - f(x0) / f1(x0);
    iter = iter + 1;
    if abs(x1-x0) < tol
           break;
    end
    x0 = x1;
end
​
x = x1;
if iter == 100
    disp('Dostignut je maksimalan broj iteracija.')
end

• Napisati M-fajl $tablica.m$ u kome su zadati vektori $Y$ i $x_0$ iste dužine $n$ , kao i vrednost $tol$. Pretpostavka je da su elementi vektora $Y$ različiti. U M-fajlu se formira tablica $(X_1,Y_1),...,(X_n,Y_n)$ gde su $X_i, i=1,...,n$ rešenja jednačina $f(X_i)=Y_i$ dobijena korišćenjem funkcije iz prethodne tačke. Vektor $x_0$ sadrži odgovarajuće početne vrednosti za iterativni proces. Tablicu štampati u komandnom prozoru u formatu $X: X(1)\ X(2) ... X(n)$ $Y: Y(1)\ Y(2) ... Y(n).$

  xxxxxxxxxx
  % vrednosti su izabrane proizvoljno
  Y = 1:5; % desna strana jednačine f(x)=y
  n = length(Y);
  tol = 10^(-5); % tolerancija
  x0 = zeros(1, n); % početne vrednosti
  X = zeros(1, n); % koreni
  ​
  for i = 1:n
      X(i) = Njutn(x0(i), Y(i), tol);
  end
  ​
  % ugrađena funkcija num2str(x,n) konvertuje broj x u string 
  % i zapisuje ga sa n značajnih cifara 
  disp(['X: ', num2str(X, 5); 'Y: ', num2str(Y, 5)])
  

• Napisati M-fajl $vredfunk.m$ sa funkcijom $y=vredfunk(x)$ koja vraća vrednost Lagranžovog interpolacionog polinoma u tački $x$ dobijenog korišćenjem svih vrednosti iz formirane tablice. (Niz čvorova $x_i, i=1,...,n$ ne mora biti rastući)

  xxxxxxxxxx
  function y = vredfunk(x)
  % Ulazni argument: x - tačka u kojoj se računa vrednost funkcije
  % Izlazni argument: y - približna vrednost funkcije
  % Test: vredfunk(-3)
  ​
    tablica;
    y = 0;
    for i = 1:n
        p = 1;
        for j = 1:n
            if j ~= i
                p = p * (x - X(j)) / (X(i) - X(j));
            end
        end
        y = y + p * Y(i);
    end
  

20. zadatak

Neka je funkcija $f(x)$ zadata tablično M-fajlom $tablica.m$ koji generiše dva niza $X=[x_1,...,x_n]$ i $F=[f_1,...,f_n]$ (od kojih je prvi strogo rastući) za tu tablično zadatu funkciju. Tablica ne mora biti ekvidistantna.

$tablica.m$

xxxxxxxxxx
% --------- I test --------------
X = linspace(0, 4, 5);
F = X - 1;
% -------- II test ---------------
% X = linspace(0,2,6);
% F = X.^3-exp(X);
% %F = exp(X).*sin(X)-cos(X);

• Napisati M-fajl $funk.m$ sa funkcijom $y=funk(x)$ koja najpre formira Njutnov interpolacioni polinom sa podeljenim razlikama na osnovu svih vrednosti vektora $X$ i $F$ iz fajla $tablica.m$ , a zatim vraća vrednost formiranog polinoma za ulazni argument $x$ .

  x
  function y = funk(x)
  ​
  tablica;
  ​
  n = length(X);
  pr = zeros(n, n-1);
  ​
  for i = 1 : n-1
      pr(i, 1) = (F(i+1) - F(i)) / (X(i+1) - X(i));
  end
  ​
  for j = 2 : n-1
      for i = 1 : n-j
          pr(i, j) = (pr(i+1, j-1) - pr(i, j-1)) / (X(i+j) - X(i));
      end
  end
  ​
  [X', F', pr]
  ​
  % Konstruišemo interpolacioni polinom 
  polinom = [F(1)];
  r = [1];
  ​
  for i = 1 : n-1
      r = conv(r, [1, -X(i)]); % množimo polinome
      % svaki sledeći sabirak je polinom stepena za jedan većeg od prethodnog
      % dodavanjem nule na početak, omogućavamo sabiranje ova dva polinoma
      polinom = [0, polinom] + r * pr(1, i);  
  end
  y = polyval(polinom, x);
  ​
  ​
  

• Napisati M-fajl $polov.m$ sa funkcijom $x=polov(tol)$ koja na intervalu $[x_1,x_n]$ računa i vraća rešenje jednačine $funk(x)=0$ metodom polovljenja intervala sa tačnošću $tol$ . Funkcija treba da proveri da li su uslovi za primenu metode polovljenja intervala ispunjeni, da prekine program i vrati poruku ukoliko nisu. Prvi i poslednji element vektora $X (x_1$ i $x_n)$ dobijaju se pozivanjem fajla $tablica.m$

  xxxxxxxxxx
  function x = polov(tol)
  ​
  tablica;
  a = X(1);
  b = X(end);
  ​
  if (funk(a) * funk(b) > 0)
      error('Funkcija ne menja znak na posmatranom segmentu')
  end
  ​
  if funk(a) == 0
      x = a;
  else
      if funk(b) == 0
          x = b;
      end
  end
  ​
  if (funk(a) ~= 0 && funk(b) ~= 0)
      iter = 0;
      while (abs(a-b) > tol)
          x = (a + b) / 2;
          iter = iter + 1;
          if funk(x) == 0
              break
          else
              if funk(a) * funk(x) < 0
                  b = x;
            else
                  a = x;
              end
          end
      end
      x=(a + b) / 2;
  end
  ​
  

  21. zadatak

  Neka je funkcija $f$ zadata tablično M-fajlom $tablica.m$ koji generiše dva niza $X=[x_1,...,x_n]$ i $F=[f_1,...,f_n]$ (od kojih je prvi strogo rastući) za tu tablično zadatu funkciju. Tablica ne mora biti ekvidistantna.

  $tablica.m$

xxxxxxxxxx
X = [0 0.3 0.55 0.60 0.98 1.07 1.27 1.33 1.6 2];
F = [0.2500    0.5455    0.7727    0.8146    1.0805    1.1272    1.2051    1.2211    1.2496    1.1593];
​
%Drugi primer
%X=[1:10];
%F=log(X)+3;
%Test: nula(1,1e-4,20)

• Napisati M-fajl $funk.m$ sa funkcijom $y=funk(x)$ koja za unetu vrednost argumenta $x$ vraća $y$, približnu vrednost funkcije u toj tački izračunatu pomoću Njutnovog interpolacionog polinoma sa podeljenim razlikama, koristeći sve vrednosti iz M-fajla $tablica.m$.

  xxxxxxxxxx
  function y=funk(x)
  ​
  tablica;
  n=length(X);
  pr=zeros(n,n-1);
  ​
  for i=1:n-1
      pr(i,1)=(F(i+1)-F(i))/(X(i+1)-X(i));
  end
  ​
  for j=2:n-1
      for i=1:n-j
          pr(i,j)=(pr(i+1,j-1)-pr(i,j-1))/(X(i+j)-X(i));
      end
  end
  ​
  L=F(1);
  p=1;
  for i=1:n-1
      p=conv(p,[1,-X(i)]);
      L=[0 L]+p*pr(1,i);
  end
  y=polyval(L,x);
  ​
  % Nije bilo neophodno konstruisati polinom (obzirom da nije naglašeno u 
  % zadatku)
  

• Napisati M-fajl $nula.m$ sa funkcijom $[x,briter]=nula(x_0,tol,iterM)$ koja računa i vraća $x$, rešenje jednačine $funk(x)=0$ metodom proste iteracije sa tačnošću $tol$ , kao i broj iteracija $briter$. Kriterijum zaustavljanja je: $|x_n-x_{n-1}|<tol.$ Broj iteracija ograničiti na $iterM$ i ispisati poruku da li je zadovoljena tačnost ili je dostignut maksimalan broj iteracija. Vrednosti funkcije u tačkama koje se javljaju kroz iteracije računati pomoću funkcije $funk().$ Pretpostavka je da je funkcija na tom intervalu kontrakcija. Za početnu tačku iterativnog niza uzeti tačku $x_0.$

  x
  function [x, briter] = nula(x0, tol, iterM)
  ​
  briter=0;
  ​
  while (briter < iterM)
      x1 = funk(x0);
      briter = briter + 1;
      if abs(x0-x1) <= tol
          break
      end
      x0 = x1;
  end
  ​
  x = x1;
  ​
  if briter == iterM
      disp('Dostignut je maksimalan broj iteracija!')
  end
  ​
  if abs(x0-x1) <= tol
      disp('Tacnost je zadovoljena!')
  end
  ​
  

• Grafički prikazati, funkcijom $grafik(x_0, iterM)$u M-fajlu $grafik.m$, zavisnost brzine konvergencije od tačnosti $tol$ ako se ona kreće od $10^{-4}$ do $10^{-3}$ sa korakom $10^{-4}$ . ( Pod brzinom konvergencije podrazumeva se broj iterativnih koraka.)

  xxxxxxxxxx
  function grafik(x0, iterM)
  ​
  tol = 1e-4:1e-4:1e-3;
  n = length(tol);
  briter = zeros(1, n);
  ​
  for i = 1 : n
      [x, briter(i)] = nula(x0, tol(i), iterM);
  end
  ​
  plot(tol,briter);
  xlabel('Tolerancija');
  ylabel('Broj iteracija');
  ​
  title('Zavisnost brzine konvergencije od tacnosti tol');
  

22. zadatak

Neka je funkcija $f$ zadata tablično M-fajlom $tablica.m$ koji generiše dva niza $X=[x_{1},x_{2}...,x_{n}]$ i $F=[f_{1},f_{2}...,f_{n}]$ (od kojih je prvi strogo rastući i ekvidistantan) za tu tablično zadatu funkciju.

$tablica.m$

xxxxxxxxxx
X = 0:0.1:0.5;
F = exp(X) - 2 * (X - 1).^2;
​
% Test: [I,Y,x]=nula(0,0.5,1e-5,15) %zadovoljena tacnost tol
% Test: [I,Y,x]=nula(0,0.5,1e-5,4)  %zadovoljen max br. iteracija 
% Test: [I,Y,x]=nula(0,0.5,1e-5,6)  %zadovoljena je tacnost i dostignut max
% br. iteracija

• Napisati M-fajl $Njutn1.m$ sa funkcijom $koef=Njutn1()$ koja vraća koeficijente I Njutnovog interpolacionog polinoma (po promenljivoj $q$ ), koristeći sve vrednosti iz tablice.

  xxxxxxxxxx
  function koef = Njutn1()
  ​
  tablica;
  n = length(X);
  ​
  KR = zeros(n, n-1); % tablica konacnih razlika
  for i = 1:n - 1
      KR(i, 1) = F(i+1) - F(i);
  end
  ​
  for j = 2:n - 1
      for i = 1:n - j
          KR(i, j) = KR(i+1, j-1) - KR(i, j-1);
      end
  end
  %disp([X' Y' KR]);
  ​
  % Konstruišemo I Njutnov interpolacioni polinom po promenljivoj q:
  yk = F(1);
  qk = [1, 0];
  for j = 1:n - 1
      yk = [0, yk] + qk * KR(1, j) / factorial(j);
      qk = conv(qk, [1, -j]);
  end
  koef = yk;
  ​
  

• Napisati M-fajl $nula.m$ sa funkcijom $[I,Y,x]=nula(x_0,x_F,tol,iterM)$ koja računa i vraća nulu $x$ tablično zadate funkcije iz $tablica.m$ metodom regula-falsi sa tačnošću $tol$ . Kriterijum zaustavljanja je $|x_n-x_{n-1}|<tol$ gde su $x_n$ i $x_{n-1}$ dve uzastopne tačke dobijene metodom regula-falsi. Broj iteracija ograničiti na $iterM$ i ispisati poruku da li je zadovoljena tačnost ili je dostignut maksimalan broj iteracija. Za fiksiranu tačku u metodi uzeti $x_F$ , a za početnu vrednost iterativnog procesa $x_0$ . Vrednosti funkcije u tački računati pomoću I Njutnovog interpolacionog polinoma dobijenog funkcijom $Njutn1()$. U vektore $I$ i $Y$ upisivati redni broj iteracije i vrednost funkcije u tačkama iz iterativnog niza, redom. (Pretpostavka je da su ispunjeni svi uslovi za primenu metode.)

  xxxxxxxxxx
  function [I, Y, x] = nula(x0, xF, tol, iterM)
  ​
  % Ulazni argumenti su:
  %     x0 - početna vrednost iterativnog niza
  %     xF - fiksirana tačka u metodi regula-falsi
  %     tol - tačnost
  %     iterM - maksimalni dozvoljeni broj iteracija
  % Izlazni argumenti su:
  %     I - vektor [1:iter] gde je iter broj iteracija
  %         dok se ne zadovolji tačnost tol ili iterM
  %     Y - vrednosti tablično zadate funkcije u tačkama iterativnog niza
  %         izračunata pomoću I Njutnovog interpolacionog polinoma
  %     x - nula tablično zadate funkcije
  ​
    tablica;
    h = X(2) - X(1);
    koef = Njutn1(); % Koeficijenti I Njutnovog int. pol. po q
  ​
    % Anonimna f-ja koja računa vrednost I Njutnovog polinoma u tački x 
    % Polinom je po q, zato nećemo direktno slati x, već odgovarajuću 
    % vrednost za q, tj. (x-X(1)/h))
    f = @(x)polyval(koef, (x - X(1))/h); 
  ​
    Y = zeros(1, iterM); % Y ima najviše iterM elemenata
    iter = 0;
  ​
    while iter < iterM
        x1 = x0 - (f(x0) / (f(xF) - f(x0))) * (xF - x0); % regula-falsi
        iter = iter + 1;
        Y(iter) = f(x1);
        greska = abs(x1-x0);
        if greska < tol
            break
        end
        x0 = x1;
    end
  ​
    I = 1:iter;
    Y = Y(1:iter);
    if iter == iterM && greska < tol
        disp('Zadovoljena je tacnost i dostignut max br iteracija.');
    elseif iter == iterM
        disp('Zadovoljen je maksimalan broj iteracija.');
    else
        disp('Zadovoljena je tacnost tol');
    end
    x = x1;
  

23. zadatak

Neka je funkcija $f$ zadata tablično M-fajlom $tablica.m$ koji generiše dva niza $X=[x_{1},x_{2}...,x_{n}]$ i $F=[f_{1},f_{2}...,f_{n}]$ (od kojih je prvi strogo rastući) za tu tablično zadatu funkciju. Tablica ne mora biti ekvidistantna.

$tablica.m$

xxxxxxxxxx
X = 1:0.1:1.5;
F = sin(X) - X.^2 + 1;
​
% Test:  x=nula(1e-4,10)

• Napisati M-fajl $Lagranz.m$ sa funkcijom $L=Lagranz()$ koja vraća koeficijente Lagranžovog interpolacionog polinoma, koristeći sve vrednosti iz tablice.

xxxxxxxxxx
function L = Lagranz()
​
tablica;
​
n = length(X);
L = zeros(1, n);
​
for i = 1:n
    p = 1;
    for j = 1:n
        if i ~= j
            p = conv(p, [1, -X(j)]) / (X(i) - X(j));
        end
    end
    L = L + p * F(i);
end

• Napisati M-fajl $nula.m$ sa funkcijom $x=nula(tol,iterM)$ koja računa i vraća nulu $x$ tablično zadate funkcije iz $tablica.m$ metodom sečice sa tačnošću $tol$ . Kriterijum zaustavljanja je $|x_n-x_{n-1}|<tol$, gde su $x_n$ i $x_{n-1}$ dve uzastopne tačke dobijene metodom sečice. Za prve dve iteracije uzeti krajeve intervala. Broj iteracija ograničiti na $iterM$ i ispisati poruku da li je zadovoljena tačnost ili je dostignut maksimalan broj iteracija. Vrednosti funkcije u tački računati pomoću Lagranžovog interpolacionog polinoma dobijenog pozivom funkcije $Lagranz()$. Funkcija treba da ispisuje redni broj iteracije i vrednost funkcije u tačkama iz iterativnog niza, redom. (Pretpostavka je da su ispunjeni svi uslovi za primenu metode).

  xxxxxxxxxx
  function x = nula(tol, iterM)
  ​
  %Ulazni argumenti su:
    %     tol - tacnost
    %     iterM - maksimalni dozvoljeni broj iteracija
  %Izlazni argumenti su:
    %     x - nula tablicno zadate funkcije
    
  ​
    tablica;
    n = length(X);
    koef = Lagranz();
    f = @(x)polyval(koef, x);
    Y = zeros(1, iterM);
    iter = 0;
  ​
    x0 = X(1); %Xn-1
    x1 = X(n); %Xn
  ​
    while iter < iterM
        x2 = x1 - (f(x1) / (f(x0) - f(x1))) * (x0 - x1); % metoda sečice
        iter = iter + 1;
        Y(iter) = f(x2);
        greska = abs(x2-x1);
        if greska < tol
          break
        end
        x0 = x1;
        x1 = x2;
    end
    I = 1:iter;
    Y = Y(1:iter);
  ​
    if iter == iterM && greska < tol
        disp('Zadovoljena je tacnost i dostignut je maksimalan broj iteracija.');
    elseif iter == iterM
        disp('Zadovoljen je maksimalan broj iteracija.');
    else
        disp('Zadovoljena je tacnost tol')
    end
    x = x2;
  ​
    disp(['I: ', num2str(I, 5)])
    disp(['Y: ', num2str(Y, 5)])