Viscous relaxation of Pluto&#8217;s clathrate-insulated ice shell below Sputnik Planitia
Authors: Kihoulou Martin , Kalousová Klára , Souček Ondřej , Čadek Ondřej 
Published: July 2021

AbstractSputnik Planitia, a 1000 km wide, Pluto-dominating feature, is located very close to Pluto-Charon tidal axis. To explain its position, a reorientation driven by a postimpact uplift of a subsurface ocean was proposed. Since pure water ice shell relaxes too quickly to maintain the reorientation up to the present, an insulating layer of high viscosity clathrates at the ice/ocean interface was proposed. Here, we solve Pluto&#8217;s ice shell evolution in a 2D spherical axisymmetric geometry with an evolving free surface, assuming a viscous rheology of both ice and clathrates. Our results show that the thermal effect of the impact and nonlinear rheology substantially decrease the relaxation timescales.&#160;IntroductionSputnik Planitia basin is situated only 400 km from the tidal axis. Since a position this close to the tidal axis is very unlikely, reorientation of the whole body was proposed [1]. Facing towards the tidal axis, Sputnik Planitia should be a&#160;positive gravity anomaly [1]. Pluto is believed to be differentiated into a rocky core and a&#160;hydrosphere [2], therefore isostatic uplift of a&#160;liquid ocean might have compensated the basin&#8217;s negative topography. Impact origin ejecta blanket and accumulated nitrogen would then provide an additional positive mass [3]. Due to the basin&#8217;s likely age of 4 Gyr [4], the gravity anomaly has to be positive up to the present. Since the warm ice close to the melting point is expected to relax quickly, Kamata et al. [5] proposed the presence of an insulating layer of high viscosity clathrates. Here, we study the effect of impact heating [6] and nonlinear ice and clathrate rheology [7] on the relaxation timescale.&#160;Numerical modelWe solve the thermal evolution of Pluto&#8217;s ice shell in a 2D spherical axisymmetric geometry with time-evolving ice/ocean interface. The initial height of the isostatic uplift is evaluated assuming Airy isostasy, a 10 km deep basin at the surface and ocean water density of 1100 kg/m3 [8]. Free surface is prescribed at the bottom boundary and free slip is prescribed elsewhere. Temperature is fixed at the top and bottom boundary (40 and 265 K, respectively) and the side boundaries are kept adiabatic. The initial condition is steady-state conductive profile assuming clathrate layer of thickness 5 or 10&#160;km. We use nonlinear viscosity of ice dependent on temperature, grain size and stress using the composite law [7], similarly&#160;for clathrates [9]. We prescribe temperature dependent thermal conductivity for ice [10] and a constant value (0.6 W/m/K) in the clathrate layer. We choose impactor velocity 4&#160;km/s, diameter 400 km and density 920 kg/m3 in order to create a 900 km wide circular basin. The resulting temperature anomaly in the isobaric core is &#8764; 80 K [6]. We use the Arbitrary Lagrangian-Eulerian formulation to address the moving boundary [11]. The problem is implemented in a freely available Finite Element Method library FEniCS [12].&#160;ResultsWe performed series of simulations combining different shell thicknesses H and clathrate layer thicknesses hc. Here, we present two combinations, which, according to Kamata et al. [5], should lead to the uplift relaxation times in the order of billion years. Models 1 and 2 represent the H/hc = 100/5 km and 200/10 km simulations, respectively, with parameters used in [5]. Models 3 and 4 are analogical to Models 1 and 2, however, with nonlinear rheology and impact heating. Note that Kamata et al. [5] define relaxation timescale as the time when the uplift volume decreases to 1/e of the initial volume, whereas we present the time when the uplift starts to subside rapidly, i.e., when it no longer supports the surface basin isostatically, which is crucial for reorientation. This timescale turns out to be substantially shorter than the timescale presented in [5].Figure 1 shows the relaxation curves, i.e., uplift height at the axis of symmetry with respect to time. We can see that Models 3 and 4 relax faster than Models 1 and 2. At 1 Gyr, the uplift of Model 4 already subsided by a few km and relaxation of Model 3 is almost over. These preliminary results suggest that none of them can compensate the negative basin topography anymore.&#160;&#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160; &#160;Figure 2 shows the moment when the uplift starts to subside rapidly for Models 3 and 4. Although the thin shell becomes quickly cold and rigid (Fig. 2a), its uplift starts to relax at t