Next Article in Journal
Mapping the ΛsCDM Scenario to f(T) Modified Gravity: Effects on Structure Growth Rate
Next Article in Special Issue
A Scenario for Origin of Global 4 mHz Oscillations in Solar Corona
Previous Article in Journal
Cosmologies with Perfect Fluids and Scalar Fields in Einstein’s Gravity: Phantom Scalars and Nonsingular Universes
Previous Article in Special Issue
X-Ray Views of Galactic Accreting Pulsars in High-Mass X-Ray Binaries
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Effects of Thermodynamics on the Concurrent Accretion and Migration of Gas Giants in Protoplanetary Disks

1
Department of Physics and Astronomy, University of Nevada, Las Vegas, 4505 S. Maryland Parkway, Las Vegas, NV 89154, USA
2
Shanghai Astronomical Observatory, Chinese Academy of Sciences, Shanghai 200030, China
*
Authors to whom correspondence should be addressed.
Universe 2025, 11(1), 1; https://doi.org/10.3390/universe11010001
Submission received: 31 October 2024 / Revised: 20 December 2024 / Accepted: 23 December 2024 / Published: 25 December 2024

Abstract

:
Accretion and migration usually proceed concurrently for giant planet formation in the natal protoplanetary disks. Recent works indicate that the concurrent accretion onto a giant planet imposes significant impact on the planetary migration dynamics in the isothermal regime. In this work, we carry out a series of 2D global hydrodynamical simulations with Athena++ to explore the effect of thermodynamics on the concurrent accretion and migration processes of the planets in a self-consistent manner. The thermodynamics effect is modeled with a thermal relaxation timescale using a β -cooling prescription. Our results indicate that radiative cooling has a substantial effect on the accretion and migration processes of the planet. As cooling timescales increase, we observe a slight decrease in the planetary accretion rate, and a transition from the outward migrating into inward migration. This transition occurs approximately when the cooling timescale is comparable to the local dynamical timescale ( β 1 ), which is closely linked to the asymmetric structures from the circumplanetary disk (CPD) region. The asymmetric structures in the CPD region which appear with an efficient cooling provide a strong positive torque driving the planet migrate outward. However, such a positive torque is strongly suppressed, when the CPD structures tend to disappear with a relatively long cooling timescale ( β 10 ). Our findings may also be relevant to the dynamical evolution of accreting stellar-mass objects embedded in disks around active galactic nuclei.

1. Introduction

With several exoplanet populations have been detected so far, the discovery of the close-in super-Earths [1] and hot Jupiters [2] has prompted the theoretical proposal that planets may migrate significantly in the protoplanetary disk (PPD) phase [3,4], although it has been suggested that the giant planet could also migrate due to a more populated asteroid belt [5]. Another compelling evidence for disk-driven migration is the discovery of multiple planets in or near mean motion resonant configurations [6,7,8,9]. The migration of giant planets has a substantial dynamic impact on the formation regions and final locations of minor bodies (e.g., [10,11]). For the migration of a low-mass protoplanet, the gravitational interaction with the PPD is usually the dominant factor which leads to the classical type I migration [12]. The torque exerted onto the protoplanet can be divided into two parts, the wave torque [13,14], and the corotation torque or horseshoe drag [15,16], which originate from the asymmetries in different regions of the PPD and usually drives the planet to migrate inward [17,18,19,20]. When the protoplanet’s mass is high enough to modify the density distribution of the PPD, the tidal interaction between the protoplanet and the disk leads to the opening of a gap around the protoplanet’s orbit, and the migration is in the realm of type II, which is believed to be locked to the viscous drift velocity of the background disk [13,21].
In the meanwhile, under the widely adopted core accretion model for planet formation, the dynamical runaway gas accretion occurs at a “critical core mass” above which the envelope fails to maintain hydrostatic equilibrium [22,23,24,25]. In the phase of the dynamical gas accretion, the gas supplied by the PPD, due to the conservation of angular momentum, could form a disk around the protoplanet, known as the circumplanetary disk (CPD) [26,27,28,29]. The accretion onto the protoplanets sensitively depend on the structures of the CPD, which finally determine the asymptotic mass of mature planetary system (e.g., [29,30,31,32,33,34,35]). The presence of the deep gap and the strong tidal barrier effect in the type II regime would prevent the gas from reaching the protoplanets’ vicinity [13,21], thus impeding their growth (e.g., [35,36]). However, further studies have shown that the gas in the PPD could still pass through the gap through the horseshoe orbit, so that the protoplanet could accrete gas anyway [29,37,38,39,40,41]. Therefore, a new paradigm for the type II migration emerges, which suggests that the type II migration rate is not at the viscous velocity of the disk, but shows a similar scaling law as the type I case after considering the steady gap structure [42].
It has been shown that the accretion process in the CPD is mainly contributed by spiral shock dissipation, which is sensitive to the disk thermodynamics [43]. While most of the studies for the CPD structures and the gas dynamical accretion adopted an isothermal equation of state (EOS), there have been some works including non-isothermal effect and radiation hydrodynamics for the CPD formation (e.g., [28,44,45]). Different treatments of the thermodynamics could also influence the circumbinary disk accretion and dynamics (e.g., [46,47,48]). Furthermore, the disk thermodynamics could alter the global disk structure [49,50,51,52] and thus the planet migration dynamics significantly [53,54,55]. However, the accretion and migration of planets in these works are typically studied in isolation. In this sense, the coupled evolution between the CPD accretion and planetary accretion are not treated in a self-consistent manner, as also highlighted in Li et al. [56].
Recently, it has been found that the concurrent accretion would imposes strong impact on the migration of the planets in PPDs when taking into account of a fully self-consistent treatment for an accreting planet [56,57]. Especially, Li et al. [56] found that the structures of CPD play an important role on the migration of the planets in the isothermal limit. As the thermodynamics of the disk could modify the CPD structures significantly, it is thus still unclear how different thermodynamics affect the concurrent accretion and migration of planet in PPDs, which is the main goal of this work. Since there have been some theoretical arguments suggesting that there could a population of stellar-mass objects (e.g., stellar-mass black holes, stars) embedded in active galactic nucleus (AGN) disks, our studies should also be relevant to the accretion and migration of stellar-mass objects embedded in AGN disks.
This rest part of the paper is organized as follows. We present the numerical methods in our hydrodynamical simulations in Section 2. The results are shown in Section 3. The conclusions and discussions are presented in Section 4.

2. Numerical Methods

2.1. Hydrodynamical Model

