erftrap

# erftrap - % erftrap.m efunc = inline( '2/sqrt(pi) * exp( -...

This preview shows pages 1–3. Sign up to view the full content.

% erftrap.m efunc = inline( '2/sqrt(pi) * exp( - x.^2 )' ); a = 0; b = 1; true = erf(1); % compute composite trapezoid and error for multiple n powers = 2:8; N = 2.^powers + 1; I = []; err = []; for n = N h = (b-a)/(n-1); x = a:h:b; f = efunc( x ); trap = h * ( f(1)/2 + sum(f(2:(n-1))) + f(n)/2 ); abserr = abs( trap - true ); I = [I, trap]; err = [err, abserr]; end figure(1); clf; subplot(2,1,1); plot( N, I, 'o-', 'LineWidth', 1.5 ); xlabel( 'n' ); title( 'Integral: erf(x) using trapezoid rule' ); subplot(2,1,2); semilogy( N, err, 'o-', 'LineWidth', 1.5 ); xlabel( 'n' ); title( 'Absolute error' );

This preview has intentionally blurred sections. Sign up to view the full version.

View Full Document
% erftrap.m efunc = inline( '2/sqrt(pi) * exp( - x.^2 )' ); a = 0; b = 1; true = erf(1); % compute composite trapezoid and error for multiple n powers = 2:8; N = 2.^powers + 1; I = []; I2 =[]; err = []; err2 = []; for n = N h = (b-a)/(n-1); x = a:h:b; f = efunc( x ); trap = h * ( f(1)/2 + sum(f(2:(n-1))) + f(n)/2 ); simp = h/3 * ( (f(1) + 4*sum(f(2:2:(n-1))) + 2*sum(f(3:2:(n-2))) + f(n) )); abserr = abs( trap - true );
This is the end of the preview. Sign up to access the rest of the document.

## This note was uploaded on 05/01/2009 for the course PSTAT 120A taught by Professor Mackgalloway during the Spring '08 term at UCSB.

### Page1 / 3

erftrap - % erftrap.m efunc = inline( '2/sqrt(pi) * exp( -...

This preview shows document pages 1 - 3. Sign up to view the full document.

View Full Document
Ask a homework question - tutors are online