The outer areas of Jupiter and Saturn have multiple zonal winds, reaching the high latitudes, that penetrate deep into the planets' interiors, as suggested by gravity measurements. These characteristics are finally replicable in numerical simulations by including both a shallow stably stratified layer, below a convecting envelope, and increasing electrical conductivity. A dipolar magnetic field, assumed to be generated by a dynamo below our model, is imposed. We find that the winds' depth into the stratified layer depends on the local product of the squared magnetic field strength and electrical conductivity. The key for the drop-off of the zonal winds is a meridional circulation which perturbs the density structure in the stable layer. In the stable region its dynamics is governed by a balance between Coriolis and electromagnetic forces. Our models suggest that a stable layer extending into weakly conducting regions could account for the observed deep zonal wind structures.

Numerical geodynamo simulations capture several features of the spatial and temporal geomagnetic field variability on historical and Holocene timescales. However, a recent analysis questioned the ability of these numerical models to comply with long-term palaeomagnetic field behaviour. Analysing a suite of 50 geodynamo models, we present here the first numerical simulations known to reproduce the salient aspects of the palaeosecular variation and time-averaged field behaviour since 10 Ma. We find that the simulated field characteristics covary with the relative dipole field strength at the core-mantle boundary (dipolarity). Only models dominantly driven by compositional convection, with an Ekman number (ratio of viscous to Coriolis forces) lower than $10^{-3}$ and a dipolarity in the range $0.34-0.56$ can capture the observed palaeomagnetic field behaviour. This dipolarity range agrees well with state-of-the-art statistical field models and represent a testable prediction for next generation global palaeomagnetic field model reconstructions.