A Pair of Planets Likely in Mean-motion Resonance From Gravitational Microlensing

TL;DR

This study combines statistical and dynamical analyses to suggest that the two massive planets in the OGLE-2012-BLG-0026L system are likely in a mean-motion resonance, revealing the potential of microlensing to probe planetary dynamics.

Contribution

It demonstrates that microlensing can be used not only to detect planets but also to infer their dynamical configurations, specifically the likelihood of mean-motion resonance in a multi-planet system.

Findings

The two planets are likely in a mean-motion resonance.

Stable configurations include resonant and nearly circular orbits.

Resonant configuration is more probable given observational and formation constraints.

Abstract

We combine statistical arguments and dynamical analysis to study the orbital configuration of the microlensing planetary system OGLE-2012-BLG-0026L. This system is known to host two massive planets, with both projected close to the Einstein ring at the time of the detection. Assuming an isotropic distribution for the orbital orientation, we find that the two planets should also be closely spaced in three-dimensional space and that the ratio of their orbital periods is almost certainly less than two. With N-body numerical integrations, we then identify two types of stable configurations: the two planets can be in first-order mean-motion resonances (MMRs) and have significant ($\gtrsim0.1$) eccentricities, or they stay out of MMRs and have nearly circular orbits. The latter is disfavored, given the absence of similar planet pairs in radial velocity (RV) observations as well as the theoretical difficulties in forming such a configuration. Therefore, the two massive planets in OGLE-2012-BLG-0026L are likely in a resonance configuration. Our work shows that the microlensing technique, which usually only measures the projected configurations, can also probe the detailed dynamical state of multi-planet systems. We also discuss theoretical implications of measuring the multiplicity and the orbital architecture of cold planets.

Figures

Tables

Table 1. Table 1: Massive planet pairs with period ratio < 2.2 absent 2.2 <2.2 .

Host name	$M_{\star }/M_{\odot }$	${(m\sin i,P,e,\omega )}_{\mathrm{in}}$aaWe use the parameters from the Keplerian fit (i.e., no planetary interactions) wherever possible. Here the minimum mass $m\sin i$ is in unit of $M_J$, the orbital period $P$ is in unit of days, and the argument of periapsis $\omega$ in unit of degrees. Under the coplanarity assumption $\varpi =\omega$.	${(m\sin i,P,e,\omega )}_{\mathrm{out}}$aaWe use the parameters from the Keplerian fit (i.e., no planetary interactions) wherever possible. Here the minimum mass $m\sin i$ is in unit of $M_J$, the orbital period $P$ is in unit of days, and the argument of periapsis $\omega$ in unit of degrees. Under the coplanarity assumption $\varpi =\omega$.	MMRbbWe consider MMR detected if the dynamical analysis (in several cases, strongly) prefers resonant solutions. Otherwise, MMR remains possible although slightly disfavored.?	Reference	Comments
HD 200964	$1.4$	$(1.9,631,0.11,223)$	$(1.3,829,0.11,301)$	Yes	Johnson et al. (2011)	(1)
HD 5319	$1.5$	$(1.8,641,0.02,97)$	$(1.2,886,0.15,252)$	Yes	Giguere et al. (2015)	(2)
HD 45364	$0.8$	$(0.19,227,0.17,163)$	$(0.66,343,0.10,7.4)$	Yes	Correia et al. (2009)	(3)
HD 33844	$1.7$	$(2.0,551,0.15,211)$	$(1.8,916,0.13,71)$	Yes	Wittenmyer et al. (2016)	(4)
HD 47366	$1.8$	$(1.8,363,0.09,100)$	$(1.9,685,0.28,132)$	Possible	Sato et al. (2016)	(5)
24 Sex	$1.5$	$(1.6,455,0.18,227)$	$(1.4,910,0.41,172)$	Yes	Johnson et al. (2011)	(6)
$\eta$ Cet	$1.7$	$(2.6,404,0.13,251)$	$(3.3,752,0.10,68)$	Yes	Trifonov et al. (2014)	(7)
HD 202696	$1.9$	$(1.9,521,0.06,259)$	$(2.0,956,0.26,129)$	Possible	Trifonov et al. (2019)	(8)
HD 73526	$1.1$	$(2.3,189,0.29,196)$	$(2.3,379,0.28,272)$	Yes	Wittenmyer et al. (2014)	(9)
HD 82943	$1.2$	$(1.9,220,0.37,117)$	$(1.7,441,0.16,300)$	Yes	Baluev & Beaugé (2014)	(10)
HD 155358	$1.2$	$(1.0,194,0.17,143)$	$(0.82,392,0.16,180)$	Yes	Robertson et al. (2012)	(11)
HD 128311	$0.8$	$(1.8,453,0.30,58)$	$(3.1,922,0.16,15)$	Possible	McArthur et al. (2014)	(12)
HD 160691	$1.1$	$(0.52,310,0.07,190)$	$(1.1,643,0.13,22)$	Possible	Pepe et al. (2007)	(13)
HD 37124	$0.9$	$(0.65,885,0.13,53)$	$(0.70,1862,0.16,0)$	Possible	Wright et al. (2011)	(14)