We simulate globally an accreting planet embedded in a thin and non-self-gravitating disk using Athena++ [58]. Most of our simulation setup are similar to Li et al. [56]. Here we only briefly describe the main physics of our models. The simulations are performed in a cylindrical coordinate system ( r , ϕ ), and the origin of the coordinate system is located at the position of the central star with mass M . We initialize a constant aspect ratio of the PPD as h H r = 0.05 over the whole disk, where H is the disk scale height. This leads to an initial temperature profile of T disk r 1.0 . The local sound speed is thus c s = h v K , where v K is the local Keplerian velocity of the disk. The disk has an axis-symmetric initial profile of Σ ( r ) = Σ 0 ( r r 0 ) p with p = 0.5 , where r 0 is the typical length scale of the disk, and Σ 0 is the surface density at r 0 . Due to the high uncertainty from the observations, the exponents of the temperature and surface density profiles exhibit considerable variations for different PPDs at different radial locations [59,60]. Therefore, we adopt specific exponents of the profiles to maintain a steady-state accretion disk, remaining within the acceptable range of observational uncertainties.
The planet is fixed at a circular orbit with a semi-major axis of a p = r 0 in a corotating frame with the planet so that the planet’s azimuthal angle fixed at ϕ p = 0 . We consider a planet with 1 Jupiter mass ( m p = M J ) which corresponds to a mass ratio between the planet and the solar-type star of q m p M = 0.001 . The gravitational potential of the planet at r p is in the form,
Φ p = G m p ( | r p r | 2 + s 2 ) 1 / 2 + q Ω p 2 r p · r ,
where s indicates the softening length of the potential which is set as 0.1 r h unless otherwise stated, and r h is the Hill radius of the planet. Other softening lengths have also been tested for the convergence, as shown in the Appendix A.
We adopt a Shakura-Sunyaev viscosity prescription ν g = α c s H with a constant α = 0.01 across the whole disk [61]. A relative large viscosity parameter is to ensure the disk reach the quasi-viscous steady state within a feasible computation time. The disk, covering a radial range of [ r in , r out ] = [ 0.3 , 3.5 ] r 0 , is resolved with a uniform base grid of 80 and 256 cells in radial and azimuthal directions, respectively. We have implemented four levels of static mesh refinement within a distance of 1 r h from the planet to effectively resolve the CPD region. With such a grid resolution, the Hill radius can be resolved by about 20 cells in each dimension. A test with different grid resolutions have been performed as shown in the Appendix A to confirm the convergence of our simulation results. For all of our simulations, we take the natural units of G = M = r 0 = 1.0 , and the simulation results are scale-free in a non-self-gravitating disk.

2.2. Cooling of the Disk

To explore the effects of thermodynamics, we conduct simulations with locally isothermal and adiabatic EOSs with β -cooling, as adopted in previous works (e.g., [50,52,62]). The EOS for the isothermal case is
P = c s 2 Σ .
The adiabatic EOS is
P = ( γ 1 ) e ,
where P is the pressure, and e is the specific internal energy of the gas. The sound speed in this case is defined as c s = γ P / Σ . In our simulations, we take γ = 1.4 for Equation (3). In addition, We perform the adiabatic EOS simulations with a simple thermal relaxation process, such that the gas internal energy is relaxed towards its initial state on a timescale controlled by the parameter β ,
d E ( t ) d t = Ω β ( E ( t ) E ( 0 ) ) ,
where Ω is the local orbital frequency of the disk. In the calculation, this is equivalent to the following equation in which E is transformed after each time step as follows,
E t + d t = E t ( 1 Ω β d t ) + E 0 Ω β d t ,
where E t + d t and E t are the internal energy of the gas in the grid before and after each time step, E 0 indicates the initial state. So the local cooling timescale t cool of the system is β / Ω . The smaller the value of β , the closer the disk is to the locally isothermal condition. In contrast, a higher value of β brings the disk closer to the adiabatic approximation. Assuming a Minimum Mass Solar Nebula, the cooling timescale is shorter (e.g., β 10 2 at 100 au) at the outer disk and much longer at the inner disk (e.g., β 10 5 at 1 au) [62]. In our numerical simulations, we examine β values across a broad range from 10 2 to 10 2 .

2.3. Planet Accretion and Boundary Conditions

We follow Li et al. [56] to implement the accretion onto the planet and the boundary condition both at the inner and outer edge of the disk. Basically, there are two parameters, which are sink hole radius r acc around the accreting planet, a removal rate f within the sink hole, to describe the dynamical accretion of the planet. The planetary accretion rates m ˙ p are recorded on-the-fly during each simulation, which are calculated as (e.g., [29,30,31,35,56])
m ˙ p = δ r < r acc f Σ d S ,
where Σ , d S , δ r are the surface density, surface area of each grid, and the distance to the planet, respectively. We adopt an accretion radius as the same of the softening length with r acc = 0.1 r h and a removal rate of f = 5 Ω , similar to Li et al. [56]. Other accretion prescriptions have also been tested for the convergence as shown in the Appendix A. The accreted mass and angular momentum are not added actively onto the planet, while their effect on the migration dynamics of the planet is analyzed by post-processing the gravitational and accretional torque components as in Li et al. [56].
For the radial boundary, we maintain a steady inflow m ˙ d from the outer boundary, and preserve a radially constant inward mass flux at the inner edge as in Li et al. [56]. A wave damping is applied at the outer and inner edges to remove the wave reflection [63]. A periodical boundary is adopted at the azimuthal direction.

3. Simulation Results

We first outline the accretion dynamics of the planet, followed by an examination of the effects of planetary migration under various cooling prescriptions.

3.1. The Accretion of the Planet

The time evolution of the planetary accretion rates m ˙ p for isothermal and different β -cooling are shown in Figure 1. All the simulations proceed to 20,000 orbits, which is roughly the viscous timescales at r out (one viscous timescale can be calculated as ( 2 π α h 0 2 ) 1 orbits). The accretion rates at 20,000 orbit for different runs are slightly below the mass supply rate m ˙ d = 3 π Σ ν at the outer boundary, similar to [56]. This suggests that most of the mass supply from the outer disk is accreted onto the planet, with a small fraction of the disk supply accreted onto the central star. The nearly constant accretion rates onto the planet around 20,000 orbits indicates the viscous quasi-steady state for the CPD accretion. We have further checked the mass flux across the Hill sphere around the planet, which is nearly constant as a function of δ r and in good agreement with the accretion rates onto the planet. The accretion rate also matches well with the global disk mass flux jump at the the planetary orbit. All of these confirm that the quasi-steady state for the CPD accretion is well established. For different cooling prescriptions, it can be seen that a higher β results in a slightly decrease in the planetary accretion rate m ˙ p , although strong variabilities appear for high β -cooling runs.
The surface density distribution and the velocity field in the frame corotating with the planet are shown in Figure 2. The gas velocity field reveals the presence of a distinct disk structure in close proximity to the planet, known as the CPD, in the isothermal case. The disk structure is particularly prominent under isothermal condition. Additionally, the horseshoe streamline is clearly visible for the isothermal and low- β cases, but it disappears for β 10 . It is evident that the accreted material is brought into the vicinity of the planet from the upper horseshoe orbits and flows outward from the spiral density waves. Notably, the rotational velocity is in the retrograde direction for β 10 , in contrast to the prograde flow observed in the low β cases.
To facilitate a more comprehensive examination of the CPD structures in isothermal and long cooling timescale cases, we plot the azimuthally-averaged rotation curves and density profiles of the CPDs in Figure 3. Based on the rotation curves presented in the left panel, it is evident that the rotation velocity within 0.2 r h under the isothermal condition closely resembles the Keplerian rotation. This is also true for the cases with β 0.1 . This indicates that the CPD is rotationally supported up to 0.2 r h and exhibits characteristics of a disk structure in these cases, where the centrifugal force effectively counterbalances the gravitational force. In the longer cooling timescale cases, the rotation curve exhibits a significant deviation from the Keplerian curve in all regions within the Hill sphere and the materials even rotate in the retrograde direction for β 10 as confirmed in Figure 2. Consequently, it is the pressure gradient due to the inefficient cooling in the CPD region that balances the gravity. This is consistent with previous 3D simulations which have shown that the rotationally supported CPD appears mostly under the isothermal EOS, while it transformed into a pressure supported envelope for an adiabatic EOS (e.g., [28]). Actually, based on the temperature (sound speed) profile within the Hill sphere shown in the right panel of Figure 2, the Bondi radius r Bondi m p / c s 2 at the expected CPD size from the planet (e.g., 0.2 r h ) is smaller than the Hill radius due to the increase of the temperature in the vicinity of the planet, which suggests that the planetary accretion proceeds indeed by the Bondi-like process.

