2D and 3D Simulation Data for: Numerical Determination of the Nonlinear Convection Threshold in Nitrogen Ice: Application to the Margins of the Cellular Plains in Sputnik Planitia, Pluto

Numerical simulations of bottom-heated convection in a layer of nitrogen ice under temperature conditions expected on Sputnik Planitia, Pluto. Numerical calculations are performed in a 12x1 box in 2D and a 12x12x1 box in 3D, where the top and bottom boundaries are maintained at constant temperatures of T_s and T_b, respectively. The top boundary is free-slip, whereas the bottom boundary is no-slip. The vertical boundaries are free-slip and thermally insulated. In all of the calculations shown here except one, the viscosity follows the power-law creep function, i.e., it depends on both temperature and stress, and the experimentally determined flow-law parameters for power-law creep in nitrogen ice are used (Yamashita et al., 2010). The 2D calculations are conducted using the Citcom code (Moresi & Solomatov, 1995) and 3D calculations using the CitcomCU code (Moresi and Gurnis, 1996; Zhong, 2006). The numerical domain is divided into cube-shaped finite elements, and the number of elements in the vertical direction is 32. For a 12x12x1 box, this resolution translates into 4,718,592 finite elements. We investigate the stability of the localized plume over a range of conditions relevant to Sputnik Planitia. These include T_s = 32 – 40 K, T_b = 47 – 63 K, and a range of experimentally constrained values for the flow-law parameters n, E, and A for power-law creep in nitrogen ice, where n is the stress exponent (n = 2.2 +/- 0.2), E is the activation energy (E = 3.5 +/- 0.5 kJ/mol), and A is the scaling coefficient. In the 3D calculations, we focus on only three sets of permissible values of these parameters, which are: n = 2.2, E = 3.5 kJ/mol, and A = 3.74x10^(-12) s^(-1) Pa^(-n) (reference values); n = 2.0, E = 3.2 kJ/mol, and A = 2.18x10^(-11) s^(-1) Pa^(-n); and n = 2.4, E = 3.8 kJ/mol, and A = 6.41x10^(-13) s^(-1) Pa^(-n). For each. set of assumptions on T_s, T_b, n, and E, the stability boundaries [in the Rayleigh number (Ra)-space] of the localized convection planform are determined using the approach described in previous studies (e.g., Solomatov, 2024). The results of our 2D and 3D simulations are provided in the attached document. The document contains 5 sheets. In the first sheet (titled '2D-Ra_loc,c'), we provide the values of the highest Rayleigh number at which localized convection decays, denoted by Ra(-)_loc,min, the lowest and highest Rayleigh numbers - Ra(+)_loc,min and Ra(-)_loc,max, respectively - at which a localized plume is stable, and the lowest Rayleigh number - Ra(+)_loc,max - at which the localized convection planform is unstable and transitions to widespread convection for each set of assumptions on n, E, T_s, and T_b in 2D simulations. The stability boundaries for localized convection in 3D are given in the sheet title '3D-Ra_loc,c'. The nondimensional values of the surface observables for the localized convection planform output by our 2D simulations, namely the width of the localized cell, w*_h, the peak amplitude of the surface topographic anomaly, A*_h, the average surface heat flux, F*_avg, and the average surface convective velocity, U*_avg, are provided in the sheet titled '2D-nondimen-params'. The same parameters output by our 3D simulations are provided in the sheet '3D-nondimen-params'. For each surface observable, the parameter values at the lowest and highest Ra for stable localized convection are given. The formulae for deriving the corresponding dimensional values are explained in the sheet "dimensionalization".针对冥王星斯普特尼克平原（Sputnik Planitia）预期的温度条件，开展氮冰层底部加热对流的数值模拟研究。本研究分别在二维12×1计算域与三维12×12×1计算域中开展数值计算，上下边界分别维持恒定温度T_s与T_b。上边界采用自由滑移条件，下边界采用无滑移条件；侧边界采用自由滑移且绝热的边界条件。除一组计算外，本次展示的所有模拟均采用幂律蠕变模型描述黏度，即黏度同时依赖于温度与应力，并使用实验测定的氮冰幂律蠕变流动定律参数（Yamashita等，2010）。二维计算采用Citcom程序（Moresi与Solomatov，1995）完成，三维计算则采用CitcomCU程序（Moresi与Gurnis，1996；Zhong，2006）。数值计算域采用六面体有限单元离散，垂直方向单元数为32；针对12×12×1的三维计算域，该离散分辨率对应4718592个有限单元。 本研究针对与斯普特尼克平原（Sputnik Planitia）相关的一系列参数条件，探究局地羽流的稳定性。所涉及的条件包括：T_s取值范围为32~40 K，T_b取值范围为47~63 K，以及一系列由实验约束得到的氮冰幂律蠕变流动定律参数n、E与A：其中n为应力指数（n=2.2±0.2），E为活化能（E=3.5±0.5 kJ/mol），A为比例系数。在三维计算中，本研究仅聚焦于三组该参数的可行取值，分别为：① n=2.2、E=3.5 kJ/mol、A=3.74×10^(-12) s^(-1)·Pa^(-n)（参考值组）；② n=2.0、E=3.2 kJ/mol、A=2.18×10^(-11) s^(-1)·Pa^(-n)；③ n=2.4、E=3.8 kJ/mol、A=6.41×10^(-13) s^(-1)·Pa^(-n)。针对每组T_s、T_b、n与E的假设条件，本研究采用已有研究提出的方法（如Solomatov，2024），确定局地对流流型的瑞利数（Rayleigh number, Ra）空间内的稳定性边界。本研究的二维与三维模拟结果均收录于随附文档中。 该文档共包含5个工作表。首个工作表（表名："2D-Ra_loc,c"）收录了二维模拟中，每组n、E、T_s与T_b假设条件下的四类瑞利数数值：局地对流衰减的最高瑞利数（记为Ra(-)_loc,min）；局地羽流稳定的最低与最高瑞利数，分别记为Ra(+)_loc,min与Ra(-)_loc,max；以及局地对流流型失稳并转变为全域对流的最低瑞利数（记为Ra(+)_loc,max）。三维局地对流的稳定性边界收录于表名"3D-Ra_loc,c"的工作表中。二维模拟输出的局地对流流型地表可观测参数的无量纲数值收录于表名"2D-nondimen-params"的工作表中，这些参数包括局地对流元的宽度w*_h、地表地形异常的峰值振幅A*_h、平均地表热通量F*_avg以及平均地表对流速度U*_avg。三维模拟输出的同类参数则收录于表名"3D-nondimen-params"的工作表中。针对每一类地表可观测参数，文档均给出了稳定局地对流对应的最低与最高瑞利数下的参数取值。各参数对应有量纲数值的推导公式详见表名为"dimensionalization"的工作表。