Equations

$$\begin{array}{rl}{q}_1=1.3\times 10^{-4},& r_{\perp ,1}=0.96\mathrm{or}1.03;\\ q_2=7.9\times 10^{-4},& r_{\perp ,2}=0.81\mathrm{or}1.25;\end{array}$$

$$\left(\begin{array}{c}{x}_j\\ y_j\end{array}\right)=r_j\left(\begin{array}{c}\cos \Omega \cos u_j-\sin \Omega \sin u_j\cos i\\ \sin \Omega \cos u_j+\cos \Omega \sin u_j\cos i\end{array}\right),$$

$$Z\equiv \frac{|e_2{e}^{i{\varpi }_2}-e_1{e}^{i{\varpi }_1}|}{\sqrt{2}},$$

$$r_{\mathrm{H}}\equiv \frac{a_{\mathrm{in}}+a_{\mathrm{out}}}{2}{\left(\frac{q_{\mathrm{in}}+q_{\mathrm{out}}}{3}\right)}^{1/3}.$$

Full text

gravitational lensing: micro — methods: numerical — planetary systems: dynamical evolution and stability

1 Introduction

Gravitational microlensing can detect planetary-mass objects around stars through their perturbations on the otherwise smooth and simple Paczynski (1986) light curves (Mao & Paczynski, 1991; Gould & Loeb, 1992). This technique is most sensitive to planets close to the Einstein ring, corresponding to $\sim$3 au separation for typical Galactic events (see recent reviews by Mao 2012 and Gaudi 2012). Although such a feature makes it efficient in probing the cold planet population, it does place challenge in obtaining orbital information of the detected planets, as their orbital periods ($\sim$10 yr) are much longer than the duration of the microlensing event ($\sim$30 days). Consequently, for the majority of microlensing planets only the instantaneously projected separations between the planets and their hosts are constrained (Gaudi & Gould, 1997).

Traditionally, there are three ways to constrain the actual three-dimensional configuration of the microlensing planetary system, but they are individually challenging and only apply to rare events. In principle, one can try to measure the radial velocities (RV) of the host star (see Yee et al. 2016 for a successful application in the binary star case). However, this requires the host star to be relatively bright and the planet to be fairly massive. If the planetary perturbation in the microlensing light curve persists for a significant fraction ($\gtrsim 10\%$) of the orbital period, then it becomes possible to constrain the full orbit (Ryu et al., 2018). If the planet imposes prominent features (i.e., caustic crossing) in the light curve, from which its locations at certain epochs can be precisely determined, then a much shorter perturbation can still provide meaningful constraints on the full orbit (Gaudi et al., 2008; Bennett et al., 2010). These latter two scenarios require special conditions, and thus they are not generally applicable.

The orbits of planets in multi-planet systems are of special interest, as they provide information on the formation and evolution of these systems. For such systems, stability analysis offers another route in constraining their orbits. A similar idea has been massively applied in the multi-planet systems found by other techniques, including radial velocity (e.g., Vogt et al., 2005; Lee et al., 2006), transit (e.g., Lissauer et al., 2011; Fabrycky et al., 2014; Tamayo et al., 2016), and direct imaging (Fabrycky & Murray-Clay, 2010; Wang et al., 2018). In particular, the direct imaging technique commonly measures the projected positions of planets, which makes it similar to microlensing. Fabrycky & Murray-Clay (2010) applied the stability analysis to the directly imaged multi-planet system HR 8799 (Marois et al., 2008), and found that, in order to keep the system long-term stable, the three planets should likely be in a double resonance configuration. Compared to direct imaging, which infers the planetary mass from measured flux, microlensing directly yields the planet-to-star mass ratio, which is more relevant for the dynamical analysis.

