function alpert_power_test ( )

%*****************************************************************************80
%
%% alpert_power_test() tests the Alpert rule on the power integrand.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    23 November 2015
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'ALPERT_POWER_TEST\n' );
  fprintf ( 1, '  Test the Alpert rule on the power integrand.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Rule    Order   J   A        N     N+2J               H        Estimate           Error\n' );
  fprintf ( 1, '\n' );

  v2 = integral_power ( );

  num_p = num_power ( );
%
%  For the righthand interval, use the regular rule of the same index.
%
  for rule = 1 : num_p

    a_p = a_power ( rule );
    j_p = j_power ( rule );
    order_p = order_power ( rule );
    [ x_p, w_p ] = rule_power ( rule, j_p );

    a_r = a_regular ( rule );
    j_r = j_regular ( rule );
    order_r = order_regular ( rule );
    [ x_r, w_r ] = rule_regular ( rule, j_r );

    n = 8;

    for nlog = 4 : 6

      n = n * 2;
      h = 1.0 / real ( n + a_p + a_r - 1 );

      x1 = h * x_p; 
      s1 = w_p(1:j_p)' * integrand_power ( j_p, x1 );

      x2 = linspace ( a_p * h, ( a_p + n - 1 ) * h, n );
      s2 = sum ( integrand_power ( n, x2 ) );

      x3 = 1.0 - h * x_r;
      s3 = w_r(1:j_r)' * integrand_power ( j_r, x3 );

      v1 = h * ( s1 + s2 + s3 );

      fprintf ( 1, '    %2d     %4.1f  %2d  %2d  %7d  %7d  %14.6g  %14.6g  %14.6g\n', ...
        rule, order_p, j_p, a_p, n, n + j_p + j_r, h, v1, abs ( v1 - v2 ) );
    end

    fprintf ( 1, '\n' );

  end

  fprintf ( 1, '\n' );
  fprintf ( 1, '                                                    Exact: %14.6g\n', v2 );

  return
end

