Scattering of Giant Planets and Implications for the Origin of the Hierarchical and Eccentric Two-planet System GJ 1148

The GJ 1148 system has two Saturn-mass planets orbiting around an M dwarf star on hierarchical and eccentric orbits, with orbital period ratio of 13 and eccentricities of both planets of 0.375. The inner planet is in the regime of eccentric warm Jupiters. We perform numerical experiments to study the planet–planet scattering scenario for the origin of this orbital architecture. We consider a third planet of 0.1M J (Jupiter's mass) in the initial GJ 1148 system with initial orbital separations of 3.5, 4, and 4.5 mutual Hill radii and initial semimajor axis of the innermost planet in the range of 0.10–0.50 au. The majority of scattering results in planet–planet collisions, followed by planet ejections, and planet–star close approaches. Among them, only planet ejections produce eccentric and widely separated two-planet systems, with some having similar orbital properties to the GJ 1148 system. We also examine the effects of general relativistic apsidal precession and a higher mass of 0.227M J for the third planet. The simulation results suggest that the GJ 1148 system may have lost a giant planet. We also perform simulations of the general problem of the origin of warm Jupiters by planet–planet scattering. As in the GJ 1148 simulations, a nontrivial number of stable two-planet systems are produced by ejection, which disagrees with the result from a previous study showing that two-planet systems arise exclusively through planet–planet collisions.

1. INTRODUCTION M dwarf stars have been a primary target in the search for exoplanets in recent years due to the large number of them in the galaxy.Their lower masses and luminosities compared to Sun-like stars make the planets in the so-called habitable zone (HZ) easier to detect with larger radial-velocity (RV) signal and shorter orbital period (Trifonov et al. 2020;Sabotta et al. 2021;González-Álvarez et al. 2022).Current observations show that most planets orbiting around M-dwarfs have masses equivalent to a few Earth to Neptune masses (Trifonov et al. 2020;González-Álvarez et al. 2022), while Jupiter-mass planets are observed more frequently around Sun-like or more massive stars (Fischer & Valenti 2005;Reffert et al. 2015).GJ 1148, GJ 876 (Rivera et al. 2010;Nelson et al. 2016;Millholland et al. 2018), and GJ 317 (Johnson et al. 2007) are three unusual M-dwarf systems that have massive planets, whose formation is difficult to explain by the core accretion theory of planet formation (Laughlin et al. 2004;Morales et al. 2019).Alternative formation mechanisms for massive planets around M-dwarfs may be needed (Schlecker et al. 2022).Left panel : Orbital evolution of the best fit of the GJ 1148 system over 1 Myr (Trifonov et al. 2020).Right panel : Phase-space diagrams of the planetary eccentricities e b and e c versus the apsidal alignment angle ∆ϖ = ϖ b − ϖ c for systems with the same total angular momentum as the GJ 1148 system.Blue crosses indicate the best-fit parameters of e b , e c and ϖ b − ϖ c and their uncertainties, while the red paths represent the trajectories starting at the best fit.
GJ 1148 is an M-type dwarf hosting two Saturn-mass planets (Trifonov et al. 2018(Trifonov et al. , 2020)).The star has mass M * = 0.354 ± 0.015 M ⊙ , effective temperature T eff = 3358 ± 51 K, and luminosity L = 0.0143 ± 0.0003 L ⊙ (Passegger et al. 2018;Schweitzer et al. 2019).The planet GJ 1148 b, with minimum mass m b sin i = 0.304M J (where M J is Jupiter mass), semimajor axis a b = 0.166 au, and orbital period of around 41.4 days, was discovered based on 37 HIRES RV measurements by Haghighipour et al. (2010).With an additional 52 RV data from CARMENES, Trifonov et al. (2018) found the second massive planet GJ 1148 c with a minimum mass of m b sin i = 0.227M J located in a relatively distant orbit with semimajor axis a c = 0.910 au and orbital period of 532.6 days.The period ratio between the two planets is around 13, and both have orbital eccentricities of around 0.375 at the current epoch.Long-term integration shows that the two planets have large eccentricity variations on a secular timescale of around 72, 000 yr, as shown in Figure 1.
The inner planet GJ 1148 b is in the regime of warm Jupiters (WJs), which are commonly defined as Jupiter-size planets with orbital periods of 10-200 days (e.g., Dawson & Johnson 2018). 1 Observations show that WJs can be classified into two populations (Sun et al. 2022) according to their eccentricity and companionship, whose origins are still hotly debated.The WJs with low eccentricity and inner super-Earth companions are suggested to form in situ (Huang et al. 2016).The WJs with moderate eccentricity and Jovian companion are proposed to form through high-eccentricity migration (Dawson & Chiang 2014;Dong et al. 2014;Petrovich & Tremaine 2016).However, Antonini et al. (2016) found that most observed WJs with characterized external giant planet companions cannot be formed in this way.Another source of eccentricity excitation is planet-planet scattering (Rasio & Ford 1996;Weidenschilling & Marzari 1996;Marzari & Weidenschilling 2002;Chatterjee et al. 2008;Jurić & Tremaine 2008;Raymond et al. 2010).Anderson et al. (2020) studied a scenario of in situ scatterings.All of their two-planet systems arise from planet-planet collision events, which produce too many low-eccentricity and closely packed two-planet systems.The dearth of stable two-planet systems from ejection events in Anderson et al. (2020) is very surprising and even doubtful, because previous studies of three-planet scatterings at several au have shown that an ejection event could produce two planets on stable and well-separated orbits (Marzari & Weidenschilling 2002;Chatterjee et al. 2008), making planet ejection a plausible way to produce systems like GJ 1148.
In this study, we aim to test whether the current GJ 1148 system could have evolved from a tightly packed system of three Saturn-mass planets by losing one of the planets through planetplanet scattering.Previously scattering studies always assumed a Sun-like star with Jupiter-mass planets.Our study explores a new parameter space: Saturn-mass planets around an M-dwarf star, although the planet-star mass ratio is comparable to Jupiter around a Sun-like star.Studying the origin of such a system could provide important implications for the formation of massive planets around M-dwarfs, i.e., it would give the possible parameter space for the configuration of the original GJ 1148 system just after the gas disk dispersal.The success in reproducing the GJ 1148 system with a planet-planet scattering scenario also gives important implications for WJs with a distant companion.
This paper is organized as follows.In Section 2, we introduce the setups for our N -body simulations.In Section 3, we show the outcomes of scattering, illustrate the properties of two-planet systems from either planet-planet collision or ejection events, and determine the optimal original configuration of the GJ 1148 system.In Section 4, we discuss the influence of general relativistic (GR) apsidal precession and the mass of the ejected planet, as well as the discrepancy with the simulations of in situ scattering of WJs by Anderson et al. (2020).We end with conclusions in Section 5.

INITIAL CONDITIONS AND PARAMETERS OF SIMULATIONS
We assume that the original GJ 1148 system includes an M-dwarf star (M ⋆ = 0.354 M ⊙ ) and three giant planets.Two planets have the same mass as GJ 1148 b (0.304 M J ) and GJ 1148 c (0.227 M J ).The third planet is considered a less-massive planet with 0.1 M J .Such a mass choice is consistent with previous scattering studies (Chatterjee et al. 2008;Anderson et al. 2020), which indicate that the planet with the lowest mass is most likely to be ejected.The planets are randomly ordered.Considering all planets are young, the planetary radius is set as 1.0 R J .The initial conditions of planet orbits are similar to Anderson et al. (2020).For each planet, we sample the argument of periapse, longitude of the ascending node, and mean anomaly in [0  ], in a uniform random manner.The planets are separated by K in units of mutual Hill radius, which is defined as To study the influence of initial separation between planets, we perform simulations in three groups with K equal to 3.5, 4, and 4.5.In each group, we firstly perform N total = 200 trial N -body simulations (Sim.trial ) with the initial innermost planet's semimajor axis a i,1 in the range of [0.1, 0.5] au to get a preliminary knowledge of the statistics of scattering outcomes (see Table 1).Subsequently, two additional sets of simulations (Sim. 1 and Sim. 2 ), are performed; they are also listed in Table 1 and will be discussed in Sections 3.3 and 3.4.We integrate the planetary orbits with the symplectic integrator SyMBA (Duncan et al. 1998), which can handle close encounters between the planets.To ensure that SyMBA can adequately integrate the close orbit of the innermost planet, we set the time step as 3.5 ×10 −4 yr.The initial three-planet system is dynamically unstable.As the system evolves, the growing instability can lead to close encounters between planets.During a close encounter between two planets, the dynamical outcome depends on the relative magnitude between the escape velocity V esc from the planet's surface and its orbital velocity V orb , which is quantified by the Safronov number (Morbidelli 2013): where a p is the semimajor axis, R p is the radius of the planet, and M p and M * are the mass of the planet and the host star, respectively.If Θ is much less than unity, planet-planet collisions will dominate the dynamical outcomes.If Θ is much greater than unity, planet ejections are favored.
Our initial planets have Θ smaller than unity (with Θ = 0.057-0.877),and thus primarily result in planet-planet collisions, followed by planet ejections and planet-star collisions.When the distance between two planets is less than the sum of their radii, they are merged together to form a new planet with the sum of the masses.A planet is ejected and removed when its distance from the host M-dwarf exceeds 1000 au.We also set an inner boundary at 0.02 au to remove planets that come too close to the star.The removed planet would probably result in a collision with the star or survive in an eccentric orbit with an extremely small periastron distance, in which additional physical processes like tidal effect or mass loss would be nonnegligible.Each three-planet system is firstly integrated with a maximum time 2π × 10 6 P i (except for systems initially separated by K = 4.5, for which we set the maximum time as 2π × 10 7 P i ), with P i the initial orbital period of the initial innermost planet.We refer to this highly chaotic phase as Phase 1.The instability often leads to orbital crossings and close encounters and results in a combination of collision, ejection, and removal at the inner boundary.Most systems end up with two planets left.For Sim. trial, we integrate systems with more than one planet for an additional 2π × 10 8 P in , with P in the inner planet's period at the end of Phase 1.This long-term integration is referred to as Phase 2. For Sim. 1 and Sim. 2, only some of the two-planet systems are integrated to the end of Phase 2 (see Table 1).
Figure 2. Dynamical evolution of a hypothetical three-planet system with K = 4.0 that evolves to a twoplanet system.Red, black, and blue colors denote GJ 1148 b, c, and the third planet, respectively.The hypothetical third planet is ejected at around 4 × 10 3 yr.

Scattering outcomes from Sim. trial
The instabilities in initially unstable three-planet systems lead to close encounters, which result in planet loss by planet-planet collision (pp), planet ejection (ej), and planet-star close approach (ps).We first made a trial simulation (Sim.trial ) to get a preliminary knowledge of the evolutionary outcomes of the unstable three-planet systems.In each group of K = 3.5, 4.0, and 4.5, we performed 200 simulations with time covering both Phase 1 and Phase 2. Almost all of our three-planet systems lose at least one planet at the end of Phase 2, (except for five systems in the simulation set with K = 4.5, as shown in Table 2) and evolve into one-planet or two-planet systems at the end of Phase 2.
Figure 2 shows the orbital evolution of a system with K = 4.0, ending with two planets in stable orbits.The planets in the initial system were GJ 1148 b, c, and the hypothetical third planet of 0.1 M J in order of semimajor axis.The innermost planet's semimajor axis a i is around 0.4 au.After about 1.3 × 10 3 yr, there is a chaotic phase in which orbital crossings repeatedly occur until the third planet is ejected.The original inner planet is scattered closer to the star, while the original middle one is left in a larger orbit.Their orbits are stable, and no orbital crossing occurs.The lower panel shows the time evolution of the eccentricities.We can see that the eccentricities change with small amplitudes before 1.3 × 10 3 yr.During the chaotic phase, the eccentricities of all planets are significantly excited.The eccentricity of the third planet increases to exceed unity and the planet ends up being ejected.In the final stable phase, we can see that the eccentricities of the two planets vary on a secular time scale due to their gravitational interactions.The final configuration with two planets that are widely separated and show moderately large eccentricity variations is qualitatively similar to the GJ 1148 system.
Figure 3 shows the time evolution of a system that loses two planets.The three planets start to evolve chaotically after a period without orbital crossing.During the chaotic phase, the planetary  eccentricities and semimajor axes change significantly.At a time of ≃ 3×10 3 yr, the original outermost planet collides with the original innermost planet, merging into a new planet in a slightly smaller orbit than the initial innermost orbit.The original middle planet is scattered to a larger eccentric orbit at the end of the chaotic phase.The orbits then evolve slowly until ≃ 2 × 10 5 yr, when the semimajor axis and eccentricity of the outer planet suddenly increase, while the orbit of the inner planet changes more mildly.The outer planet is finally ejected at a time ≃ 7 × 10 5 yr, making the inner planet's eccentricity larger and its semimajor axis smaller.
Figure 4 shows the branching ratio of how initially unstable three-planet systems evolve into twoplanet and one-planet systems for K = 4.0 simulations.The most common events that result in Note-The first two columns are the initial separation between the planets in units of mutual Hill radii (K) and the number of simulations (N total ).The quantities N pp,2p , N ps,2p , N ej,2p are the numbers of two-planet (2p) systems from planet-planet collision (pp), planet-star close approach (ps), and planet ejection (ej), respectively.The quantities N (1) pp,1p , N (1) ps,1p , and ps,1p , and N (2) ej,1p are similar, but for the loss of the first (superscript 1) and second (superscript 2) planets of one-planet (1p) systems.The last three columns are the numbers of final one-planet (N F, 1P ), two-planet (N F, 2P ), and three-planet (N F, 3P ) systems.
the first planet loss are planet-planet collisions (N pp,2p + N (1) pp,1p = 160), followed by planet ejections (N ej,2p + N (1) ej,1p = 27) and planet-star close approach (N ps,2p + N (1) ps,1p = 13).The scattering outcomes dominated by collisions are consistent with the small value of the Safronov number Θ discussed in Section 2. For the second-planet loss, the most common events are planet ejections (N  2. As expected, the K = 3.5 simulations produce the highest number of planet-planet collisions since the three planets' mean Θ value is the smallest.We also notice that both K = 3.5 and 4 simulations do not have three-planet systems left at the end of integration, but five three-plant systems are left in K = 4.5 simulations.This is also not surprising because dynamical instabilities take a long time to grow for a large initial separation of K = 4.5.In addition, we observe that the fraction of final one-planet systems seems to decrease with K.It indicates that systems with larger K are less likely to experience second-planet loss. The two-planet systems that follow a planet-star close approach (ps) are not studied here because the planet is removed artificially, and the properties of the remaining two planets are not realistic.We focus on the remaining two-planet systems following an ejection (post-ejection two-planet system) or a collision event (post-collision two-planet system).We will show their properties and compare them with the observation of the GJ 1148 system in the following.

Properties of post-collision two-planet systems from Sim. trial
Figure 5 shows another view of the properties of the post-collision two-planet systems for simulations with initial spacing K = 3.5, 4.0, and 4.5.The upper panels show the semimajor axes of the remaining inner planets a f, i (red) and outer planets a f, o (black) versus a i, 1 .The semimajor axes of the initial innermost a i, 1 (pink) and outermost a i, 3 (grey) planets are included for comparison.We can see that a linear relationship can fit the semimajor axis of the remaining inner planet (a f, i ) versus the semimajor axis of the initial innermost planet (a i, 1 ) for all three sets of simulations very well.The semimajor axis a f, i (red circles) shows a small scatter around a i,1 (pink circles), and the slopes of the three linear models from fitting a f, i versus a i,1 are close to unity, showing that the position of the inner planet almost remains the same.The semimajor axis a f, o seems to have a wider range of variations, but it has a tendency to increase with a i,1 and is generally close to the semimajor axis of the initial outermost planet a i,3 .So the positions of the remaining planets do not change much compared to the innermost and outermost planets of the initial three-planet systems.Compared to the observation, we find that the post-collision two-planet systems with the inner planet's semimajor axis around 0.166 au could reproduce the semimajor axis of GJ 1148 b.However, the orbit of the produced outer planet is too small to match the semimajor axis of GJ 1148 c.
The middle panels of Figure 5 show the eccentricities (e f, i and e f, o ) of the remaining two planets versus a i,1 .The horizontal grey line shows the observed eccentricity of 0.375 for both planets (Trifonov et al. 2020).The red and black dashed lines represent the 90th percentiles of e f, i and e f, o , respectively, which are lower than the grey line, showing that post-collision two-planet systems have lower eccentricities compared to the observation.The property of low eccentricities for the planets after a collision is consistent with previous studies (Ford et al. 2001;Anderson et al. 2020).Never-theless, the eccentricities of GJ 1148 b and c are not constant but vary with large amplitude due to their secular interaction with each other.So, low eccentricity is allowed, but only for one of the planets.However, only few post-collision two-planet systems have either e f, i or e f, o exceeding 0.375.In addition, we observe that, unlike the semimajor axis, there is no evident dependence between the final eccentricities and initial a i,1 .
The lower panels of Figure 5 show the inner and outer planets' inclinations.For K = 3.5 simulations, 90% of inner (outer) planets have orbital inclination below 10 • (6 • ).The orbital inclination of the inner (outer) planets in K = 4.0 simulations has 90th percentiles around 6 • (5 • ).In K = 4.5 simulations, the orbital inclinations of both inner and outer planets have 90th percentiles near 5 • , which are slightly lower than those of K = 4.0 simulations.We can see that post-collision two-planet systems generally have low orbital inclinations.
In summary, Figure 5 shows that post-collision two-planet systems have smaller orbital spacing and eccentricities, so they cannot produce systems like GJ 1148.
3.3.Properties of post-ejection two-planet systems from Sim. trial and Sim. 1 The fraction of post-ejection two-planet systems in Sim.trial only constitutes a small percentage of the remaining two-planet systems.For example, in simulations with K = 4.0, it accounts for approximately 11% (refer to Table 2).Therefore, the fraction with the same planetary mass order as the GJ 1148 system is even smaller.To get more statistics about post-ejection two-planet systems, we performed an additional 1000 simulations (Sim. 1 ) and only fully integrated the two-planet systems arising from the ejection of the hypothetical third planet (Table 1).
Figure 6 shows the orbital properties of the planets in the post-ejection two-planet systems with the same mass order as the GJ 1148 system.The upper panels show that, in each group, a f, i versus a i,1 can be well fit by a linear relationship, similar to the post-collision cases, but the slope of the linear model is smaller than the post-collision cases.The red circles are lower than the pink circles, showing that a f, i is smaller than a i,1 and that ejections push the remaining inner planet to a smaller orbit than the initial innermost planet.In addition, we can see the slope of the linear model slightly increases with K, indicating that the broader initial systems produce inner planets in orbit closer to the initial inner planets.The semimajor axis a f, o also seems to have a linear relationship with a i,1 but with large scatter.The semimajor axis a f, o (black circles) is generally above a i,3 (grey circles), suggesting that ejection results in the outer planet being in a larger orbit than the initial outermost planet.So, in general, ejection would cause one planet to orbit closer to the star and the other to move to a larger orbit.Comparing simulations with different K, we see no major differences in orbital properties.However, at a i,1 ≈ 0.21 au, which is where a f,i can match the semimajor axis of GJ 1148 b, only K = 4.0 and 4.5 simulations have several cases in which a f,o can simultaneously match the semimajor axis of GJ 1148 c.
The middle panels show the eccentricities of the remaining two planets.For all cases, we can see that the distributions of eccentricities do not have an evident correlation with a i,1 , similar to the post-collision cases.The eccentricities are generally higher than the post-collision cases (with 90th percentiles of eccentricities for both planets exceeding 0.375), which indicates that ejections have significantly excited the eccentricities of the planets.In addition, the 90th percentile of the eccentricity of the inner planet is larger than that of the outer planet in the simulations of K = 3.5, while they are nearly equal for K = 4 and 4.5.The bottom panels show the inclinations of the two planets following ejections.The inclinations for both planets are higher than those of the post- collision cases in Figure 5, suggesting that ejection also excites inclination.The inclination of the inner planet is higher than that of the outer planet for all K.
Planet-planet collisions produce two-planet systems with small separation and low eccentricities, thus unlikely to reproduce the GJ 1148 system.In contrast, planet ejection tends to produce eccentric and widely separated two-planet systems.The semimajor axis a f,i correlates with a i,1 , and it can reproduce a b when a i,1 is around 0.21 au.Nevertheless, matching a c is not easy because a f,o varies over a wide range.We expect systems like GJ 1148 could be discovered by searching through more post-ejection two-planet systems with a i,1 close to 0.21 au.  and 6).Right panel: the period ratio P f, o /P f, i vs. a i,1 .The shaded area represents the range in which the difference between the ratio of planetary periods to that of GJ 1148 does not exceed 1.The systems in this area are selected and marked yellow.

Searching for GJ 1148-like two-planet systems from Sim. 1 and Sim. 2
To search for initial configurations of three-planet systems that can produce stable two-planet systems comparable to the GJ 1148 system, we perform Sim. 2, which includes 800 simulations with K = 4.5 and a i,1 in a small range of [0.15, 0.3] au and keeps other setups the same.We only fully integrate the systems in which one planet is ejected, the remaining planet with the mass of GJ 1148 b is in the inner orbit, and the planet with the mass of GJ 1148 c is in the outer orbit.We select systems (including those with a i,1 in the range of [0.15, 0.3] au from Sim. 1 ) close to the period ratio (or spacing) of the GJ 1148 system.Then we analyze each system's secular eccentricity variation and apsidal alignment angle ∆ϖ to select the optimal systems.
Figure 7 shows the properties of the post-ejection two-planet systems with a i,1 in the range of [0.15, 0.3] au.The left panel shows a f, o and a f, i compared to the observation, in which we can see several systems with semimajor axes close to both GJ 1148 planets.The right panel shows the period ratio P f, o /P f, i versus a i,1 .We find five cases (yellow squares) whose period ratio is within ±1 of GJ 1148's period ratio.Three of them have a i,1 > 0.25 au, which produces two-planet systems with an overall distance farther from the central star than the GJ 1148 planets.We observe that about 24% of the systems shown in Figure 7 can reach period ratio ≳ 10.If the initial system has a more significant separation, such as 5R m, H , more widely spaced systems close to the GJ 1148 system probably could arise, but the instability would take longer to grow, making the numerical simulations impractical.
We analyze the five selected planetary systems by comparing their eccentricity and apsidal alignment angle variations to the GJ 1148 system in Figure 8.All of the systems are qualitatively similar to GJ 1148 (right panel of Figure 1) in showing significant eccentricity variations and ∆ϖ librating around 0 • or circulating near the boundary between libration and circulation.We observe that there are three systems with ∆ϖ circulating and two with ∆ϖ librating around 0 • .The green one has  semi-amplitude of ∆ϖ around 60 • and eccentricity variations close to the best-fit GJ 1148 system.We regard this system as the most similar one to the observation found in our scattering experiments.
Figure 9 shows the evolution of this system over a million years, which is similar to the evolution of the best-fit GJ 1148 system shown in the left panel of Figure 1.The semimajor axes of the two planets in this hypothetical system do not have significant variation and are generally close to the GJ 1148 system.The semimajor axis of the inner planet is about 0.201 au, which is 0.035 au larger than GJ 1148 b.The mean semimajor axis of the outer planet is about 1.085 au, which is 0.175 au larger than that of the GJ 1148 c.Therefore, the orbital spacing between the two planets in this system is slightly larger than that of the GJ 1148 system.In post-ejection two-planet systems, planets typically exhibit larger inclinations.For instance, the 90th percentile of the inner planet's inclination is approximately 24 • , while that of the outer planet's inclination is roughly 15 • .Nevertheless, the planets in this system have inclinations smaller than 13 • .The eccentricities of the two planets vary with large amplitudes, with e f, i from 0.06 to 0.45 and e f, o from 0.3 to 0.44, similar to the GJ 1148 planets.We can see that the eccentricities of the two planets are anticorrelated.When e f, i reaches its maximum value, e f, o , reaches its minimum value.The timescale of the eccentricity variations is about 95, 000 yr, which is about 23, 000 yr larger than GJ 1148.The longer timescale is due to the larger orbit semimajor axes of this system.As suggested by Equation ( 40) in Lee & Peale (2003), the timescale for secular evolution of eccentricity is proportional to where α = a 1 /a 2 and n 1 is the mean motion of the inner planet.The α of this system is close to the GJ 1148 system, but n 1 is smaller than that of GJ 1148 b, making the secular timescale t e longer.

The influence of GR apsidal precession
GR apsidal precession has been incorporated in N -body codes such as REBOUND (Tamayo et al. 2020), which has been used to study the origin of warm Jupiters (Anderson et al. 2020).However, the effects of GR apsidal precession on planet-planet scattering outcomes have been rarely explored.
Here we want to explore whether incorporating GR apsidal precession can produce post-ejection two-planet systems with different properties than without it.We have modified the SyMBA code to include GR apsidal precession by using the potential of Nobili & Roxburgh (1986). 2 The kick to the linear momentum of each planet due to this potential is implemented in barycentric coordinates to conserve total angular momentum.We use the same initial conditions as the K = 4.5 simulations (denoted as K4.5) but with GR apsidal precession.We perform a total of 800 simulations (denoted as K4.5-GR), the same number as the K4.5 simulations in Sim. 2. As in Section 3.4, we only analyze those two-planet systems produced by planet ejection, with GJ 1148 b in an inner orbit and GJ 1148 c in an outer orbit.At the end of Phase 2, only 25 systems remain stable, which is consistent with 31 out of 800 in Sim. 2 of K4.5 within statistical uncertainty.
Figure 10 (a) shows the cumulative distributions of the properties of the inner planet for the simulations K4.5 and K4.5-GR.The distributions of a f, i (upper left panel), e f, i (lower left panel), pericenter distance a f, i (1 − e f, i ) (upper right panel), and inclination I f, i (lower right panel) are generally close to each other, but the inclination distribution for K4.5 is higher than that for K4.5-GR when I f, i ≲ 20 • .Figure 10    the outer planet for the simulations K4.5 and K4.5-GR.The distributions are again close to each other.However, the lower right panel shows that the inclination distribution for K4.5 is lower than that for K4.5-GR when I f, o ≳ 3 • .To quantitatively determine whether K4.5 and K4.5-GR produce significantly different results, we have performed the two-sample Kolmogorov-Smirnov (KS) test on the orbital elements of both cases.Suppose a sample has size m with cumulative distribution F (x), and the other sample has size n with cumulative distribution G(x).D m,n is the maximum difference between F (x) and G(x).The null hypothesis H 0 is that the two samples have the same distribution.H 0 should be rejected if D m,n > D crit , where D crit is the critical value defined as where α is the significance level and c(α) is the inverse of the Kolmogorov distribution at α.When we select α as 0.05, the corresponding D crit is 0.365.The calculated D m,n values for the orbital elements of the inner and outer planets are listed in Table 3, in which we find no D m,n larger than D crit .Thus, the KS test suggests no significant differences between K4.5 and K4.5-GR at the significance level of 0.05.

The influence of the mass of the hypothetical third planet
We assume that the original GJ 1148 system has lost a hypothetical third planet via planet ejection in its history.Since the planet with the lowest mass is most likely the ejected one during close encounters (Chatterjee et al. 2008), we assume that the mass of the third planet is not larger than the mass of the less-massive planet (GJ 1148 c with a mass of 0.227M J ) in the current GJ 1148 system.We have so far performed scattering experiments where we assume a third planet of mass 0.1M J .
Here we explore the situation where the lost planet has a mass of 0.227M J , the same as GJ 1148 c.The other initial conditions and the integration schemes remain the same.Since the properties of the post-collision two-planet systems cannot account for the GJ 1148 system, we only analyze the post-ejection two-planet systems.We only fully integrate the post-ejection two-planet systems with the same mass order as the GJ 1148 system.To explore the influence of initial separation on scattering results, we perform three groups of simulations with K = 3.5, 4.0, and 4.5, with 1400 simulations in each group.
Figure 11 shows the properties of two-planet systems arising from three-planet systems with the initial separation of 3.5, 4.0, and 4.5 R m, H by ejecting a 0.227M J planet.The upper panels show a f, i and a f, o versus a i, 1 .We can see that a f, i versus a i,1 can be well-fitted by a linear relationship, and the slopes are similar for simulations with different K.However, the slopes are smaller than the simulations assuming a hypothetical 0.1M J planet in the original GJ 1148 system.This suggests that the inner planet can move to a smaller orbit if the ejected planet is more massive.Such a result is expected because the system's energy is conserved during the whole process (Marzari & Weidenschilling 2002).The total energy of the initial system is If the hypothetical third planet is more massive, then the absolute value of total energy in the initial system is greater.Once the third planet is ejected, its orbital energy is negligible, and this energy is redistributed to the remaining two planets.The outer planet always moves to a distant orbit, and its energy gets smaller.Thus, the inner planet can move to a smaller orbit and dominate the energy of the final two-planet system.The upper panels of Figure 11 show that systems with a i,1 close to 0.29 au can reproduce the observed semimajor axis (0.166 au) of the inner planet of GJ 1148, which is larger than a i,1 ≈ 0.21 au in the simulations that lost a 0.1 M J planet.The semimajor axis of the outer planet in the final system also seems to increase as a i,1 increases, but their range of values is extensive, similar to what we see in the 0.1 M J cases.In addition, we note that fewer systems have the outer planet's semimajor axis less than 0.910 au compared to the 0.1 M J cases, which indicates that the outer planet tends to scatter to a larger orbit when the ejected planet is more massive.
The distributions of eccentricities are shown in the middle panels of Figure 11.We can see that the 90th percentiles of eccentricities of both planets are larger than 0.6 for all three groups of simulations.We also notice a more significant portion of systems with eccentricities larger than 0.375 compared to the 0.1 M J cases, suggesting that ejecting a more massive planet can increase the excitation of eccentricities.In the bottom panels, the 90th percentiles of the inclination of the inner planet are larger than those of the outer planet, similar to the 0.1 M J cases.In all three groups of simulations, the 90th percentiles of the inner planet are about 40 • , while those of 0.1 M J cases are less than 30 • , suggesting that ejecting a more massive planet can also increase the excitation of inclinations.
We did not study the scenario in which the hypothetical third planet has a mass greater than 0.227 M J for the following two reasons.First, less-massive planets are prone to being ejected during their dynamical evolution, implying that GJ 1148 c is likely to be ejected.Second, if a more massive planet were to be ejected, it would result in two-planet systems with even larger separations, further reducing the likelihood of finding a two-planet system similar to the GJ 1148 system.2020) explored the origin of warm Jupiters via the in situ planet-planet scattering scenario by using the N -body code REBOUND (Rein & Liu 2012).They performed N -body simulations starting with a range of initial separations and planetary masses, including general relativistic (GR) apsidal precession.It is surprising that there are no post-ejection two-planet systems at all in their results.However, in our study, the search for GJ 1148-like planets is based on the existence of stable post-ejection two-planet systems.We have tried to check some results of the fiducial K = 4.0 simulations in Anderson et al. (2020) by repeating their simulations with the N -body codes SyMBA including the GR effect (SyMBA-GR) and REBOUND including the GR effect.The results of our simulations are presented in Appendix A. In particular, the simulations using the two N -body codes give similar results, and they produce a significant number (about 20%) of stable post-ejection two-planet systems.

