83 views (last 30 days)

Hello everyone,

basically the problem is exactly what is stated in the title. I want to create an unitary (2D) matrix of random complex numbers, such that the elements of each column of this matrix sum up to exactly 1. Is there any way to do this? As long as it's still reasonable, computation speed and/or memory usage are not that important.

Thank you to everyone trying to help!

David Goodmanson
on 22 Nov 2020

Edited: David Goodmanson
on 4 Dec 2020

MODIFIED to replace previous random function

Hi Michael,

here is one way. It's based on the idea that if the unitary matrix U is nxn, and onz = [1 1 1 1 1 1... ] (length n), then the sum-of-each-column condition is

[1 1 1 1 1 1... ]*U = [1 1 1 1 1 1... ]

so

n = 5;

onz = ones(1,n);

onzc = onz'; % column vector

na = null(onzc');

% construct an (n-1)x(n-1) unitary matrix by employing random numbers

% uniformly distributed on {-1, 1} x {-i, i}

h = 2*(rand(n-1,n-1) + i*rand(n-1,n-1)) -(1+i);

% method 1

[u, ~] = qr(h);

% method 2 (slower)

% [u, ~] = eig(h+h');

% construct the result

U = na*u*na' + (1/n)*onzc*onzc';

% checks, all of these should be small

U'*U -eye(size(U))

U*U' -eye(size(U))

onz*U - onz

onz*U' - onz % U' works too

David Goodmanson
on 4 Dec 2020

Hi Bruno,

That's a good point about eig, every eigenvector having one real component, And in your example

r = rand(10)+1i*rand(10);

[u,~] = eig(r);

it's the component with largest absolute value. However for eig(r+r') the real component is usualy not the one with largest absolute value. But just the fact that there is always a real component means that eig needs a phase shuffle.

And I agree that when m is complex gaussian, the plateau value for qr is 2. How did I mess that up?

I probably didn't explain the plateau calculation all that well, The idea is that if you let r = abs(u(:)) and sort the r's, then for every r, all the r's with smaller array index are contained inside a circle of radius r. So for every r, it's a calculation of

rho_average = [number of points < r]/(pi*r^2).

If the density is constant, then rho_average does not change. From the plot there is a large region where this is true (after getting past large fluctuations at small r due to small sample size). Eventually you get to wider scattered points, not uniformly distributed in radius, and rho_average drops down.

I don't agree that [ uniform random variable and qr ] followed by a phase shuffle will raise the plateau value from 1.6 to 2, because phase shuffle has no effect on the r values.

Bruno Luong
on 5 Dec 2020

Yes I was wrong on the plateau value of qr since I though the RAND/RANDN do not have the effect, so I run my RandUnitary code and the only thing that add is phase shuffle.

Indeed the drop of plateau value is due to RAND as we both observe now. I knew RAND is not a good input at least for QR.

I'm still stumped that EIG can be somewhat pass few uniformity tests when using with RAND. I admit that I still don't fully understand it.

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!
## 0 Comments

Sign in to comment.