function v=fdm(primer, h)
close
%% ucitavamo podatke
run(primer)
%% formiramo mrezu
x=a:h:b;
n=length(x);
%% matrica sistema
%poddijagonala
a=-1/h^2-p(x)/(2*h);
a(n)=0; %prvi element je fiktivni, necemo ga koristiti u racunu
%dijagonala
b=2/h^2+q(x);
b(1)=1;b(n)=1;
%naddiajgonala
c=-1/h^2+p(x)/(2*h);
c(1)=0;%poslednji element je fiktivni, ne koristi se u racunu
%M=diag(b)+diag(a(2:n),-1)+diag(c(1:n-1),1);
%disp(M);
F=f(x');
F(1)=A;
F(n)=B;
%v=M\F;
v=tomas(a,b,c,F);
if nargout==0
plot(x,ue(x),'b')
hold on
plot(x,v,'--or')
hold off
end
end
function x=tomas(a,b,c,f)
n=length(b);
%sistem svodimo na gornje trougaounu matricu
u=b; d=f;
for k=2:n-1
pom=a(k)/u(k-1);
u(k)=b(k)-c(k-1)*pom;
d(k)=f(k)-d(k-1)*pom;
end
%resavamo sistem
x=d; %inicijalizujemo sa d
for k=n-1:-1:1
x(k)=(d(k)-c(k)*x(k+1))/u(k);
end
end