In this paper, we perform the dynamical analysis of a two-planet system detected via microlensing, OGLE-2012-BLG-0026L (hereafter OB120026L; Han et al. 2013). We first summarize the known properties of this system in Section 2. Then, in Section 3, we take into account the projection effect to constrain the orbital spacing between two planets. In Section 4 we perform N-body simulations to determine the long-term stable configurations and, based on theoretical reasonings and supporting evidence from known systems, we argue that the two planets in OB120026L are likely in one of the mean-motion resonances (MMRs). Finally, in Section 5, we discuss the implications from this particular system to future microlensing detections and the general demographics of exoplanets.

2 Known Properties of OB120026L

OB120026L is known to host two planets with planet-to-star mass ratios and projected separations from the host star (in units of Einstein ring radius $r_{\rm E}$)

[TABLE]

respectively. Furthermore, the opening angle between the two planets is $\psi=137^{\circ}$ (Han et al., 2013). The ambiguities in $r_{\perp}$ values come from the close-wide degeneracy (Griest & Safizadeh, 1998), as both planets were detected through their central caustic perturbations. Even so, the ratio of projected separations differs by $<30\%$ for all four solutions, corresponding to a period ratio $<1.5$ if the actual semi-major axes are proportional to the projected separations. The fact that the two planets are not aligned (or anti-aligned) with respect to the central star also suggests that the orbital plane cannot be highly inclined, a property we will use to constrain the orbital orientation in Section 3.

For demonstration purposes, in the paper we only address one of the four microlensing models, namely model (B) of Han et al. (2013). This model has $r_{\perp,1}=1.03$ and $r_{\perp,2}=0.81$, and thus the separation ratio ($r_{\perp,\rm out}/r_{\perp,\rm in}=1.26$) is the second largest among the four models. We provide in Appendix A the results from using the other three solutions.

Although the dynamical analysis is mostly determined by the well constrained mass ratios, we notice that this system also has a well constrained stellar mass. The mass of OB120026L host star was measured to be $1.06\pm 0.05~{}M_{\odot}$, by combining the lens flux and finite source effect measurements (Beaulieu et al., 2016). This makes it a planetary system with a Sun-like host, which is not common for Galactic microlenses (e.g., Zhu et al., 2017). The masses of the two planets are therefore $M_{1}=0.15~{}M_{\rm J}$ and $M_{2}=0.86~{}M_{\rm J}$, respectively, and the Einstein ring radius $r_{\rm E}=4~{}$au. The well determined physical properties of this system make it possible to connect with known exoplanet demographics around Sun-like stars, as we will discuss in Section 5.

3 Orbital Orientation of OB120026L

The microlensing light curve provides stringent constraints ($\lesssim 3\%$) on the mass ratio $q$ and the projected positions (i.e., dimensionless separation $s$ and azimuthal separation $\psi$), 111The uncertainties on masses and absolute projected separations are relatively large ($\sim 10\%$). These are dominated by the uncertainties on the lens mass measurement (Beaulieu et al., 2016). However, the dynamics of the system can be well described by the dimensionless parameters $q$, $s$ and $\psi$, and the absolute scale only has marginal effect. but it has nearly zero constraint on the projected velocities as well as positions and velocities along the line-of-sight direction (Han et al., 2013). In this section, we explore the impact of these unknown parameters on the relative positions between two planets in the orbital plane.

The projected position of the planet $j$ ($j=1$ or 2) relative to the host is given by

[TABLE]

where $\Omega$ is the longitude of ascending node, $i$ is the orbital inclination, $u$ is the argument of latitude (i.e., the sum of the argument of periapsis $\omega$ and the true anomaly $f$), and $r$ is the radial separation between the planet and the star. We define the $z$-axis to be the direction toward the observer, and the $x$-axis to be the direction toward the projected position of the lower-mass planet (see Figure 2 of Han et al. 2013). We assume that the two planets are coplanar. The deviation from this assumption will be addressed at the end of this section.