3.2. Mechanism of the Accretion

To explore the accretion process, we examine the CPD’s angular momentum conservation equation in the vicinity of a planet following Zhu et al. [43],
t Σ δ r v ϕ ϕ = 1 δ r r ( δ r 2 Σ v r v ϕ ϕ ) + δ r × F ϕ ,
where δ r is the radial distance relative to the planet and X ϕ represents π π X d ϕ . The quantities with ′ indicate that they are measured in the frame relative to the planet. Here we have neglected the viscous stress associated with the explicit alpha viscosity in the CPD since it is found that the explicit α viscosity is much smaller than the effective viscosity due to shock dissipation. The sum of all external forces, denoted as F , encompasses various components, including stellar gravity, centrifugal forces and Coriolis forces within the reference frame that orbits the star with the planet’s orbital period. The planet’s gravity does not contribute to the external forces because it is in parallel to the radial vector δ r .
By decomposing the velocity in ϕ direction into two components, the averaged velocity v ¯ (it is worth noting that we use the azimuthally averaged velocity in the ϕ direction in instead of the local Keplerian velocity due to the distinct deviation from Keplerian motion under long cooling timescales) and the disturbance velocity δ v ϕ , and afterwards multiplying both sides of Equation (7) by δ r 2 to ensure angular momentum dimensions, the equation below can be obtained,
δ r 2 t Σ δ r v ϕ ϕ = M ˙ δ r r ( δ r v ¯ ) δ r r ( δ r 2 Σ v r v ϕ ϕ ) + δ r 2 δ r × F ϕ ,
where M ˙ CPD = Σ δ r v r ϕ . It turns out the mass flux in the CPD M ˙ CPD well matches with the planetary accretion rate m ˙ p , which again suggests a steady state for the CPD accretion. The term on the left side represents the time-dependent evolution term of angular momentum, the first term on the right side denotes the angular momentum transfer resulting from accretion, the second term corresponds to the gradient of Reynolds stress or the gradient of angular momentum transfer flow, and the last term represents the external torque.
When the state of accretion equilibrium for the CPD has been established as in our cases, the time-dependent evolution term in Equation (8) is negligible. Consequently, Equation (8) can be expressed as follows,
M ˙ CPD δ r r ( δ r v ¯ ) δ r r ( δ r 2 Σ v r v ϕ ϕ ) + δ r 2 δ r × F ϕ .
By performing integration on both sides of Equation (9), the components of the aforementioned equation is written as,
α eff = α rey + α con + α T ,
where
α eff = M ˙ CPD 3 π Σ c s H ,
α rey = Σ v r δ v ϕ ϕ 3 π Σ c s 2 ,
α con = C 3 π Σ c s H δ r v ¯ ,
α T = δ r δ r × F ϕ d r 3 π Σ c s H δ r v ¯ ,
where α eff represents the effective accretion coefficient of the CPD, while α rey denotes the accretion coefficient resulting from Reynolds stress, α T is the accretion coefficient induced by external force. α con is associated with the accretion coefficient determined by a constant free parameter C associated with the accretion boundary of the CPD.
As shown in Figure 4, it is evident that the total accretion coefficient ( α rey + α T + α con ) by combing different components can be in reasonable agreement with the effective α eff based on measured planetary accretion rate m ˙ p for both isothermal and β = 100 cases. It should be noted that due to the softening lengths and the inner accreting boundary in the vicinity of the planet, the two curves do not align perfectly. Therefore, the consistency of Equation (10) is another signature of a steady state within the CPD region. Furthermore, we can see that α eff is quite close to α rey in most parts of the disk, and both α eff and α rey is much higher than the explicit α viscosity implemented in the global disk. These suggest that the angular momentum transport and thus the accretion within the CPD could be mainly induced by spiral shocks in the CPD regions instead of the global disk viscosity [43].
The effective accretion coefficients α eff for different β -cooling models are shown in the left panel of Figure 5, illustrating its variation across various thermodynamic factors. It reveals that the effective accretion coefficient is highest under the isothermal condition and progressively decreases as the cooling timescale increases, ultimately reaching its minimum under adiabatic conditions. As the net accretion rates onto the planet vary slightly for different β -cooling, the decreasing of α eff is partly compensated by a higher CPD temperature resulting from the adiabatic heating with a longer cooling timescale, as shown in the right panel of Figure 3. Nevertheless, the mass flux across the CPD and the accretion rate onto the planet under the isothermal limit are still the highest one after considering the lowest disk temperature for allowing the formation of a rotationally supported disk.
Furthermore, the Reynolds stress coefficient α rey and the external torque coefficient α T under different cooling timescales are represented in the right panel of Figure 5. It is evident that the α T is anti-correlated with the effective accretion coefficient α eff among different thermodynamical models. The Reynolds stress coefficient α rey , however, exhibits noticeable positive correlation with α eff for different cooling timescales, with the highest values observed under isothermal conditions. It, therefore, again suggests that the CPD accretion for different thermodynamic models is mostly attribute to the Reynolds stress.
In the long cooling timescale case, the gas is unable to efficiently dissipate the compressional heat, leading to an elevated temperature in the vicinity of the planet, as shown in the right panel of Figure 3. The sound speed appeared in Equation (12) is larger compared to the isothermal limit. In addition, as the azimuthal motion relative to the planet is significantly suppressed with a longer cooling timescale, the Reynolds stress in Equation (12) also becomes weaker. As a result, the value of α rey is reduced with the increasing of the cooling timescale, and the shock-driven accretion becomes less efficient.

3.3. The Migration of the Planet

To study the effects of thermodynamics on the migration of the planets, we calculate the torques exerted onto the planet. The torques can be divided into two parts, the first part is from the gravitational force from the PPD, the other one is from the mass and angular momentum transfer in the sink hole area, even though we do not actively add them onto the planet in our simulations [35,64,65]. Considering the small mass ratio between the planet and star ( q = 10 3 ), the gravitational torque and accretional torque associated with migration can be calculated as follows [56],
Γ gra = r p × Σ Φ p d S ,
and
Γ acc = r p × δ r < r acc ( v v p ) d m ˙ p ,
where Φ p is the gravitational potential of the planet, v is the velocity of the gas within the sink hole radius in the inertial frame. The migration rate for the planet’s semi-major axis a p can be expressed as
a ˙ p a p = 2 p p = 2 ( Γ gra + Γ acc ) m p p ,
where p = G M a p is the specific angular momentum of the planet. Here we have neglected the accretion term onto the central star M ˙ / M as this is usually negligible.
The migration rates of the semi-major axis a p for two typical models are shown in Figure 6. It is evident that both the isothermal case and the one with β = 100 achieve a quasi-steady state at approximately 20,000 orbits. In the isothermal case, the planet migrates outward, whereas it migrates inward in the case of β = 100 .
We present the results for the migration torques for all different models in Figure 7. The migration torque is usually dominated by the gravitational component except for the cases with β 10 , where the accretional component becomes comparable to the gravitational torque; however, it is still not the primary factor determining the direction of migration. It can be seen that as the cooling timescale increases, the gravitational torque Γ gra and also the total torque Γ tot drop monotonically from positive to negative, with Γ tot > 0 for β 1 and Γ tot < 0 for β 1 . This suggests that when the cooling timescale exceeds the local orbital timescale (i.e., β 1 ), migration transitions to inward. The positive torque for an accreting Jupiter-mass planet in the isothermal case with a high disk viscosity is consistent with the results from Li et al. [56], Laune et al. [57].
In the type I/II regime for non-accreting planets, the inward migration rate is a ˙ p / a p 4 q / h 0 2 Σ p r 0 2 Ω 0 / M 500 m ˙ d / M [18], where the Σ p is the surface density at the planetary orbit after considering the gap opening effect [42]. We can see that the inward migration rate in the high β -cooling cases is a factor of ∼3 slower than that of the non-accreting planet. In the adiabatic limit for the non-accreting planet, the migration torque can be positive for our density and temperature profiles due to the non-linear corotation torque [54]. With a long β -cooling timescale, which is close to the adiabatic case, the negative gravitational torque on the accreting planet is due to the suppression of the positive corotation torque associated with accretion.

