Optimize function on $x^2 + y^2 + z^2 \leq 1$ Optimize $f(x,y,z) = xyz + xy$ on $\mathbb{D} = \{ (x,y,z) \in \mathbb{R^3} : x,y,z \geq 0 \wedge x^2 + y^2 + z^2 \leq 1  \}$. The equation $\nabla f(x,y,z) = (0,0,0)$ yields $x = 0, y = 0, z \geq 0 $ and we can evaluate $f(0,0,z) = 0$.
Now studying the function on the boundary $x^2 + y^2 + z^2 = 1$ gets really hairy. I tried replacing $x$ with $\sqrt{1 - y^2 - z^2}$ in order to transform $f(x,y,z)$ into a two-variable function $g(y,z)$ and optimize it on $y^2 + z^2 \leq 1$ but $g(y,z)$ is a pain to differentiate. I then tried spherical coordinates which really did not make it any much easier.
Got any suggestions on how to tackle it?
 A: Now $f(x,y,z) = xyz + xy=xy(1+z)$
Note that on $x^2 + y^2 + z^2 = 1$, $|z| \le 1$ so $1+z \ge 0$.
Hence, for a given $z$, you maximise $f(x,y,z)$ by maximising $xy$ on $x^2 + y^2 = 1-z^2$
$x^2 + y^2 = 1-z^2$ is a circle of radius $\sqrt{1-z^2}$, so we can parametrise $(x,y)$ by $(\sqrt{1-z^2} \;sin(\theta), \sqrt{1-z^2} \;cos(\theta))$, for $0 \le \theta \lt 2\pi$
Hence we want to maximise $xy = \sqrt{1-z^2}*sin(\theta) * \sqrt{1-z^2} *cos(\theta) = (1-z^2) \dfrac{sin(2\theta)}{2}$.
$sin(2\theta)$ has a obvious maximum of $1$, so the maximum of $xy$ is $\dfrac{(1-z^2)}{2}$
Now we choose $z$ to maximise $\dfrac{(1-z^2)}{2}(1+z) = \dfrac{(1+z-z^2-z^3)}{2}$
Setting the derivative to zero gives $\dfrac{(1-2z-3z^2)}{2} \equiv 0$
$(3z-1)(z+1) \equiv 0$
i.e. $z=\dfrac{1}{3}$ or $z=-1$
$z=-1$ corresponds to $f(x,y,z)=0$ so this must be a minimum.
$z=\dfrac{1}{3}$ corresponds to $f(x,y,z)=\dfrac{(1 - \dfrac{1}{9})}{2}(1 + \dfrac{1}{3}) = \dfrac{16}{27}$
A: $f=xy(z+1) \le \dfrac{x^2+y^2}{2}(z+1)=\dfrac{(z+1)(k-z^2)}{2},k\le1$
$f'_z=\dfrac{k-z^2-2z^2-2z}{2}=0,z=\dfrac{-2 \pm \sqrt{4+12k}}{6}=\dfrac{-1 \pm \sqrt{1+3k}}{3}$
for max: $z=\dfrac{-1 + \sqrt{1+3k}}{3},f_{max}=\dfrac{(2+\sqrt{1+3k})(6k-1+2\sqrt{1+3k})}{2*27}\le \dfrac{32}{2*27}$ ,
when $k=1$ get max becasue $6k-1+2\sqrt{1+3k}$ is mono increasing funtion 
$z= \dfrac{1}{3}, x=y=\dfrac{2}{3}$
A: Another way:
$f$ is convex in each variable $x, y, z$ so the maximum can be obtained only on the boundary.  This means $x^2+y^2+z^2=1$.  Further,
$$4 \,f = 4xy(1+z) \leqslant 2(x^2+y^2)(1+z) = 2(1-z^2)(1+z) = (1+z)(2-2z)(1+z)$$
The RHS now is a product of three positive terms with a constant sum, so it gets maximised when all terms are the same, viz. when $1+z = 2- 2z \implies z=\frac13$.  As the first inequality becomes an equality iff $x=y$, the maximum is indeed achieved when $z = \frac13, \; x=y=\frac23$ and $\displaystyle f_{max} = \frac{16}{27}$.
