# Formation of the terrestrial planets in the solar system around 1 au via radial concentration of planetesimals

Masahiro Ogihara, Eiichiro Kokubo, Takeru Suzuki, Alessandro Morbidelli

► To cite this version:

Masahiro Ogihara, Eiichiro Kokubo, Takeru Suzuki, Alessandro Morbidelli. Formation of the terrestrial planets in the solar system around 1 au via radial concentration of planetesimals. Astronomy & Astrophysics - A&A, 2018, 612, pp.L5. 10.1051/0004-6361/201832654. hal-02307855

HAL Id: hal-02307855

https://hal.science/hal-02307855v1

Submitted on 16 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
LETTER TO THE EDITOR

# Formation of the terrestrial planets in the solar system around 1 au via radial concentration of planetesimals

Masahiro Ogihara<fnref n="1" />, Eiichiro Kokubo<fnref n="1" />, Takeru K. Suzuki<fnref n="2" />, and Alessandro Morbidelli<fnref n="3" />

Received 17 January 2018 / Accepted 29 March 2018

## ABSTRACT

*Context.* No planets exist inside the orbit of Mercury and the terrestrial planets of the solar system exhibit a localized configuration. According to thermal structure calculation of protoplanetary disks, a silicate condensation line (~1300 K) is located around 0.1 au from the Sun except for the early phase of disk evolution, and planetesimals could have formed inside the orbit of Mercury. A recent study of disk evolution that includes magnetically driven disk winds showed that the gas disk obtains a positive surface density slope inside ~1 au from the central star. In a region with positive midplane pressure gradient, planetesimals undergo outward radial drift.

*Aims.* We investigate the radial drift of planetesimals and type I migration of planetary embryos in a disk that viscously evolves with magnetically driven disk winds. We show a case in which no planets remain in the close-in region.

*Methods.* Radial drifts of planetesimals are simulated using a recent disk evolution model that includes effects of disk winds. The late stage of planet formation is also examined by performing *N*-body simulations of planetary embryos.

*Results.* We demonstrate that in the middle stage of disk evolution, planetesimals can undergo convergent radial drift in a magnetorotational instability (MRI)-inactive disk, in which the pressure maximum is created, and accumulate in a narrow ring-like region with an inner edge at ~0.7 au from the Sun. We also show that planetary embryos that may grow from the narrow planetesimal ring do not exhibit significant type I migration in the late stage of disk evolution.

*Conclusions.* The origin of the localized configuration of the terrestrial planets of the solar system, in particular the deficit of close-in planets, can be explained by the convergent radial drift of planetesimals in disks with a positive pressure gradient in the close-in region.

**Key words.** planets and satellites: formation – protoplanetary disks – planet–disk interactions – methods: numerical

## 1. Introduction

The properties of the terrestrial planets of our solar system are well understood compared to other systems. A number of *N*-body simulations of planetary accretion were performed to reproduce these properties, such as small eccentricities and late moon formation (e.g., Nagasawa et al. 2005; Ogihara et al. 2007; Raymond et al. 2009; Morishima et al. 2010).

One of the major issues facing the theory of planet formation is reproducing the current orbital configuration of the solar system's terrestrial planets. Their orbital locations are confined to a region between 0.38 and 1.5 au from the Sun. Several models have been proposed to reproduce the outer boundary of the localized distribution which includes the grand tack model (e.g., Walsh et al. 2011). In the meantime, it is quite hard to explain why there is no large planets in close-in orbits. If planetesimals form when the condensation line for silicates is located inside *r* = 0.7 au, protoplanets/planets should form in close-in orbits. However, no planets exist inside the orbit of Mercury in our solar system. Regarding this problem, Ogihara et al. (2015) demonstrated that protoplanets larger than Mars that form in close-in orbits can undergo outward type I migration in disks with positive surface density profiles, which would help to reproduce the deficit of close-in planets.