We are mostly interested in the probability distribution of the period ratio between the two planets. To reach this goal, we first randomly choose $i$ from an isotropic distribution ($\propto\sin{i}$) and $\Omega$ from a uniform distribution between [math] and $2\pi$. We then solve the above set of equations for $r_{j}$ and $u_{j}$ ($j=1,~{}2$). For circular orbits, the semi-major axis, $a$, is given by $r$ and then the period ratio can be derived from the semi-major axis ratio. We always report the outer-to-inner period ratio. To account for the commonly adopted log-flat distribution in $a$, we weight the derived period ratio by $1/(a_{1}a_{2})$. From $50,000$ random sets of $i$ and $\Omega$ we construct the probability distribution of the period ratio as well as the distribution of the azimuthal separation between the two planets in the orbital plane (i.e., $u_{2}-u_{1}$), both of which are illustrated in Figure 1. The two planets likely have small period ratios even after the correction of the projection effect. Although the period ratio distribution has a long tail toward larger values, it is unlikely ($<14\%$) for the period ratio to exceed 2. There is also a non-negligible probability that the two planets are in a super compact (period ratio close to unity) configuration, the majority of which is dynamically unstable. This is further explored in Section 4. The observed azimuthal separation $\psi$ also serves as a good proxy for the azimuthal separation between the two planets in the orbital plane.

Introducing eccentric orbits only changes both probability distributions at a quantitative level. To demonstrate this, we carry out the following exercise. For each orientation of the orbital plane, we also introduce to each planet an eccentricity vector with a fixed amplitude $e=0.3$ but a random direction $\omega$. The semi-major axis is then given by $a=r[1+e\cos(u-\omega)]/(1-e^{2})$, and the rest procedure is the same as in the circular case. The resulting probability distributions are also shown in Figure 1. Even with such extreme eccentricities, the key feature remains: the period ratio is most likely smaller than 2.

We now discuss the impact of deviations from coplanar orbits. If the orbit of the apparently inner planet (planet 2) is more inclined with respected to the other one, the period ratio becomes even smaller than that derived under the coplanarity assumption, unless the mutual inclination becomes larger than $\cos^{-1}{[r_{\perp,2}/r_{\perp,1}]}=38^{\circ}$. After that the apparently inner planet becomes the outer one. Such a configuration becomes a hierarchical system, which we defer to Section 5 for further discussion. Inclining the orbit of the apparently outer planet (planet 1) would increase the period ratio, but unless the outer planet is significantly inclined, the period ratio derived from the projected separations serves as a reasonably good approximation of the actual period ratio. 222Retrograde orbits can in principle produce a more stable system than the prograde orbits (Nesvorný et al., 2003). However, such a configuration is considered unlikely and therefore not discussed in this paper.

To summarize, the closely spaced configuration of the two-planet system OB120026L in the projected plane strongly suggests that the two planets are also closely spaced in 3D space.

4 Two Planets Likely in Mean-Motion Resonances

Here we investigate the impact of dynamical instability on the system configuration. Because the two planets are closely spaced, mean motion resonances (MMRs), which are not captured by the stability criteria either derived analytically (e.g., Deck et al., 2013; Hadden & Lithwick, 2018) or empirically (e.g., Petrovich, 2015), make a non-negligible contribution to all stable configurations. Therefore, we use N-body numerical integrations to determine the stability of any given configuration. These are done with the REBOUND (Rein & Liu, 2012) code with the WHFast integrator (Rein & Tamayo, 2015), which is a fast and accurate implementation of the Wisdom & Holman (1991) symplectic algorithm.

We further assume that the planetary system is face-on ($i=0$). Because the spacing between two planets, which is most relevant for this dynamical analysis, is very small even when different orientations of the orbital plane are taken into account (see Section 3), our result below is not affected by the choice of this initial condition.

