function [sum] = trap(alpha,n)

%  Computes an approximation to the integral from 0 to 1 of
%  cos(alpha x^2) dx  using the trapezoidal rule with n points.

dx = 1/n;

sum = 0;
for i=0:n,              % Loop over nodes.
  x = i*dx;
  cosaxsq = cos(alpha*x^2);  % Evaluate integrand at x.
  if i==0 | i==n, 
    sum = sum + .5*cosaxsq;  %  Add .5 * Value of integrand at endpoints
  else
    sum = sum + cosaxsq;     %  Add 1 * Value of integrand at interior points.
  end;
end;

sum = sum*dx;        %  Multiply final result by dx.