  Tobasco et al. [Physics Letters A, 382:382-386, 2018; see
https://doi.org/10.1016/j.physleta.2017.12.023] recently suggested that
trajectories of ODE systems that optimize the infinite-time average of a
certain observable can be localized using sublevel sets of a function that
arise when bounding such averages using so-called auxiliary functions. In this
paper we demonstrate that this idea is viable and allows for the computation of
extremal unstable periodic orbits (UPOs) for polynomial ODE systems. First, we
prove that polynomial optimization is guaranteed to produce auxiliary functions
that yield near-sharp bounds on time averages, which is required in order to
localize the extremal orbit accurately. Second, we show that points inside the
relevant sublevel sets can be computed efficiently through direct nonlinear
optimization. Such points provide good initial conditions for UPO computations.
As a proof of concept, we then combine these methods with a single-shooting
Newton-Raphson algorithm to study extremal UPOs for a nine-dimensional model of
sinusoidally forced shear flow. We discover three previously unknown families
of UPOs, one of which simultaneously minimizes the mean energy dissipation rate
and maximizes the mean perturbation energy relative to the laminar state for
Reynolds numbers approximately between 81.24 and 125.
