[SciPy-User] Conditional bivariate normal

[Ted To]

Hi All,

I have a puzzle that I'm having trouble figuring out.  Suppose U=X+Y
where X and Y are independent normals.  Sampling X and Y conditional on
U>=\bar U takes an inordinate amount of time since at times \bar U can
be fairly large so I've been thinking about how to shorten the time.  I
thought I came upon a solution by using scipy.stats.truncnorm.rvs to
first draw U=u and then draw an X=x conditional on U=u.  Using U=u and
X=x, Y=y=u-x.

I'm getting the correct means and the correct standard deviation for U
but the standard deviations for X and Y are too small.  Is there
something wrong with my logic or have I incorrectly derived the sd for
X|U=u?

Thanks,
Ted