CONCLUSIONS
We have explored a scenario in which the GJ 1148 system evolves from a three-planet system.Dynamical instabilities can lead to close encounters between the planets, resulting in planet-planet collisions, planet-star close approaches, or planet ejections.We first assumed that there is a third planet with a mass of 0.1 M J in the original GJ 1148 system.By performing scattering experiments with initial orbital separations of 3.5, 4, and 4.5 R m,H , we found that most scattering outcomes are planet-planet collisions, followed by planet ejections, and then a planet getting too close to the star.The properties of the resulting two-planet systems show that the post-collision two-planet systems are unlikely to reproduce the GJ 1148 system because they have small orbital spacings and small eccentricities.However, the post-ejection two-planet systems have more similar properties to the GJ 1148 system.
The semimajor axis of the final inner planet has a strong linear relationship with the semimajor axis of the initial innermost planet.We found that systems with a i,1 close to 0.21 au can reproduce the semimajor axis of the inner planet of GJ 1148.However, the outer planet's semimajor axis varies in a wide range and requires more simulations to search for ones that can match well the observed system.By performing simulations with K = 4.5 in a narrower range of the semimajor axis of the initial innermost planet, we found several systems that are similar to GJ 1148, including one that can reproduce the best-fit properties of the GJ 1148 planets (but with both planets on slightly larger orbits than the GJ 1148 planets).
We also checked the influence of GR precession on ejection outcomes.The Kolmogorov-Smirnov (KS) test shows that GR precession does not significantly affect the properties of post-ejection twoplanet systems.To explore the influence of the mass of the ejected planet, we performed simulations with a hypothetical third planet having the same mass as GJ 1148 c (0.227 M J ).We found that the optimal initial a i,1 moves to around 0.29 au, and the excitation of eccentricities and inclinations is increased.
Figure 12.The fractions of one-, two-, and three-planet systems as a function of time for the fiducial simulations using SyMBA-GR (solid lines) and REBOUND (dash lines).For comparison, it is plotted to look similar to Figure 1 of Anderson et al. (2020).Left: Phase 1, in which we integrate the initial threeplanet systems for a time span of 10 6 periods of the initial inner planet.Right: Phase 2, in which the remaining two-and three-planet systems are integrated for a time span of 10 8 orbital periods of the inner planet.
Figure 1.Left panel : Orbital evolution of the best fit of the GJ 1148 system over 1 Myr(Trifonov et al. 2020).Right panel : Phase-space diagrams of the planetary eccentricities e b and e c versus the apsidal alignment angle ∆ϖ = ϖ b − ϖ c for systems with the same total angular momentum as the GJ 1148 system.Blue crosses indicate the best-fit parameters of e b , e c and ϖ b − ϖ c and their uncertainties, while the red paths represent the trajectories starting at the best fit.

