The C2,1 and S2,1 spherical harmonic coefficients of glacial isostatic adjustment (GIA) models are dominated by true polar wander (TPW), the secular drift of the position of the Earth’s rotational axis with respect to the Earth’s surface. The Earth’s rotational state and its interior are closely linked. The deformation of the body induced by a surface ice load depends on the local viscosity, while the adjustment of the equatorial bulge to the position perpendicular to the rotational axis depends on the average global viscosity of the body. Therefore, lateral viscosity variations could play a role in GIA-induced TPW. Here our goal is to extend an existing GIA model based on a finite element method (FEM) with a simulation for TPW (Hu et al., JGR 2017).

The model is based on FE software ABAQUS coupled to the solution of the Laplace equation. The perturbed gravitational potential is a function of the radial displacements. Therefore, an iterative process is required to solve for the displacements in the body. For the TPW algorithm we assume that during the process of TPW the equatorial bulge readjustment is fast enough such that the equatorial bulge is always nearly perpendicular to the rotational axis. In that case the linearized Liouville equation can be used by using coordinate transformations for each time step. This method allows for large-angle TPW and takes into account non-stationary surface loads with respect to the rotational axis. In the spin-up phase the flattening of the Earth is simulated by applying a centrifugal force for a long-enough duration.

We find that when the perturbed gravitational potential in the first time step is not fully converged, it affects the perturbed gravitational potential in future time steps and thus TPW. Furthermore, when a surface ice load is applied to the model, TPW is triggered. The centrifugal potential changes based on the new position of the rotational axis, and this also affects the perturbed gravitational potential for the following time steps. As a result, the perturbed gravitational potential and centrifugal potential cannot simultaneously for all time steps be iterated for using our TPW approach. Therefore, to be able to study the effect of lateral viscosity variations on GIA-induced TPW in the numerical model, the iteration between radial displacements and centrifugal potential needs to occur per time step.