You are just asking to compute $p=P(X>Y)$$p=P(X>Y)=P(Z<0)$, where $(X,Y)$ has the bivariate normal distribution with given $EX=\mu_1$, $EY=\mu_2$, $Var\,X=\sigma_1^2:=\Sigma_{1,1}$, $Var\,Y=\sigma_2^2:=\Sigma_{2,2}$, and $\rho:=corr(X,Y)=\Sigma_{1,2}/(\sigma_1 \sigma_2)$, and $Z:=Y-X\sim N(\mu,\sigma^2)$, where $\mu:=\mu_2-\mu_1$ and $\sigma^2:=\sigma _1^2+\sigma _2^2-2 \rho \sigma _1 \sigma_2$. It is easy to see that this probability isSo, $$p=\frac{1}{2} \text{erfc}\Big(\frac{\mu _2-\mu _1}{\sqrt{2} \sqrt{\sigma _1^2+\sigma _2^2-2 \rho \sigma _1 \sigma_2}}\Big) =\Phi\Big(\frac{\mu _1-\mu _2}{ \sqrt{\sigma _1^2+\sigma _2^2-2 \rho \sigma _1 \sigma_2}}\Big),$$$$p=\Phi\Big(\frac{\mu _1-\mu _2}{ \sqrt{\sigma _1^2+\sigma _2^2-2 \rho \sigma _1 \sigma_2}}\Big)$$ or, equivalently, $$\frac{\mu _1-\mu _2}{ \sqrt{\sigma _1^2+\sigma _2^2-2 \rho \sigma _1 \sigma_2}}=\Phi^{-1}(p), \tag{1}$$ where $\Phi$ is the standard normal cdf and $\Phi^{-1}$ is the function inverse to $\Phi$. We have $p\in(\frac12,1)$ iff $\mu_1>\mu_2$.