  There is an ongoing quest to improve on the spectroscopic quality of nuclear
energy density functionals (EDFs) of the Skyrme type through extensions of its
traditional form. One direction for such activities is the inclusion of terms
of higher order in gradients in the EDF. We report on exploratory
symmetry-breaking calculations performed for an extension of the Skyrme EDF
that includes central terms with four gradients at next-to-next-to-leading
order (N2LO) and for which the high-quality parametrization SN2LO1 has been
constructed recently [P. Becker et al, Phys. Rev. C 96, 044330 (2017)]. Up to
now, the investigation of such functionals with higher-order terms was limited
to infinite matter and spherically symmetric configurations of singly- and
doubly-magic nuclei. We address here nuclei and phenomena that require us to
consider axial and non-axial deformation, both for reflection-symmetric and
also reflection-asymmetric shapes, as well as the breaking of time-reversal
invariance. Achieving these calculations demanded a number of formal
developments. These all resulted from the formulation of the N2LO EDF requiring
the introduction of new local densities with additional gradients that are not
present in the EDF at NLO. Their choice is not unique, but can differ in the
way the gradients are coupled. While designing a numerical implementation of
N2LO EDFs in Cartesian 3d coordinate-space representation, we have developed a
novel definition and a new unifying notation for normal and pair densities that
contain gradients at arbitrary order. The resulting scheme resolves several
issues with some of the choices that have been made for local densities in the
past, in particular when breaking time-reversal symmetry. Guided by general
practical considerations, we propose an alternative form of the N2LO
contribution to the Skyrme EDF that is built from a different set of densities.
