function nint_exactness ( quad_filename, degree_max )

%% MAIN is the main program for NINT_EXACTNESS.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license. 
%
%  Modified:
%
%    03 July 2008
%
%  Author:
%
%    John Burkardt
%
  timestamp;
  fprintf ( 1, '\n' );
  fprintf ( 1, 'NINT_EXACTNESS\n' );
  fprintf ( 1, '  MATLAB version\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Investigate the polynomial exactness of a quadrature\n' );
  fprintf ( 1, '  rule by integrating all monomials of a given degree\n' );
  fprintf ( 1, '  over the [0,1] hypercube.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  The rule will be adjusted to the [0,1] hypercube.\n' );
%
%  Get the quadrature file root name:
%
  if ( 1 <= nargin )
%
%  A commandline argument can be blank.  But to let a user
%  specify the shortcut names, also let "D" indicate them.
%
    if ( quad_filename == 'D' )
      quad_filename = '';
    end

  else

    fprintf ( 1, '\n' );
    fprintf ( 1, 'NINT_EXACTNESS:\n' );

    quad_filename = input ( '  Enter the "root" name of the quadrature files.' );

  end
%
%  Create the names of:
%    the quadrature X file;
%    the quadrature W file;
%    the quadrature R file;
%    the output "exactness" file.
%
  if ( 0 < s_len_trim ( quad_filename ) )
    quad_x_filename = strcat ( quad_filename, '_x.txt' );
    quad_w_filename = strcat ( quad_filename, '_w.txt' );
    quad_r_filename = strcat ( quad_filename, '_r.txt' );
    quad_exact_filename = strcat ( quad_filename, '_exact.txt' );
  else
    quad_x_filename = 'x.txt';
    quad_w_filename = 'w.txt';
    quad_r_filename = 'r.txt';
    quad_exact_filename = 'exact.txt';
  end
%
%  The second command line argument is the maximum degree.
%
  if ( 2 <= nargin )

  else

    fprintf ( 1, '\n' );
    fprintf ( 1, 'NINT_EXACTNESS:\n' );

    degree_max = input ( '  Please enter the maximum total degree to check.' );

  end

  exact_unit = fopen ( quad_exact_filename, 'wt' );

  if ( exact_unit == -1 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'NINT_EXACTNESS - Fatal error!\n' );
    fprintf ( 1, '  Could not open "%s".\n', quad_exact_filename );
    error ( 'NINT_EXACTNESS - Fatal error!' );
  end
%
%  Summarize the input.
%
  string = timestring;
  fprintf ( exact_unit, '%s\n', string );
  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, 'NINT_EXACTNESS\n' );
  fprintf ( exact_unit, '  MATLAB version\n' );
  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, '  Investigate the polynomial exactness of a quadrature\n' );
  fprintf ( exact_unit, '  rule by integrating all monomials of a given degree\n' );
  fprintf ( exact_unit, '  over the [0,1] hypercube.\n' );
  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, '  The rule will be adjusted to the [0,1] hypercube.\n' );
  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, 'NINT_EXACTNESS: User input:\n' );
  fprintf ( exact_unit, '  Quadrature rule X file = "%s".\n', quad_x_filename );
  fprintf ( exact_unit, '  Quadrature rule W file = "%s".\n', quad_w_filename );
  fprintf ( exact_unit, '  Quadrature rule R file = "%s".\n', quad_r_filename );
  fprintf ( exact_unit, '  Maximum total degree to check = %d', degree_max );
%
%  Read the X file.
%
  [ dim_num, point_num ] = dtable_header_read ( quad_x_filename );

  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, '  Spatial dimension = %d\n', dim_num );
  fprintf ( exact_unit, '  Number of points  = %d\n', point_num );

  x = dtable_data_read ( quad_x_filename, dim_num, point_num );
