function U = random_orthogonal_matrix(integerN)
  %
  % Generate a random orthogonal matrix.
  %
  % INPUT:
  %
  %   - ``integerN`` -- The dimension of the resulting square matrix.
  %
  % OUTPUT:
  %
  %   - ``U`` -- An orthogonal matrix of dimension integerN.
  %
  % REFERENCES:
  %
  %   1. G.W. Stewart, Efficient Generation of Random Orthogonal Matrices,
  %      SIAM J. Numer. Analysis, 1980, pp. 403--409, Section 3.
  %

  % We begin by computing a random matrix A.
  A = rand(integerN);

  % The Q-R decomposition of A will give us an orthogonal factor of A.
  % See the reference for why this doesn't suck.
  [U, R] = qr(A);
end
