Consider "random" x, y, z such that x + y + z = 1.

A number of different definitions of "random" might make sense here, but I'm interested in those for which expected value of all three random variables is 1/3 and where the possible range includes everything from 0 to 1.

I think choosing a random point in an equilateral triangle and finding its barycentric coordinates fits this bill.

Question 1: is this the case?

But i'm also interesting in finding a simple algorithm for actually generating these values. And, in particular, one that scales easily to any number of dimensions:

x1 + x2 + ... + xn = 1

Generate random x_i where each has expected value 1/n but can take on any value between 0 and 1.

Question 2: What is the simplest algorithm to accomplish generating such numbers?

I think so.

My first reflex is to try something like this: Pick three random integers A, B, and C from 0 to N (uniform) with the condition that not all three are zero (redraw all three if all are zero) and construct the following value:

(A/(A + B +C), B/(A + B +C), C/(A + B +C))

This is symmetric and it satisfies that the sum of the coordinates is 1. This should mean that the average of the individual coordinates is 1/3. It's also easy to generalize. So it satisfies your conditions.

However, I don't know whether this is a "good" distribution for your application. I don't think it reflects the idea of a randomly chosen point on the triangle with uniform distribution. (It might for "large" N. But I didn't give this too much thought.)

Generate two random numbers between 0 and 1 using a uniform distribution. Let's call them a and b. Discard if a+b > 1. Add b/2 to a. Plot these points on your triangle defined by (1,0,0) (0,1,0) (0,0,1) making the point a=0, b=0 correspond to (1,0,0), a=1, b=0 correspond to (0,1,0), and a=1/2, b=1 correspond to (0,0,1). Recover the x,y,z values of the point plotted. That last bit is the hard part, I know, but shouldn't be too hard.

Generate 3 iid uniform(0,1)s. Let their sum be S.

Return (x1/S, x2/s, x3/s). The mean of each component is 1/3, as desired; and it sums to 1.

This generalises to n dimensions easily.

thanks all. aaron / rivertilt's suggestion should work fine for what I need.

As noted by everyone else, taking y_i i.i.d. with any distribution supported on some [0,a] (a possibly infinity) and dividing by the sum will satisfy your requirements, but it will typically not give you something which is uniformly distributed on the region of interest.

Useful keyword here is "simplex" (and "uniform"), google has a lot of material...

According to this math stackexchange link, taking the y_i exponential (-log(uniform)) then normalizing gives you something uniformly distributed on the simplex.

checktheriver, thanks for the link.