I'm aware that $Y = \sum X_i$ where $X_i \sim U(0,1)$ is the Irwin-Hall distribution. I'm interested in the distribution $$ Z_i = \frac{X_i}{\sum X_i}. $$ I'm aware this maybe difficult for $n$ random variables. So I tried it for $n = 2$ and obtained $$ f_Z(z) = \begin{cases} \frac{1}{2(1-z)^2}, & 0 < z < 1/2 \\ \frac{1}{2z^2} & 1 \ge z \ge 1/2. \end{cases} $$ The case $n = 3$ is more difficult to visualize how the limits of integration should be obtained.
Attempt: \begin{align} F_Z(z) & = P(Z \le z) \\ & = P \left(\frac{X_1}{X_1 + X_2 + X_3} \le z \right) \\ & = P(X_1(1-z) \le z(X_2 + X_3)). \end{align} At this point I'm unsure how to set up the triple integral... $$ F_Z(z) = \int_0^1 \int_?^? \int_?^? dx_3 dx_2 dx_1. $$
Find the distribution of $1/Z_1$, which equals
$$\frac{1}{Z_1}=1+\frac{X_2+X_3}{X_1}.$$
You can find the distribution of $X_2+X_3$, which is a triangular distribution, and then since $X_2+X_3$ is independent of the uniformly distributed $X_1$, you can use the product of their distributions to find the distribution of their ratio, shift it by $1$ and do the reciprocal transform. I think this method generalizes to $n$ variables. For the $n=3$ case, the probability density of $Z_1$ is given by \begin{align} p(z)&=\int_0^1dx\int_0^2dy\,(1-|y-1|)\,\delta\left(\frac{y}{x}+1-\frac{1}{z}\right)\left|\frac{d}{dz}\frac{1}{z}\right|, \end{align} where $\delta(\cdot)$ is the Dirac-delta function, $x=X_1$ is uniformly distributed on the support $[0,1]$, and $y=X_2+X_3$ is triangularly distributed with the density $1-|y-1|$ on the support $[0,2]$. The last absolute value factor is the Jacobian of transforming from $1/z$ to $z$, which equals $1/z^2$, since we have $z\in(0,1)$. So we have \begin{align} p(z)&=\frac{1}{z^2}\int_0^2dy\,(1-|y-1|)\int_0^1dx\,\delta\left(\frac{y}{x}+1-\frac{1}{z}\right)\\ &=\frac{1}{z^2}\int_0^2dy\,(1-|y-1|)\,\left.\Theta(0\leq x\leq 1)\,\frac{x^2}{y}\right|_{x=\frac{yz}{1-z}}\\ &=\frac{1}{(1-z)^2}\int_0^2ydy\,(1-|y-1|)\,\Theta\!\left(y\leq\frac{1-z}{z}\right), \end{align} where $\Theta(\cdot)$ is the Heaviside step function which equals $1$ if the condition in the bracket is satisfied and $0$ otherwise. The rest is for Wolfram alpha. The result splits into $z\in[0,1/3], (1/3,1/2], (1/2,1]$ and is continuous. To do calculations for $n>3$, just replace the probability density of $y$ by the general $(n-1)$-variable Irwin-Hall distribution. The calculation is harder but still doable.