The simplest method is based on the QR decomposition of a random complex matrix with independent elements following a Normal distribution.

In Mathematica one can generate a random complex number with Normal distribution for both the real and imaginary parts with the following function

NormalRandomComplex[] := (#[[1]] + #[[2]]*I) &@

RandomReal[MultinormalDistribution[{0, 0}, {{1, 0}, {0, 1}}]]

where {0, 0} specifies the origin and {{1, 0}, {0, 1}} the covariance matrix (standard deviation 1)

The function that generates the random n x n unitary matrix is

RandomUnitaryMatrix[n_] := Module[{z, q, r},

z = Table[NormalRandomComplex[], {n}, {n}];

{q, r} = QRDecomposition[z];

q.DiagonalMatrix[Sign[Tr[r, List]]]

] z = Table[NormalRandomComplex[], {n}, {n}];

{q, r} = QRDecomposition[z];

q.DiagonalMatrix[Sign[Tr[r, List]]]

In this code, the QR decomposition of the complex random matrix with normal distribution outputs a unitary matrix "q" and a triangular matrix "r". The sought random unitary matrix is a correction of "q", which must be multiplied by a diagonal matrix filled with 1 and -1 according to the sign of the corresponding diagonal element of "r".

## No comments:

## Post a Comment