You can simulate Brownian motion on $[0, 1]$ for instance by splitting it into $K$ intervals and at each time $k \Delta t$ add $N(0, \Delta t)$ to your running total, where $\Delta t = 1/K$.
If you only have a fair coin and can't sample from $N(0, \Delta t)$ you can instead flip the coin and add $\sqrt{\Delta t}$ if heads, $-\sqrt{\Delta t}$ if tails.
In some rigorous sense this approximation converges to Brownian motion.
My question is if instead of a fair coin you have a die, how would you use it to sample from $N(0, \Delta t)$?
I thought of using a discrete version of the probability integral transform, but that doesn't seem to give the right variance.
Edit: the kind of solution I'm looking for would be to find the "best" vector $v \in \mathbb{R}^n$ to map the die's outcomes with.
For $n=2$ (the two-sided die) the solution is $v = [-\sqrt{\Delta t}, \sqrt{\Delta t}]'$.
For $n=3$, it seems like $v = [-\sqrt{(3/2) \Delta t}, 0, \sqrt{(3/2) \Delta t} ]'$ is the right pick: it's symmetric and has the right variance.
For $n=4$, we can use the tetrahedron as two coins (like a binomial tree) and get the vector $v = [-\sqrt{2 \Delta t}, 0, 0, \sqrt{2 \Delta t}]'$. There are many other vectors that have the same variance and are symmetric, for example $[- \sqrt{\Delta t}, -\sqrt{\Delta t}, \sqrt{\Delta t}, \sqrt{\Delta t}]'$. This is what we get if we take the parity of the outcome.
Arguably both these vectors are "suboptimal" in the sense that it intuitively makes sense that there should be one value from each quartile.
I guess the question is, what am I trying to optimise?
If X represents the outcome of the roll of a fair die, the X is a random variable with mean (7/2) and standard deviation $\sqrt{(35/12)}$.
$(X-7/2)\sqrt{(12/35)}$ then gives you a random variable with mean 0 and standard deviation 1.
For a large number of rolls $\sum (X_i-7/2)\sqrt{(12/35)\Delta t}$ is approximately a standard normal random variable.