This is the set of Hermitian \( N \times N \) matrices \( H \) with probability density \begin{equation*} P(H) dH = \mathcal{N} \prod_{ i < j} d \mathrm{Re} H_{ij} \prod_{ i < j} d \mathrm{Im} H_{ij}\ \prod_{ i} d H_{ii}\exp \left( -\frac{1}{2}\ \mathrm{tr}(H^2)\right) \end{equation*} where \( \mathcal{N} \) is the normalization constant. The spectral density \( \rho (\lambda) \) is defined by \begin{equation*} \rho(x) dx = \sum_{k=1}^N P(x \le \lambda_k \le x + dx) \end{equation*} with \( \lambda_k \) the \(N \) eigenvalues of the matrix \( H \). In chapter 6 of his book, Mehta calculates that \( \rho(x) = K_N(x,x) \) with \begin{equation*} K_N(x,y) = \frac{H_N(x) H_{N-1}(y) - H_N(y) H_{N-1}(x)}{x-y}e^{-\displaystyle\frac{x^2}{4} -\displaystyle\frac{y^2}{4}} \end{equation*} Thus \begin{equation}\label{eq:20151114c} \rho(x) = \left( H_{N-1}(x) H_N'(x) - H_{N}(x) H_{N-1}'(x) \right) e^{-\displaystyle\frac{x^2}{2}} \end{equation} Here \( H_N(x) \) are the Hermite polynomials. I now check formula \eqref{eq:20151114c} for \( N = 1, 2, 3 \) with a Monte Carlo simulation in Mathematica. I generate matrices from the GUE as follows: I draw matrices \( A \) and \( B \) whose elements \( A_{ij} \) and \( B_{ij} \) are independent normal random variables and take \begin{equation*} H = \frac{1}{2} \left( A + A^T \right) + \frac{i}{2} \left( B -B^T \right) \end{equation*} Then I calculate the \( N \) eigenvalues of \( H \). This is the code that I use for this.

```
Clear[randomsamples]
(* This code calculates the eigenvalues of NSamples matrices from GUE.
The size of the matrices is NN times NN *)
```

```
randomsamples[NN_, NSamples_] :=
Table[A = RandomVariate[NormalDistribution[], {NN, NN}];
```

`B = RandomVariate[NormalDistribution[], {NN, NN}]; `

`H = 1/2 (A + Transpose[A]) + I / 2 (B - Transpose[B]); `

```
Eigenvalues[H], {k, NSamples}] // Flatten
```

__N = 1__

In this case \begin{equation*} \rho(x) = \frac{1}{\sqrt{2 \pi }} \exp \left(-\frac{x^2}{2}\right) \end{equation*} I produce a graph with the following commands

` `

`NSamples = 10000; `

`NN = 1; `

`randomeigs = randomsamples[NN, NSamples]; `

`fig1 = Histogram[randomeigs, {-3, 3, 0.1}, "PDF"]; `

`fig2 = Plot[1/Sqrt[2 Pi] Exp[-x^2/2], {x, -3, 3}, PlotStyle -> Red]; `

```
Show[fig1, fig2]
```

__N = 2__

In this case \begin{equation*} \rho(x) = \frac{e^{-\frac{x^2}{2}} \left(x^2+1\right)}{2 \sqrt{2 \pi }} \end{equation*} I produce a graph with the following commands

```
NSamples = 50000;
NN = 2;
```

`randomeigs = randomsamples[NN, NSamples]; `

`fig1 = Histogram[randomeigs, {-4, 4, 0.1}, "PDF"]; `

```
fig2 = Plot[ (E^(-(x^2/2)) (1 + x^2))/(2 Sqrt[2 \[Pi]]), {x, -4, 4},
PlotStyle -> Red];
```

```
Show[fig1, fig2]
```

__N = 3__

In this case \begin{equation*} \rho(x) = \frac{e^{-\frac{x^2}{2}} \left(x^4+3\right)}{6 \sqrt{2 \pi }} \end{equation*} I produce a graph with the following commands

```
NSamples = 50000;
NN = 3;
```

`randomeigs = randomsamples[NN, NSamples]; `

`fig1 = Histogram[randomeigs, {-5, 5, 0.1}, "PDF"]; `

```
fig2 = Plot[ (E^(-(x^2/2)) (3 + x^4))/(6 Sqrt[2 \[Pi]]), {x, -5, 5},
PlotStyle -> Red];
Show[fig1, fig2]
```

Warning: In this post, I have not paid attention to overall constants in expressions, I only made sure the normalizations were consistent in the Mathematica code.

## No comments:

## Post a Comment