Member
offline
3 posts

Hi community,

I'm working on incorporating a model of a scattering medium into a simulation and I need to check, whether my model is correct. It is assumed that the conditions for random matrix theory to be a valid description are met. I apply SVD to compare to the results of Popoff et al. (doi.org/10.1088/1367-2630/13/12/123021). For SVD implementation I relied on the paper by Edelman&Wang (web.eecs.umich.edu/~rajnrao/Acta05rmt.pdf). My problem lies in the normaliziation of the singular values.

I assume normal distribution N(0,1) for both real and imaginary part of the transmission matrix. The matrix is normalized to the biggest absolute value, as I assume a passive transmission medium with no amplification possible.
Now, the quarter circle law distribution for the singular values is to be expected for this case of gamma=1. The probability distribution follows this curve, but I obtain an incorrect range for the singular values (quarter circle law predicts [0,2]).
After the scaling of K, the variance its elements will not be 1 anymore. So, I found the previously applied scaling factor, divided by sqrt(2), to be a working conversion to the [0,2] interval predicted by the quarter circle law.

However, the physical interpretation is not clear to me. By visualizing the singular value I want to verify, that I use a description previously shown to be valid for certain scattering media.

Would someone know, how to correctly scale the singular values in order to show, that the example follows the quarter circle law?
I'd be grateful for any hint.
Cheers, Matthias

Here's the numeric Python code I used:

from numpy import sqrt, angle, amax, sum
from numpy.random import randn
from scipy.linalg import svd
import matplotlib.pyplot as plt
m = 1024
gamma = 1
n = gamma * m

K = randn(m,n)+ 1j*randn(m,n) # transmission matrix
normFactor = amax(amax(abs(K))) # biggest amplitude element in K
K = K/normFactor # normalization
U,s,V = svd(K)

# different normalizations
sNorm1 = s/sqrt(m) # as in Edelman & Wang
sNorm2 = s/sqrt(sum(s**2)) # as in Popoff et al.
sNorm3 = normFactor/sqrt(2)*s/(sqrt(m)) # results in [0,2] interval (as expected)

#%% plotting
plt.figure()
plt.hist(sNorm1, bins = 50, normed = True)
plt.title('histogram of singular values, normalization 1/sqrt(m)')
plt.xlabel('normalized singular value')
plt.ylabel('probability density')

plt.figure()
plt.hist(sNorm2, bins = 50, normed = True)
plt.title('histogram of singular values, normalization 1/sqrt(sum(s**2))')
plt.xlabel('normalized singular value')
plt.ylabel('probability density')

plt.figure()
plt.hist(sNorm3, bins = 50, normed = True)
plt.title('histogram of singular values, normalization normFactor/sqrt(2)*/(sqrt(m))')
plt.xlabel('normalized singular value')
plt.ylabel('probability density')

plt.show()