It is well known that Gaussian elimination without pivoting is numerically unstable, and in practice Gaussian elimination is done with row pivoting (partial pivoting). A theorem of Wilkinson states that if $\tilde{x}$ is the computed solution to $Ax = b$ using Gaussian elimination with partial pivoting, with $A \in \mathbb{C}^{n \times n}$, then $(A + \Delta A) \tilde{x} = b$ with $\Delta A$ satisfying the norm inequality
$$ \| \Delta A\|_{\infty} \le n^2 \frac{ 3n \epsilon}{1 - 3n \epsilon} \rho_n \|A\|_{\infty}, $$
where $\epsilon$ is the floating point error and $\rho_n$ is the *growth factor* given by
$$ \rho_n = \max_{i,j,k} |A_{ij}^{(k)}| / \max_{ij} |A_{ij}|$$
where $A^{(k)}$ is the manipulated matrix after $k$ steps of Gaussian elimination. The growth factor achieves its upper bound of $\rho_n = 2^{n-1}$, but it is also observed in practice that no typically matrix comes even close to this upper bound in terms of instability, and that in fact it seems that almost every matrix has a small growth factor. This leads to the following unresolved problem:

What is the probability distribution of $\rho_n$ over the space of $n \times n$ random matrices?

We can take the space of random matrices to be those matrices with entries drawn from i.i.d. normal distributions, typically, though one can also pose the question when entries are drawn from other distributions.

This is all pretty well-known as a mystery in numerical linear algebra, but my question is mainly this:

What are the primary hurdles that random matrix theory has yet to overcome in the analysis of the stability of Gaussian elimination with partial pivoting?