Order Statistics-Expected Value of Random Length Let $Y_1<Y_2 $ denote the order statistics of a random sample of size 2 from a distribution that is $N\left( \mu,\sigma^2 \right) $, where $\sigma^2$ is known. Compute the expected value of the random length $Y_2-Y_1$.
I can see that the answer is $\frac{2\sigma}{\sqrt{\pi}}$ but I do not know how to get there since I cannot evaluate the double integral:
$$ \int_{-\infty}^{\infty} \int_{-\infty}^{y_2} \left( y_2-y_1 \right) \frac{1}{2\pi \sigma^2} exp \left\{ -\frac{1}{2\sigma^2}\left[ \left( y_1-\mu \right)^2 +\left( y_2-\mu \right)^2 \right]\right\} \mathrm {dy_1 dy_2}$$
Any ideas on how to compute this are greatly appreciated, thank you!
 A: Here is a quick check using a computer algebra system. I am using the mathStatica package for Mathematica (I am one of the developers of the former) to automate the nitty gritties for me ... 
Given: The parent pdf is $N(\mu, \sigma^2)$ with pdf $f(y)$:

Then, the joint pdf of the 1st and 2nd order statistics $(Y_1, Y_2)$, in a sample of size 2, denoted say $g(y_1,y_2)$ can be easily obtained using the OrderStat function:

Note that the constant multiplier here is $\frac{1}{\pi \sigma^2}$, not $\frac{1}{2\pi \sigma^2}$ in your equation.
Because the constraint $Y_1 < Y_2$ is already entered into the pdf definition, we can enter the domain of support on the real line as:

Finally, the expectation you seek is: 

which agrees with your stated solution. 

@Ioannis wrote: 

Is there a way to see how the integral is computed though?

One can activate VerboseMode, which shows all the integrands being sent off for calculation. With VerboseMode[On], one can see the intermediary integrands ...

You might need to open the pic manually in a separate window to see the detail ...
A: As was mentioned in the comment section, the sample range $Y_2-Y_1$ can be expressed as $Y_2-Y_1=X_{(2)}-X_{(1)}=|X_1-X_2|$ where $(X_1,X_2)$ is the random sample under consideration.
Since $(X_1,X_2)$ is i.i.d $\mathcal{N}(\mu,\sigma^2)$, by the reproductive property of normal distribution, we have $X_1-X_2\sim\mathcal{N}(0,2\sigma^2)$. So, $E(|X_1-X_2|)$ is simply the mean absolute deviation of $X_1-X_2$ about its mean. 
It is well-known (and can be easily shown) that, for some $X\sim\mathcal{N}(\mu,\sigma^2)$, one has mean absolute deviation about mean $E(|X-\mu|)=\sqrt{\frac{2}{\pi}}\sigma$.
So, for $X\sim\mathcal{N}(0,\sigma^2)$, we readily get $E(|X|)=\sqrt{\frac{2}{\pi}}\sigma$
Hence, $E(Y_2-Y_1)=E(|X_1-X_2|)=\sqrt{\frac{2}{\pi}}\cdot\sqrt{2}\sigma=\frac{2\sigma}{\sqrt{\pi}}$
This, in my opinion, is easier to show than to find the expectation from the distribution of $Y_2-Y_1$ or to find $E(Y_2)-E(Y_1)$ from the distributions of $Y_1$ and $Y_2$.

Interestingly, by the same logic, it can be shown that expected value of the sample range $R$ from the same population when the sample size is three is $E(R)=\frac{3\sigma}{\sqrt{\pi}}$.