3.4. How Thermodynamics Affects Migration

The gravitational torque could be the dominant factor that determines the planet migration. To understand what causes the transition of the migration, we explore the influence of thermodynamics on the gravitational torque. As the regions where the Lindblad torque dominates and the corotation torque dominates overlap with each other, we do not distinguish their contributions in this paper. Following Li et al. [56], we plot the spatial torque distribution in the vicinity of the planet and the results are shown in Figure 8. The CPD structure is clear in the isothermal case in which the leading horseshoe region provides positive torque (red color region) while the trailing region provides negative torque (blue color). In the β = 100 case, the disk-like structures break down due to thermal pressure and the torque within the Hill sphere is highly symmetrical with respect to the ϕ = 0 plane. We then calculate the asymmetrical part of the torque defined as Γ asy = Γ ( r , ϕ ) + Γ ( r , ϕ ) as shown in the right panels in the figure. It can be seen that in the CPD region for the isothermal case, the positive torque from the ϕ > 0 hemisphere is stronger than the negative torque from the ϕ < 0 hemisphere so the net torque is positive, as expected from the asymmetric spiral arms feeding from the global disk into the Hill sphere [56]. While for the high β -cooling case, the net positive torque in the vicinity of the planet almost vanishes, due to that the asymmetric disk-like structure is replaced by a spherical-like envelope within the Hill sphere.
The asymmetry of the gravitational torque comes from the asymmetry of the density wave structure in the vicinity of the planet. We plot the surface density distribution around the planet in Figure 2. The asymmetry in the spiral arms feeding the planet occurs in the isothermal case, but vanishes in the β = 100 case. In the accretion process of the planet, the materials from the outer regions of the disk approach the planet through outer horseshoe tracks and accumulate at the CPD region. Part of the materials are accreted onto the planet with the left materials moving to the inner horseshoe tracks and will enter the CPD region again after one orbit. Therefore, the asymmetries in density and torque originate from the accretion process of the planet. In contrast, for the β = 100 case, the CPD structure is destroyed and supported by the thermal pressure gradient, which tends to be nearly isotropic, such that the accretion tends to be Bondi-like.
The cumulative radial torque profiles by azimuthally averaging the 2D spatial gravitational torque distribution > r ϕ Γ ( r , ϕ ) d r d ϕ as a function of radius for two different models are shown in Figure 9. We can see that the torque outside the Hill sphere ( r < r p r h , r > r p + r h ), which is dominated by the Lindblad resonance, is negative both in isothermal and long cooling timescales. These differential Lindblad torque usually drives the planet migrate inward. However, as shown in Li et al. [56], the gravitational torque within the planet’s Hill sphere is significantly positive in cases close to the isothermal EOS, allowing it to counterbalance the negative differential Lindblad torque. This is found to be associated with the asymmetric structures in the CPD region as discussed above. As the cooling timescale becomes longer (e.g., β = 100 ), a rotationally supported CPD disappears, the positive gravitational torque from the CPD region thus becomes too weak to balance the negative differential Lindblad torque, in line with Figure 8. As a result, the planet migrates inward.

4. Conclusions and Discussion

In this work, we perform a series of 2D global hydrodynamical simulations using Athena++ to study the impact of thermodynamics on gas giants’ concurrent accretion and migration embedded within their natal disks. The disk is evolved over the viscous timescale with a desired mass flux from the outer boundary, allowing the system to reach a quasi-steady state in a self-consistent manner. A parameterized thermal relaxation process with β -cooling is implemented to simulated the critical role of disk thermodynamic on the accretion and migration dynamics of the embedded object. The planetary accretion is modeled as a sink particle, with a high resolution around the planet to resolve the CPD region. Our findings are summarized as follows:
  • The planetary accretion rates for varying cooling timescales are primarily constrained by the disk supply rate from the outer disk, exhibiting a slight decrease with increasing cooling timescales, even though different thermodynamical states have a significant impact on the CPD structures.
  • The CPD is much smaller in size and less rotationally supported with longer cooling timescales. When the cooling timescale is an order of magnitude longer than the local dynamical timescale (e.g., β 10 ), the rotation of the CPD can even shift to a retrograde orbit, which could strongly influence satellite formation within the CPD.
  • The planetary accretion is chiefly driven by spiral shock dissipation in the CPD region. We observe that the Reynolds stress coefficient varies significantly with different cooling timescales, following a trend similar to that of the effective accretion coefficients, making it a dominant factor in controlling the accretion process onto the planet.
  • We confirm that in an isothermal disk, the planet migrates outward for a high disk viscosity, consistent with recent findings by Li et al. [56], Laune et al. [57]. However, there is a tendency for inward migration as the cooling timescale increases ( β 1 ), different from the general trend where non-accreting planets migrate outward in the adiabatic limit.
  • The transition in migration under the longer cooling timescales is closely linked to the gravitational torque from the co-orbital region of the planet. In longer cooling timescales, the positive gravitational torque from this region is significantly suppressed due to the absence of a rotationally supported CPD. These results highlight the critical role of thermodynamics in the accretion and migration of planets in PPDs.
The physical process in this work could be also relevant to the dynamics of the accreting stellar-mass black holes/compact objects embedded in AGN disks. Their dynamical evolution in AGN disks could provide be an important channel for binary black hole mergers detected by LIGO/Virgo/KAGRA gravitational wave observations (e.g., [46,64,65]). In this case, the mechanical and radiative feedback from the embedded compact objects could be significant to modify the circumsingle disk structures, and thus their dynamical evolution. In addition, the inner solar system are subjected to extreme radiation environments from young stars [66]. Massive asymptotic giant branch stars could also play a significant role in delivering short lived radionuclides to heat rock-forming materials in the disk [67]. All of these would necessitate a more self-consistent approach to cooling and heating with radiative hydrodynamics simulations for the CPD, rather than relying on the β -cooling prescription. Second, a laminar viscous flow with a α -viscosity is assumed in our work. Magnetic turbulence in disks, which could modify the coherent spiral waves lunched by the embedded objects, may also play an important role in the dynamical evolution of embedded objects. Third, our focus has been primarily on a parameter space involving a Jupiter mass planet embedded in a highly viscous disk fixed in a circular orbit. A systematic exploration of different planet masses, disk viscosities, disk scale heights, and planetary eccentricities is still needed, which could be crucial for future population synthesis studies. Lastly, our simulations are limited to 2D cases. While our previous studies have shown that 3D simulations yield similar accretion and migration dynamics as 2D cases with an isothermal EOS [56], the extent of the 3D effects for our adiabatic equation of state with varying β -cooling remains unclear. All these aspects will be addressed in future work.