Based on three-dimensional (3D) magnetohydrodynamic (MHD) simulations, Suzuki et al. (2016) derived the long-term global disk evolution including effects of magnetically driven disk winds. They considered the mass loss due to disk winds (Suzuki & Inutsuka 2009) and the wind-driven accretion (Bai & Stone 2013; Gressel et al. 2015; magnetic braking) in addition to the standard viscous evolution. They found that the disk profile in the close-in region (r ≲ 1 au) can be significantly altered from the minimum-mass solar nebular (MMSN) model (Weidenschilling 1977; Hayashi 1981). In particular, the slope of the gas surface density is almost flat inside r ≈ 1 au for MRI-active disk, and the slope can even be positive for MRI-inactive disk. In a separate paper (Ogihara et al. 2018, hereafter Paper II), we investigated the orbital evolution of planetary embryos larger than Mars in such disks under various conditions, and found that type I migration of protoplanets can be significantly suppressed in many cases compared to the type I migration rate estimated for the MMSN model. On the other hand, it is expected that smaller particles (dust/planetesimals) can undergo significant radial drift in MRI-inactive disks of Suzuki et al. (2016).

In this work, we examine another path to form the localized configuration of the solar system's terrestrial planets using the recent disk evolution model of Suzuki et al. (2016). We examine whether or not the following formation scenario can take place

<footnote n="1">¹ Division of Theoretical Astronomy, National Astronomical Observatory of Japan, 2-21-1, Osawa, Mitaka, 181-8588 Tokyo, Japan
e-mail: masahiro.ogihara@nao.ac.jp</footnote>

<footnote n="2">² School of Arts & Sciences, University of Tokyo, 3-8-1, Komaba, Meguro, 153-8902 Tokyo, Japan</footnote>

<footnote n="3">³ Laboratoire Lagrange, Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Blvd de l'Observatoire, CS 34229, 06304 Nice Cedex 4, France</footnote>
by numerical simulations: planetesimals undergo convergent radial drift to form localized planetesimal disks, then planetary embryos grow to terrestrial planets without undergoing significant type I migration.

The paper is organized as follows. In Sect. 2, we present the numerical model of disk evolution. In Sect. 3, we show the outcomes of simulations that consider the radial drift of planetesimals. In Sect. 4, we present results of N-body simulations from planetary embryos. In Sect. 5, we make concluding remarks.

## 2. Disk model

For evolution of disk surface density, we numerically solve the diffusion equation that includes effects of disk winds based on MHD simulations (the same as used by Suzuki et al. 2016; see also Simon et al. 2015; Hasegawa et al. 2017):

$$ \frac{\partial \Sigma_g}{\partial t} = \frac{1}{r} \frac{\partial}{\partial r} \left[ \frac{2}{r\Omega} \left\{ \frac{\partial}{\partial r} \left( r^2 \Sigma_g \overline{\alpha_{r,\phi} c_s^2} \right) + r^2 \overline{\alpha_{\phi,z}} \frac{\Sigma_g H \Omega^2}{2\sqrt{\pi}} \right\} \right] - C_w \frac{\Sigma_g \Omega}{\sqrt{2\pi}}, \quad (1) $$

where $\Sigma_g$, $r$, $\Omega$, $c_s$, and $H$ are the gas surface density, the radial distance, the Keplerian frequency, the sound speed, and the disk scale height, respectively. An effective turbulent viscosity $\overline{\alpha_{r,\phi}}$ is based on a Shakura & Sunyaev (1973) $\alpha$ viscosity parameterization. Parameters $\overline{\alpha_{\phi,z}}$ and $C_w$ describe measures of the angular momentum loss due to the wind torque and mass loss due to disk winds, respectively. According to results of MHD simulations, we chose the parameters as follows. We adopt two values of the effective turbulent viscosity; namely, $\overline{\alpha_{r,\phi}} = 8 \times 10^{-3}$ (MRI-active case) and $\overline{\alpha_{r,\phi}} = 8 \times 10^{-5}$ (MRI-inactive case). Other parameters are fixed as $\overline{\alpha_{\phi,z}} = 10^{-5}(\Sigma_g/\Sigma_{g,ini})^{-0.66}$ (Bai 2013), $C_w = 2 \times 10^{-5}$ (MRI-active case; Suzuki & Inutsuka 2009) and $C_w = 10^{-5}$ (MRI-inactive case; Suzuki et al. 2010). As described in Sect. 2.3 of Suzuki et al. (2016), the mass-loss rate due to disk winds $C_w$ is also constrained by energetics.

