Multiple Integration non rectangular domain I'm trying to do the following question from an advanced calculus past paper at university:
Suppose that $f(x)$ is continuous in $0\le x  \le 1$ show that
$$\int_0^1 \int_x^1 \int_x^y f(x)f(y)f(z) dzdydx = \frac{1}{6} \left( \int_0^1 f(x) dx \right)^3$$ 
I've found an example in two dimensions (so answer is area square /2) but can't see how to apply the same method to 3 dimensions. Not sure it's correct, but I get down to.
$$I=\frac{1}{2}\int_0^1 \int_0^1 \int_x^y f(x)f(y)f(z) dzdydx$$
Any pointers would be appreciated. I have a feeling the $\int_0^x$ and $\int_x^y$ and $\int_y^1$ for z will be equal giving the factor of 3?
 A: I knew there is a simpler way:
$$\eqalign{\left(\int_0^1 f(x)\>dx\right)^3&=\int_{[0,1]^3}f(x)f(y)f(z)\>{\rm d}(x,y,z)\cr &=6\int_{0\leq x\leq z\leq y\leq1}f(x)f(y)f(z)\>{\rm d}(x,y,z)\cr &=6\int_0^1\int_x^1\int_x^y f(x)f(y)f(z)\>dz\,dy\,dx\ .\cr}$$
A: It's just a repeated routine calculation, so I'll just get you started.
$f(x)$ being continuous on $[0,1]$ $\displaystyle \Rightarrow \exists \ F(x) = \int_0^x f(s) \, \mathrm{d}s$ for $x \in [0,1]$.
This implies that
$$
\int_0^1 \int_x^1 \int_x^y f(x) f(y) f(z) \, \mathrm{d}z \, \mathrm{d}y \, \mathrm{d}x = \int_0^1 \int_x^1 f(x) f(y)[F(y) - F(x)] \, \mathrm{d}y \, \mathrm{d}x
$$
You can then let $u = F(y)$ to get
$$
\int_0^1 \left [f(x) \int_x^1 f(y)F(y) \, \mathrm{d}y - F(x)f(x)\int_x^1 f(y) \, \mathrm{d}y \right ] \, \mathrm{d}x \\ = \int_0^1 \frac{F(1)^2}{2}f(x) + \frac{F(x)^2f(x)}{2} - F(1)F(x)f(x) \, \mathrm{d}x
$$
Making the same substitution and using the fact that $\displaystyle F(1) = \int_0^1 f(x) \, \mathrm{d}x $, you will reach your desired result.
A: For your two-dimensional example, I suppose you have in mind the fact that
$$
 \int_0^1 \!\!\int_x^1 f(x)f(y)\,dy\,dx
 = \int_0^1 \!\!\int_0^x f(x)f(y)\,dy\,dx
= \frac12 \left( \int_0^1 f(x)\,dx \right)^2.
$$
This is not a carte-blanche license to exchange 
$\int_x^1 \cdots dy$ with $\int_0^x \cdots dy$
in an arbitrary multiple integral.
Rather, this exploits a symmetry of the function $g(x,y) = f(x)f(y)$ on the
square $[0,1]\times[0,1],$ namely, you can swap the two parameters without changing the value: $g(x,y) = g(y,x).$
Moreover, the region of integration for
$$
\int_0^1 \!\!\int_0^x \cdots\,dy\,dx
$$
is exactly the same as the region of integration for
$$
\int_0^1 \!\!\int_y^1 \cdots\,dx\,dy
$$
(namely, it is the triangle bounded by the lines $y=0,$ $x=1,$ and $y=x$).
It follows that
\begin{align}
\int_0^1 \!\!\int_0^x g(x,y)\,dy\,dx
&= \int_0^1 \!\!\int_y^1 g(x,y)\,dx\,dy
 && \text{same region of integration} \\
&= \int_0^1 \!\!\int_x^1 g(y,x)\,dy\,dx 
 &&\text{substitute $(y,x)$ for $(x,y)$}\\
&= \int_0^1 \!\!\int_x^1 g(x,y)\,dy\,dx && g(x,y) = g(y,x),
\end{align}
and of course
\begin{align}
\int_0^1 \!\!\int_0^x f(x)f(y)\,dy\,dx + \int_0^1 \!\!\int_x^1 f(x)f(y)\,dy\,dx
&= \int_0^1 \!\!\int_0^1 f(x)f(y)\,dy\,dx \\
&= \int_0^1 f(x)\left(\int_0^1 f(y)\,dy\right) \,dx \\
&= \left( \int_0^1 f(x)\,dx \right)^2.
\end{align}
This symmetry about the line $y=x$ extends to a symmetry about the plane $y=x$ in your triple integral.
If we let $h(x,y,z)=f(x)f(y)f(z),$ then
\begin{align}
\int_0^1 \!\!\int_0^x \!\!\int_x^y &\, h(x,y,z)\,dz\,dy\,dx \\
&= \int_0^1 \!\!\int_y^1 \!\!\int_x^y h(x,y,z)\,dz\,dx\,dy
 && \text{same region} \\
&= \int_0^1 \!\!\int_x^1 \!\!\int_y^x h(y,x,z)\,dz\,dy\,dx 
 &&\text{$(y,x) \mapsto (x,y)$}\\
&= \int_0^1 \!\!\int_x^1 \!\!\int_y^x h(x,y,z)\,dz\,dy\,dx
 && h(x,y,z) = h(y,x,z)\\
&= -\int_0^1 \!\!\int_x^1 \!\!\int_x^y h(x,y,z)\,dz\,dy\,dx
 && \text{swap integration limits for $z$}.
\end{align}
Therefore 
\begin{align}
\int_0^1 \!\!\int_0^1 \!\!\int_x^y &\, h(x,y,z)\,dz\,dy\,dx \\
&=  \int_0^1 \!\!\int_x^1 \!\!\int_x^y  h(x,y,z)\,dz\,dy\,dx 
  + \int_0^1 \!\!\int_0^x \!\!\int_x^y  h(x,y,z)\,dz\,dy\,dx \\
&=  \int_0^1 \!\!\int_x^1 \!\!\int_x^y  h(x,y,z)\,dz\,dy\,dx 
  - \int_0^1 \!\!\int_x^1 \!\!\int_x^y  h(x,y,z)\,dz\,dy\,dx \\
&= 0,
\end{align}
so whatever you did to conclude that
$I = \frac12 \int_0^1 \int_0^1 \int_x^y f(x)f(y)f(z)\,dz\,dy\,dx,$
there is an error somewhere in those calculations.
The full symmetry to use in the case of the triple integral is to
realize that it is possible not only to swap the variables $x$ and $y,$
but indeed any permutation of the three variables $x,$ $y,$ and $z$
will produce an integral with the same value.
Moreover, except on a set of measure zero,
every point in the cube $[0,1]\times[0,1]\times[0,1]$
satisfies $\sigma(x) < \sigma(y) < \sigma(z)$ for exactly one permutation
$\sigma$ of the variables $x,$ $y,$ and $z.$
Therefore the union of the six regions of integration we get by permuting the 
variables is the entire cube,
and since any pair of these regions overlaps only on a set of measure zero,
the sum of the integrals over these regions is 
$\int_0^1 \int_0^1 \int_0^1 f(x)f(y)f(z)\,dz\,dy\,dz.$
Christian Blatter's answer states the same conclusion much more succinctly.