In principle, the two detected planets can both have eccentric orbits, but we consider eccentric orbits for two planets separately and assume initially circular orbit for the other planet. In each case, we run $10^{5}$ simulations. The initial state of the planet in a circular orbit is fully determined by the known $(x,~{}y)$ positions. For the eccentric orbit, we first randomly choose eccentricity, $e$, uniformly between 0 and 0.5 and mean anomaly, $l$, uniformly between 0 and $2\pi$, and then solve Kepler’s equation for eccentric anomaly, $E$. Together with the known $(x,~{}y)$ coordinates, these parameters yield the semi-major axis, $a$, and the argument of periapsis, $\omega$, and thus the initial state of this eccentric planet is also fully determined. For each simulation, we run up to 3000 orbits of the longer-period planet, with a step size of one 100th of the smaller orbital period. To distinguish between regular and chaotic orbits, we use the Mean Exponential Growth of Nearby Orbits (MEGNO) chaos indicator (Cincotta et al., 2003) built into REBOUND and use MEGNO=2 as the threshold. 333This is a conservative threshold since the MEGNO value of a regular orbit can be slightly larger than 2 for a finite integration time. However, our result is not sensitive to the choice of the MEGNO threshold, as similar simulations have shown (e.g., Figure 6 of Hadden & Lithwick 2018). The system is also determined to be unstable if either of the planets is too far ($>20$ au) away from or too close ($<0.01$ au) to the central star. This last criterion is introduced in order to make the computations more efficient. For a subset (139) of the stable simulations, we integrate up to $10^{5}$ orbits. A large fraction of these integrations lead to unstable orbits, but there is no systematic preference over any specific initial condition (low-$e$, non-resonant or high-$e$, resonant). Therefore, the results we report below are not affected by the integration runtime.

We show the relative eccentricities and period ratios of the regular configurations (MEGNO$\leq 2$) in Figure 2. Here the relative eccentricity is defined as

[TABLE]

where $\varpi_{1}$ and $\varpi_{2}$ are the longitudes of periapsis of the two planets, respectively (Hadden & Lithwick, 2018). Given our choices of the initial condition, this is approximately the initial eccentricity of the eccentric planet. As Figure 2 shows, there are two different types of regular configurations regardless of which planet is set eccentric initially: planets have nearly circular ($e\lesssim 0.1$) orbits, or planets have significant eccentricities and can be in one of the first-order MMRs. In fact, as seen in Figure 2, there are regular configurations under which the two planets have crossing orbits, and it is MMR that protects the system from instability.

Although both types of configurations are allowed from the stability point of view, the MMR configurations are preferred once the formation history of such a planet pair is taken into account. Figure 3 shows the selected planet pairs from NASA Exoplanet Archive (Akeson et al., 2013). 444We also included $\eta$ Cet and HD 202696. Both systems have planet pairs that meet our criteria. Parameters of individual systems are provided in Table 1. Similar to the two planets in OB120026L, these planets are relatively massive ($>0.15~{}M_{\rm J}$) and well-separated ($P>100~{}$d) from their hosts. In the period ratio range $1<P_{\rm out}/P_{\rm in}<2.2$, the majority of the planet pairs have period ratios close to commensurabilities, with half close to the 2:1 ratio. Dynamical analyses have suggested the existence of MMRs in most cases, in particular if the period ratio is $<2$.

The lack of giant planet pairs in closely spaced but non-MMR configurations in Figure 3 is also well understood theoretically. When protoplanets grow, orbital repulsion keeps the separation between neighboring planets wider than 5–10 times the mutual Hill radius of the protoplanets (Kokubo & Ida, 1998). Here the mutual Hill radius is given by

[TABLE]

For OB120026L, this means that the two planets should start with at least $P_{\rm out}/P_{\rm in}>1.6$ and more likely $P_{\rm out}/P_{\rm in}>2.2$. The initial orbits are also nearly circular ($e<0.1$), due to the interactions between protoplanets and the embedded disk. To produce planet pairs in smaller period ratios, convergent migration is required, and together with eccentricity damping this usually leads to planet pairs trapped in MMRs (Lee & Peale, 2002). It is difficult to form closely spaced massive planet pairs that are not in MMRs without fine-tuning the conditions.

The above theoretical arguments and the lack of planet pairs outside MMRs in small period ratios in Figure 3 both suggest that the MMR configurations are preferred for the two-planet system OB120026L.

5 Discussion

