function dot_forplus_graph ( log_max )

%*****************************************************************************80
%
%% DOT_FORPLUS_GRAPH generates and plots timings for the dot product operation.
%
%  Discussion:
%
%    A FOR loop is used for the calculation.
%
%    The FOR loop is "unrolled" to a depth of 8.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    17 December 2011
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer LOG_MAX, the log base 2 of the maximum value of N
%    to try.  A value of 25 is roughly as big as you can go without
%    hitting memory restrictions.
%
  if ( nargin < 1 )
    log_max = 24;
  end

  speed_limit = 2800000000;

  fprintf ( 1, '\n' );
  fprintf ( 1, 'DOT_FORPLUS_GRAPH\n' );
  fprintf ( 1, '  Compute the dot product of two N-vectors,\n' );
  fprintf ( 1, '  using a FOR loop.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Do this for vectors of size 2^0, 2^1, ..., 2^%d\n', log_max );

  n = 1;
  k = 0;

  for logn = 0 : log_max

    x = rand ( n, 1 );
    y = rand ( n, 1 );

    if ( n < 8 )
      tic;
      z = 0.0;
      for i = 1 : n
        z = z + x(i) * y(i);
      end
      t = toc;
    else
      tic;
      z = 0.0;
      for i = 1 : 8 : n
        z = z + x(i)   * y(i);
              + x(i+1) * y(i+1);
              + x(i+2) * y(i+2);
              + x(i+3) * y(i+3);
              + x(i+4) * y(i+4);
              + x(i+5) * y(i+5);
              + x(i+6) * y(i+6);
              + x(i+7) * y(i+7);
      end
      t = toc;
    end

    k = k + 1;
    xplot(k) = 4*n+2;
    yplot(k) = t;
    n = n * 2;

  end

  plot ( r8_log_2 ( xplot ), r8_log_2 ( yplot ), 'r-*', 'LineWidth', 2, 'MarkerSize', 10 );
  grid on
  xlabel ( 'Log2(Number of vector entries)' );
  ylabel ( 'Log2(TicToc Time in Seconds)' );
  title ( 'Dot product time versus vector size' );
  axis equal
  print -dpng 'dot_forplus_time_graph.png'

  clf
  hold on
  plot ( r8_log_2 ( xplot ), r8_log_2 ( xplot ./ yplot ), 'b-*', 'LineWidth', 2, 'MarkerSize', 10 );
  zplot(1:log_max+1) = r8_log_2 ( speed_limit );
  plot ( r8_log_2 ( xplot ), zplot, 'r-', 'LineWidth', 2 );
  hold off
  grid on
  xlabel ( 'Log2(Number of vector entries)' );
  ylabel ( 'Log2(Computational Rate)' );
  title ( 'Computational rate versus vector size' );
  axis equal
  print -dpng 'dot_forplus_rate_graph.png'

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Log2(N)  Log2(W)  Log2(T)  Log2(N/T)            N    T                N/T\n' );
  fprintf ( 1, '\n' );

  for logn = 0 : log_max
    rate = xplot(logn+1) / yplot(logn+1);
    fprintf ( 1, '  %6d  %8.2f  %8.2f  %8.2f  %12d  %f  %14.2g\n', ...
      logn, r8_log_2 ( xplot(logn+1) ), r8_log_2 ( yplot(logn+1) ), r8_log_2 ( rate ), ...
      xplot(logn+1), yplot(logn+1), rate );
  end

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Maximum possible rate W/T assumed to be                       %14.2g\n', speed_limit );

  return
end
function value = r8_log_2 ( x )

%*****************************************************************************80
%
%% R8_LOG_2 returns the logarithm base 2 of |X|.
%
%  Discussion:
%
%    R8_LOG_2 ( X ) = Log ( |X| ) / Log ( 2.0 )
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license. 
%
%  Modified:
%
%    21 March 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real X, the number whose base 2 logarithm is desired.
%    X should not be 0.
%
%    Output, real VALUE, the logarithm base 2 of the absolute
%    value of X.  It should be true that |X| = 2**R8_LOG_2.
%
  if ( x == 0.0 )
    value = -inf;
  else
    value = log ( abs ( x ) ) / log ( 2.0 );
  end

  return
end

