function pisim ( trial_num )

%*****************************************************************************80
%
%% pisim() estimates PI by sampling points in a square.
%
%  Discussion:
%
%    Sample N points (x,y) in a square.
%    The number of points M with x^2+y^2 <= 1 lie inside the quarter circle.
%    M/N is an estimate of pi/4. 
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    29 January 2018
%
%  Author:
%
%    Paul Nahin
%
%  Reference:
%
%    Paul Nahin,
%    Dueling Idiots and Other Probability Puzzlers,
%    Princeton, 2012,
%    ISBN: 978-0691155005.
%
%  Input:
%
%    integer trial_num: the number of trials.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'pisim():\n' );
  fprintf ( 1, '  Sample N points (x,y) in a square.\n' );
  fprintf ( 1, '  The number of points M with x^2+y^2 <= 1 lie inside the quarter circle.\n' );
  fprintf ( 1, '  M/N is an estimate of pi/4. \n' );

  darts = 0;
 
  for i = 1 : trial_num
    x = rand ( );
    y = rand ( );
    r = x * x + y * y;
    if ( r < 1 )
      darts = darts + 1;
    end
  end

  pi_estimate = 4 * darts / trial_num;

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Estimate for PI is %g\n', pi_estimate );

  return
end