We combine statistical arguments and dynamical analysis to study the orbital configuration of the two-planet system OB120026L that was found through gravitational microlensing (Han et al., 2013). The two massive planets, one with $0.86~{}M_{\rm J}$ and the other with $0.15~{}M_{\rm J}$ (Beaulieu et al., 2016), were found to be closely spaced in the sky-projected plane, with the ratio of the projected separations $<1.3$. Under the reasonable assumptions that the orbital orientation follows an isotropic distribution and that the mutual inclinations are typically small, this observed configuration strongly suggests that the two planets should also be closely spaced in the 3D space, and that the period ratio is unlikely to exceed two (Figure 1).

We then study the eccentric orbits with numerical integrations and find two types of stable configurations. The two planets can be in MMRs and have significant ($e\gtrsim 0.1)$ eccentricities, or they have nearly circular ($e\lesssim 0.1$) orbits and stay out of nominal MMRs. The out-of-MMR configurations are disfavored, given the absence of similar planet pairs from RV observations and the theoretical difficulties in forming such configurations. Therefore, the two massive planets in OB120026L are probably in MMR.

The above conclusion builds on the assumption that the two planets in OB120026L have nearly coplanar ($\lesssim 38^{\circ}$, see Section 3) orbits. Given our understanding of the architecture of planetary systems from both observations and theories, this is a fair assumption, although other scenarios cannot be entirely ruled out. For example, a hierarchical, highly mutually inclined two-planet system can also produce the observed configuration by chance. We consider such a scenario unlikely, given the small a priori probability for the two planets in such a hierarchical system to have similar projected separations. 555Assuming isotropic orientation, the chance for a planet with semi-major axis $a$ to have a projected separation $r_{\perp}$ is simply $r_{\perp}/a$. With $r_{\perp}\approx 4$ au, this is $10\%(a/{\rm 40~{}au})^{-1}$. The two planets in such hierarchical systems will likely have dynamical interactions (i.e., Lidov-Kozai libration Lidov 1962; Kozai 1962), the inclusion of which will further decrease the probability. In principle, this low by-chance projection probability can be compensated by a much higher occurrence rate. This, if true, will lead to an enhanced rate of microlensing detections in which a normal planetary event is followed (or preceded) by a short single-lens event due to the distant planet. 666See OGLE-2008-BLG-092 (Poleski et al., 2014) for an example, although in this case the distant companion is a star. Future microlensing observations will be able to test whether this is true or not.

Our work shows that the microlensing technique, which usually only measures the projected configuration, can also probe the detailed dynamical state of the multi-planet system. The two features of OB120026L that are key to our conclusion are that 1) the two planets are well separated azimuthally, and 2) they are both close to the Einstein ring. For a comparison, the first two-planet system, OGLE-2006-BLG-109L, has two planets almost aligned ($\sim 13^{\circ}$) with the central star and the projected separation ratio $\sim 1.7$ (Gaudi et al., 2008; Bennett et al., 2010), making our method inapplicable. So is OB120026L-type configuration rare in microlensing detections? The detections of multiple planets in the same system are largely independent (Gaudi et al., 1998; Zhu et al., 2014a, b; Shin et al., 2015). Hence, it is more likely to detect two-planet systems with both planets close to the Einstein ring, given that microlensing sensitivity is a steep function of the distance from the Einstein ring. Furthermore, there is no obvious preference to any particular azimuthal separation in two-planet detections, as shown in previous simulations (Zhu et al., 2014a). Therefore, we expect that a significant fraction of future two-planet systems (and likely higher multiples) should have configurations similar to OB120026L, and thus the method we developed here will be generally applicable.

The orbital architecture and occurrence rate of planetary systems like OB120026L, which hosts both a cold Jupiter ($q\sim 10^{-3}$) and a cold Neptune ($q\sim 10^{-4}$), are important for our better understanding of the planet formation process. Capturing planets into MMRs requires some dissipative migration through interactions with the gas disk or planetesimals (e.g., Lee & Peale, 2002), so the fraction of planet pairs in MMRs is an important proxy for understanding how planetary architectures are set during formation. Ground-based RV observations find that about 1/3 of well-characterized multi-planet systems contain planet pairs with apparent low-order period commensurability (Wright et al., 2011). Such a high concentration cannot be explained by random fluctuations, but is consistent with the expectation from migration theory (Goldreich & Tremaine, 1980). By contrast, the Kepler mission found that the close-in small planets do not show a strong preference for MMRs (Lissauer et al., 2011), suggesting that migration, or at least large scale migration, might not have happened for these planets (Petrovich et al., 2013). These two types of planets differ in both mass and orbital separation, and therefore the cold Neptunes that microlensing is sensitive to (Gould et al., 2006), which have similar separations to the RV cold giants but similar masses to the Kepler super Earths, convey important clues on what drives the different architectures.

