Random triple with uniform marginals Output a random triple \$(x,y,z)\$ of real numbers such that


\$x+y+z=\frac{3}{2}\$
Each of \$x,y,z\$, taken alone, is uniformly distributed from \$0\$ to \$1\$.


Please explain in your answer what distribution your code produces and why it meets these conditions. In particular, please include a demonstration that with many samples produced by running your code many times, the distribution of each variable alone is approximately uniform.

Allowances for precision

Since computers can't produce actual real numbers, you can approximate them by outputting any of:


Floats
Decimals
Rationals


I've tried to make the rules here as friendly as possible. 

Your outputs need to be accurate within \$\pm 10^{-3}\$, that is 3 decimal digits. That means you can imagine that each real-number output from your true distribution has been fudged by at most \$\pm 10^{-3}\$ in any way. It's fine if this makes \$x+y+z\$ be a bit off from \$\frac{3}{2}\$.

You may output rationals as (numerator, denominator) pairs of integers or similar. You may also choose to output just one denominator that applies to all three values, that is, output \$a,b,c\$ and \$d\$ to represent \$\frac{a}{d},\frac{b}{d},\frac{c}{d}\$. Rationals don't have to be in reduced form.

Randomness

You can assume any standard PRNG is in fact random. For instance, if it outputs random floats from 0 to 1, you can assume these are uniformly random reals that are independent from each other, and not worry about any endpoint issues. It's also fine if your code fails in a probability-zero event.

To decide whether a pseudo-random source such as time or memory contents is sufficiently random, consider an empirical test where your code is run many times in sequence to check if the randomness varies sufficiently between runs. You may have these runs be done within your language such as in a loop, or by repeatedly calling the code from outside like with a shell command. You may not assume any intentional time delay or other actions taken between runs.
APL (Dyalog Unicode), 15 bytesSBCS



1|(-,.5∘-,+⍨)?0


Try it online!

A full program that generates three floating-point numbers \$ X,Y,Z \$.

The distribution

Uses the same method as Luis Mendo's MATL answer, which in turn is from this Math.SE answer. The Math.SE answer presents the following random variables:

$$
X = \text{uniform random value in } (0,1) \\
Y = \begin{cases} X + \frac12, \quad X &lt; \frac12 \\
X - \frac12, \quad X \geq \frac12 \end{cases} \\
Z = \frac32 - (X + Y)
$$

I found that this is equivalent to

$$
X = \text{uniform random value in } (0,1) \\
Y = \{ X + 1/2 \}, Z = \{-2X\}
$$

which, in turn, has the same distribution with

$$
r = \text{uniform random value in } (0,1) \\
X = \{ -r \}, Y = \{ 1/2 - r \}, Z = \{ 2r \}
$$

where the notation \$ \{ x \} \$ means the positive fractional part of \$ x \$.

How the code works

1|(-,.5∘-,+⍨)?0
             ?0  ⍝ r←Generate a random number from (0,1)
  ( ,    ,  )    ⍝ Create a 3-element array from r:
   -             ⍝ -r
     .5∘-        ⍝ 0.5 - r
          +⍨     ⍝ r + r (= 2r)
1|               ⍝ Extract positive fractional part from each number


Test script

⍝ Take 10000 samples of X, Y, Z's
data←f¨⍳10000

⍝ Tests that X+Y+Z = 1.5 in all cases
∧/1.5=+/¨data

⍝ Shows the distribution of X,Y,Z in 10 bins
{(⊂∘⍋⌷⊢)↓,∘≢⌸⍵}¨⌊10×↓⍉↑data
                    ↓⍉↑data  ⍝ Reform the data into X's, Y's, and Z's
                ⌊10×         ⍝ Classify into 10 bins
{             }¨     ⍝ For each of X's, Y's, or Z's,
        ↓,∘≢⌸⍵       ⍝ Gather the list of (bin, count)
 (⊂∘⍋⌷⊢)             ⍝ Sort the bins

⍝ Example result
X:  0 1012  1 1016  2 993   3 1025  4 1002  5 969   6 1002  7 977  8 1049  9 955  
Y:  0 969   1 1002  2 977   3 1049  4 955   5 1012  6 1016  7 993  8 1025  9 1002 
Z:  0 1013  1 944   2 1022  3 1052  4 990   5 980   6 1039  7 979  8 959   9 1022 




Alternative method, 21 bytesSBCS



(3*9),3⊥¨3|(⍳3)+⊂?9/3


Try it online!

This one uses the OP's formulation from the same Math.SE question:


Generate a uniform random number \$ X \$ from \$ (0, 1) \$.
Take its ternary digits \$ X = 0.X_1 X_2 X_3 \cdots \$ and generate ternary digits of \$ Y \$ and \$ Z \$:


$$
Y_i = (X_i + 1) \% 3 \\
Z_i = (X_i + 2) \% 3
$$


Then the resulting numbers \$ X,Y,Z \$ are uniform in range \$ (0, 1) \$, and their sum is


$$
\frac3{3^1} + \frac3{3^2} + \frac3{3^3} + \cdots = \frac32
$$

Theoretically it requires infinite number of digits; the code approximates the scheme with just 9 digits, which makes \$ X+Y+Z=1.499923792 \$. I'm not exactly sure if it meets the OP's spec though.

How the code works

(3*9),3⊥¨3|(⍳3)+⊂?9/3
                 ?9/3  ⍝ Generate 9 ternary digits
         3|(⍳3)+⊂      ⍝ Ternary digits for X,Y,Z
      3⊥¨  ⍝ Convert each list of digits to numerator for the fraction
(3*9),     ⍝ Prepend the denominator