Fast top-n from samples of many different normal distributions Let's say I have 1 million normal distributions, each with a different mean and stdev. I want to sample from each distribution and take the top 10 samples (for a Thompson sampling application.)
Is there a faster approximate way to sample than generating all million samples and taking the top n?
This process gets repeated many times, so I am happy to precompute something if it helps with the amortized cost.

As whuber points out in his answer, for non-IID variables it is not possible to use the standard results for order statistics for IID values to directly generate the order statistics of interest, without generating the underlying sample.  That makes this type of problem quite tricky, and it is likely that full simulation and sorting is an appropriate method.
The only thing I can think of here that might give a reasonable approximation is to  see if you can make a preliminary identification of a subset of distributions where the mean and variance are low enough that the top value from these distributions is highly unlikely to be in the top-10 values overall.  If you could successfully make a preliminary identification of this kind, it would allow you to remove those random variables from consideration while still giving an approximately correct method.  Of course, you would need to be careful here --- this would only work if identifying a set of random variables with this property is less computationally expensive than just including them in the sampling.

How might this work: As an example, suppose that you have $N$ random variables in your analysis but you can identify a subset of $K$ random variables (call them $X_1,...,X_K$) with maximum mean and variance bounded as follows:
$$\begin{align}
\max \{ \mu_1,..., \mu_K \} \leqslant \mu_0, \\[6pt]
\max \{ \sigma_1,..., \sigma_K \} \leqslant \sigma_0, \\[6pt]
\end{align}$$
for some low values of $\mu_0$ and $\sigma_0$.  Then for any $r > \mu_0$ you have:
$$\begin{align}
\mathbb{P}(\max X_1,...,X_K \geqslant r)
&= 1 - \mathbb{P}(\max X_1,...,X_K < r) \\[12pt]
&= 1 - \prod_{i=1}^K \mathbb{P}(X_i < r) \\[6pt]
&= 1 - \prod_{i=1}^K \Phi \bigg( \frac{r-\mu_i}{\sigma_i} \bigg) \\[6pt]
&\leqslant 1 - \Phi \bigg( \frac{r-\mu_0}{\sigma_0} \bigg)^K. \\[6pt]
\end{align}$$
Similarly, suppose you can identify a separate set of $R$ random variables (call them $Y_1,...,Y_R$) with minimum mean and maximum variance bounded as follows:
$$\begin{align}
\min \{ \mu_{1*},..., \mu_{R*} \} \geqslant \mu_*, \\[6pt]
\max \{ \sigma_{1*},..., \sigma_{R*} \} \leqslant \sigma_*, \\[6pt]
\end{align}$$
for some high value of $\mu_*$ and low value of $\sigma_*$.  Then for any $r < \mu_*$ you have:
$$\begin{align}
\mathbb{P}(\min Y_1,...,Y_R \leqslant r)
&= 1 - \mathbb{P}(\min Y_1,...,Y_R > r) \\[12pt]
&= 1 - \prod_{i=1}^R \mathbb{P}(Y_i > r) \\[6pt]
&= 1 - \prod_{i=1}^R \bigg( 1- \Phi \bigg( \frac{r-\mu_{i*}}{\sigma_{i*}} \bigg) \bigg) \\[6pt]
&\leqslant 1 - \bigg( 1- \Phi \bigg( \frac{r-\mu_*}{\sigma_*} \bigg)^R \bigg). \\[6pt]
\end{align}$$
Putting these results together, and choosing any $\mu_0 < r < \mu_*$ you get:
$$\begin{align}
\mathbb{P}(\max X_1,...,X_K \geqslant \min Y_1,...,Y_R) 
&\leqslant \mathbb{P}(\max X_1,...,X_K \geqslant r \geqslant \min Y_1,...,Y_R) \\[12pt]
&= \mathbb{P}(\max X_1,...,X_K \geqslant r) \cdot \mathbb{P}(\min Y_1,...,Y_R \leqslant r)) \\[12pt]
&\leqslant \Bigg[ 1 - \Phi \bigg( \frac{r-\mu_0}{\sigma_0} \bigg)^K \Bigg] \Bigg[ 1 - \bigg( 1- \Phi \bigg( \frac{r-\mu_*}{\sigma_*} \bigg)^R \bigg) \Bigg]. \\[6pt]
\end{align}$$
Now, suppose that you can identify subsets with these requirements for a value of $\mu_0 \ll r \ll \mu_*$ (relative to the values $\sigma_0$ and $\sigma_*$) such that:
$$\Phi \bigg( \frac{r-\mu_0}{\sigma_0} \bigg)^K \approx 1
\quad \quad \quad \quad \quad 
\Phi \bigg( \frac{r-\mu_*}{\sigma_*} \bigg)^R \approx 0.$$
This would be enough to guarantee that:
$$\begin{align}
\mathbb{P}(\max X_1,...,X_K \geqslant \min Y_1,...,Y_R) 
\approx 0,
\end{align}$$
which means that we can be highly confident that the maximum value of the first subset of $K$ values will not be as large as the minimum of the second subset of $R$ values.  If we can make this identification for $R = 10$ in this case, this means that we can be highly confident that none of the values in the first subset will be in the top-10 values of the sample overall.  This would then allow us to remove consideration of the first subsample from analysis and conduct the standard simulation and sorting for the remaining sample of $N-K$ values.
