Radiative scale-height and shadows in protoplanetary disks
Matías Montesinos, Nicolás Cuello, Johan Olofsson, Jorge Cuadra, Amelia Bayo, Gesa H.-M. Bertrang, Clément Perrot
DD raft version F ebruary
8, 2021Typeset using L A TEX twocolumn style in AASTeX63
Radiative scale-height and shadows in protoplanetary disks M at ´ ıas M ontesinos , N icol ´ as C uello , J ohan O lofsson , J orge C uadra , A melia B ayo , G esa H.-M. B ertrang , and C l ´ ement P errot Instituto de F´ısica y Astronom´ıa, Universidad de Valpara´ıso, Chile Chinese Academy of Sciences South America Center for Astronomy, National Astronomical Observatories, CAS, Beijing 100012, China N´ucleo Milenio de Formaci´on Planetaria (NPF), Chile Univ. Grenoble Alpes, CNRS, IPAG, F-38000 Grenoble, France Instituto de Astrof´ısica, Pontificia Universidad Cat´olica de Chile, Santiago, Chile Departamento de Ciencias, Facultad de Artes Liberales, Universidad Adolfo Ib´a˜nez, Av. Padre Hurtado 750, Vi˜na del Mar, Chile Max Planck Institute for Astronomy, K¨onigstuhl 17, 69117 Heidelberg, Germany LESIA, Observatoire de Paris, Universit´e PSL, CNRS, Sorbonne Universit´e, Univ. Paris Diderot, Sorbonne Paris Cit´e, 5 place Jules Janssen, 92195 Meudon,France (Received February 8, 2021; Revised; Accepted)
Submitted to AJABSTRACTPlanets form in young circumstellar disks called protoplanetary disks. However, it is still di ffi cult to catchplanet formation in-situ. Nevertheless, from recent ALMA / SPHERE data, encouraging evidence of the directand indirect presence of embedded planets has been identified in disks around young stars: co-moving pointsources, gravitational perturbations, rings, cavities, and emission dips or shadows cast on disks.The interpretation of these observations needs a robust physical framework to deduce the complex disk geom-etry. In particular, protoplanetary disk models usually assume the gas pressure scale-height given by the ratioof the sound speed over the azimuthal velocity H / r = c s / v k . By doing so, radiative pressure fields are oftenignored, which could lead to a misinterpretation of the real vertical structure of such disks.We follow the evolution of a gaseous disk with an embedded Jupiter mass planet through hydrodynamicalsimulations, computing the disk scale-height including radiative pressure, which was derived from a general-ization of the stellar atmosphere theory. We focus on the vertical impact of the radiative pressure in the vicinityof circumplanetary disks, where temperatures can reach (cid:38) ff ects create a vertical optically thick column of gas and dust at the proto-planetlocation, casting a shadow in scattered light. This mechanism could explain the peculiar illumination patternsobserved in some disks around young stars such as HD 169142 where a moving shadow has been detected, orthe extremely high aspect-ratio H / r ∼ . Keywords:
Planetary-disk interactions – Hydrodynamical simulations – Protoplanetary disks – Radiative trans-fer simulations INTRODUCTIONDuring the last years, high angular resolution imaging ofdisks around young stars has revealed di ff erent patterns inthe gas and dust structures. These structures are usually ob-served in gas-rich disks with an inner cavity around low and Corresponding author: Mat´ıas [email protected] intermediate-mass stars, classified as transition disks. Spa-tially resolved observations reveal astonishing images of cir-cumstellar disks. For instance, in the disk around HL TauALMA unveiled for the first time a series of gaps at millime-ter wavelengths (ALMA Partnership et al. 2015). Several re-markable other examples from high angular resolution im-ages (at near-infrared and sub-mm wavelengths) show mul-tiple gaps, spirals, and rings. It is thought that embeddedplanets cause most of these features (e.g., Long et al. 2018a; a r X i v : . [ a s t r o - ph . E P ] F e b M ontesinos et al .Keppler et al. 2018). A complete review of disk substruc-tures’ observations at high angular resolution can be foundin the DSHARP project paper series (Andrews et al. 2018;Huang et al. 2018a,b).Besides these astonishing gas / dust structures, intriguing il-lumination features — both in scattered light and thermal(sub)millimeter emission — have been observed. For in-stance, Avenhaus et al. (2014) reported two intensity nullsseen in the circumbinary disk of HD 142527. These werelater identified as shadows cast by a tilted inner circumpri-mary disk (Marino et al. 2015; Casassus et al. 2015a). More-over, multiple shadow features had been discovered throughhigh-contrast polarimetric di ff erential imaging (PDI) withVLT / SPHERE in the outer disk of HD 135344B (Stolkeret al. 2016) and the disk around TW Hya (Debes et al. 2017).Interestingly, Quanz et al. (2013) reported VLT / NACOobservations of the disk around the isolated Herbig Ae / Bestar HD 169142 showing a narrow emission dip located at ∼ ◦ . More recently, Bertrang et al. (2018, 2020) presentedVLT / SPHERE / ZIMPOL observations of the same system athigher spatial resolution. The latter reveals an inner ring at24 au from the star and a narrow emission dip in surfacebrightness located now at ∼ ◦ . Remarkably, the movingshadow could be caused by a rotating optically thick bumppossibly located at ∼
12 au from the star, as suggested byBertrang et al. (2018).Most of these features seem to be related to planet forma-tion processes, in which planet-disk interactions sculpt thedisk surface and mid-plane. However, direct detection ofplanets embedded in disks remains elusive. To date, onlyone planetary-mass companion has been imaged coexistingwith a protoplanetary disk (Keppler et al. 2018; M¨uller et al.2018), corresponding to the discovery of a forming planetwithin the gap in PDS 70, while some planet candidateshave been proposed through kinematic signatures (Pinte et al.2019; Pinte et al. 2020).An important quantity that characterizes the circumstellardisk geometry is its pressure scale-height, which is neededto explain several observations (especially the non-resolvedones). For instance, transitional disks exhibit an excess ofNIR emission related to unusually high aspect-ratios ( H / r ∼ .
2) at the inner rim, where a pu ff ed-up wall is created bydust evaporation (Natta et al. 2001; Dullemond et al. 2001;Dullemond & Monnier 2010; Olofsson et al. 2013). There isalso some inferred di ff erence in the aspect-ratio between thedisks of groups I and II stars (Meeus et al. 2001). Variationson the vertical scale-height in flat or only moderately flareddisks are also supposed to be responsible for self-shadowinge ff ects (Garufi et al. 2014; Stolker et al. 2016)The vertical disk extension is directly related to the fluidpressure field and gravity. In general, circumstellar gaseousdisks are expected to be gas pressure-dominated. This con- dition is used to assume equilibrium between gravity and gaspressure, leading to the well-known relation for the aspectratio H / R = c s / v φ , where c s and v φ are the sound speed andthe azimuthal velocity of the gas, respectively.However, it is worth asking if this condition holds every-where in the disk — especially at high temperatures. In par-ticular, upon which circumstances (if any) does the scale-height formula break down? In such a case, what kind ofequation should be used instead?One important modification should be implemented in thevicinity of a planet, where its gravity dominates and mustbe taken into account (M¨uller et al. 2012). Also, very hightemperatures are expected to develop in some regions of pro-toplanetary disks. One of these hot regions arises at theinner boundaries of transitional disks where the dust subli-mates, reaching temperatures in the range of 1000 − ff ed-up inner rim (Natta et al. 2001).Also, highly luminous (10 − − − L (cid:12) ) forming planets areexpected to form during early stages ( ∼ ∼ ≥ ff ects through radiative transfer calcu-lations. This could explain, for instance, the moving shadowobserved in HD 169142 by associating it to the presence ofan inner planet (Bertrang et al. 2018, 2020).In Section §2, we discuss the physical background of theatmosphere model to compute the disk scale-height. In §3,we run numerical hydro-simulations following the evolutionof a protoplanetary disk with the adopted formulation. In §4,we use our hydro-models to interpret observations by per- adiative scale - height SCALE-HEIGHTThe vertical structure of a disk can be computed underthe assumption of hydro-static equilibrium between pres-sure forces and gravity. i.e., the vertical gradient pressure isbalanced by the vertical tidal acceleration towards the mid-plane. Then, we have:d P d z = − (cid:90) ρ ∇ Ψ dz , (1)where z corresponds to the distance from the mid-plane, P and ρ are the total pressure and density, respectively. Ψ is thegravitational potential given by contributions from the star,and the embedded planet, namely, Ψ = Ψ ∗ + Ψ planet .2.1. Scale-height near the planet
Taking into account the gravitational potential from the starand the planet, the vertical hydro-static equation 1 can bewritten in cylindrical coordinates as:1 ρ ∂ P ∂ z = − GM ∗ ( r + z ) (3 / − GM p ( s + z ) (3 / , (2)where M ∗ , and M p are the mass of the star and planet, respec-tively, and s is the distance from the projected position of theplanet ( z = P g = ρ c s , the scale-height H in the absence of a planet is defined as: H = c s Ω , (3)where c s is the vertically isothermal sound speed, and Ω theangular velocity of the disk.The embedded planet changes the disk structure, leadingto a reduced thickness near it. We follow the procedure byM¨uller et al. (2012) to compute the scale-height H planet , andthe density ρ planet of the disk in the vicinity of the planet. Fora vertically isothermal disk, equation 2 can be integrated: ρ planet = ρ exp (cid:40) − (cid:32) z H + | z | H planet (cid:33)(cid:41) , (4)where the scale-height near the planet is given by, H planet = s H qr , (5)where q = M p / M ∗ , and H defines the scale-height in theabsence of the gravitational influence of the planet (Eq. 3).The planet influence is characterized by the transition dis-tance (M¨uller et al. 2012): s t = (cid:32) r qH (cid:33) / . (6) When s < s t , the planet dominates and therefore H = H planet (Eq. 5), otherwise H is given by equation 3.2.2. Radiative scale-height
In some cases, radiation pressure may be relevant, suchas in circumplanetary disks around accreting planets wheretemperatures could reach high temperatures above 1000 K.We use an analytical model of the scale-height derived byHubeny (1990) to compute the vertical structure, obtainedfrom a generalization of the classical stellar atmospheric the-ory applied to a plane parallel geometry compatible with adisk structure. We combine this procedure with the calcula-tions of the scale-height in the vicinity of the planet (Eq. 5).The total pressure is assumed to be the contribution of thegas, and radiation pressure, i.e., P = P g + P r , respectively.The gas pressure can be written as: P g = ρ c . (7)It is important to note that the isothermal sound-speed c s isrelated to the gas pressure P g only, and not to the total pres-sure P . The radiation pressure P r , can be written as: P r = π c K ν , (8)where c is the speed of light, and K ν is the second momentumradiation flux (proportional to the radiation pressure tensor)and defined as: K ν ≡ π (cid:90) I ν µ d Ω , (9)where I ν is the radiation intensity field in the solid angle Ω and its direction specified by µ ≡ cos( θ ) (being θ the eleva-tion angle respect to the surface normal). For a semi-isotropicradiation field we have µ (cid:39) .Assuming an energy balance in the vertical direction inwhich the energy dissipated per unit volume Q + equals thenet radiation loss per unit volume Q − , from the total pressure P of the fluid, and the definition of the moment equation K ν from the stellar atmospheres theory (Mihalas 1978), the hy-drostatic equilibrium condition in Equation 1 can be rewrittenas a density equation (Hubeny 1990):d ρ d z = − H ρ z + H r H ρ [1 − θ ( z )] , (10)where H g is defined as the pure gas pressure scale-height,given by: H g ≡ H planet (Eq. 5) , if s < s t H g (Eq. 3) , if s ≥ s t , (11) For a detailed discussion on Stellar atmospheres, and the radiative transferequations consult Mihalas (1978); Mihalas & Weibel Mihalas (1984). M ontesinos et al .and H r represents a pure radiative pressure scale-height,given by: H r ≡ ( σ/ c ) T ff κ/ Ω , (12)where T e ff is the disk e ff ective temperature, Ω K the Kep-lerian velocity, κ the flux mean opacity, and σ the Stefan-Boltzmann constant. The e ff ective temperature is relatedwith the mid-plane temperature through: T ff = T τ e ff , (13)where we use for the e ff ective optical depth τ e ff (valid foroptically thick case, Hubeny 1990): τ e ff = √ + τ + τ , (14)and the optical depth is obtained from the mean opacity κ : τ = κρ. (15)To solve Equation 10 one can use a monotonically increas-ing function between 0 and 1 for θ ( z ): θ ( z ) = − ( z / H ) , if z < H , if z ≥ H , (16)where H may be called the density scale-height, given by: H = Σ / ( √ πρ ) , (17)where Σ = (cid:82) ∞ ρ ( z )d z corresponds to the integrated surfacedensity, and ρ = ρ ( z = ff erential density Eq. 10 canbe solved to obtain ρ ( z ): ρ ( z ) = ρ exp (cid:26) − (cid:20) − z H (cid:16) − H r H (cid:17) + | z | H palanet (cid:21)(cid:27) , if z < H ρ exp (cid:40) − (cid:34) (cid:18) z − H r H g (cid:19) + | z | H planet (cid:35)(cid:41) exp (cid:26) − ( H − H r ) H g H r H g (cid:27) , if z ≥ H .(18)The term | z | / H planet in the above equation is only consid-ered in regions where the planet gravity dominates, i.e., when s < s t .We have now all the ingredients to derive the density scale-height H of the disk. We follow the treatment from Hubeny(1990) by introducing the dimensionless parameters: h ≡ H / H g , h r ≡ H r / H g , (19)which allows us to construct an algebraic equation to com-pute H in the form of the dimensionless parameter h : h = √ π (cid:16) hh − h r (cid:17) / [1 − erf( { h ( h − h r ) } / )] + erf( h − h r ) exp [ − ( h − h r ) h r ] , (20) where erf( x ) = (cid:82) ∞ x exp ( − t )d t is the error function.Once a numerical solution for h is found, the final scale-height is simply given by: H ≡ H g + r = hH g , (21)where H g is defined in Equations 11. We call H g + r the densityscale-height obtained from Equation 21 and 20, which cor-responds to the scale-height, containing contributions fromboth gas and radiation pressure. It is worth to mention that H r in Equation 12 represents the semi-thickness of the diskonly when radiative forces are dominant. When the radiationpressure is negligible, i.e. H r (cid:28) H g , the scale-height is givenby: H (cid:39) H g , (22)recovering the standard calculation defined in Eq. 11In this limit, the vertical density profile in Equation 18 hasthe standard form: ρ ( z ) (cid:39) ρ exp {− ( z / H g ) } , if s ≥ s t ρ exp {− [ ( z / H g ) + | z | / H planet ] } , if s < s t , (23)where ρ = Σ / ( √ π H g ).Summarizing, H g + r (Eq. 21) represents the circumstellarscale-height when gas and radiation pressures are taken intoaccount, and H g (Eq. 11) when only gas pressure is included. NUMERICAL HYDRO-SIMULATIONSWe ran a set of hydro-simulations to follow the evolutionof a gaseous disk with an embedded Jupiter mass planet, inwhich the planet radiates away a fraction of its internal en-ergy at a constant rate. We modify the the public FARGO-ADSG code (Baruteau & Masset 2008) following the sameprocedure as in Montesinos et al. (2015, 2016). Assum-ing hydrostatic equilibrium, we computed the vertical scale-height following the expression for H g + r (Eq. 21) to recreatea 3D structure, with an initial aspect-ratio H / r set to 0.05 foreach model. The gravitational force of the planet uses a soft-ening length (cid:15) = . α -viscosity prescription (Shakura & Sunyaev 1973), with α = − . For the radiative cooling mechanism we assumeblack-body emission Q − = σ T ff . We also take into accounta radiative heating source (feedback) from the planet as inMontesinos et al. (2015).The disk extends from 4 to 25 au, with an initial surfacedensity given by: Σ ( r ) = Σ r p r , (24) adiative scale - height Figure 1.
Temperature profile of the disk after 10 years. The left panel shows a model without planet feedback L p = ∼
200 K), while the right panel includes a L p = × − L (cid:12) luminous planet (temperature peak ∼ Figure 2.
Slice of the temperature field along the radial direction— passing through the planet located at r p =
10 au — for modelswith and without feedback in red and blue (respectively). where Σ =
30 g cm − corresponds to the density at the planetlocation r p . The total disk mass is therefore M disk ≈ × − M (cid:12) . The embedded luminous planet ( L p = − L (cid:12) ) islocated at 10 au.The grid resolution for all the simulations was set to n r = n θ = π . We initially put aplanet located at r p =
10 au, φ p = ◦ on a circular orbit( φ = ◦ is the north). The planet is not allowed to migrate.We follow the evolution of the system in the reference frameof the planet. We present here two cases: one in which theplanet feedback is turned-o ff (i.e. L p = L p = × − L (cid:12) , which we have shown is achievable even whenconsidering accretion feedback (G´arate et al. 2020).Figure 1 shows the mid-plane temperature field after 10 years of disk evolution (or ∼
316 planetary orbits) at whichthe system has reached a quasi-steady state. The left panelshows a model without feedback ( L p = r passing through the planet atfixed φ p . As in Figure 1, we observe that when the feedbackis turned-on, the temperature peaks at 1060 K, while it peaksat 200 K when the planet is not emitting.In Fig. 3 we show how the disk aspect-ratio H / r is a ff ectedby temperature changes due to the planet feedback. We plot H / r for a slice passing over the planet, using the pure gasscale-height H g model (Eq. 11) and the radiative scale-height H g + r (Eq. 21). The left panel of Figure 3 corresponds to amodel with no feedback. In this case, both calculations ( H g and H g + r ) matched, indicating that radiative pressure is neg-ligible when the feedback is turned o ff . In the right panelof Fig. 3, we plot the model with a luminous planet. Inthis case, the radiative pressure locally overcomes the planetgravity, enhancing the scale-height above the planet. TheCPD region increased its height — with respect to H g — by( H r − H g ) / H g ∼ (0 . − . / . =
19. Radiative forces aredominant above the CPD when the planet feedback is acti-vated. OBSERVATIONAL EFFECTS OF THE RADIATIVESCALE-HEIGHT4.1.
Vertical optical depth
To create a 3D distribution, we vertically extend the diskusing Equation 21, where the volume density ρ ( r , φ, z ) iscomputed from Equation 18. A scale-height enhancementas the one observed in Figure 3 (right panel) could act as abump blocking or scattering a fraction of the stellar radia-tion. Therefore, it is expected to cast a shadow into the outerdisk in the radial direction away from the planet. However, M ontesinos et al . Figure 3.
Radial slices of the aspect ratio h = H / r — passing through the planet location — computed with the pure gas formula h g = H g / r (blue) from equation 11, and the radiative one h g + r = H g + r / r (red) from Equation 21. The Left panel corresponds to a model without feedback.In this case both curves ( h g and h g + r ) matches. The middle panel corresponds to a model with planet feedback. The radiative component in theexpression h g + r overcomes gravitational forces from the planet, locally enhancing the aspect ratio. The right panel is a zoom of the vicinity ofthe planet (around ∼
10 au) for the feedback model. The observed fluctuation around ∼ .
95 au responds to a perturbation in the temperaturefield due to the planetary feedback.
Figure 4.
Optical depth τ ( r , z ) of the disk as a function of the ver-tical and radial direction for a model with L p = − L (cid:12) . The toppanel was computed using the pure gas formulation H g (Eq. 11),the bottom panel uses H g + r (Eq. 21). At the planet location ( r p = bump is enhanced by a factor of 19. this phenomenon only occurs if, and only if, the “bump” isoptically thick in the vertical and radial directions. We compute the optical depth in Equation 15 from the vol-ume density given by Eq. 18, using for κ the Rosseland meanopacity appropriate for protoplanetary disks (Semenov et al.2003). Figure 4 shows a slice cut of the vertical optical depth τ ( z , r ) = (cid:82) ∞ z κ ( r ) ρ ( r , z ) dz , passing through the planet for themodel with planetary feedback. We note that, when usingthe radiative formulation H g + r (Eq. 21), an optically thick“bump” of about ∼ r p =
10 au), compared with a no-bump situation in this case “radial” refers to a spherical coordinate, rather than the cylin-drical r coordinate used in the simulations. Figure 5.
Grain size (10 − −
10 cm) distribution in r , z for a modelwith L p = − L (cid:12) . The slice passes through the planet located at r p =
10 au. The top figure shows a model computed with H g , thebottom one includes radiative pressure H g + r . The upper layers ofthe disk are populated with micron-size particles as expected, whilebiggest particles ( ∼ cm) settle to the mid-plane. if the pure gas formulation H g (Eq. 11) is used, meaning thatplanet gravity dominates over gas pressure alone.As mentioned before, if the “bump” is also optically thickin the radial direction, then a shadow is expected to be cast.We explore this possibility in the next section through radia-tive transfer calculations.4.2. Monte Carlo radiative transfer: Dust vertical structure
We perform radiative transfer calculations using the radia-tive transfer code radmc d (Dullemond et al. 2012). Wefeed radmc d with the gas distribution obtained from our adiative scale - height Figure 6.
Images for λ = µ m obtained from radmc d . Panel a): model without feedback ( L p = H g . Panel b):model with L p =
0, but using the new formula H g + r . Panel c): model with feedback L p = × − L (cid:12) , using the pure gas formula H g . Panel d):same feedback as c), but using H g + r . A shadow cast from the planet location is observed in this model. hydro-simulations (see section §3), corresponding to an evo-lutionary step after ≈ yrs. The dust of the disk is as-sumed to be composed by astrosilicates of intrinsic den-sity ρ intr = − following a power-law size distributiond n ( a ) ∝ a − . d a , with a ranging from 0.1 µ m to 1 cm. Thedust density is normalized by imposing a gas-to-dust ratioof 100. We assume local thermodynamic equilibrium, i.e. T gas = T dust . The absorption e ffi ciencies were computed us-ing Mie theory (Bohren & Hu ff man 1983).We decompose the dust size range into 12 logarithmicallyspaced bins representing 12 dust i -species. We obtain a dustdensity distribution Σ i for each i -bin in such a way that thesum of individual species surface density gives the total dustdensity i.e., Σ dust ( r , φ ) = (cid:88) i Σ i ( r , φ ), where the sum is per- M ontesinos et al .formed from a min to a max . The Stokes number is computedfrom St = ρ intr a Ω K /ρ c s , where ρ intr is the intrinsic density ofparticles, a the particle radius, Ω K the Keplerian velocity, ρ the volume density, and c s the isothermal sound speed.The dust scale-height is assumed to be H d = H g √ α/ ( α + St i ) , where α is the turbulent viscosity, St i theaverage Stokes number of the i th -species, and H g the densityscale-height computed from H g (Eq. 11) or H g + r (Eq. 21).This method allows us to take into account di ff erent verticaldistributions for di ff erent dust species i (Youdin & Lithwick2007).In Figure 5, we plot the dust scale-height H d ( r , z ) for thecase with planetary feedback. The top panel corresponds tocomputations with H g (Eq. 11), while the bottom one usedthe generalized H g + r expression (Eq. 21). When using H g alone, the bump over the CPD disappears.For the 3D radiative transfer calculation, we assume thatthe central star is of solar-type with an e ff ective temperatureof 6000K. We incline the disk by 13 ◦ , and we set the systemat 140 pc from the Earth.In Figure 6 we show the relative intensity ( I / I max ) for λ = µ m. We compare the e ff ect of the feedback (switchedon / o ff ) when using the di ff erent scale-height formulations( H g and H g + r ) to compute the 3D distribution as an input in radmc d . We remark that:1. A shadow is cast in scattered light from the CPD pro-jected to the outer disk only when the feedback is acti-vated, and the formulation H g + r is used.2. No shadows are observed when L p =
0, whether H g or H g + r is used to compute the synthetic image.The features described above can be explained as fol-lows. At the CPD location, the gravity from the planet dom-inates over gas pressure forces. When the feedback is acti-vated (producing a local temperature increment from 200 to1060 K), radiative forces play a relevant role. They can over-come the gravity from the planet, and the material is pu ff ed-up in the vertical direction (Fig. 3, left panel). This materialis optically thick, producing a narrow shadow cast from theCPD to the outer disk regions (Fig. 6).We repeat the same calculation for λ = µ m. In thiscase, the shadows disappear since large dust tends to settle,making the bump optically thin at those frequencies. Theshadow is observed only in reflected light and not in the dusttemperature. This suggests that the “optimal” wavelength forshadow detection should be dominated by scattered light ( ∼ − µ m) rather than thermal emission ( > µ m). However, tovalidate these results, proper dust modeling is required with from this we have that small grains with small Stokes number (St (cid:28) α )have scale-heights ∼ gas scale-height. large grain sizes ( ≥ µ m), where the dust may not be wellcoupled to the gas as assumed in this work.It is worth noting in Figure 6 that in models with L p = − L (cid:12) (bottom panels) some azimuthal features appear,such as wakes within the gap and the disk. Their origincomes from the planetary feedback, which heats-up the gap,producing noticeable turbulence (see Montesinos et al. 2015;G´arate et al. 2020 for a discussion about the feedback).In Figure 7, we plot the same disk emissivity reported for λ = µ m (Figure 6), but following a slice arc between 260 ◦ to 290 ◦ for di ff erent radii ranging from 16 to 20 au. The toppanels of Figure 7 correspond to models without feedback,where the emissivity was calculated using H g (top left), and H g + r (top right). The bottom panels include the L p = − L (cid:12) feedback model, with H g (left) and H g + r (right). Looking atthe feedback model with H g + r (bottom, right panel) an emis-sion dip is observed up to 20 au. DISCUSSIONWe model the evolution of a protoplanetary disk with anembedded Jupiter mass planet, in which the planet presentsan intrinsic luminosity of 10 − L (cid:12) . When this feedback isactivated, the region around the planet reaches a mid-planetemperature of about ∼ ∼ Myr old planet) and can last for about ∼ L p = adiative scale - height Figure 7.
Slice cuts between 260 ◦ to 290 ◦ degrees for di ff erent radial cuts — ranging from 16.1 to 20.5 au — of the emissivity reported for λ = µ m in Fig. 9. The top panels correspond to models without feedback, where the emissivity was calculated using H g (top left) and H g + r (top right), respectively. The bottom panels include the 10 − L (cid:12) feedback models, with H g (left) and H g + r (right). The feedback model report astrong brightness dip when H g + r is used. proaches, H g (Eq. 11) and H g + r (Eq. 21), produce the samescale-height (see Figure 3), and no shadows are observed.In protoplanetary disks, hot regions can be reached at leastin two di ff erent locations: at the inner rim of a transitionaldisk and a planet-forming region. gThe inner rim of a tran-sitional disk is expected to be at temperatures of about ∼ H g + r (Eq. 21)to explain the extremely high aspect-ratios H / r ∼ . ∼ ff ects the vertical disk geometry.More specifically, by surpassing the planet’s gravity, radia-tion pressure locally enhances its height by 19% comparedto a pure gas pressure model (left vs. right panels in Fig. 3). From radiative transfer calculations, we showed that the“bump” created above the CPD is optically thick, casting ashadow in scattered light, which extends to the outer regionsof the disk in about ∼
20 au from the planet. The shadowsare only observed at wavelengths of ∼ − µ m, mostly be-cause small particles (micron-size) get lifted to the upper disklayers making the bump optically thick at those wavelengths(Fig. 4). At longer wavelengths, the shadows disappear,as thermal emission tends to dominate over scattered light.For close-in planets, even though they are relatively brightat 1 micron (Fig. 6), their direct detection can be harmed byspeckles and other observational artifacts. On the other hand,since shadows extend much further away in the disk, they aremore favorable for detection.To conciliate our findings with observations, recent zim - pol/sphere images obtained by Bertrang et al. (2018, 2020)revealed small-scale structures could be due to planet-diskinteractions. The observations also show a moving narrowsurface brightness dip with an azimuthal width of ∼ ◦ .The moving shadow seems to be cast by a large amount ofoptically thick material that blocks a fraction of the stel-lar emission, which could be interpreted as an undetected0 M ontesinos et al .CPD (Bertrang et al. 2018). Such dips are precisely the kindof predictions obtained with our prescription for the scale-height.A compelling case is the first direct image of a forming-planet within the gap of PDS 70 (Keppler et al. 2018; M¨ulleret al. 2018), where an accreting 10 M J mass planet sur-rounded by a protoplanetary disk was suggested (Keppleret al. 2019). Also, Christiaens et al. (2019) recently presentedobservational evidence of the presence of a CPD around theprotoplanet PDS 70 b. Under these circumstances, a shadowshould be cast by the putative CPD. In fact, some shadow-ing has been observed in the outer disk of PDS 70 (Longet al. 2018b). Its origin remains however unclear. High-contrast observations of micron-sized dust from PDS 70bcould, in principle, reveal new shadows arising from theplanet-forming region. This system is an exceptional labo-ratory to test our ideas, where a dedicated model would beneeded to predict possible shadowing e ff ects.Other interesting illumination patterns — non-related toaccreting planets — have also been reported. For instance, inHD 142527 two diametrically opposed shadows have beendiscovered in polarized scattered light by Avenhaus et al.(2014). These intensity nulls were later explained by thepresence of a misaligned inner disk able to block a frac-tion of the stellar emission (Marino et al. 2015; Casassuset al. 2015b; Price et al. 2018). The main di ff erence with theshadows obtained in our simulations is that the latter causeboth a deeper intensity null and a broader azimuthal exten-sion (e.g. ∼ ◦ for HD 142527 as opposed to ∼ ◦ inour models). Moreover, it has been shown that these strongshadows can trigger spirals in the gas and dust distribution(Montesinos et al. 2015; Cuello et al. 2019), respectively.Also, if the inclined inner disk is precessing, then the pro-jected shadow rotates as well (Facchini et al. 2018). Remark-ably, this could create planetary-like spirals in the gaseousdisk at the co-rotating region with the shadow (Montesinos& Cuello 2018). Further observations of intriguing mov-ing shadows have been reported by Stolker et al. (2017) andDebes et al. (2017). However, there are still doubts abouttheir origin. They could be caused by a precessing innerwarped disks (Nealon et al. 2019), or by moving opticallythick clouds / bumps at the inner disk rims. In this work, theshadows from the CPD move at Keplerian rotation, and wedo not expect any impact on the gas dynamics due to suchillumination e ff ects. We emphasize that the obtained dipsare not deep enough to cool down the gas and trigger az- imuthal perturbations in the density field, as in Montesinoset al. 2016. SUMMARYIn this work we showed that for any region in a circum-stellar disk with temperatures larger than ∼ ff ect should be taken into account tocompute the disk scale-height, which can be approximated byan analytical expression (Eq. 21) derived from stellar atmo-spheres theory adapted to a disk geometry by Hubeny (1990).In the vicinity of an accreting proto-planet, radiative forcessurpass the planet’s gravitational force, producing an opti-cally thick bump that may cast a narrow shadow into the outerdisk. Such shadow should move at Keplerian speed with theplanet. For a non-accreting planet, radiative pressure can beneglected, and no shadows are observed. Due to the narrowwidth of the predicted shadows ( ∼ ◦ ), observations at highangular resolution are required to reach high signal-to-noiselevels. These dips should be better observed in scattered lightat wavelengths close to ∼ − µ m. The shadows are expectedto disappear at longer wavelengths since the bump above theCPD becomes optically thin for lower frequencies, and thedisk thermal emission contaminates any residual shadowingpattern. ACKNOWLEDGMENTSThe authors thank the referee for constructive commentsand recommendations, which help to improve this pa-per. MM acknowledges financial support from the Chi-nese Academy of Sciences (CAS) through a CAS-CONICYTPostdoctoral Fellowship administered by the CAS SouthAmerica Center for Astronomy (CASSACA) in Santiago,Chile. NC acknowledges financial support provided byFONDECYT grant 3170680. GHMB acknowledges fund-ing from the European Research Council (ERC) under theEuropean Union Horizon 2020 research and innovation pro-gramme (grant agreement No. 757957). JO acknowledgesfinancial support from Fondecyt (grant 1180395). AB ac-knowledges support from FONDECYT grant 1190748. CPacknowledges financial support from FONDECYT (grant3190691). MM, NC, JC, JO, AB and CP acknowledge sup-port from Iniciativa Cient´ıfica Milenio via the N´ucleo Mile-nio de Formaci´on Planetaria. This project has received fund-ing from the European Union’s Horizon 2020 research andinnovation programme under the Marie Skłodowska-Curiegrant agreement No 210021.REFERENCES ALMA Partnership, Brogan, C. L., P´erez, L. M., et al. 2015, ApJL,808, L3, doi: 10.1088 / / / / L3 Andrews, S. M., Huang, J., P´erez, L. M., et al. 2018, ApJ, 869,L41, doi: 10.3847 / / aaf741 adiative scale - height Armitage, P. J. 2015, arXiv e-prints, arXiv:1509.06382.https: // arxiv.org / abs / / / / / / / mnras / stx3052Bertrang, G. H. M., Flock, M., Keppler, M., et al. 2020, arXive-prints, arXiv:2007.11565. https: // arxiv.org / abs / ff man, D. R. 1983, Absorption and scatteringof light by small particlesCasassus, S., Marino, S., P´erez, S., et al. 2015a, ApJ, 811, 92,doi: 10.1088 / / / / / / / / // arxiv.org / abs / / / / / / / / // ascl.net / / annurev-astro-081309-130932Facchini, S., Juh´asz, A., & Lodato, G. 2018, MNRAS, 473, 4459,doi: 10.1093 / mnras / stx2523G´arate, M., Cuadra, J., Montesinos, M., & Ar´evalo, P. 2020,MNRAS, doi: 10.1093 / mnras / staa3860Garufi, A., Quanz, S. P., Schmid, H. M., et al. 2014, A&A, 568,A40, doi: 10.1051 / / / / aaf740Huang, J., Andrews, S. M., P´erez, L. M., et al. 2018b, ApJL, 869,L43, doi: 10.3847 / / aaf7a0Hubeny, I. 1990, ApJ, 351, 632, doi: 10.1086 / / / // arxiv.org / abs / / / aae8e1Long, Z. C., Akiyama, E., Sitko, M., et al. 2018b, ApJ, 858, 112,doi: 10.3847 / / aaba7c Marino, S., Perez, S., & Casassus, S. 2015, ApJL, 798, L44,doi: 10.1088 / / / / L44Masset, F. S. 2002, A&A, 387, 605,doi: 10.1051 / / / / / / / mnrasl / sly001Montesinos, M., & de Freitas Pacheco, J. A. 2011, A&A, 526,A146, doi: 10.1051 / / / / / / L8Mordasini, C., Alibert, Y., Klahr, H., & Henning, T. 2012, A&A,547, A111, doi: 10.1051 / / / / / / / / mnras / stz346Olofsson, J., Benisty, M., Le Bouquin, J.-B., et al. 2013, A&A,552, A4, doi: 10.1051 / / / s41550-019-0852-6Pinte, C., Price, D. J., M´enard, F., et al. 2020, ApJL, 890, L9,doi: 10.3847 / / ab6ddaPrice, D. J., Cuello, N., Pinte, C., et al. 2018, MNRAS, 477, 1270,doi: 10.1093 / mnras / sty647Quanz, S. P., Avenhaus, H., Buenzli, E., et al. 2013, ApJL, 766, L2,doi: 10.1088 / / / / L2Semenov, D., Henning, T., Helling, C., Ilgner, M., & Sedlmayr, E.2003, A&A, 410, 611, doi: 10.1051 / / / / / / / ff , B., et al. 2017, ApJ, 849, 143,doi: 10.3847 / / aa886aVinkovi´c, D. 2014, A&A, 566, A117,doi: 10.1051 / / / j.icarus.2007.07.012 ontesinos et al . Zhu, Z., Stone, J. M., Dong, R., Rafikov, R., & Bai, X.-N. 2015, inAAS / Division for Extreme Solar Systems Abstracts, Vol. 3,AAS //