  Invariant solutions of the Navier-Stokes equations play an important role in
the spatiotemporally chaotic dynamics of turbulent shear flows. Despite the
significance of these solutions, their identification remains a computational
challenge, rendering many solutions inaccessible and thus hindering progress
towards a dynamical description of turbulence in terms of invariant solutions.
We compute equilibria of three-dimensional wall-bounded shear flows using an
adjoint-based matrix-free variational approach. To address the challenge of
computing pressure in the presence of solid walls, we develop a formulation
that circumvents the explicit construction of pressure and instead employs the
influence matrix method. Together with a data-driven convergence acceleration
technique based on dynamic mode decomposition, this yields a practically
feasible alternative to state-of-the-art Newton methods for converging
equilibrium solutions. We compute multiple equilibria of plane Couette flow
starting from inaccurate guesses extracted from a turbulent time series. The
variational method outperforms Newton(-hookstep) iterations in successfully
converging from poor initial guesses, suggesting a larger convergence radius.