The initial gas profile is given by $\Sigma_{g,ini} = 1.7 \times 10^4(r/1 \text{ au})^{-3/2} \exp(-r/30 \text{ au}) \text{ g cm}^{-2}$. In order to consider the early stage of disk evolution, which is connected to star formation, the initial gas disk is about ten times more massive than the MMSN model. We note that planetesimals and planetary embryos take some time to form, therefore we start simulations some time after the disk formation. In the same way as in our previous works, we start simulations at $t = 10^5 \text{ yr}$ in the disk evolution of Suzuki et al. (2016), and refer to this time as “initial”.

Figure 1 shows the time evolution of the gas surface density for MRI-active disks and MRI-inactive disks. The evolution is the same as in Suzuki et al. (2016). We see that disk profiles are altered in the close-in region ($r \le 1 \text{ au}$). An important point is that the disk obtains a positive surface density slope in a broad close-in region in MRI-inactive disks. As discussed in Suzuki et al. (2016), the angular momentum loss due to the wind torque is inversely correlated with the gas density. In addition, the mass loss rate is proportional to the Keplerian time. Therefore, the gas density is carved out from inside based on the model of Suzuki et al. (2016). It is likely that in the late stage of disk evolution, the disk is MRI active due to smaller gas surface density. For the evolution of the disk temperature, the stellar irradiation and viscous heating are considered. We refer to Sect. 2.4 of Suzuki et al. (2016) for a detailed description.

<figure><img src="image_2.png" /><figcaption>Fig. 1. Time evolution of gas surface density based on Suzuki et al. (2016). Blue and red lines represent MRI-active disk and MRI-inactive disk, respectively.</figcaption></figure>

## 3. Radial drift of planetesimals

We first examine the radial drift of planetesimals before the formation of planetary embryos. Particles with the size of planetesimals or smaller undergo radial migration due to gas drag. As stated in Sect. 2, it is considered that the MRI activity increases with decreasing gas surface density. Therefore, here we use the MRI-inactive disk for the orbital evolution phase of planetesimals. As discussed in Suzuki et al. (2016), disk winds affect the radial drift of planetesimals. As the radial pressure gradient force and hence the rotation velocity of gas are altered, planetesimals can feel a tail wind in the region of a positive pressure gradient, leading to outward drift.

Here we perform simulations of the orbital evolution of planetesimals in MRI-inactive disks with mass $M = 1 \times 10^{16}$ g ($R \approx 1$ km assuming the internal density of 3 g cm$^{-3}$) that are placed between $r = 0.3$ and $3 \text{ au}$. Mutual interactions between planetesimals (i.e., gravity, collision) are ignored. For aerodynamical gas drag force, we use the formula of Adachi et al. (1976). The drift timescale (see Appendix A) is inversely proportional to a scaling factor of pressure gradient force $\eta = -1/2(c_s/v_K)^2(d\ln P/d\ln r)$. See Fig. A.1 for the evolution of $\eta$. Figure 2 shows the result of our simulation, in which the evolution of solid surface density is indicated. We find that planetesimals undergo convergent radial drift. As the disk evolves, the location of the pressure maximum moves outward. After $t = 0.3$ Myr, planetesimals accumulate in a narrow ring-like region because all planetesimals inside $\sim 0.7 \text{ au}$ undergo outward drift. Therefore, it is possible that a localized configuration is attained by radial drift of planetesimals before formation of planetary embryos. We note that typical sizes of planetesimals that form via streaming instability would be larger than 1 km (e.g., Simon et al. 2016). Larger planetesimals suffer less efficient radial drift. Thus, if planetesimals larger than 10 km are distributed between 0.3 and 3 au, the strong convergent drift seen in Fig. 2 does not occur.