It is now known that cold Jupiters are almost certainly accompanied by inner small planets (Zhu & Wu, 2018), but how often do they also have cold Neptune companions? While this question can serve a test of planet formation theories, it is also relevant for understanding the prevalence of the solar system-like architecture. By the end of 2018, there were 35 microlensing planetary systems with $q>4\times 10^{-4}$ listed on NASA Exoplanet Archive (Akeson et al., 2013), out of which two were two-planet systems with cold Neptunes. 777There is also a candidate two-planet system OGLE-2014-BLG-1722, but both planets have $q>4\times 10^{-4}$ (Suzuki et al., 2018). If one simply takes the detection efficiency of an additional planet from the simulation ($\sim 5\%$, Zhu et al. 2014a), these numbers seem to suggest that all cold Jupiters can be accompanied by cold Neptunes. However, the above result has to be taken with caution because of the mismatch between observations and the simulation. Both two-planet systems were detected in high-magnification events that received intensive follow-up observations (e.g., Gould et al., 2010), whereas the simulation done by Zhu et al. (2014a) assumed an observing strategy with a consistent cadence (10 min) similar to the Korean Microlensing Telescope Network (KMTNet, Kim et al., 2016). Therefore, more detailed works and larger samples are needed in order to settle this issue.

It is also interesting to notice that the first two multi-planet systems from microlensing both contain one planet with intermediate mass ($30-100~{}M_{\oplus}$): OGLE-2006-BLG-109c has $86~{}M_{\oplus}$ (Bennett et al., 2010), and OGLE-2012-BLG-0026c has $46~{}M_{\oplus}$ (Beaulieu et al., 2016). 888Both events have finite-source effect, microlensing parallax, and lens flux measurement, and thus the mass measurements are fairly secure. In the standard core accretion scenario (Pollack et al., 1996), hydrostatic gas accretion can no longer be maintained once the gaseous envelop doubles the planetary core mass ($\sim 10~{}M_{\oplus}$), and subsequent run-away gas accretion pushes the total mass rapidly to the gas giant regime ($\gtrsim 100~{}M_{\oplus}$). Therefore, the population synthesis models based on the core accretion theory predicted that planets with masses in the intermediate mass range should be very rare (Ida & Lin, 2004; Mordasini et al., 2009). This is not supported by microlensing observations, which have found that the transition from Neptune mass to Jupiter mass is rather smooth (Suzuki et al., 2016). One possibility to reconcile the theory and observations, as suggested by the first two multi-planet systems, is that the multiplicity rate may be higher than model predicts, and that competitions for gas material starve some of the embryos from fully growing into giants. Such a scenario is particularly possible for OB120026L, since in this system the two planets are so closely spaced. Future microlensing studies on the multiplicity rate can confirm whether this explanation stands or not.

We would like to thank Cristobal Petrovich, Yanqin Wu, Dan Tamayo, and Mohamad Ali-dib for useful discussions, and Andy Gould for comments on an earlier version of the paper. We also thank the anonymous referee for comments and suggestions that helped improve the manuscript. W.Z. also thanks Dan Fabrycky for some early discussions. This work is the result of the 2018 Summer Undergraduate Research Program (SURP) in astronomy & astrophysics at the University of Toronto. W.Z. was supported by the Beatrice and Vincent Tremaine Fellowship at CITA.

Appendix A Results from other Solutions

Our main results do not change even when other allowed microlensing solutions are used. Here we show the results from the other three allowed solutions given in Han et al. (2013). Figure 4 shows the probability distributions of the period ratio and the azimuthal separation in the orbital plane (i.e., before projection effect). Figures 5, 6, and 7 show the stable configurations out of the N-body integrations. For each case we have $10^{4}$ realizations.

Appendix B Packed Massive Planet Pairs from RV

We list in Table 1 the information of RV planet pairs used in Figure 3. Some comments on the dynamical analysis are also provided.