Author Contributions

Conceptualization, methodology, supervision, funding acquisition, Y.-P.L.; writing—original draft preparation, H.W.; software, validation, formal analysis, visualization, investigation, resources, writing—review and editing, H.W. and Y.-P.L. All authors have read and agreed to the published version of the manuscript.

Funding

This work is supported in part by the Natural Science Foundation of China (grants Nos. 12373070, and 12192223), the Natural Science Foundation of Shanghai (grant No. 23ZR1473700).

Data Availability Statement

All the data that are used in this article are available from the corresponding author upon reasonable request.

Acknowledgments

We thank the referees for many constructive suggestions. We also thank beneficial discussions with Zhaohuan Zhu. The calculations have made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Shanghai Astronomical Observatory.

Conflicts of Interest

The authors declare no conflicts of interest.

Appendix A

Here we test the convergence of the accretion prescription and numerical resolution for a model in the transition region (i.e., β = 1 ), mainly focused on their impact on the accretion dynamics. There are two parameters for the planetary accretion, the removal rate f, and the accretion radius (or sink hole radius) r acc .
Figure A1. The accretion rates for different parameters. The number following the letter ‘b’ indicates the value of β , ‘l’ refers to the mesh refinement level, ‘r’ denotes the sink hole radius, which is the same as the softening radius, and ‘f’ represents the removal rate of the sink hole region, e.g., ‘b1l4r01f5’ indicates β = 1 , 4 levels of mesh refinement, r acc = 0.1 r h , and f = 5 . It is clear that all of these parameters have a negligible effect on planetary accretion rates.
Figure A1. The accretion rates for different parameters. The number following the letter ‘b’ indicates the value of β , ‘l’ refers to the mesh refinement level, ‘r’ denotes the sink hole radius, which is the same as the softening radius, and ‘f’ represents the removal rate of the sink hole region, e.g., ‘b1l4r01f5’ indicates β = 1 , 4 levels of mesh refinement, r acc = 0.1 r h , and f = 5 . It is clear that all of these parameters have a negligible effect on planetary accretion rates.
Universe 11 00001 g0a1
The sink hole radius r acc should be small enough to represent the size of the accreting object accurately. However, using a sink hole radius comparable to that of a Jupiter-mass planet necessitates an extremely high resolution around the planet in our global simulations, which is impractical even with mesh refinement. Therefore, it is advisable to select an appropriate value for r acc that does not significantly affect the dynamical accretion of the embedded object. By comparing the simulation results with r acc = 0.1 r h , 0.05 r h , and different removal rates f, we observe no significant differences in the accretion rate, as illustrated in Figure A1. This indicates good convergence with the accretion prescription.
In addition, we have conducted convergence tests at various numerical resolutions in the vicinity of the planet. With 4 and 5 levels of mesh refinement around the planet, we found that the differences in accretion rates among these cases are negligible. Similar tests were performed for different accretion prescriptions and global numerical resolutions (up to n r = 1800 , and n ϕ = 2500 ) without any mesh refinement using FARGO3D [68]. All of these indicate satisfactory convergence for our simulation results.

