function value = beta ( x, y )

%% BETA returns the value of the Beta function.
%
%  Formula:
%
%    BETA(X,Y) = ( GAMMA(X) * GAMMA(Y) ) / GAMMA(X+Y)
%
%  Restrictions:
%
%    Both X and Y must be greater than 0.
%
%  Properties:
%
%    BETA(X,Y) = BETA(Y,X).
%    BETA(X,Y) = Integral ( 0 <= T <= 1 ) T**(X-1) (1-T)**(Y-1) dT.
%    BETA(X,Y) = GAMMA(X) * GAMMA(Y) / GAMMA(X+Y)
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    22 July 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real X, Y, the two parameters that define the Beta function.
%    X and Y must be greater than 0.
%
%    Output, real VALUE, the value of the Beta function.
%
  if ( x <= 0.0 | y <= 0.0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'BETA - Fatal error!\n' );
    fprintf ( 1, '  Both X and Y must be greater than 0.\n' );
    error ( 'BETA - Fatal error!' );
  end

  value = exp ( gamma_log ( x ) + gamma_log ( y ) - gamma_log ( x + y ) );

  return
end