We develop a computational method for extremal Steklov eigenvalue problems and apply it to study the problem of maximizing the pth Steklov eigenvalue as a function of the domain with a volume constraint. In contrast to the optimal domains for several other extremal Dirichlet- and Neumann-Laplacian eigenvalue problems, computational results suggest that the optimal domains for this problem are very structured. We reach the conjecture that the domain maximizing the pth Steklov eigenvalue is unique (up to dilations and rigid transformations), has p-fold symmetry, and has at least one axis of symmetry. The pth Steklov eigenvalue has multiplicity 2 if p is even and multiplicity 3 if p >= 3 is odd.