Although it is reasonable to assume that the disk is MRI inactive in the early/middle phase of planet formation, we also comment on radial drift of planetesimals in MRI active disks. When the surface density distribution exhibits flat profiles (e.g., MRI-active disk), planetesimals undergo slow inward drift instead of outward drift. In this case, narrow ring-like regions would not form; however, as the drift velocity is slowed down significantly, it would be useful to overcome the drift barrier of dust/planetesimals.

We note that we adopt several simplifications that affect the realism of the simulation. First, a single size population of planetesimals is considered, and their growth is also ignored. In reality, planetesimals have size distribution, which changes with time. Second, mutual interactions between planetesimals and the
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Time evolution of solid surface density. Planetesimals with $R \approx 1$ km are placed between $r = 0.3$ and $3$ au. The disk is considered to be MRI-inactive.</figcaption></figure>

back reaction of planetesimals on the gas disk are neglected. When solid particles accumulate in a narrow ring, where the dust-to-gas ratio is close to unity, the back reaction on the gas disk should be important. The evolution in Fig. 2 is therefore an idealized result. We will investigate these effects in future works.

## 4. Orbital evolution of planetary embryos via type I migration

Now we examine the late stage of planet formation from planetary embryos. As we showed in the previous subsection, planetesimals can accumulate in a narrow-ring region. We therefore consider the subsequent evolution for the initial embryo distribution. Hansen (2009) adopted a localized configuration for the initial embryo distribution and performed *N*-body simulations. We use the same initial conditions as Hansen (2009). As an initial condition of *N*-body simulations, 400 planetary embryos with mass $M = 0.005 M_⊕$ are randomly placed between 0.7 and 1 au according to a uniform surface density. Their initial eccentricities and inclinations are set to be small ($\approx 10^{-2}$). Although Hansen (2009) considered the presence of Jupiter, we do not include any giant planets in our simulations.

For the disk evolution, as we stated in the previous section, it is reasonable to consider that the disk is MRI-inactive in the early phase, which evolves to the MRI-active disk. However, the MRI activity is uncertain when embryos form; therefore we perform *N*-body simulations of embryos both in MRI-active and -inactive disks. The disk evolution model is already shown in Sect. 2.

For the formulae of type I migration, we use those described in Paardekooper et al. (2011). The unsaturated torque is given by $\Gamma = (-0.6 + 2.3p - 2.8q)\Gamma_0$, where $p = \text{dln } \Sigma_g/\text{dln } r$ is the local surface density slope and $q = \text{dln } T/\text{dln } r$ is the local temperature slope (Paardekooper et al. 2010). In our disk model, the temperature slope is roughly given by $q = -1/2$ (see Fig. 5 of Suzuki et al. 2016), therefore the rate and direction of type I migration is mainly determined by the surface density slope $q$. The actual torque is calculated by considering the saturation effect. Here the saturation is determined by the viscous diffusion coefficient $\nu$, which is proportional to $\bar{\alpha}_{r,\phi}$ (see Sect. 4.1 of Ogihara et al. 2017).

Figure 3a shows snapshots of the system for MRI inactive disks (see also Paper II for the migration map). No significant inward migration is seen. Although some protoplanets move

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Snapshots of the system for MRI-inactive disk (a) and MRI-active disk (b). The filled circles represent the size of particles. When particles grow to more massive than $0.1 M_⊕$, they are indicated by filled blue circles. The red lines indicate the gas surface density (right axis). The red circles represent the initial and final orbital configurations of simulations from isolation-mass protoplanets.</figcaption></figure>

outward to $r \approx 2$ au, the innermost larger protoplanets ($M \gtrsim 0.1 M_⊕$) are located at $r = 0.56$ au at $t = 100$ Myr.

Figure 3b shows results of the *N*-body simulation for MRI-active disks. We find no clear tendency for planets to migrate inward or outward. In the final state ($t = 100$ Myr), four larger protoplanets remain between $r = 0.46$ and $1.1$ au. The region where planets exist spreads compared to the initial ring-like region. This is not because planets exhibit migration but because they undergo scattering between protoplanets.

