function cuda_loop_test ( )

%*****************************************************************************80
%
%% cuda_loop_test() tests cuda_loop().
%
%  Discussion:
%
%    A CUDA kernel "kernel()" is invoked by a command of the form
%      
%      kernel << blocks, threads >> ( args )
%
%    where blocks and threads are each vectors of up to 3 values,
%    listing the number of blocks and number of threads to be used.
%
%    If a problem involves N tasks, then tasks are allotted to 
%    specific CUDA processes in an organized fashion.  Some processes
%    may get no tasks, one task, or multiple tasks.  
%
%    Each process is given variables that can be used to determine
%    the tasks to be performed:
%
%      gridDim.x, gridDim.y, gridDim.z: the block dimensions as
%      given by the user in "blocks"
%
%      blockDim.x, blockDim.y, blockDim.z: the thread dimensions as
%      given by the user in "threads"
%
%      blockIdx.x, blockIdx.y, blockId.z: the block indices for this process.
%
%      threadIdx.x, threadIdx.y, threadIdx.z: the thread indices for this process.
%
%    Essentially, a process can determine its linear index K by:
%
%      K = threadIdx.x
%        +  blockdim.x  * threadIdx.y
%        +  blockDim.x  *  blockDim.y  * threadIdx.z
%        +  blockDim.x  *  blockDim.y  *  blockDim.z  * blockIdx.x
%        +  blockDim.x  *  blockDim.y  *  blockDim.z  *  gridDim.x  * blockIdx.y
%        +  blockDim.x  *  blockDim.y  *  blockDim.z  *  gridDim.x  *  gridDim.y  * blockIdx.z
%
%    Set task T = K.
%
%    while ( T < N )
%      carry out task T
%      T = T + blockDim.x * blockDim.y * blockDim.z * gridDim.x * gridDim.y * gridDim.z.
%
%    This program suggests how a specific set of block and thread parameters 
%    would determine the assignment of individual tasks to CUDA processes.
%
%  Licensing:
%
%    This code is distributed under the MIT license. 
%
%  Modified:
%
%    13 December 2018
%
%  Author:
%
%    John Burkardt
%
%  Local:
%
%    integer BLOCKS(3), the CUDA block values.  These 
%    should be nonnegative.  Typically, the third entry is 1.  Generally, the 
%    first two values cannot be greater than 35,535.
%
%    integer THREADS(3), the CUDA thread values.  These should be 
%    nonnegative.  Typically, there is a maximum value imposed on these 
%    quantities, which depends on the GPU model.
%
%    integer N, the number of tasks to be carried out.
%
  addpath ( '../cuda_loop' );

  timestamp ( );
  fprintf ( 1, '\n' );
  fprintf ( 1, 'CUDA_LOOP_TEST():\n' );
  fprintf ( 1, '  MATLAB/Octave version %s.\n', version ( ) );
  fprintf ( 1, '  Test cuda_loop()\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Simulate the way CUDA breaks into iterative task, using\n' );
  fprintf ( 1, '  blocks and threads.\n' );
%
%  Linear array of blocks and threads.
%  Essentially, blocks = your hands and threads = your fingers.
%  Now count up to 23..
%
  blocks(1) = 2;
  blocks(2) = 1;
  blocks(3) = 1;
  threads(1) = 5;
  threads(2) = 1;
  threads(3) = 1;
  n = 23;
 
  cuda_loop ( blocks, threads, n );
%
%  Unit arrays of blocks and threads.
%  Waste your GPU by having a single block and thread do everything.
%
  blocks(1) = 1;
  blocks(2) = 1;
  blocks(3) = 1;
  threads(1) = 1;
  threads(2) = 1;
  threads(3) = 1;
  n = 23;

  cuda_loop ( blocks, threads, n );
%
%  2D block array, 3D thread array.
%  More processes than tasks.
%
  blocks(1) = 2;
  blocks(2) = 3;
  blocks(3) = 1;
  threads(1) = 2;
  threads(2) = 1;
  threads(3) = 4;
  n = 40;
 
  cuda_loop ( blocks, threads, n );
%
%  One block, 8 threads.
% 
  blocks(1) = 1;
  blocks(2) = 1;
  blocks(3) = 1;
  threads(1) = 2;
  threads(2) = 2;
  threads(3) = 2;
  n = 23;
 
  cuda_loop ( blocks, threads, n );
%
%  Terminate.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'CUDA_LOOP_TEST():\n' );
  fprintf ( 1, '  Normal end of execution.\n' );
  fprintf ( 1, '\n' );
  timestamp ( );

  rmpath ( '../cuda_loop' );

  return
end
function timestamp ( )

%*****************************************************************************80
%
%% timestamp() prints the current YMDHMS date as a timestamp.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    14 February 2003
%
%  Author:
%
%    John Burkardt
%
  t = now;
  c = datevec ( t );
  s = datestr ( c, 0 );
  fprintf ( 1, '%s\n', s );

  return
end

