Bivariate Normal Conditional Variance
Rather than embarking on some pretty involved computations of conditional distributions, one should rely on one of the main assets of Gaussian families, namely, the...
Key feature: In Gaussian families, conditioning acts as a linear projection.
Hence, as the OP suggested, one could do worse than to start from a representation of $(X,Y)$ by standard i.i.d. Gaussian random variables $U$ and $V$, for example, $$ X=\mu_x+\sigma_xU\qquad Y=\mu_y+\sigma_y(\rho U+\tau V)$$ where the parameter $\tau$ is $$\tau=\sqrt{1-\rho^2} $$ Since $\sigma_x\ne0$, the sigma-algebra generated by $X$ is also the sigma-algebra generated by $U$ hence conditioning by $X$ or by $U$ is the same. Furthermore, constants and functions of $X$ or $U$ are all $U$-measurable while functions of $V$ are independent on $U$, thus, $$ \mathrm E(Y\mid X)=\mu_y+\sigma_y(\rho U+\tau \mathrm E(V))=\mu_y+\sigma_y \rho U $$ which is equivalent to $$ \color{red}{\mathrm E(Y\mid X)=\mu_y+\rho\frac{\sigma_y}{\sigma_x}(X-\mu_x)} $$ Likewise, when computing conditional variances conditionally on $X$, deterministic functions of $X$ or $U$ should be considered as constants, hence their conditional variance is zero, and functions of $V$ are independent on $X$, hence their conditional variance is their variance. Thus, $$ \mbox{Var}(Y\mid X)=\mbox{Var}(\sigma_y\tau V\mid X)=\sigma_y^2\tau^2\mbox{Var}(V\mid X)=\sigma_y^2\tau^2\mbox{Var}(V) $$ that is, $$ \color{red}{\mbox{Var}(Y\mid X)=\sigma_y^2(1-\rho^2)} $$ Finally, the event $$A=[X>\mu_x,Y>\mu_y]$$ is also $$ A=[U>0,\rho U+\tau V>0]. $$ To evaluate $\mathrm P(A)$, one can turn to the planar representation of couples of independent standard Gaussian random variables, which says in particular that the distribution of $(U,V)$ is invariant by rotations. The event $A$ means that the direction of the vector $(U,V)$ is between the angle $\vartheta$ in $(-\pi/2,\pi/2)$ such that $$\tan(\vartheta)=-\rho/\tau$$ and the angle $\pi/2$. Thus, $$\mathrm P(A)=\frac{\pi/2-\vartheta}{2\pi}$$ that is, $$ \color{red}{\mathrm P(X>\mu_x,Y>\mu_y)=\frac14+\frac1{2\pi}\arcsin\rho} $$ Numerical application: If $\mu_x=2$, $\mu_y=-1$, $\sigma_x=2$, $\sigma_y=1$ and $\rho=-\sqrt3/2$, then $$ \mathrm E(Y\mid X)=-1+\sqrt3/2-(\sqrt3/4)X\qquad \mbox{Var}(Y\mid X)=1/4 $$ and $\tau=1/2$, hence $\vartheta=\pi/3$ and $$\mathrm P(A)=1/12$$
First, the joint PDF $f(x,y)$ is obvious, just plug in your parameters. Bivariate Normal. Then you can find the marginal density for $X$, which gives you the conditional density of $Y$ given $X=x$: $$f_{Y|X}(y|x)=\frac{f(x,y)}{f_X(x)}.$$ Now use the conditional density you can evaluate both conditional expectation and conditional variance : $$\mathbb{E} (Y|X=x)=\int_{-\infty}^\infty y f_{Y|X}(y|x)dy,$$ and $$\text{Var} (Y|X=x)=\int_{-\infty}^\infty (y-h(x))^2 f_{Y|X}(y|x)dy=\frac14,$$ where $h(x)=\mathbb{E} (Y|X=x)=-\frac{\sqrt 3}4(x-2)-1$.
And with the joint PDF, $P(X>\mu_x, Y > \mu_y)$ is just an integration: $$P(X>\mu_x, Y > \mu_y)=\int_{\mu_x}^\infty\int_{\mu_y}^\infty f(x,y)dydx=\frac1{12},$$ though I guess there's an easier way to compute.