For comparison, we perform *N*-body simulations for the MMSN model. We find that after planetary embryos grow to Mars-mass protoplanets ($M \approx 0.1 M_⊕$), they rapidly migrate inward and are lost to the central star with a migration timescale of $\approx 0.1$ Myr. At $t = 10$ Myr, the total mass of remaining objects is only 0.445 $M_⊕$. Therefore we can confirm that disk winds
play a role in suppressing type I migration, which leads to localized orbital configurations. We note that from geochemical constraints, Mars grew to half of its present size in about 1.8 Myr (Dauphas & Pourmand 2011). If Mars grew on this timescale and the MMSN decays exponentially with a timescale of 1 Myr, Mars would migrate marginally and not fall onto the star.

We have also performed additional simulations that start from protoplanets with an isolation mass instead of embryos with $M = 0.005 \ M_⊕$. For the calculation of the isolation mass, see Eq. (13) of Kokubo & Ida (2002) for example. We chose a relatively small initial orbital separation of 5 $r_H$, where $r_H$ represents the mutual Hill radius of neighboring protoplanets. The initial and final orbital configuration is overplotted in Fig. 3 by circles. We find that our conclusion does not change. No planets remain in the close-in region ($a < 0.4$ au) and planetary orbits are confined to a localized region. Detailed results of N-body simulations of isolation-mass protoplanets are presented in Paper II.

## 5. Conclusions

We investigate the formation of the solar system's terrestrial planets in disks evolving with disk winds. We demonstrate that planetesimals formed in close-in orbits can undergo outward radial drift in MRI-inactive disks and the inner edge of the planetesimal ring is created at $r \sim 0.7$ au. Subsequently, we examine the orbital evolution and growth of embryos from a narrow ring by N-body simulation, and find that migration of protoplanets can be significantly suppressed. As a result, we confirm that no large planets remain inside $\sim 0.5$ au. Hansen (2009) showed that the current orbital configuration of the solar system's terrestrial planets can be reproduced if a narrow solid ring is formed in a gas-free environment. Here, in this work, even if we start simulations from an earlier stage when gas is present, we show that the localized configuration can be maintained. In a previous work, Ogihara et al. (2015) showed that protoplanets can undergo outward migration to $\sim 1$ au in a region with a positive surface density slope. In the present work, we showed that even if planetary embryos do not undergo convergent migration, the localized configuration can be created due to convergent radial drift of planetesimals. The reason why all planetesimals that formed in the close-in region undergo outward radial drift, which creates the inner edge of planetesimal ring, is that the gas disk obtains a positive pressure gradient in a broad close-in region, which expands outwards in time. This results in outward drift of all planetesimals in the close-in region, leaving no planetesimals distributed in the broad close-in region. This is qualitatively different from other narrower pressure traps (e.g., the snow line). Such traps can stop inward drift, but cannot drag planetesimals outwards over the full close-in region. If there was no clear disk inner edge, the narrow-trap case would give a ring-like planetesimal distribution because planetesimals inside the pressure trap would fall onto the star. However, if there existed an inner edge, inward drifting planetesimals would be accumulated there; the planetesimal distribution could be very different from our model. We note also that the outer boundary is determined by the initial distribution of planetesimals/embryos. Other models are required (e.g., the grand tack model) to explain the origin of the outer boundary in the solar system's terrestrial planets.

In Sect. 3, planetesimals are initially distributed between 0.3 and 3 au. Although we do not specify their origin, they would be formed via radial drift of dust and pebbles, which undergo rapid drift. We note here that no pressure bump exists in the early

phase of disk evolution of Suzuki et al. (2016), and it takes some time ($\sim 0.1$ Myr) to form positive slope in the MRI-inactive disk. Therefore, in the earlier stage, pebbles can drift inwards with no limit and planetesimals can form everywhere in the close-in region, which justifies our choice of initial planetesimal distribution. On the other hand, if a pressure bump exists before the planetesimal formation, dust can accumulate at the pressure maximum and planetesimals only form in a local region at the pressure bump. In this case, the planetesimal distribution would differ from that used in Sect. 3. However, our conclusion would still be valid because this situation is consistent with our assumption in Sect. 4 that embryos form in a narrow ring-like region. In this case, the radial drift of planetesimals is not needed and big ($>10$ km) planetesimals could form. In either case, as we stated in Sect. 3, a simulation of growth and radial drift of dust particles in a consistent manner is required. We note also that an alternative mechanism to create the narrow planetesimal ring by drifting pebbles and their backreaction to the disk has been proposed (Drążkowska et al. 2016).

