function p = ellipse_uniform ( n )

%% ELLIPSE_UNIFORM returns sample points from a particular ellipse.
%
%  Discussion:
%
%    This routine returns N points sampled uniformly at random
%    from within the ellipse defined by
%
%      x' * A * x <= R*R
%
%    with A = ( 9  6 ) and R = 6.
%            ( 6 20 )
%
%    The routine requires the Cholesky factor of A.  For this
%    problem, the formulas are worked out explicitly.  The approach
%    generalizes to other ellipses, and with a little more work,
%    to higher dimensions as well.
%
%  Modified:
%
%    19 May 2006
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of points to generate.
%
%    Output, real P(2,N), the sample points.
%

%
%  Step 1: Generate N points in the unit circle.
%  (Generate N normally distributed points.
%   Normalize them (now they lie on the circumference).
%   Assign a radius.)
%
  p(1:2,1:n) = randn(2,n);

  p_norm(1:n) = sqrt ( dot ( p, p, 1 ) );

  r(1:n) = rand(1,n);

  for i = 1 : 2
    p(i,1:n) = sqrt ( r(1:n) ) .* p(i,1:n) ./ p_norm(1:n);
  end
%
%  Step 2.  The ellipse is x' * A * x <= R*R.
%  A has the Cholesky factorization A = U' * U.
%  Take each point from the unit ball, and:
%  *) multiply by R;
%  *) multiply by inverse(U).
%
%  For our problem, we have A = ( 9  6 )
%                               ( 6, 20 )
%
%  and R = 6
%
%  Hence
%
%  U = ( 3 2 )
%      ( 0 4 )
%
%  inverse ( U ) = ( 1/3 -1/6 )
%                  (  0   1/4 )
%
%  
%  so we can do this in one step.
%
  p(1,1:n) = 2.0 * p(1,1:n) - p(2,1:n);
  p(2,1:n) = 1.5 * p(2,1:n);