% Eksempel på bruk:
% GaussInt2(@(x) sin(x)/2, 0, pi, 2)
% integrerer sin(x)/2 fra x=0 til x=pi med 2 integrasjonspunkt
function f = GaussInt2(funk,x1,x2,n)
% Bestemmer integrasjonspunktene og tilhørende vekter
beta   = (1:n-1)./sqrt(4*(1:n-1).^2 - 1);
J      = diag(beta,-1) + diag(beta,1);    
[V,D]  = eig(J);
[x,ix] = sort(diag(D));  
w      = 2*V(1,ix)'.^2;
% Gjør om til integrasjon fra -1 til +1 og summerer alle leddene
a      = 0.5*(x2-x1);
b      = 0.5*(x1+x2);
f=0;
for i=1:n
    f=f+w(i)*funk(x(i)*a+b)*a;
end
end

