The fastest way to generate Gaussian Orthogonal Ensemble (GOE) matrix and find its eigenvalues

I want to generate a Gaussian Orthogonal Ensemble (GOE) matrix (definition please refer to https://www.lpthe.jussieu.fr/~leticia/TEACHING/Master2019/GOE-cuentas.pdf) and find its largest and the second largest eigenvalues. Is there any faster way?

My code is like below now:

def simulate_f(n):
    main = np.sqrt(2) * np.random.normal(size=(1, n))
    off = np.random.normal(size=(n, n))
    tril = np.tril(off, -1)
    W_n = tril + tril.T
    np.fill_diagonal(W_n, main)

    eigenvalues = LA.eigvals(W_n)
    w = np.flip(np.sort(eigenvalues))
    # GOE_L12_dist[:,i] = w[0:2]
    return w[0:2]