In addition, disk winds would contribute to a solution of the “snow line problem”. According to radiative transfer calculations, the ice condensation line (snow line) at the midplane moves inside 1 au when the accretion rate of the disk drops below $\sim 10^{-9} \ M_⊕\ yr^{-1}$ (e.g., Oka et al. 2011). If planetesimals form at this time and the Earth forms in situ, the Earth should accrete a large amount of icy material, and it is hard to explain the origin of the amount of water (0.023 wt%) on the Earth. Here we apply the disk-wind model to this problem. When the snow line moves inside 1 au, we can assume that the pressure bump is located outside 1 au. If planetary embryos or planetesimals that grow to form the Earth were to form at 1 au by the time of the passage of the snow line through 1 au, the pressure maximum outside 1 au would possibly prevent the inward drift of icy pebbles in the vicinity of the Earth. This could be another mechanism of a fertilization of the snow line proposed by Morbidelli et al. (2016). This point should be explored in future works.

*Acknowledgements.* We thank Cornelis Dullemond for constructive comments. Numerical computations were conducted on the general-purpose PC farm at the Center for Computational Astrophysics, CfCA, of the National Astronomical Observatory of Japan. This work was supported by JSPS KAKENHI Grant Numbers 16H07415 and 17H01105.

## References

Adachi, I., Hayashi, C., & Nakazawa, K. 1976, Prog. Theor. Phys., 56, 1756

Bai, X.-N. 2013, ApJ, 772, 96

Bai, X.-N., & Stone, J. M. 2013, ApJ, 769, 76

Dauphas, N., & Pourmand, A. 2011, Nature, 473, 489

Drążkowska, J., Alibert, Y., & Moore, B., 2016, A&A, 591, A105

Gressel, O., Turner, N. J.; Nelson, R. P., & McNally, C. P. 2015, ApJ, 801, 84

Hansen, B. M. S. 2009, ApJ, 703, 1131

Hasegawa, Y., Okuzumi, S. Flock, M., & Turner, N. J. 2017, ApJ, 845, 31

Hayashi, C. 1981, Prog. Theor. Phys. Suppl., 70, 35

Kokubo, E., & Ida, S. 2002, ApJ, 581, 666

Morbidelli, A., Bitsch, B., Crida, A., et al. 2016, Icarus, 267, 368

Morishima, R., Stadel, J., & Moore, B. 2010, Icarus, 207, 517

Nagasawa, M., Lin, D. N. C., & Thommes, E. 2005, ApJ, 635, 578

Ogihara, M., Ida, S., & Morbidelli, A. 2007, Icarus, 188, 522

Ogihara, M., Kobayashi, H., Inutsuka, S., & Suzuki, T. K. 2015, A&A, 579, A65

Ogihara, M., Kokubo, E., Suzuki, T. K., Morbidelli, A., & Crida, A. 2017, A&A, 608, A74

Ogihara, M., Kokubo, E., Suzuki, T. K., & Morbidelli, A. 2018, A&A, in press
DOI:10.1051/0004-6361/201832720

Oka, A., Nakamoto, T., & Ida, S. 2011, ApJ, 738, 141

Paardekooper, S.-J., Baruteau, C., Crida, A., & Kley, W. 2010, MNRAS, 401, 1950
Paardekooper, S.-J., Baruteau, C., & Kley, W. 2011, MNRAS, 410, 293

Raymond, S., O'Brien D. P., Morbidelli, A., & Kaib, N. A. 2009, Icarus, 203, 644

Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337

Simon, J. B., Lesur, G., Kunz, M. W., & Armitage, P. J. 2015, MNRAS, 454, 1117