References

  1. Borucki, W.J.; Koch, D.G.; Basri, G.; Batalha, N.; Brown, T.M.; Bryson, S.T.; Caldwell, D.; Christensen-Dalsgaard, J.; Cochran, W.D.; DeVore, E.; et al. Characteristics of Planetary Candidates Observed by Kepler. II. Analysis of the First Four Months of Data. Astrophys. J. 2011, 736, 19. [Google Scholar] [CrossRef]
  2. Mayor, M.; Queloz, D. A Jupiter-mass companion to a solar-type star. Nature 1995, 378, 355–359. [Google Scholar] [CrossRef]
  3. Lin, D.N.C.; Papaloizou, J. On the Tidal Interaction between Protoplanets and the Protoplanetary Disk. III. Orbital Migration of Protoplanets. Astrophys. J. 1986, 309, 846. [Google Scholar] [CrossRef]
  4. Levison, H.F.; Morbidelli, A.; Gomes, R.; Backman, D. Planet Migration in Planetesimal Disks. In Protostars and Planets V; Reipurth, B., Jewitt, D., Keil, K., Eds.; University of Arizona Press: Tucson, AZ, USA, 2007; p. 669. [Google Scholar]
  5. Clement, M.S.; Morbidelli, A.; Raymond, S.N.; Kaib, N.A. A record of the final phase of giant planet migration fossilized in the asteroid belt’s orbital structure. Mon. Not. R. Astron. Soc. 2020, 492, L56–L60. [Google Scholar] [CrossRef]
  6. Lissauer, J.J.; Fabrycky, D.C.; Ford, E.B.; Borucki, W.J.; Fressin, F.; Marcy, G.W.; Orosz, J.A.; Rowe, J.F.; Torres, G.; Welsh, W.F.; et al. A closely packed system of low-mass, low-density planets transiting Kepler-11. Nature 2011, 470, 53–58. [Google Scholar] [CrossRef] [PubMed]
  7. Goldreich, P.; Schlichting, H.E. Overstable Librations can Account for the Paucity of Mean Motion Resonances among Exoplanet Pairs. Astron. J. 2014, 147, 32. [Google Scholar] [CrossRef]
  8. Xu, W.; Lai, D.; Morbidelli, A. Migration of planets into and out of mean motion resonances in protoplanetary discs: Overstability of capture and non-linear eccentricity damping. Mon. Not. R. Astron. Soc. 2018, 481, 1538–1549. [Google Scholar] [CrossRef]
  9. Yang, H.; Li, Y.P. Mean-motion resonances with interfering density waves. Mon. Not. R. Astron. Soc. 2024, 534, 485–501. [Google Scholar] [CrossRef]
  10. Brownlee, D.; Tsou, P.; Aléon, J.; Alexander, C.M.O.D.; Araki, T.; Bajt, S.; Baratta, G.A.; Bastien, R.; Bland, P.; Bleuet, P.; et al. Comet 81P/Wild 2 Under a Microscope. Science 2006, 314, 1711. [Google Scholar] [CrossRef] [PubMed]
  11. Huang, H.; Wu, Y. Review of Jupiter-Trojan asteroids research. Rev. Geophys. Planet. Phys. 2024, 55, 175–183. [Google Scholar] [CrossRef]
  12. Ward, W.R. Protoplanet Migration by Nebula Tides. Icarus 1997, 126, 261–281. [Google Scholar] [CrossRef]
  13. Goldreich, P.; Tremaine, S. Disk-satellite interactions. Astrophys. J. 1980, 241, 425–441. [Google Scholar] [CrossRef]
  14. Artymowicz, P. On the Wave Excitation and a Generalized Torque Formula for Lindblad Resonances Excited by External Potential. Astrophys. J. 1993, 419, 155. [Google Scholar] [CrossRef]
  15. Ward, W.R. Horsehoe Orbit Drag. Lunar Planet. Sci. Conf. 1991, 22, 1463. [Google Scholar]
  16. Masset, F.S. On the Co-orbital Corotation Torque in a Viscous Disk and Its Impact on Planetary Migration. Astrophys. J. 2001, 558, 453–462. [Google Scholar] [CrossRef]
  17. Papaloizou, J.; Lin, D.N.C. On the tidal interaction between protoplanets and the primordial solar nebula. I—Linear calculation of the role of angular momentum exchange. Astrophys. J. 1984, 285, 818–834. [Google Scholar] [CrossRef]
  18. Tanaka, H.; Takeuchi, T.; Ward, W.R. Three-Dimensional Interaction between a Planet and an Isothermal Gaseous Disk. I. Corotation and Lindblad Torques and Planet Migration. Astrophys. J. 2002, 565, 1257–1274. [Google Scholar] [CrossRef]
  19. Kley, W.; Nelson, R.P. Planet-Disk Interaction and Orbital Evolution. Ann. Rev. Astron. Astrophys. 2012, 50, 211–249. [Google Scholar] [CrossRef]
  20. Paardekooper, S.; Dong, R.; Duffell, P.; Fung, J.; Masset, F.S.; Ogilvie, G.; Tanaka, H. Planet-Disk Interactions and Orbital Evolution. In Protostars and Planets VII; Inutsuka, S., Aikawa, Y., Muto, T., Tomida, K., Tamura, M., Eds.; Astronomical Society of the Pacific Conference Series; Astronomical Society of the Pacific: San Francisco, CA, USA, 2023; Volume 534, p. 685. [Google Scholar] [CrossRef]
  21. Lin, D.N.C.; Papaloizou, J. On the Tidal Interaction between Protoplanets and the Primordial Solar Nebula. II. Self-Consistent Nonlinear Interaction. Astrophys. J. 1986, 307, 395. [Google Scholar] [CrossRef]
  22. Bodenheimer, P.; Pollack, J.B. Calculations of the accretion and evolution of giant planets: The effects of solid cores. Icarus 1986, 67, 391–408. [Google Scholar] [CrossRef]
  23. Pollack, J.B.; Hubickyj, O.; Bodenheimer, P.; Lissauer, J.J.; Podolak, M.; Greenzweig, Y. Formation of the Giant Planets by Concurrent Accretion of Solids and Gas. Icarus 1996, 124, 62–85. [Google Scholar] [CrossRef]
  24. Ida, S.; Lin, D.N.C. Toward a Deterministic Model of Planetary Formation. I. A Desert in the Mass and Semimajor Axis Distributions of Extrasolar Planets. Astrophys. J. 2004, 604, 388–413. [Google Scholar] [CrossRef]
  25. Chen, Y.X.; Li, Y.P.; Li, H.; Lin, D.N.C. The Preservation of Super-Earths and the Emergence of Gas Giants after Their Progenitor Cores Have Entered the Pebble-isolation Phase. Astrophys. J. 2020, 896, 135. [Google Scholar] [CrossRef]
  26. Machida, M.N.; Kokubo, E.; Inutsuka, S.i.; Matsumoto, T. Angular Momentum Accretion onto a Gas Giant Planet. Astrophys. J. 2008, 685, 1220–1236. [Google Scholar] [CrossRef]
  27. Tanigawa, T.; Ohtsuki, K.; Machida, M.N. Distribution of Accreting Gas and Angular Momentum onto Circumplanetary Disks. Astrophys. J. 2012, 747, 47. [Google Scholar] [CrossRef]
  28. Fung, J.; Zhu, Z.; Chiang, E. Circumplanetary Disk Dynamics in the Isothermal and Adiabatic Limits. Astrophys. J. 2019, 887, 152. [Google Scholar] [CrossRef]
  29. Li, Y.P.; Chen, Y.X.; Lin, D.N.C. 3D global simulations of accretion onto gap-opening planets: Implications for circumplanetary disc structures and accretion rates. Mon. Not. R. Astron. Soc. 2023, 526, 5346–5364. [Google Scholar] [CrossRef]
  30. Kley, W.; D’Angelo, G.; Henning, T. Three-dimensional Simulations of a Planet Embedded in a Protoplanetary Disk. Astrophys. J. 2001, 547, 457–464. [Google Scholar] [CrossRef]
  31. D’Angelo, G.; Kley, W.; Henning, T. Orbital Migration and Mass Accretion of Protoplanets in Three-dimensional Global Computations with Nested Grids. Astrophys. J. 2003, 586, 540–561. [Google Scholar] [CrossRef]
  32. Machida, M.N.; Kokubo, E.; Inutsuka, S.I.; Matsumoto, T. Gas accretion onto a protoplanet and formation of a gas giant planet. Mon. Not. R. Astron. Soc. 2010, 405, 1227–1243. [Google Scholar] [CrossRef]
  33. Bodenheimer, P.; D’Angelo, G.; Lissauer, J.J.; Fortney, J.J.; Saumon, D. Deuterium Burning in Massive Giant Planets and Low-mass Brown Dwarfs Formed by Core-nucleated Accretion. Astrophys. J. 2013, 770, 120. [Google Scholar] [CrossRef]
  34. Choksi, N.; Chiang, E.; Fung, J.; Zhu, Z. The maximum accretion rate of a protoplanet: How fast can runaway be? Mon. Not. R. Astron. Soc. 2023, 525, 2806–2819. [Google Scholar] [CrossRef]
  35. Li, Y.P.; Chen, Y.X.; Lin, D.N.C.; Zhang, X. Accretion of Gas Giants Constrained by the Tidal Barrier. Astrophys. J. 2021, 906, 52. [Google Scholar] [CrossRef]
  36. Dobbs-Dixon, I.; Li, S.L.; Lin, D.N.C. Tidal Barrier and the Asymptotic Mass of Proto-Gas Giant Planets. Astrophys. J. 2007, 660, 791–806. [Google Scholar] [CrossRef]
  37. Duffell, P.C.; Haiman, Z.; MacFadyen, A.I.; D’Orazio, D.J.; Farris, B.D. The Migration of Gap-opening Planets is Not Locked to Viscous Disk Evolution. Astrophys. J. Lett. 2014, 792, L10. [Google Scholar] [CrossRef]
  38. Fung, J.; Shi, J.M.; Chiang, E. How Empty are Disk Gaps Opened by Giant Planets? Astrophys. J. 2014, 782, 88. [Google Scholar] [CrossRef]
  39. Dürmann, C.; Kley, W. Migration of massive planets in accreting disks. Astron. Astrophys. 2015, 574, A52. [Google Scholar] [CrossRef]
  40. Robert, C.M.T.; Crida, A.; Lega, E.; Méheut, H.; Morbidelli, A. Toward a new paradigm for Type II migration. Astron. Astrophys. 2018, 617, A98. [Google Scholar] [CrossRef]
  41. Chen, Y.X.; Zhang, X.; Li, Y.P.; Li, H.; Lin, D.N.C. Retention of Long-period Gas Giant Planets: Type II Migration Revisited. Astrophys. J. 2020, 900, 44. [Google Scholar] [CrossRef]
  42. Kanagawa, K.D.; Tanaka, H.; Szuszkiewicz, E. Radial Migration of Gap-opening Planets in Protoplanetary Disks. I. The Case of a Single Planet. Astrophys. J. 2018, 861, 140. [Google Scholar] [CrossRef]
  43. Zhu, Z.; Ju, W.; Stone, J.M. Shock-driven Accretion in Circumplanetary Disks: Observables and Satellite Formation. Astrophys. J. 2016, 832, 193. [Google Scholar] [CrossRef]
  44. Szulágyi, J.; Masset, F.; Lega, E.; Crida, A.; Morbidelli, A.; Guillot, T. Circumplanetary disc or circumplanetary envelope? Mon. Not. R. Astron. Soc. 2016, 460, 2853–2861. [Google Scholar] [CrossRef]
  45. Krapp, L.; Kratter, K.M.; Youdin, A.N.; Benítez-Llambay, P.; Masset, F.; Armitage, P.J. A Thermodynamic Criterion for the Formation of Circumplanetary Disks. Astrophys. J. 2024, 973, 153. [Google Scholar] [CrossRef]
  46. Li, Y.P.; Dempsey, A.M.; Li, H.; Li, S.; Li, J. Hot Circumsingle Disks Drive Binary Black Hole Mergers in Active Galactic Nucleus Disks. Astrophys. J. Lett. 2022, 928, L19. [Google Scholar] [CrossRef]
  47. Sudarshan, P.; Penzlin, A.B.T.; Ziampras, A.; Kley, W.; Nelson, R.P. How cooling influences circumbinary discs. Astron. Astrophys. 2022, 664, A157. [Google Scholar] [CrossRef]
  48. Wang, H.Y.; Bai, X.N.; Lai, D. On the Role of Dynamical Cooling in the Dynamics of Circumbinary Disks. Astrophys. J. 2023, 943, 175. [Google Scholar] [CrossRef]
  49. Miranda, R.; Rafikov, R.R. Planet-Disk Interaction in Disks with Cooling: Basic Theory. Astrophys. J. 2020, 892, 65. [Google Scholar] [CrossRef]
  50. Zhang, S.; Zhu, Z. The effects of disc self-gravity and radiative cooling on the formation of gaps and spirals by young planets. Mon. Not. R. Astron. Soc. 2020, 493, 2287–2305. [Google Scholar] [CrossRef]
  51. Ziampras, A.; Nelson, R.P.; Rafikov, R.R. Modelling planet-induced gaps and rings in ALMA discs: The role of in-plane radiative diffusion. Mon. Not. R. Astron. Soc. 2023, 524, 3930–3947. [Google Scholar] [CrossRef]
  52. Zhang, M.; Huang, P.; Dong, R. The Dependence of the Structure of Planet-opened Gaps in Protoplanetary Disks on Radiative Cooling. Astrophys. J. 2024, 961, 86. [Google Scholar] [CrossRef]
  53. Masset, F.S.; Casoli, J. On the Horseshoe Drag of a Low-Mass Planet. II. Migration in Adiabatic Disks. Astrophys. J. 2009, 703, 857–876. [Google Scholar] [CrossRef]
  54. Paardekooper, S.J.; Baruteau, C.; Crida, A.; Kley, W. A torque formula for non-isothermal type I planetary migration—I. Unsaturated horseshoe drag. Mon. Not. R. Astron. Soc. 2010, 401, 1950–1964. [Google Scholar] [CrossRef]
  55. Lega, E.; Crida, A.; Bitsch, B.; Morbidelli, A. Migration of Earth-sized planets in 3D radiative discs. Mon. Not. R. Astron. Soc. 2014, 440, 683–695. [Google Scholar] [CrossRef]
  56. Li, Y.P.; Chen, Y.X.; Lin, D.N.C. Concurrent Accretion and Migration of Giant Planets in Their Natal Disks with Consistent Accretion Torque. Astrophys. J. 2024, 971, 130. [Google Scholar] [CrossRef]
  57. Laune, J.; Li, R.; Lai, D. Migration of Accreting Planets and Black Holes in Disks. arXiv 2024, arXiv:2405.00296. [Google Scholar] [CrossRef]
  58. Stone, J.M.; Tomida, K.; White, C.J.; Felker, K.G. The Athena++ Adaptive Mesh Refinement Framework: Design and Magnetohydrodynamic Solvers. Astrophys. J. Suppl. Ser. 2020, 249, 4. [Google Scholar] [CrossRef]
  59. Andrews, S.M. Observations of Protoplanetary Disk Structures. Ann. Rev. Astron. Astrophys. 2020, 58, 483–528. [Google Scholar] [CrossRef]
  60. Miotello, A.; Kamp, I.; Birnstiel, T.; Cleeves, L.C.; Kataoka, A. Setting the Stage for Planet Formation: Measurements and Implications of the Fundamental Disk Properties. In Protostars and Planets VII; Inutsuka, S., Aikawa, Y., Muto, T., Tomida, K., Tamura, M., Eds.; Astronomical Society of the Pacific Conference Series; Astronomical Society of the Pacific: San Francisco, CA, USA, 2023; Volume 534, p. 501. [Google Scholar] [CrossRef]
  61. Shakura, N.I.; Sunyaev, R.A. Black holes in binary systems. Observational appearance. Astron. Astrophys. 1973, 24, 337–355. [Google Scholar]
  62. Zhu, Z.; Dong, R.; Stone, J.M.; Rafikov, R.R. The Structure of Spiral Shocks Excited by Planetary-mass Companions. Astrophys. J. 2015, 813, 88. [Google Scholar] [CrossRef]
  63. de Val-Borro, M.; Edgar, R.G.; Artymowicz, P.; Ciecielag, P.; Cresswell, P.; D’Angelo, G.; Delgado-Donate, E.J.; Dirksen, G.; Fromang, S.; Gawryszczak, A.; et al. A comparative study of disc-planet interaction. Mon. Not. R. Astron. Soc. 2006, 370, 529–558. [Google Scholar] [CrossRef]
  64. Li, Y.P.; Dempsey, A.M.; Li, S.; Li, H.; Li, J. Orbital Evolution of Binary Black Holes in Active Galactic Nucleus Disks: A Disk Channel for Binary Black Hole Mergers? Astrophys. J. 2021, 911, 124. [Google Scholar] [CrossRef]
  65. Li, R.; Lai, D. Hydrodynamical evolution of black-hole binaries embedded in AGN discs. Mon. Not. R. Astron. Soc. 2022, 517, 1602–1624. [Google Scholar] [CrossRef]
  66. Ribas, I.; Guinan, E.F.; Güdel, M.; Audard, M. Evolution of the Solar Activity over Time and Effects on Planetary Atmospheres. I. High-Energy Irradiances (1-1700 Å). Astrophys. J. 2005, 622, 680–694. [Google Scholar] [CrossRef]
  67. Trigo-Rodríguez, J.M.; García-Hernández, D.A.; Lugaro, M.; Karakas, A.I.; van Raai, M.; García Lario, P.; Manchado, A. The role of massive AGB stars in the early solar system composition. M&PS 2009, 44, 627–641. [Google Scholar] [CrossRef]
  68. Benítez-Llambay, P.; Masset, F.S. FARGO3D: A New GPU-oriented MHD Code. Astrophys. J. Suppl. Ser. 2016, 223, 11. [Google Scholar] [CrossRef]