%
%  Read the W file.
%
  [ dim_num2, point_num2 ] = dtable_header_read ( quad_w_filename );

  if ( dim_num2 ~= 1 )
    fprintf ( exact_unit, '\n' );
    fprintf ( exact_unit, 'NINT_EXACTNESS - Fatal error!\n' );
    fprintf ( exact_unit, '  The quadrature weight file should have exactly\n');
    fprintf ( exact_unit, '  one value on each line.\n' );
    error ( 'NINT_EXACTNESS - Fatal error!' );
  end

  if ( point_num2 ~= point_num )
    fprintf ( exact_unit, '\n' );
    fprintf ( exact_unit, 'NINT_EXACTNESS - Fatal error!\n' );
    fprintf ( exact_unit, '  The quadrature weight file should have exactly\n' );
    fprintf ( exact_unit, '  the same number of lines as the abscissa file.\n' );
    error ( 'NINT_EXACTNESS - Fatal error!' );
  end

  weight = dtable_data_read ( quad_w_filename, 1, point_num );
%
%  Read the R file.
%
  [ dim_num2, point_num2 ] = dtable_header_read ( quad_r_filename );

  if ( dim_num2 ~= dim_num )
    fprintf ( exact_unit, '\n' );
    fprintf ( exact_unit, 'NINT_EXACTNESS - Fatal error!\n' );
    fprintf ( exact_unit, '  The quadrature region file should have the same\n' );
    fprintf ( exact_unit, '  number of values on each line as the abscissa file\n' );
    fprintf ( exact_unit, '  does.\n' );
    error ( 'NINT_EXACTNESS - Fatal error!' );
  end

  if ( point_num2 ~= 2 )
    fprintf ( exact_unit, '\n' );
    fprintf ( exact_unit, 'NINT_EXACTNESS - Fatal error!\n' );
    fprintf ( exact_unit, '  The quadrature region file should have two lines.\n' );
    error ( 'NINT_EXACTNESS - Fatal error!' );
  end

  x_range = dtable_data_read ( quad_r_filename, dim_num, 2 );
%
%  Rescale the weights, and translate the abscissas.
%
  volume = abs ( prod ( x_range(1:dim_num,2) - x_range(1:dim_num,1) ) );
  weight(1:point_num) = weight(1:point_num) / volume;

  for dim = 1 : dim_num
    x(dim,1:point_num) = ( x(dim,1:point_num) - x_range(dim,1) ) ...
                       / ( x_range(dim,2) - x_range(dim,1) );
  end
%
%  Explore the monomials.
%
  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, '      Error    Degree  Exponents\n' );
  fprintf ( exact_unit, '\n' );

  for degree = 0 : degree_max

    expon = [];
    more = 0;
    h = 0;
    t = 0;

    while ( 1 )

      [ expon, more, h, t ] = comp_next ( degree, dim_num, expon, more, h, t );

      quad_error = monomial_quadrature ( dim_num, expon, point_num, weight, x );

      fprintf ( exact_unit, '  %24.16f   %2d  ', quad_error, degree );
      for dim = 1 : dim_num
        fprintf ( exact_unit, '%2d', expon(dim) );
      end
      fprintf ( exact_unit, '\n' );

      if ( ~more )
        break
      end

    end

    fprintf ( exact_unit, '\n' );

  end

  fprintf ( exact_unit, '\n' );
  fprintf ( exact_unit, 'NINT_EXACTNESS:\n' );
  fprintf ( exact_unit, '  Normal end of execution.\n' );
  fprintf ( exact_unit, '\n' );
  string = timestring;
  fprintf ( exact_unit, '%s\n', string );

  fclose ( exact_unit );

  fprintf ( 1, '\n' );
  fprintf ( 1, 'NINT_EXACTNESS:\n' );
  fprintf ( 1, '  The exactness results were written to "%s".\n', ...
    quad_exact_filename );

  fprintf ( 1, '\n' );
  fprintf ( 1, 'NINT_EXACTNESS:\n' );
  fprintf ( 1, '  Normal end of execution.\n' );
  fprintf ( 1, '\n' );
  timestamp;