function value = sexpo ( )

%*****************************************************************************80
%
%% SEXPO samples the standard exponential distribution.
%
%  Discussion:
%
%    This procedure corresponds to algorithm SA in the reference.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    19 September 2014
%
%  Author:
%
%    Original FORTRAN77 version by Barry Brown, James Lovato.
%    MATLAB version by John Burkardt.
%
%  Reference:
%
%    Joachim Ahrens, Ulrich Dieter,
%    Computer Methods for Sampling From the
%    Exponential and Normal Distributions,
%    Communications of the ACM,
%    Volume 15, Number 10, October 1972, pages 873-882.
%
%  Parameters:
%
%    Output, real VALUE, a random deviate from the standard
%    exponential distribution.
%
  q = [ ...
       0.6931472, ...
       0.9333737, ...
       0.9888778, ...
       0.9984959, ...
       0.9998293, ...
       0.9999833, ...
       0.9999986, ...
       0.9999999 ];

  a = 0.0;
  u = r4_uni_01 ( );

  while ( true )

    u = u + u;

    if ( 1.0 < u )
      break
    end

    a = a + q(1);

  end

  u = u - 1.0;

  if ( u <= q(1) )
    value = a + u;
    return
  end

  i = 1;
  ustar = r4_uni_01 ( );
  umin = ustar;

  while ( true )

    ustar = r4_uni_01 ( );
    umin = min ( umin, ustar );
    i = i + 1;

    if ( u <= q(i) )
      break;
    end

  end

  value = a + umin * q(1);

  return
end
