The problem of solving `4^x = x^4` was just posed on the Wikipedia reference
desk. I posted there some code for finding all the solutions with mpmath.

With SymPy, one can find one solution like this:
```
solve(4**x-x**4, x)
[-4_LambertW(-1/4_log(4))/log(4)]
```
Four more solutions are:
```
p = -log(4)/4
(LambertW(p)/p).evalf()
2.00000000000000
(LambertW(-p)/p).evalf()
-0.766664695962123
(LambertW(I_p)/p).evalf()
-0.270279229207913 + 0.869544816180665_I
(LambertW(-I_p)/p).evalf()
-0.270279229207913 – 0.869544816180665_I
```
Then, one obtains infinitely many solutions by replacing each `LambertW(x)`
above with `LambertW(x,k) for k = …-2,-1,0,1,2,…` (SymPy needs to be
fixed to support the two-arg LambertW).

What solve needs to do is, to begin with, to include all four branches when
computing `x^(1/4)` in tsolve.

Secondly, it would be nice it could return a solution containing the free
parameter when one exists; something like `Set(LambertW(-p,C1)/p,
is_integer(C1))`.