Solutions of a set of trigonometric equations I am trying to find maxima of a linear combination of functions $\text{trig}(x)\text{trig}(y)\text{trig}(z)$, where $\text{trig}$ is $\sin$ or $\cos$ (which of course gives 8 such terms):
$$
f(x,y,z) = 
A\cos{x}\cos{y}\cos{z}
+B\cos{x}\cos{y}\sin{z}
+C\cos{x}\sin{y}\cos{z}
+D\cos{x}\sin{y}\sin{z}\\
+E\sin{x}\cos{y}\cos{z}
+F\sin{x}\cos{y}\sin{z}
+G\sin{x}\sin{y}\cos{z}
+H\sin{x}\sin{y}\sin{z}
$$
Coefficients of this linear combination are arbitrary. After a series of substitutions, I was able to reduce it to maximizing the following function:
$$
g(u,v,w)=K\cos(u+v+w+e) + L\cos u + M\cos v + N\cos w.
$$
Equating its gradient to 0, we are presented with a set of equations:
$$
\begin{cases}
K\sin(u+v+w+e) + L\sin u = 0, \\
K\sin(u+v+w+e) + M\sin v = 0, \\
K\sin(u+v+w+e) + N\sin w = 0.
\end{cases}
$$
Unfortunately I am not able to proceed further. Do you have any ideas how to find a solution?
Analytical solution would be great, however in the end it is used in computer calculations, so numerical solutions are also accepted, however they should be stable and work for arbitrary coefficients. For some specific values of free parameters maxima can be lines or planes, so these cases have to be detected (I am probably able to figure it out myself).
Thanks.
 A: $\def\c{{\rm c}}\def\s{{\rm s}}$
Hint ($c_x=\cos(x),s_x=\sin(x)$ and so on):
One can approach the problem in the following way. Write:
$$
\begin{align}
f(x,y,z)
&=[A\,c_xc_y+C\,c_xs_y+E\,s_xc_y+G\,s_xs_y]\cos(z)\\
&+[B\,c_xc_y+D\,c_xs_y+F\,s_xc_y+H\,s_xs_y]\sin(z)\\
&=g(x,y) \cos(z-\phi(x,y))\tag1
\end{align}
$$
with
$$
\begin{align}
g^2(x,y)
&=[A\,c_xc_y+C\,c_xs_y+E\,s_xc_y+G\,s_xs_y]^2\\
&+[B\,c_xc_y+D\,c_xs_y+F\,s_xc_y+H\,s_xs_y]^2\tag2
\end{align}
$$
Clearly the function (1) has maximum value $g(x,y)$ at $z=\phi(x,y)$. Therefore we need to maximize $g^2(x,y)$. This can be dealt exactly as we have done with $f(x,y,z)$ expanding in (2) the squares and using
$$\cos^2y=\frac{1+\cos2y}2,\quad
\sin^2y=\frac{1-\cos2y}2,\quad \cos y\sin y=\frac{\sin2y}2.\tag3$$
After doing this you will obtain an expression
$$
\lambda(x)+h(x)\cos(2y-2\psi(x)),
$$
which after setting $y=\psi(x)$ reduces to
a simple problem of maximizing $\lambda(x)+h(x)$.
PS. By my calculation
$$\begin{align}
\lambda(x)&=\frac12\Big[(A^2+B^2+C^2+D^2)c_x^2+(E^2+F^2+G^2+H^2)s_x^2
+2(AE+BF+CG+DH)c_xs_x\Big]\\
h^2(x)&=\frac14\Big[(A^2+B^2-C^2-D^2)c_x^2+(E^2+F^2-G^2-H^2)s_x^2
+2(AE+BF-CG-DH)c_xs_x\Big]^2\\
&+\quad\Big[(AC+BD)c_x^2+(EG+FH)s_x^2+2(AG+BH+CE+DF)c_xs_x\Big]^2
\end{align}$$
but I could make a mistake.
A: Assuming that $e$ is known, rewrite the equations as
$$Z + L\sin (u) = 0 \tag 1$$
$$Z + M\sin (v) = 0 \tag 2$$
$$Z + N\sin (w) = 0 \tag 3$$
$$Z =K\sin(u+v+w+e) \tag 4$$
Hoping that there is modulo problems, using $(1),(2),(3)$ we have
$$u=-\sin ^{-1}\left(\frac{Z}{L}\right)\qquad\qquad v=-\sin ^{-1}\left(\frac{Z}{M}\right)\qquad\qquad w=-\sin ^{-1}\left(\frac{Z}{N}\right)\tag 5$$
Plug $(5)$ into $(1)$ to obtain
$$K \sin \left(e-\sin ^{-1}\left(\frac{Z}{L}\right)-\sin ^{-1}\left(\frac{Z}{M}\right)-\sin
   ^{-1}\left(\frac{Z}{N}\right)\right)-Z=0 \tag 6$$ Expanding the sine, this gives an algebraic equation in $Z$.
When $Z$ is known, go back to $(5)$ to get $(u,v,w)$.
May be, this could work
A: This accepted answer was most helpful in devising a solution, however I made some modifications, so I am posting my own additional answer to describe them.
In the cited answer the last step was to find solutions of the equation $\frac{d}{dx}[h(x) + \lambda(x)] = 0$. One could write it as $\frac{dh(x)}{dx} = -\frac{d\lambda(x)}{dx}$, compute derivatives, square both sides to eliminate square roots, multiply by denominators and divide by $\cos^8x$ ending up with 8-th degree polynomial in $\tan x$. It can then be solved numerically (using specialized polynomial root finding algorithms, such as Aberth or Jenkins-Traub methods). One then examines all solutions and checks which one gives the actual maximum.
I was able to improve numerical stability by reducing the polynomial's degree to 6. First we solve the problem in 2D:
$$
f_{2D}(x,y) = A'\cos x \cos y + B'\cos x \sin y + C'\sin x \cos y + D'\sin x \sin y
$$
which we rewrite as
$$
\begin{align}
f_{2D}(x,y) &= \underbrace{\frac{1}{2}\sqrt{(A'-D')^2 + (C'+B')^2}}_{\mathcal{N}_1} \cos(x+y-a) + \underbrace{\frac{1}{2}\sqrt{(A'+D')^2 + (C'-B')^2}}_{\mathcal{N}_2} \cos(x-y-b), \\
\tan(a) &= \frac{C'+B'}{A'-D'}, \\
\tan(b) &= \frac{C'-B'}{A'+D'}.
\end{align}
$$
The maximum is obviously $\mathcal{N} = \mathcal{N}_1 + \mathcal{N}_2$ and its position is given by $x+y = a, x-y = b$.
In 3D we restrict the maximalization to a plane with a constant $z$ reducing the problem to 2D and we find that
$$
\begin{align}
A'(z) &= A\cos z + B\sin z, \\
B'(z) &= C\cos z + D\sin z, \\
C'(z) &= E\cos z + F\sin z, \\
D'(z) &= G\cos z + H\sin z.
\end{align}
$$
We know how that the maximum in this fixed plane is $\mathcal{N}(z)$, so we just have to maximize this function over $z$. This is done in an analogous way as for $h(x) + \lambda(x)$, but we end up with 6-th degree polynomial instead of 8-th.
The whole solution is implemented here - I left a link if someone wants to see it in detail.