Figure 1. The left panel displays the running time-averaged planetary accretion rates in a locally isothermal disk and disks with varying β -cooling timescales, measured in unit of Σ 0 r 0 2 Ω 0 . The black dashed line indicates the material supply rate from the outer boundary. The right panel illustrates the short-term fluctuations in accretion rates for several representative simulations, without any time averaging.
Figure 1. The left panel displays the running time-averaged planetary accretion rates in a locally isothermal disk and disks with varying β -cooling timescales, measured in unit of Σ 0 r 0 2 Ω 0 . The black dashed line indicates the material supply rate from the outer boundary. The right panel illustrates the short-term fluctuations in accretion rates for several representative simulations, without any time averaging.
Universe 11 00001 g001
Figure 2. The density distribution δ Σ in the vicinity of the planet with different cooling timescales. δ Σ = Σ Σ where Σ is the azimuthally averaged surface density from the perturbed disk at the same orbital time. The arrows in the plots represent the velocity vectors relative to the planet. All plots are presented at 20,000 orbits. For the cases with β = 10 and β = 100 , time averages over 100 and 500 orbits, respectively, were performed to achieve steady flow structures. The red dotted lines enclose the classical horseshoe region around the planetary orbit defined as x s = a p ( q / h ) 1 / 2 (e.g., [19]).
Figure 2. The density distribution δ Σ in the vicinity of the planet with different cooling timescales. δ Σ = Σ Σ where Σ is the azimuthally averaged surface density from the perturbed disk at the same orbital time. The arrows in the plots represent the velocity vectors relative to the planet. All plots are presented at 20,000 orbits. For the cases with β = 10 and β = 100 , time averages over 100 and 500 orbits, respectively, were performed to achieve steady flow structures. The red dotted lines enclose the classical horseshoe region around the planetary orbit defined as x s = a p ( q / h ) 1 / 2 (e.g., [19]).
Universe 11 00001 g002
Figure 3. Left panel: the azimuthally-averaged rotation velocity relative to the planet for different models in unit of the local Keplerian velocity orbiting around the star. The black line represents the Keplerian rotation velocity around the planet. Middle panel: the surface density profile around the planet in unit of Σ 0 . Right panel: the azimuthally-averaged sound speed c s in the vicinity of the planet in code unit.
Figure 3. Left panel: the azimuthally-averaged rotation velocity relative to the planet for different models in unit of the local Keplerian velocity orbiting around the star. The black line represents the Keplerian rotation velocity around the planet. Middle panel: the surface density profile around the planet in unit of Σ 0 . Right panel: the azimuthally-averaged sound speed c s in the vicinity of the planet in code unit.
Universe 11 00001 g003
Figure 4. Different accretion coefficient calculated by Equation (10) under isothermal EOS (left panel) and β = 100 (right panel) cases. The the summation of α rey , α con and α T is shown as red lines.
Figure 4. Different accretion coefficient calculated by Equation (10) under isothermal EOS (left panel) and β = 100 (right panel) cases. The the summation of α rey , α con and α T is shown as red lines.
Universe 11 00001 g004
Figure 5. Left panel: the effective accretion coefficients in Equation (11) under different cooling timescales. Right panel: the accretion coefficients resulting from Reynolds stress ( α rey ; solid lines) and external forces ( α T ; dashed lines) for different cooling timescales. The same color represents the same cooling timescale.
Figure 5. Left panel: the effective accretion coefficients in Equation (11) under different cooling timescales. Right panel: the accretion coefficients resulting from Reynolds stress ( α rey ; solid lines) and external forces ( α T ; dashed lines) for different cooling timescales. The same color represents the same cooling timescale.
Universe 11 00001 g005
Figure 6. The semi-major axis evolution due to the gravitational (blue lines), accretional (red lines) and total (black lines) torques in isothermal (left panel) and β = 100 (right panel) cases. For the case of β = 100 , a running-time averaged over 500 orbits is carried out to smooth out the short timescale variability.
Figure 6. The semi-major axis evolution due to the gravitational (blue lines), accretional (red lines) and total (black lines) torques in isothermal (left panel) and β = 100 (right panel) cases. For the case of β = 100 , a running-time averaged over 500 orbits is carried out to smooth out the short timescale variability.
Universe 11 00001 g006
Figure 7. The torque exerted on the planet under different cooling timescales in code unit, where 0 = r 0 2 Ω 0 is the specific angular momentum of the disk at r 0 . The blue and red circle points represent the gravitational torque and accretional torque respectively. The black square points represent the sum of the two torques. The positive (negative) torque means that it drives the planet to migrate outward (inward).
Figure 7. The torque exerted on the planet under different cooling timescales in code unit, where 0 = r 0 2 Ω 0 is the specific angular momentum of the disk at r 0 . The blue and red circle points represent the gravitational torque and accretional torque respectively. The black square points represent the sum of the two torques. The positive (negative) torque means that it drives the planet to migrate outward (inward).
Universe 11 00001 g007
Figure 8. Spacial distribution of gravitational torques in the vicinity of the planet around 20,000 orbits. The upper panel represents the isothermal condition and the lower panel shows the β = 100 case. A time-averaging over 500 orbits is done for the β = 100 case. The contour in the left panels represented by solid and dashed lines being the same magnitude of the torques but with different signs. The vertical dotted lines correspond to the region of 1 r h from the planet. The right panels show the net torque by summing the leading (upper) horseshoe and trailing (lower) horseshoe region, i.e., Γ ( r , ϕ ) + Γ ( r , ϕ ) .
Figure 8. Spacial distribution of gravitational torques in the vicinity of the planet around 20,000 orbits. The upper panel represents the isothermal condition and the lower panel shows the β = 100 case. A time-averaging over 500 orbits is done for the β = 100 case. The contour in the left panels represented by solid and dashed lines being the same magnitude of the torques but with different signs. The vertical dotted lines correspond to the region of 1 r h from the planet. The right panels show the net torque by summing the leading (upper) horseshoe and trailing (lower) horseshoe region, i.e., Γ ( r , ϕ ) + Γ ( r , ϕ ) .
Universe 11 00001 g008
Figure 9. The cumulative radial distribution (integration from r in to r) of the gravitational torque exerted onto the planet, where the blue (red) line represents the isothermal ( β = 100 ) case. The torque for β = 100 is time-averaged over 500 orbits around 20,000 orbits. The two vertical dotted lines show the Hill radius. The torque is in code unit with 0 = r 0 2 Ω 0 being the specific angular momentum of the disk at r 0 .
Figure 9. The cumulative radial distribution (integration from r in to r) of the gravitational torque exerted onto the planet, where the blue (red) line represents the isothermal ( β = 100 ) case. The torque for β = 100 is time-averaged over 500 orbits around 20,000 orbits. The two vertical dotted lines show the Hill radius. The torque is in code unit with 0 = r 0 2 Ω 0 being the specific angular momentum of the disk at r 0 .
Universe 11 00001 g009
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Wu, H.; Li, Y.-P. Effects of Thermodynamics on the Concurrent Accretion and Migration of Gas Giants in Protoplanetary Disks. Universe 2025, 11, 1. https://doi.org/10.3390/universe11010001

AMA Style

Wu H, Li Y-P. Effects of Thermodynamics on the Concurrent Accretion and Migration of Gas Giants in Protoplanetary Disks. Universe. 2025; 11(1):1. https://doi.org/10.3390/universe11010001

Chicago/Turabian Style

Wu, Hening, and Ya-Ping Li. 2025. "Effects of Thermodynamics on the Concurrent Accretion and Migration of Gas Giants in Protoplanetary Disks" Universe 11, no. 1: 1. https://doi.org/10.3390/universe11010001

APA Style

Wu, H., & Li, Y.-P. (2025). Effects of Thermodynamics on the Concurrent Accretion and Migration of Gas Giants in Protoplanetary Disks. Universe, 11(1), 1. https://doi.org/10.3390/universe11010001

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop