function clenshaw_curtis_test03 ( )

%% TEST03 calls CC_GRIDS_CONSTRAINED to collect constrained grids.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    13 March 2007
%
%  Author:
%
%    John Burkardt
%
  dim_num = 2;

  alpha = [ 2.0, 3.0 ];
  q_max = 13.0;

  fprintf ( 1, '\n' );
  fprintf ( 1, 'TEST03:\n' );
  fprintf ( 1, '  CC_GRIDS_CONSTRAINED returns all Clenshaw Curtis grids\n' );
  fprintf ( 1, '  satisfying a set of constraints.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  ORDER(I), the order of the 1D rule in dimension I,\n' );
  fprintf ( 1, '  is constrained by\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '    ORDER_MIN(I) <= ORDER(I) <= ORDER_MAX(I)\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  We also define the total weighted order Q\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '    Q = ALPHA(1) * ORDER(1) + ... + ALPHA(N) * ORDER(N)\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  and further constrain our grids to satisfy\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '    Q <= Q_MAX = %f\n', q_max );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Spatial dimension of grids = %d\n', dim_num );

  order_min(1:dim_num) = 1;
  order_max(1:dim_num) = 5;

  fprintf ( 1, '\n' );
  fprintf ( 1, ' Dimension Order_min Order_max     Alpha\n' );
  fprintf ( 1, '\n' );

  for dim = 1 : dim_num
    fprintf ( 1, '  %8d  %8d  %8d  %12f\n', ...
      dim, order_min(dim), order_max(dim), alpha(dim) );
  end

  [ grid_num, point_num ] = cc_grids_constrained_size ( dim_num, q_max, ...
    alpha, order_min, order_max );

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Number of grids = %d\n', grid_num );
  fprintf ( 1, '  Number of points in the grids = %d\n', point_num );

  [ grid_order, grid_point ] = cc_grids_constrained ( dim_num, q_max, ...
    alpha, order_min, order_max, grid_num, point_num );
%
%  Now we're done.  Print the merged grid data.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, '         I         Q          Grid orders:\n' );
  fprintf ( 1, '\n' );
  for j = 1 : grid_num

    q = alpha(1:dim_num) * grid_order(1:dim_num,j);

    fprintf ( 1, '  %8d  %14f', j, q );
    for dim = 1 : dim_num
      fprintf ( 1, '  %12d', grid_order(dim,j) );
    end
    fprintf ( 1, '\n' );
  end

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Grid points:\n' );
  fprintf ( 1, '\n' );
  for j = 1 : point_num
    fprintf ( 1, '  %8d', j );
    for dim = 1 : dim_num
      fprintf ( 1, '  %12f', grid_point(dim,j) );
    end
    fprintf ( 1, '\n' );
  end

  return
end