Figure 3 .
Figure 3. Similar to Figure 2, but for a hypothetical three-planet system with K = 4.0 that ends up with only one planet left.Note that after the merger of GJ 1148 b and the third planet at ≃ 3 × 10 3 yr, the blue curve denotes the new and more massive planet.

Figure 4 .
Figure 4. Statistics of initially unstable three-planet systems evolving into two-planet systems and oneplanet systems in trial simulations (Sim.trial ) with K = 4.0 (see also Table2).

Figure 5 .
Figure 5. Properties of post-collision two-planet systems for simulations with initial separation K = 3.5, 4, and 4.5.Only systems from the first 100 simulations are shown for clarity.

Figure 6 .
Figure 6.Properties of post-ejection two-planet systems (with the same planet mass order as the GJ 1148 system) for simulations (Sim.trial and Sim. 1 ) with initial separation K = 3.5, 4, and 4.5.

Figure 7 .
Figure 7. Properties of post-ejection two-planet systems for simulations with initial a i,1 in range of [0.15, 0.3] au and separation K = 4.5.Left panel: a f, i and a f, o vs. a i,1 (similar to the upper panels of Figures 5and 6).Right panel: the period ratio P f, o /P f, i vs. a i,1 .The shaded area represents the range in which the difference between the ratio of planetary periods to that of GJ 1148 does not exceed 1.The systems in this area are selected and marked yellow.