Simon, J. B., Armitage, P., Li, R., & Youdin, A. N. 2016, ApJ, 822, 55

Suzuki, T. K., & Inutsuka, S. 2009, ApJ, 691, L49

Suzuki, T. K., Muto, T., & Inutsuka, S. 2010, ApJ, 718, 1289

Suzuki, T. K., Ogihara, M., Morbidelli, A., Crida, A., & Guillot, T. 2016, A&A, 596, A74

Walsh, K. J., Morbidelli, A., Raymond, S. N., O'Brien D. P., & Mandell, A. M. 2011, Nature, 475, 206

Weidenschilling, S. J. 1977, Ap&SS, 51, 153

## Appendix A: Radial drift due to aerodynamical gas drag

<figure><img src="image_5.png" /><figcaption>Fig. A.1. Time evolution of the absolute value of the scaling factor of pressure gradient force $\eta$ for an MRI-inactive disk based on Suzuki et al. (2016). Thick lines correspond to the region where $\eta$ is negative.</figcaption></figure>

We use a force formula of aerodynamical gas drag on planetesimals developed by Adachi et al. (1976);

$$F_{\text{aero}} = -\frac{1}{2} C_D \pi R^2 \rho_g \Delta u \Delta u, \qquad (A.1)$$

where $C_D$, $\rho_g$, and $\Delta u$ are the gas drag coefficient, the gas density, and the relative velocity of the planetesimal to the disk gas. The

relative velocity is expressed by

$$\Delta u \approx \left(\frac{5}{8}e^2 + \frac{1}{2}i^2 + \eta^2\right)^{1/2} v_{\text{K}}. \qquad (A.2)$$

Here we use the random velocity relative to the local Keplerian motion and the gas velocity ($v_{\text{gas}} = [1 - \eta]v_{\text{K}}$). The stopping time is given by $t_{\text{stop}} = M\Delta u/|F_{\text{aero}}|$. Then the drift timescale is

$$t_{\text{drift}} \approx \frac{t_{\text{stop}}}{2|\eta|} \qquad (A.3)$$

$$\begin{aligned} & \approx 2 \times 10^5 \left(\frac{\eta}{10^{-3}}\right)^{-1} \left(\frac{C_D}{2}\right)^{-1} \left(\frac{M}{10^{16} \text{ g}}\right)^{1/3} \left(\frac{\rho}{3 \text{ g cm}^{-3}}\right)^{2/3} \\ & \times \left(\frac{\Sigma_g}{10^3 \text{ g cm}^{-2}}\right)^{-1} \left(\frac{H/r}{0.05}\right) \left(\frac{r}{1 \text{ au}}\right)^{3/2} \\ & \times \left(\frac{\sqrt{\frac{5}{8}e^2 + \frac{1}{2}i^2 + \eta^2}}{0.01}\right)^{-1} \text{yr}, \end{aligned} \qquad (A.4)$$

where $\rho$ is the internal density of planetesimals. In Sect. 3, $e = 0.01$ and $i^2 = e^2/4$ are used. To obtain a more accurate orbital evolution of planetesimals, we need to calculate the evolution of $e$ and $i$, which should be considered in a future work.

Figure A.1 shows the time evolution of $\eta$ for an MRI-inactive disk. A similar plot is shown in Fig. 10 of Suzuki et al. (2016). The direction of radial drift is determined by the sign of $\eta$. Although Eq. (A.3) can be a little more complicated for $|\eta| < e^2$ (see Eq. (4.21) of Adachi et al. 1976), outward drift is clearly seen for $\eta \lesssim -1 \times 10^{-4}$. Thick lines represent the region of negative $\eta$. There exists a region with positive pressure gradient in the inner region ($<1$ au), which move outwards in time. The positive pressure gradient is seen in region between 0.05 and 0.4 au at the initial time, while the region is between 0.2 and 0.6 au at $t = 0.1$ Myr. Therefore, all planetesimals inside ~1 au can undergo outward radial drift before $t = 1$ Myr in results of Sect. 3.