Romberg-algoritmen (Ma4n: Oppgave1, øving8)

Kopier teksten nedenfor til en teksteditor, og lagre som "romberg.m" (f.eks).

 function R = romberg(f,a,b,n)
    h = b-a;
    R = zeros(n,n);
    R(1,1) = (feval(f,a) + feval(f,b))*h/2;
    for i=2:n 
      h = 0.5*h;
      sum = 0;
      for k = 1:2:2^(i-1)-1
        sum = sum + feval(f,a+k*h);
      end
      R(i,1) = 0.5*R(i-1,1) + sum*h;
      for j = 2:i 
        R(i,j) = R(i,j-1) + (R(i,j-1)-R(i-1,j-1))/(4^(j-1)-1);
      end
    end