Figure 8 .
Figure 8. Trajectories of eccentricities versus apsidal alignment angle ∆ϖ = ϖ b − ϖ c for five selected two-planet systems with comparable period ratios to the GJ 1148 system.The trajectory in green color is the one most similar to the best-fit GJ 1148 system.

Figure 9 .
Figure9.Million-year orbital evolution of a post-ejection two-planet system, which is very similar to the best-fit configuration of the GJ 1148 system (see left panel of Figure1).The time-averaged semimajor axis of its inner planet is 0.201 au (red solid line), and that of the outer planet is about 1.085 au (black solid line).The eccentricity variations of the two planets are very similar to the best-fit model for the GJ 1148 system but on a longer time scale.The apsidal alignment angle ∆ϖ = ϖ b − ϖ c varies from −60 • to 60 • .
(b)  shows the cumulative distributions of the orbital properties of

Figure 11 .
Figure 11.Properties of post-ejection two-planet systems with the ejected planet having a mass of 0.227M J .

Table 1 .
Parameters of Different Sets of Simulations

Table 2 .
Scattering Outcomes of Trial Simulations (Sim.trial ) Starting with Different Initial Separations K N total N pp,2p N ps,2p N ej,2p N

Table 3 .
Two-sample Kolmogorov-Smirnov (KS) Test for the Cumulative Distributions of the Orbital Elements of the Inner and Outer Planet for K4.5 and K4.5-GR