a) og b) og c)
Code Block | ||
---|---|---|
| ||
function result = simpsons(a,b,Nn,f) N = N+isEven(N)-1; h = calcDelta(a,b,N); result = h/3*f(h*0); mult = [4,2]; for i = 1:N-1 result = result + h/3*f(h*i)*mult(mod(i,2)+1); end result = result + h/3*f(h*N); end function even = isEven(n) even = mod(n,2)+1; end function h = calcDelta(a,b,n) h = (b-a)/n; end |
c)
Code Block |
---|
1.9832 = simpsons(0,pi,14,@(x)sin(x)) |
d)
if mod(n, 2) ~= 0
error('Parameteren n må være et partall, var %d\n', n);
end
dx = (b-a)/n;
y= f(a:dx:b);
yodd = 2*sum(y(3:2:n-1));
yeven = 4*sum(y(2:2:n));
ysum = y(1) + yodd + yeven + y(n+1);
result = dx/3 * ysum;
end |
b)
Code Block | ||
---|---|---|
| ||
function result = simpsons_error(start, stop, error, fn)
N = 2; % N må være et partall
Si = simpsons(start, stop, N, fn);
S2i = simpsons(start, stop, 2*N, fn);
while abs(Si-S2i) >= error
N = 2 * N;
Si = S2i;
S2i = simpsons(start, stop, 2*N, fn);
end
fprintf('Antal ledd: %d\n', N);
result = Si | ||
Code Block | ||
function err = deviation(a,b,N,f,corr)
err= corr-simpsons(a,b,N,f);
end |