clear all
close all
clc

x=linspace(-1,1,101);
y =sqrt(1-x.^2); %funkcja
figure;
plot(x,y(x))
hold on;

xn = [-1 -0.5 0 0.5 1]; 
y = [y(xn)];

ylgr = []; 
ynwt = [];

for i=1:length(x) 
ylgr = [ylgr lagrange_int(xn,y,x(i))]; 
end 



 function xv = lagrange_int(xn,yn,x)
licznik = 1;
 mianownik = 1;
 sum = 0;
  N = length(xn);
     for j = 1:N
         for k = 1:N
             if (j~=k)
               mianownik = mianownik*(xn(j)-xn(k));
                 licznik = licznik*(x-xn(k));
             end
         end
         sum = sum+yn(j)*licznik/mianownik;
         licznik = 1;
         mianownik = 1;
     end
     xv = sum;
   end
   
   ____________________________________________________________________________________________
   
   clear all
close all
clc
  % policzyć błędy bezwzględne i względne dla poszczególnych funkcji
  % wyrysować błędy (w razie potrzeby zastosować skalę 
  % logarytmiczną w osi rzędnych)
  
  % efekt rungego - oscylacje na krańcach przedziału
y = @(x) 1./(1+25*x.^2);
  
 % moduł x, brak ciągłej pochodnej
y1 = @(x) abs(x);
 
 % nieograniczona pochodna na krańcach przedziału
y2 = @(x) sqrt(1-x.^2);
 
  % nieróżniczkowalna dla x=0
y3 = @(x) x.*exp(x);
x = linspace(-1,1,100);
figure;
plot(x,y(x),'r')
hold on;
 
xn = [-1 -0.5 0 0.5 1];
yn = [y(xn)];
  
ylgr = [];
ynwt = [];
blad_bzwzg=0;
blad_wzg=0;
for i=1:length(x)
      ylgr = [ylgr lagrange_int(xn,yn,x(i))];
%       ynwt = [ynwt newton_int(xn,yn,x(i))];
end
%  plot(x,ylgr,'b',x,ynwt,'g');
plot(x,ylgr,'b');
% plot(x,ynwt,'g')

%blad bezwzgl
figure
for i=1:length(x)
          blad_bzwzg=abs(ylgr-y(x))
end
plot(x,blad_bzwzg,'-.g');


%blad wzgledny
figure
for i=1:length(x)
    
          blad_wzg=[abs((ylgr-y(x))./y(x))]
          if y(x)==0
              blad_wzg=0;
          end
end
plot(x,blad_wzg,'r');


%  Lagrange

 function xv = lagrange_int(xn,yn,x)
licznik = 1;
 mianownik = 1;
 sum = 0;
  N = length(xn);
     for j = 1:N
         for k = 1:N
             if (j~=k)
               mianownik = mianownik*(xn(j)-xn(k));
                 licznik = licznik*(x-xn(k));
             end
         end
         sum = sum+yn(j)*licznik/mianownik;
         licznik = 1;
         mianownik = 1;
     end
     xv = sum;
   end
  
% Newton
%    function xv = newton_int(xn,yn,x)
%        xes = 1;
%        suma = yn(1);
%        N = length(xn);
%   
%        for j=1:N-1
%            ynew = [];
%            xes = xes*(x-xn(j));
%           for i=1:N-j
%               ynew=[ynew(yn(i+1)-yn(i))/(xn(i+j)-xn(i))];
%           end
%           yn = ynew;
%           suma = suma+ynew(1)*xes;
%       end
%       xv = suma;
%    end

-------------------------------------------------------------------------------------------------------
   clear all
   clc
   
   % policzyć błędy bezwzględne i względne dla poszczególnych funkcji
   % wyrysować błędy (w razie potrzeby zastosować skalę 
   % logarytmiczną w osi rzędnych)
   
   % efekt rungego - oscylacje na krańcach przedziału
   y = @(x) 1./(1+25*x.^2);
   
  % moduł x, brak ciągłej pochodnej
  %y = @(x) abs(x);
  
  % nieograniczona pochodna na krańcach przedziału
  %y = @(x) sqrt(1-x.^2);
  
  % nieróżniczkowalna dla x=0
  %y = @(x) x.*exp(x);
  
  
  x = linspace(-1,1,100);
  
  figure(1);
  plot(x,y(x),'r')
  grid on
  
  hold on;
  
  xn = [-1 -0.5 0 0.5 1];
  yn = [y(xn)];
  
  ylgr = [];
  ynwt = [];
  for j=1:length(x)
      ylgr = [ylgr lagrange_int(xn,yn,x(j))];
     % ynwt = [ynwt newton_int(xn,yn,x(j))];
  end
  
  %figure()
 % plot(x,ylgr,'b',x,ynwt,'g');
   plot(x,ylgr,'black',x,y(x),'red');
   title('Funkcja interpolujaca(czarna) i interpolowana(czerwona)')
   grid on
   hold off
   
   blad_b = 0;
   blad_w = 0;
   
  %%blad bezwzgledny
  figure(3)
  
  for j=1:length(x)
      
      blad_b = abs(ylgr - y(x));
      
  end
  plot(x,blad_b,'black');
  title('Blad bezwzgledny')
  grid on

  %%blad wzgledny
  figure(4)
  
  for j=1:length(x)
      
      if y(x) == 0 % zabezpieczam przed dzieleniem przez 0      
      blad_w = 0 % wartosc dla zachowania ciaglosci przebiegu
      %(chociaz matlab chyba sam przyjmuje wtedy f(x)=0)
      else
      blad_w = (abs(ylgr - y(x))./y(x));
      end
  end
  plot(x,blad_w,'black');
  title('Blad wzgledny')
  grid on
  
  %%funkcje
  
    function xv = lagrange_int(xn,yn,x)
     licznik = 1;
     mianownik = 1;
     sum = 0;
     N = length(xn);
     for j = 1:N,
         for k = 1:N,
             if (j~=k)                              
                 licznik = licznik * (x-xn(k));
                 mianownik = mianownik * (xn(j)-xn(k));        
            end
         end
         sum = sum + yn(j) * (licznik/mianownik);
         licznik = 1;
         mianownik = 1;
     end
     xv = sum;
  end


%    function xv = newton_int(xn,yn,x)
%        xes = ...;
%       suma = yn(1);
%       N = length(xn);
%   
%        for j=1:N-1,
%           ynew = [];
%            xes = ...;
%           for i=1:N-j,
%               ynew=[...];
%           end
%           yn = ...;
%           suma = ...;
%       end
%       xv = suma;
%    end

 --------------------------------------------------------------------