Solutions of variational PDEs can often be interpreted as local minima of an associated energy functional (local minima being special instances of critical points) arising in applications such as physics and mechanics. Classical approaches typically focus solely on discretizations of the underlying PDEs. By contrast, we present a perspective that explicitly incorporates the underlying energy structure into the numerical approximation, thereby naturally exploiting more of the available structure of the problem. Specifically, we study linearized iterative Galerkin schemes that guarantee energy reduction at every iteration and give rise to an energy-based paradigm for adaptive finite element discretizations, called variational adaptivity. Applications to semilinear elliptic PDEs from ecology and quantum chemistry illustrate the proposed approach.