function v = quadint_p(f, a, b, n)
  v = 0;
  h = (b - a) / (n - 1);
  parfor i = 1:n
    x  = ((n - i)*a + (i - 1)*b) / (n - 1);
    fx = f(x);
    % normally, the first and last point would be 
    % outside of the loop...
    if (i == 1) || (i == n)
        v = v + 0.5*h*fx;
    else
        v = v + h*fx;
    end
  end
return
