JWST/MIRI Detection of a Carbon-rich Chemistry in the Disk of a Solar Nebula Analog

JWST

/

MIRI Detection of a Carbon-rich Chemistry in the Disk of a Solar Nebula Analog

María José Colmenares

, Edwin A. Bergin

, Colette Salyk

, Klaus M. Pontoppidan

, Nicole Arulanantham

Jenny Calahan

, Andrea Banzatti

, Sean Andrews

, Geoffrey A. Blake

, Fred Ciesla

, Joel Green

, Feng Long

(

龙

凤

)

, Michiel Lambrechts

, Joan Najita

, Ilaria Pascucci

, Paola Pinilla

, Sebastiaan Krijt

Leon Trapman

and

the JDISCS Collaboration

Department of Astronomy, University of Michigan, 1085 South University Avenue, Ann Arbor, MI 48109, USA;

mjcolmen@umich.edu

Vassar College, 124 Raymond Avenue, Poughkeepsie, NY 12604, USA

Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, Pasadena, CA 91109, USA

Division of Geological and Planetary Sciences, California Institute of Technology, MC 150-21, 1200 E California Boulevard, Pasadena, CA 91125, USA

Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA

Center for Astrophysics, Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA

Department of Physics, Texas State University, 749 North Comanche Street, San Marcos, TX 78666, USA

Dept. of the Geophysical Sciences, The University of Chicago, Chicago, IL 60637, USA

Lunar and Planetary Laboratory, University of Arizona, Tucson, AZ 85721, USA

Center for Star and Planet Formation and Natural History Museum of Denmark, Globe Institute, University of Copenhagen, Oster Voldgade 5-7, 1350 Cope

nhagen,

Denmark

NSFs NOIRLab, 950 N. Cherry Avenue, Tucson, AZ 85719, USA

Mullard Space Science Laboratory, University College London, Holmbury St Mary, Dorking, Surrey RH5 6NT, UK

School of Physics and Astronomy, University of Exeter, Stocker Road, Exeter, EX4 4QL, UK

Department of Astronomy, University of Wisconsin-Madison, 475 N Charter Street, Madison, WI 53706, USA

Received 2024 July 2; revised 2024 October 22; accepted 2024 October 23; published 2024 December 11

Abstract

It has been proposed, and con

fi

rmed by multiple observations, that disks around low-mass stars display a molecule-

rich emission and carbon-rich disk chemistry as compared to their hotter, more massive solar counterparts. In this

work, we present JWST Disk Infrared Spectral Chemistry Survey MIRI-MRS observations of the solar-mass star

DoAr 33, a low-accretion rate T Tauri star showing an exceptional carbon-rich inner disk. We report detections of

O, OH, and CO

, as well as the more complex hydrocarbons, C

and C

. Through the use of

thermochemical models, we explore different spatial distributions of carbon and oxygen across the inner disk and

compare the column densities and temperatures obtained from LTE slab model retrievals. We

fi

nd the best match to

the observed column densities with models that have carbon enrichment, and the retrieved emitting temperature

and area of C

with models that have C

–

4 inside the 500 K carbon-rich dust sublimation line. This

suggests that the origin of the carbon-rich chemistry is likely due to the sublimation of carbon-rich grains near the

soot line. This would be consistent with the presence of dust processing as indicated by the detection of crystalline

silicates. We propose that this long-lived hydrocarbon-rich chemistry observed around a solar-mass star is a

consequence of the unusually low M-star-like accretion rate of the central star, which lengthens the radial mixing

timescale of the inner disk, allowing the chemistry powered by carbon grain destruction to linger.

Uni

fi

ed Astronomy Thesaurus concepts:

Protoplanetary disks

(

1300

)

;

Astrochemistry

(

)

;

Planet formation

(

1241

)

1. Introduction

The predominant approach for understanding the connec-

tions between disk composition and the atmospheres of giant

planets relies on the elemental carbon-to-oxygen

(

)

ratio

(

K. I. Öberg et al.

2011b

; N. Madhusudhan

2012

)

.Ata

baseline level, the temperature-dependent condensation

sequence of volatile carriers of oxygen and carbon, such as CO,

, and H

O, can lead to spatially varying C

O ratios in both

the gas and solid components

(

K. I. Öberg et al.

2011b

)

Beyond the CO snowline, the expectation is that all of the C

and O carriers reside mostly in the solids, frozen as icy mantles

atop dust grain cores. As dust grains move closer to the star and

cross the corresponding local snowlines, the icy mantles will

sublimate and alter the abundance of C and O in the gas. This

implies that the material available for planet formation will

change across the disk, potentially yielding a variety of planet

compositions.

However, a variety of mechanisms can in

fl

uence the gaseous

O ratio. For instance, extensive inward pebble drift can

supply abundant water ice to the inner disk. Once this water ice

sublimates it would substantially lower the gaseous inner disk

O ratio

(

R. A. Booth & J. D. Ilee

2019

; A. Banzatti et al.

2020

2023

)

. Extensive planetesimal formation and water

sequestration beyond the snowline should have the opposite

effect

(

raising inner disk gas-phase C

O; J. R. Najita et al.

2013

)

. The trapping of ices can also strongly in

fl

uence the C

ratio

(

N. F. W. Ligterink et al.

2024

)

. These concepts focus on

the most volatile and abundant components that are found as

molecular ices coating the more refractory solid grain cores,

which are themselves composed of solid-state organics and

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

https:

doi.org

10.3847

1538-4357

ad8b4f

(

)

. Published by the American Astronomical Society.

NASA Hubble Fellowship Program Sagan Fellow.

Corresponding author.

Original content from this work may be used under the terms

of the

Creative Commons Attribution 4.0 licence

. Any further

distribution of this work must maintain attribution to the author

(

)

and the title

of the work, journal citation and DOI.

silicate minerals

(

B. T. Draine

2011

)

. Closer to the star, within

1 au, there will exist additional dust sublimation fronts that

return this less volatile material to the gas. Of particular

importance for the composition of planets forming within this

region is the soot line. This is the theorized sublimation front of

carbonaceous organic material

(

M. E. Kress et al.

2010

)

that is

believed to activate near

∼

500 K

(

J. Li et al.

2021

)

. This

sublimation will

fl

ood the gas with carbon, elevating the local

gaseous C

O ratio while at the same time lowering C

O in the

solids. Interior to the soot line, this would lead to the formation

of rocky planets with reduced carbon content, such as the Earth

(

J. Li et al.

2021

)

. Beyond the soot line, it would facilitate the

formation of worlds with signi

fi

cant carbon inventories in their

mantles, which, through outgassing, would impact the nascent

atmospheric composition

(

E. A. Bergin et al.

2023

)

Within this context, spectroscopic surveys of protoplanetary

disks by the Spitzer Space Telescope

(

M. W. Werner et al.

2004

)

isolated a dichotomy where disks surrounding very low-

mass stars

(

0.2

e

)

have elevated C

HCN and

HCN

fl

ux ratios

(

I. Pascucci et al.

2009

2013

)

comparison to solar-mass T Tauri disks. This is suggested to be

the result of an increase in the gaseous C

O ratio in the

innermost regions close to the star

(

J. R. Najita et al.

2013

;

I. Pascucci et al.

2013

)

. Additional support for this dichotomy

has also been seen in recent spectra of disks surrounding very

low-mass stars obtained with the James Webb Space Telescope

(

JWST

)

. The disk around 2MASS-J160532-1933159

(

0.14

e

)

exhibits broad optically thick spectral features of

, as well as detections of benzene

(

)

, and

(

B. Tabone et al.

2023

)

, also found in ISO-ChaI 147

(

A. M. Arabhavi et al.

2024

)

and Sz 28

(

J. Kanwar et al.

2024a

)

, along with emission from HC

N, C

, and C

Other than C

, these species were not isolated in Spitzer

spectra, and the chemistry surrounding these very low-mass

stars appears to be hydrocarbon dominated.

The origin of the hydrocarbon chemistry is uncertain. One

possibility is the destruction of carbon-rich organics near the

soot line. However, this points toward a conundrum: if this

chemistry is a result of carbon-rich grain destruction near

500 K, then these emission features should be prevalent

surrounding more luminous early-type stars, but they are not.

J. Mah et al.

(

2023

)

provide an alternative explanation. They

suggest that the supply of carbon-rich gas from the outer disk

may provide the fuel to generate hydrocarbons in the inner disk

gas. This is more prevalent in lower-mass stars, which have ice

lines that lie close to their stars where the viscous timescales

are shorter, in addition to the radial drift of particles being more

fi

cient

(

P. Pinilla et al.

2013

)

. This scenario was proposed to

explain the spectrum of Sz 114, which showed strong water

emission lines despite having the same stellar mass as J160532-

1933159, likely due to the difference in their pebble disks

(

C. Xie et al.

2023

)

. Alternatively, C. Walsh et al.

(

2015

)

propose that a combination of weak far-ultraviolet

(

FUV

)

and

X-ray-driven chemistry in M dwarf disks is the cause of

abundant C

and HCN observed in the atmosphere of

the disk.

With the idea that carbon-rich chemistry so far has been

expected to be dominant in low-mass stars, in this work, we

present JWST MIRI observations of the hydrocarbon-rich

source, DoAr 33, a 1.1

e

star. We detect the emission of H

, OH, and C

, including its isotopologue

, and tentatively CH

and HC

N, for all of which we

retrieve column densities, temperatures, and emitting areas, as

described in Section

. With the objective of disentangling the

origin of a rich carbon and oxygen chemistry in the inner disk

of solar analogs, in Section

we introduce thermochemical

models with varying abundances of the main carbon and

oxygen carriers. We calculate abundances, column densities,

and

fl

ux-weighted temperatures for comparison with the

detected molecules in Section

. We then discuss

(

Section

)

how DoAr 33 may shed important light on the origin of the

hydrocarbon dichotomy found within disks surrounding solar-

mass and very low-mass stars.

2. Observations and Retrievals

2.1. Observations and Data Reduction

We present JWST

(

J. Rigby et al.

2023

)

, MIRI

(

G. S. Wright

et al.

2023

)

, MRS

(

M. Wells et al.

2015

; I. Argyriou et al.

2023

)

observations of DoAr 33, observed as a part of the JWST

Cycle 1 GO1584 program

(

PIs C. Salyk and K. Pontoppidan

)

and the JWST Disk Infrared Spectral Chemistry Survey

(

JDISCS; K. M. Pontoppidan et al.

2024

)

. The properties of

DoAr 33 are listed in Table

. This source was acquired by

JWST using the MRS target acquisition procedure along with

the neutral density

fi

lter. We used the four-point dither pattern

in the negative direction and observed in all three MIRI-MRS

grating settings such that we had full spectral coverage

(

4.9

–

)

. The exposure time per module was set to

1687.224 s in order to achieve a signal-to-noise ratio

∼

100 at

the longest wavelengths.

The spectrum was reduced with the standard JDISCS

reduction pipeline, explained in detail in K. M. Pontoppidan

et al.

(

2024

)

. This custom pipeline produces a higher spectral

contrast than the standard JWST pipeline by using observations

of asteroids as spectral response function calibrators, which

produces a higher signal-to-noise ratio for Channels 3 and 4

when compared to a traditional calibration using stars. In the

case of DoAr 33, the 1.15.0 version of the JWST Calibration

Pipeline

(

H. Bushouse et al.

2024

)

and the jwst_1253.pmap

Calibration Reference Data System context was used. The

reduced spectrum is shown in Figure

. The process to compute

the continuum involves an iterative procedure as described in

K. M. Pontoppidan et al.

(

2024

)

. This method begins by

applying a median

fi

lter to the spectrum to smooth it out and

then comparing it with past iterations. During each iteration,

wavelength channels exhibiting a lower

fl

ux density than the

previous iteration are kept, while rejected ones are replaced

through linear interpolation; a second-order Savitzky

–

Golay

fi

lter is utilized to further smooth the result. The regions

13.4

–

14.1

m and 14.95

–

m are excluded from the

continuum calculation to avoid overestimation because broad

molecular features from the

branch of organic molecules

produce a pseudocontinuum.

2.2. LTE Slab Model Retrievals

To benchmark the thermochemical models of DoAr 33, we

fi

rst retrieve best-

fi

t column densities, temperatures, and

emitting areas for each detected molecule using the local

thermal equilibrium

(

LTE

)

slab models in

spectools

ir

(

C. Salyk

2022

)

. These models assume an isothermal slab of

gas that can be characterized by only its temperature, line-of-

sight column density, and emitting area. The emitting area can

be interpreted as a circle of radius

slab

, in which case an

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

emitting radius can be obtained by assuming

slab

However, the emission may instead arise from one or more

rings with total area

. The code uses the line parameters from

the HITRAN database

(

I. E. Gordon et al.

2022

)

to model the

emission spectra of different molecules, allowing us to directly

compare to the observations. The code takes into account

mutual line overlap, which has been shown to cause saturation

in regions of dense line clustering for the organics

(

B. Tabone

et al.

2023

)

and in ortho

–

para line pairs for water

(

A. Banzatti

et al.

2024

)

. Despite their simplicity, LTE models have been

able to successfully reproduce infrared molecular spectra

observed with Spitzer

(

e.g., J. S. Carr & J. R. Najita

2011

;

C. Salyk et al.

2011

)

and more recently JWST

(

e.g., D. Gasman

et al.

2023

; C. Xie et al.

2023

; C. E. Muñoz-Romero et al.

2024

; M. Temmink et al.

2024a

)

We follow a procedure similar to previous works

(

C. Salyk

et al.

2011

; S. L. Grant et al.

2023

)

to obtain the parameters

from the slab models. We run a grid of slab models with

col

–

−

in steps of

log

0.166

col

and

slab

100

–

1500 K in steps of 25 K. These models are

generated within the spectral window of 12

–

16.5

m and then

compared with the MRS spectrum through a

minimization.

The only exception is the OH models, which are run until

slab

3000 K and within a window of 12

–

m to maximize

the amount of features included in the

fi

t. We convolve the

models at each wavelength to the corresponding MIRI-MRS

resolving power by adopting the spectral resolution measured

by K. M. Pontoppidan et al.

(

2024

)

, which, in the spectral range

considered, is in agreement with previous measurements

(

A. Labiano et al.

2021

; O. C. Jones et al.

2023

)

. Then, we

select a spectral window for each molecule to calculate the

reduced

and

fi

nd the best-

fi

t parameters. The details of the

calculation are explained in Appendix

, and the spectral

ranges used to calculate it are shown in Figures

and

. For

each column density and temperature pair in the grid, we adjust

the emitting area to minimize

for that pair. Then, we

fi

nd the

lowest of those

values among the pairs in the grid. The

resulting

maps for each molecule can be found in

Appendix

Molecules like C

, HCN, H

O and OH have a signi

fi

cant

amount of overlap in their emission; therefore, we implement

an additional step to further isolate their emission. First, we

fi

the best-

fi

t model for H

O. Then, prior to simulating CO

,we

subtract the best-

fi

O model to eliminate interfering features

that could in

fl

uence the

fi

tting. Similarly, we subtract the best

model from the residual spectrum before

fi

tting the rest of

the molecules. This ensures that we are only taking into

account the most signi

fi

cant feature of each species. In the case

of C

, we simultaneously

fi

t its spectrum with

and

HCN due to the overlap in wavelength of their emission

features. We test all possible combinations of column densities

and temperatures within the speci

fi

ed grid, allowing the values

to vary independently. We report the parameters of the best-

fi

models for each molecule in Table

, and we show the

fi

t of the

–

16.5

m region in Figure

. The individual models, along

with the residuals, are shown in Figures

and

In the case of H

O, we only focus on the rotational

transitions detected around 13

–

16.5

m, for which we found a

single component with

600 K,

= ́

310cm

and

slab

0.27 auis able to reproduce the features. More detailed

Figure 1.

JWST MIRI-MRS spectrum of DoAr 33. The insets show isolated emission features of some of the molecules detected, and the calculated continuum is

shown as a pink line. Most of the unmarked emission lines correspond to water vapor. We discuss the presence of dust emission features in Appendix

Table 1

Stellar and Disk Properties of DoAr 33

Property

Value

Reference

Spectral type

J. Bouvier & I. Appenzeller

(

1992

)

Distance

139 pc

Gaia Collaboration et al.

(

2018

)

eff

4467 K

S. M. Andrews et al.

(

2010

)

Luminosity

1.5

e

S. M. Andrews et al.

(

2010

)

Mass

1.1

e

S. M. Andrews et al.

(

2018

)

Inclination

J. Huang et al.

(

2018

)

Accretion rate 2.52

−

e

−

L. A. Cieza et al.

(

2010

)

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

modeling, including more than one temperature component or a

temperature gradient, could be implemented to simultaneously

fi

t for the rotational transitions across the MIRI bandpass

(

A. Banzatti et al.

2023

; C. E. Muñoz-Romero et al.

2024

;

K. R. Schwarz et al.

2024

; M. Temmink et al.

2024b

)

;

however, this lies beyond the scope of this work. The water

rotational emission is believed to arise from larger disk radii

than for the vibrational mode emission

(

D. Gasman et al.

2023

)

; therefore, this approach allows us to isolate water

emission that is cospatial with that of the organics in the

inner disk.

fi

t the

branch of C

and

∼

13.7

while limiting the emitting area and temperature of both species

to be the same. We

fi

nd that models with

500

slab

0.17 au,

−

, and C

4.64 are able to reproduce the full acetylene feature. Given the

local interstellar medium

C ratio of 68, this suggests that

emission is optically thick, although not optically thick

enough to blend into a pseudocontinuum as observed in other

carbon-rich sources

(

B. Tabone et al.

2023

)

. Because of this

and the large degeneracies that are present when simulta-

neously

fi

tting three species, we determine the actual C

column from

with an assumed isotopic ratio. We

discuss this further in Section

4.2

We also detect the presence of rotational lines of OH, which

fi

t across a larger wavelength range to obtain a better

constraint of the parameters. The regions used for the

fi

t are

shown shaded in Figure

. The best-

fi

t temperature,

950 K, is higher than that of water

(

600

)

with a low column density of

−

. Despite its

relatively high

fl

ux at longer wavelengths, the column density

and temperature have large uncertainties. We also checked for

the presence of lines due to prompt emission between 10 and

m caused by UV photodissociation of water

(

J. R. Najita

et al.

2010

)

. For reference, we created a high-temperature

(

6000 K

)

OH slab model, which we show in the top panel of

Figure

. The insets show that some of the lines between 10

and 12

m show the

-doublet asymmetry previously observed

with Spitzer

(

J. S. Carr & J. R. Najita

2014

)

and predicted by

B. Tabone et al.

(

2021

2024

)

, indicative of prompt emission.

These features are very weak compared to the thermal OH

observed at longer wavelengths, suggesting that this source has

a low UV

fi

eld.

In the case of diacetylene

(

)

, we only detect emission

from the

branch at

∼

15.9

m, which we are able to

fi

with a model with

= ́

210

−

150 K, and

slab

1.36 au. We also explored models where we

fi

xed the

temperature and emitting area to the one found for C

. This

is discussed further in Section

4.3.1

. We also report a tentative

detection of methane

(

)

∼

7.65

m and HC

Nat

∼

15.07

m, which we further discuss in Section

5.5

. We also

detect the presence of CO emission around 5

m; however,

because these features are from optically thin transitions at high

energy levels

(

e.g., M. Temmink et al.

2024a

)

, they cannot be

used to determine the CO column density.

3. Methods

3.1. Disk Modeling

For modeling the disk, we used the Dust and Lines

(

DALI

)

thermochemical code

(

S. Bruderer et al.

2012

; S. Bruderer

2013

)

, with the extensions implemented by A. D. Bosman et al.

(

2022b

)

to reproduce the inner disk regions more accurately.

Figure 2.

Continuum-subtracted spectrum and best-

fi

t slab models for the molecules detected. The data were offset vertically for clarity. The pink-shaded region

shows the composite model of C

, HCN, CO

O, C

, and HC

N. We show a closer look at the HC

N feature in Figure

Table 2

Best-

fi

t Parameters from the Slab Models

Molecule

log

col

slab

(

−

)(

)(

)

18.5

0.3

600

0.27

17.8

0.5

1.0

200

0.77

17.0

0.3

500

0.17

16.3

(

500

)(

0.17

)

HCN

14.0

600

3.99

15.3

150

1.36

14.5

950

1.14

Note.

The uncertainties are reported only for cases with a closed 1

contour,

based on Figure

(

)

This column is uncertain; see Section

for a discussion.

(

)

The emitting area and temperature are

fi

xed to that of C

during the

fi

tting.

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

These extensions include the increased formation of H

, water

UV shielding, and dissociation heating

(

A. E. Glassgold &

J. R. Najita

2015

)

. As DoAr 33 shares a similar spectral type

classi

fi

cation

(

)

with AS 209

(

)

, we use the AS 209 stellar

spectrum from K. Zhang et al.

(

2021

)

. The UV emission,

generated by accretion

(

E. Gullbring et al.

2000

)

,isan

important driver of the surface chemistry in low-mass T Tauri

systems

(

P. Woitke et al.

2009

; I. Kamp et al.

2024

)

. Based on

the lower accretion rate of DoAr 33

(

2.5

−

e

−

)

compared to AS 209

(

−

e

−

)

, we decreased the

FUV emission by a factor of

∼

100 to account for the difference

in accretion rates between the two sources. The low UV

fl

ux is

consistent with the lack of prompt OH emission due to H

photodissociation, as discussed in Section

2.2

. We show the

spectra of both sources in Figure

In terms of the physical structure of the disk, we consider a

viscous accretion disk

(

D. Lynden-Bell & J. E. Pringle

1974

)

that follows the form

()

⎜⎟⎜⎟

⎛

⎝

⎞

⎠

⎡

⎣

⎢

⎛

⎝

⎞

⎠

⎤

⎦

⎥

S=S

exp

gas

where

gas

is the surface density of the gas,

is the critical

radius,

is the surface density at

, and

is the power-law

slope. The scale height angle is de

fi

ned by

()

⎜⎟

⎛

⎝

⎞

⎠

where

is the

fl

aring angle, and the physical scale height is

(

)

∼

. We include two different dust populations made up

of small

(

0.005

–

)

and large

(

0.005

–

1000

)

dust grains,

distributed in the disk according to

()

()

()

()

⎜⎟

⎡

⎣

⎢

⎛

⎝

⎞

⎠

⎤

⎦

⎥

-S

exp

dust,small

dust

()

()

()

⎜⎟

⎡

⎣

⎢

⎛

⎝

⎞

⎠

⎤

⎦

⎥

exp

dust,large

dust

Here,

dust

is the surface density of small and large dust,

the large dust fraction, and

describes the settling of large

dust. The values for all of the parameters used can be found in

Table

. The opacities of both populations follow the DSHARP

opacities

(

T. Birnstiel et al.

2018

)

used by K. Zhang et al.

(

2021

)

. The photometry of DoAr 33 shows an emission bump

around

∼

100

(

Y. Liu et al.

2022

)

. In order to reproduce this

feature, we increased the scale height of the large dust outside

10 au. This effectively mimics a wall of dust, consistent with

the gap substructure at

∼

9 au inferred from millimeter

observations of the source

(

J. Huang et al.

2018

)

. Our spectral

energy distribution

(

SED

)

fi

tting procedure is further explained

in Appendix

. The resulting 2D gas and dust density pro

fi

les

are shown in Figure

In order to properly model the full hydrocarbon chemistry

and include the relevant species, we use an expanded chemical

network as presented in S. E. Duval et al.

(

2022

)

. The original

network follows that implemented by C. Walsh et al.

(

2015

)

with the modi

fi

cations added by A. D. Bosman et al.

(

2022b

)

We also vary the initial abundances in the disk in order to

mimic the environment of an elevated C

O ratio. We have two

different scenarios that give rise to a change in C

O: an excess

of CH

or an excess of CH

accompanied by H

O depletion.

At the same time, in order to explore the origin of excess

carbon, we apply the changes in abundances to two regions:

inside a dust temperature of

dust

250 K or inside of

dust

500 K. The latter temperature has been argued to be

associated with the sublimation of carbon-rich grains at the

“

soot line

”

by J. Li et al.

(

2021

)

. This is a process that would

locally elevate the C

O ratio as carbon grains carry

∼

50% of

elemental carbon in the interstellar medium

(

ISM; A. Mishra &

A. Li

2015

; E. A. Bergin et al.

2015

)

. These two setups will

allow us to differentiate between a scenario where the excess

carbon is produced due to the sublimation of carbon grains at

the soot line

(

500 K model

)

, which would happen in the inner

disk, therefore resulting in high temperatures and a scenario

where the carbon is sublimated off icy grains in the outer disk

and then transported inward via advection resulting in a more

Table 3

Thermochemical Model Parameters

Parameter

Symbol

Value

Critical radius

20 au

Gas surf. dens. at

13 g cm

−

Surf. dens. power-law slope

0.9

Disk opening angle

0.045

Disk

fl

aring angle

0.08

Large dust fraction

0.999

Large dust settling

0.1 for



10 au

0.8 for

10 au

Figure 3.

Gas

(

top

)

and dust

(

bottom

)

density pro

fi

les from the thermochemical

models. The orange line in the bottom plot shows the continuum

1 line at

m. Inside 10 au, the large dust population is mostly settled in the midplane.

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

spatially extended distribution with overall lower temperatures

(

250 K models

)

. This would be roughly consistent with the

model of the origin of hydrocarbon-rich chemistry suggested

by J. Mah et al.

(

2023

)

. Outside of the 250 K and 500 K

regions, the abundances of carbon and oxygen carriers are set

so that the C

O ratio is solar. As a standard for comparison, we

also include a model with solar abundance everywhere in the

disk. This results in a grid of 19 models, for which we indicate

the names and abundances considered in Table

In this work, we refrain from comparing the observations to

simulated spectra, in part due to the uncertainty that can arise

from accounting for mutual line overlap, which can be

signi

fi

cant for molecules like C

, as demonstrated in

B. Tabone et al.

(

2023

)

. But, more importantly, we aim to

focus on the retrieved outputs from the data, rather than on the

details of the radiative transfer. Therefore, we will compare to

column densities, emitting areas, and

fl

ux-weighted tempera-

tures

(

discussed below

)

. This reduces the time required to

forward model this source using complex thermochemical

models. This solution should approximate the key effects and

enable a direct exploration of why DoAr 33 is a rare K star

exhibiting rich hydrocarbon emission.

3.2. Flux-weighted Temperatures

Column density estimates derived from observations can

intrinsically reduce the uncertainty arising from not knowing

how the gas is distributed vertically in the disk since they are

integrated along the line of sight; therefore, they can be directly

compared to the column densities calculated in the thermo-

chemical models, without having to match a 2D gas

distribution. This is not the case for the emitting temperature,

which requires that we consider where in the disk this emission

is coming from. From the thermochemical models, we are able

to obtain a 2D temperature distribution for both the gas and the

dust components; we show these distributions in Figure

However, in order to compare to the estimates from the slab

models, we must calculate a representative value for each

molecule in the thermochemical models. Here, we describe the

procedure we applied to condense the 2D temperature

structures into a single emitting temperature.

In the LTE slab models, the temperature retrieved for each

species is a representative temperature at which the gas we can

observe is emitting. In order to draw a similar parallel, we make

use of contribution functions, which trace where in the disk an

emission line is formed. Once the necessary information about

energy levels, statistical weights, and Einstein A coef

fi

cients is

provided to DALI in the form of LAMDA

fi

les

(

F. L. Schöier

et al.

2005

)

, the code is able to perform excitation calculations

for any given molecule included in the chemical network. With

this information, we can calculate the 2D contribution function

of a speci

fi

ed transition

(

S. Bruderer et al.

2012

)

to understand

where it is originating physically. It is important to note that the

region where a speci

fi

c line is generated does not necessarily

correspond to the region where that species is most abundant

due to, e.g., the effects of optical depth. Therefore, this method

allows us to mimic an emitted temperature instead of a

temperature simply weighted by abundance. In the case of

, the spectroscopic information for the LAMDA

fi

le was

compiled based on the ExoMol database

(

K. L. Chubb et al.

2020

)

In practice, we select a speci

fi

c transition for each of the

species observed in the MIRI spectrum and calculate its

contribution function in DALI. To select the transition, we use

the spectral line analysis tool

iSLAT

(

E. Jellison et al.

2024

)

which takes line transition information from HITRAN

(

I. E. Gordon et al.

2022

)

. This tool allows us to input the

parameters obtained from the

fi

tting of slab models and

generate a synthetic spectrum. Once the spectrum is generated,

a spectral range can be selected for which the strongest line

transition is shown. We then select the strongest transition of

each molecule in the spectral range between 12 and 17

m for

which we calculate the contribution function using DALI. The

selected transitions are indicated in Table

and plotted for the

Table 4

Input Abundances of the Thermochemical Models Inside the 250 K and 500 K Contours

Name

Total C

Total O

OCOCO

Solar

5.00

−

1.44

−

0.47

1.20

−

8.50

−

4.15

−

8.50

−

ch4_excess_co_1_t500

2.03

−

4.50

−

1.20

−

8.50

−

4.15

−

1.62

−

ch4_excess_co_2_t500

4.91

−

1.02

−

1.20

−

8.50

−

4.15

−

4.49

−

ch4_excess_co_1_t250

2.03

−

4.50

−

1.20

−

8.50

−

4.15

−

1.62

−

ch4_excess_co_2_t250

4.91

−

1.02

−

1.20

−

8.50

−

4.15

−

4.49

−

co_dep_co_1_t500

1.20

−

2.49

−

1.20

−

8.50

−

1.20

−

co_dep_co_2_t500

2.49

−

5.06

−

1.20

−

8.50

−

2.49

−

co_dep_co_4_t500

5.06

−

1.02

−

1.20

−

8.50

−

5.06

−

co_dep_co_1_t250

1.20

−

2.49

−

1.20

−

8.50

−

1.20

−

co_dep_co_2_t250

2.49

−

5.06

−

1.20

−

8.50

−

2.49

−

co_dep_co_4_t250

5.06

−

1.02

−

1.20

−

8.50

−

5.06

−

co_water_dep_co_1_t500

1.20

−

3.25

−

1.20

−

8.50

−

1.20

−

co_water_dep_co_2_t500

3.25

−

7.35

−

1.20

−

8.50

−

3.25

−

co_water_dep_co_4_t500

7.35

−

1.56

−

1.20

−

8.50

−

7.35

−

co_water_dep_co_1_t250

1.20

−

3.25

−

1.20

−

8.50

−

1.20

−

co_water_dep_co_2_t250

3.25

−

7.35

−

1.20

−

8.50

−

3.25

−

co_water_dep_co_4_t250

7.35

−

1.56

−

1.20

−

8.50

−

7.35

−

water_dep_no_co

8.50

−

2.55

−

1.75

1.20

−

8.50

−

8.50

−

water_dep_co

1.13

−

2.55

−

1.36

1.20

−

5.67

−

2.83

−

8.50

−

Note.

The abundances are reported with respect to H.

https:

home.strw.leidenuniv.nl

moldata

https:

github.com

spexod

iSLAT

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

solar model in Figure

. We map these contribution functions

to the gas temperature in the disk model and calculate a single

contribution function-weighted temperature for each model.

These temperatures are discussed in Section

, where we

compare them to the temperatures derived from slab models.

4. Results

The thermochemical modeling process returns 2D distribu-

tions of abundances in the disk for all of the species included in

the chemical network. Here, we only focus on the oxygen and

carbon carriers that are detected in the MIRI-MRS spectrum of

DoAr 33. In order to directly compare the modeling products to

the parameters retrieved from the LTE slab models, we

compute column densities and

fl

ux-weighted temperatures for

all of the models described in Table

based on the abundance

maps. We calculated the radial distributions of column

densities by vertically integrating the abundance from the

surface of the disk down to the

1 surface at 15

(

shown

as an orange line in Figure

)

. For the contribution function-

weighted temperature, we followed the procedure outlined in

Section

3.2

In Sections

4.1

–

4.3

, we individually discuss the results for

each molecule. We include all of the abundance maps for H

, and C

in Appendix

. To make the analysis more

tractable, we did not include HCN and HC

N and did not

attempt to match the retrieved columns for these species.

4.1. H

O, CO

, and OH

In the top left panel of Figure

, we show the column density

of water as a function of radius for all of the models, colored

according to the C

O ratio assumed. The line styles indicate

whether we assume the increased C

O is a product of only

carbon enhancement or the combination of carbon enhance-

ment and oxygen depletion, and we include a model with solar

abundances for comparison.

The solar abundance model predicts the overall highest

column density for water, which then decreases beyond 1 au.

The high column density of the solar abundance model appears

to be inconsistent with the reference values retrieved through

the LTE slab models

(

Table

)

. In fact, all of the models that

consider the default water abundance of

∼

−

, shown in solid

lines, overpredict the column density of water. Only the models

where we include water depletion are able to match the water

column. However, there are several important caveats to

determining an overall water abundance through column

densities. For one, we do not model the entire water spectrum,

and therefore, the

fi

t results reported for water in Table

are

only approximately valid

(

see A. Banzatti et al.

2023

;

D. Gasman et al.

2023

; C. E. Muñoz-Romero et al.

2024

)

We note that the retrieved water column is consistent with

values obtained for other full disks like GW Lup

(

S. L. Grant

et al.

2023

)

and FZ Tau

(

K. M. Pontoppidan et al.

2024

)

, but in

this latter disk, the retrieved water column of

∼

−

was

found to be consistent with solar abundances. In our case, the

models with the canonical water abundance of

∼

−

predict

columns closer to 10

−

. This discrepancy could be, in

part, because some of the columns we modeled might not have

the right temperature for IR emission, as has been noted by

R. Meijerink et al.

(

2009

)

and A. D. Bosman et al.

(

2022a

)

. The

modeled columns depend on how far into the disk we assume

we can see, which is highly dependent on the gas-to-dust ratio.

In our models, the water emission is arising from material with

a gas

dust ratio of

∼

, which is consistent with previous

fi

ndings of a well-settled dust disk

(

Y. Liu et al.

2012

)

.Itis

possible that a potential solution to the water column issue

(

e.g., R. Meijerink et al.

2009

)

is a lower hydrogen column,

which would also have implications for the overall carbon

chemistry. This is explored in Section

5.2

. Stronger observa-

tional constraints on the total gas and dust content are needed to

be able to model the water column more accurately.

We then examined the column densities of the second most

signi

fi

cant oxygen tracer in the MIRI spectrum, CO

. In the top

right panel of Figure

, we show the column densities

calculated for this molecule. In contrast to water, the column

density of this molecule is highest around 1 au. Both H

O and

share the same precursor, OH. A. D. Bosman et al.

(

2022b

)

and S. E. Duval et al.

(

2022

)

suggested that CO

formation in the warm

(

400 K

)

gas is curtailed as OH is more

Figure 4.

Gas

(

top

)

and dust

(

bottom

)

temperature pro

fi

les from the

thermochemical models. The colored contours in the top

fi

gure represent the

contribution functions for the

water

dep

co

model

(

see Table

)

, where

90% of the emission is produced. The transitions for each molecule are

indicated in Table

. The lines in the bottom plot represent the 250 K and

500 K contours, inside of which we change the abundances to increase the C

ratio.

Table 5

Strongest Line Transitions in the 12

–

m Region for Best-

fi

t Slab Model

Parameters

Species

Einstein A Coeff.

Wavelength

(

−

)(

)(

)

3.48

3950.8

14.2096

1.54

1113.3

14.9777

6.06

1358.3

13.6992

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

likely to react with H

than CO. The models are only able to

reproduce the column density retrieved from slab models in the

zone around 1 au. For reference, the emitting area estimated

from the slab models is shown as a shaded annulus of the same

area. This indicates that the CO

emission is likely coming

from a ring of gas around 1 au. One implication of this result is

that it would suggest that CO

emission would trace cooler gas

than regions with abundant water vapor, which has been found

in some instances via analysis of JWST spectra

(

C. Xie et al.

2023

; D. Gasman et al.

2023

)

In the bottom left panel of Figure

, we compared the

retrieved OH column density. As expected, the models with a

lower C

O ratio produce the highest OH columns. More

speci

fi

cally, the model with a solar C

O ratio predicts a column

of almost 10

−

, much higher than the retrieved best

fi

tof

−

. Nonetheless, the uncertainty of the column

density does not allow us to rule out C

1, pointing more

toward a C

1. We note that despite the weak prompt

emission at shorter wavelengths, OH emission may not be in

non-LTE

(

see discussion for H

O by B. Tabone et al.

2021

)

;

therefore, the best-

fi

t parameters from the LTE models might

not fully represent the local gas properties.

In Figure

, we show the model emitting temperatures for

,CO

, and H

O, along with the estimate from the slab

model

fi

tting. We show only the water depletion scenarios

since they are consistent with the water column densities from

the slab models. Overall, the thermochemical models predict

water has the highest temperature, and CO

has the lowest,

which is consistent with the expectation from the slab model

retrievals from observations. Nonetheless, the CO

and H

temperatures do not allow us to differentiate between the two

different carbon enhancement scenarios proposed

(

see

Section

3.1

)

4.2. C

and

As mentioned in Section

2.2

, in order to derive the actual

column, we assume an isotopic ratio. This ratio is

unknown, and B. Tabone et al.

(

2023

)

assumed 35 based on the

models of P. M. Woods & K. Willacy

(

2009

)

. We adopt this

value and report the C

column density to be

CCH

13 12

which results in

= ́

7.5 10 cm

. We note this value

could be substantially higher as the

C ratio of

hydrocarbons and nitriles in the outer disk of TW Hya is

closer to the interstellar value

(

P. Hily-Blant et al.

2019

;

Figure 5.

Column density for H

O, CO

, and C

as a function of radius for the models in Table

. The color bar re

fl

ects the value of C

O for each model, and the

linestyle indicates the method for increasing the C

O ratio. The dashed

–

dotted horizontal line shows the best-

fi

t value found from LTE slab models, along with the

shaded 1

uncertainty. The shaded gray vertical region corresponds to the emitting area estimated from slab models. In the case of CO

, we assumed this value to be

an annulus of the same area instead of a continuous distribution, and we placed the ring where it best matched the thermochemical models.

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

E. A. Bergin et al.

2024

; T. C. Yoshida et al.

2024

)

. Therefore,

the uncertainties of the C

column shown in Figure

stem

from assuming

C can go from

≈

5, as directly retrieved

from slab models to 68, the interstellar value.

From the thermochemical models, it is clear to see that,

similarly to water, the C

emitting column peaks in the inner

disk and then decreases between 0.4 and 1 au. The difference

in the location of the column density drop is due to where the

O ratio is increased. In the models where C

O is increased

within

500 K, the carbon-rich chemistry is restricted to a

smaller area, whereas in the other models, the C

column is

more spatially extended. It is clear that most models that are

able to match the expectation from the slab model

fi

t have

Contrary to CO

and H

O, in the case of C

, the selection

of the temperature where we transition from solar to a high

O ratio

(

i.e., 250 K versus 500 K

)

plays an important role in

determining the

fi

nal weighted temperature, as we can see in

Figure

, most of the models cluster around a particular value.

In the models where the carbon-rich chemistry is allowed to be

more spatially extended, the C

has a lower temperature.

This implies that a high C

O ratio extended over larger radii

would produce an overall lower emitting temperature for the

hydrocarbons. Similarly, by restricting the extent of the high

O ratio

(

i.e., the 500 K transition models

)

, we limit the

presence of C

to the innermost disk, which in turn results in

a higher emitting temperature since the material is closer to the

star. Comparing to the temperature retrieved from slab models,

we are not able to strictly discriminate between the two

scenarios since most of the

fl

ux-weighted temperatures are

within the uncertainty range. This is also true for water, where

there is an overlap in the

fl

ux-weighted temperatures.

In our approach, the key to distinguishing between a radially

extended or limited carbon-rich chemistry is the emitting area.

Slab models retrieve a maximum emitting area of

()

0.25 au

for C

, corresponding to the 3

contour shown in Figure

The 500 K contour for the gas temperature extends to 0.4 au in

the IR-emitting layer

(

Figure

)

, whereas it goes past 1 au for

the 250 K contour. This small emitting area, in combination

with the high emitting temperature of C

, indicates that the

carbon chemistry is only released in the innermost disk, where

the sublimation of carbonaceous grains is possible. If the

carbon-rich material present in this system traced by C

caused by an advection of material from further out in the disk

(

J. Mah et al.

2023

)

, then it should produce more extended

emission than inferred here.

4.3. Other Hydrocarbons

4.3.1. C

We compare the retrieved C

column to that predicted by

our thermochemical models in Figure

. Based on the

maps

(

Figure

)

from the slab model

fi

tting, this column likely

represents an upper limit. In turn, this makes it consistent with

the models where we assumed C

2. The C

emission is

extended over a larger area than C

and is therefore

compatible with a lower-emitting temperature

(

150 K

)

.To

provide a more de

fi

nite constraint, we explored

fi

tting C

with the same emitting temperature as C

in Appendix

and

found that the higher temperature is incompatible with the

branch feature. This indicates the presence of a warm

(

400 K

)

hydrocarbon-rich reservoir in addition to the hot

(

400 K

)

. We discuss the potential origin of these two reservoirs in

Section

5.3

4.3.2. Nondetections

The chemical network used in this work is able to predict

chemical abundances for over 400 species. Therefore, we

checked the predictions for the abundance of hydrocarbons that

have been detected in hydrocarbon-rich low-mass sources

(

B. Tabone et al.

2023

; A. M. Arabhavi et al.

2024

; J. Kanwar

et al.

2024a

)

, in order to see if we can obtain any additional

constraints on the C

O ratio. In the top row of Figure

,we

show the predicted column densities for C

, and

. In the case of C

and C

, the predicted columns are

lower than

∼

−

for all of the models. Based on this, we

assumed an average column for C

and C

and produced

LTE slab models with two different sets of temperatures and

emitting areas, mimicking the ones found for C

and C

Figure 6.

Contribution function-weighted temperatures of carbon and oxygen-

bearing species for the water depletion models. The symbols represent where in

the model the C

O ratio is changed, inside

250 K or

500 K. The

markers are shifted horizontally for clarity. For comparison, we include the

temperature retrieved from slab models as a black square, with errors

corresponding to the 1

uncertainty.

Figure 7.

Same as Figure

, but for C

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

The lower panel of Figure

shows the produced spectra, as

well as the parameters assumed for each one. It is clear that the

emission from these molecules would be too faint to be

observable in the DoAr 33 spectrum. This indicates that our

model prediction is consistent with the nondetection of C

and C

. The columns for these two molecules do not have

any discernible trend in terms of the C

O ratio, making them

poor C

O tracers for this source. Despite this prediction being

consistent with our observations, we note that the predicted

columns for these species are lower than the ones being

detected in low-mass sources

(

A. M. Arabhavi et al.

2024

;

J. Kanwar et al.

2024a

)

. This could be because the formation of

hydrocarbons larger than C

might be sensitive to the

explicit X-ray luminosity and the presence or absence of

cosmic rays

(

E. Raul et al.

2024

)

or due to the completeness of

the hydrocarbon reactions in the chemical network

(

J. Kanwar

et al.

2024b

)

In contrast, C

follows a similar trend to C

, where high

column densities are predicted for the high C

O ratio models,

as shown in the top right panel of Figure

. Because of this, we

attempted a retrieval in the 14.8

–

14.9

m range, where the

main feature is located, using the spectroscopic

fi

les and

partition functions provided by A. M. Arabhavi et al.

(

2024

)

fi

rst subtracted the contribution from all of the other

detected molecules, adopting the parameters from the best-

fi

slab models Table

. Then, we estimated the column density of

, assuming a given emitting area and temperature. Similar

to the other hydrocarbons, we tested two cases, one where the

is cospatial with C

, i.e., a small emitting area and high

temperature, and one where the C

is cospatial with C

having a lower temperature and bigger emitting area. The

resulting column density retrieved for each case is overplotted

with the model predictions as horizontal lines in the top right

panel of Figure

. The bottom right panel shows the spectra of

both sets of best-

fi

t parameters, overplotted with the observed

spectrum and the composite model of the detected molecules.

Both sets of solutions point toward C

–

4, but we are not

able to con

fi

rm the presence of C

in this spectrum.

5. Discussion

fi

nd that C

1 is needed to explain the hydrocarbon-

rich environment around DoAr 33. The two methods to

increase the C

O ratio are to deplete oxygen and

or enrich

carbon. Here, we

fi

rst discuss how C

O can be increased

through water sequestration. Then, we argue that C

cannot be achieved through oxygen depletion alone and

propose that additional carbon must be supplied to the gas

phase. We explore how much excess carbon is inferred to be

present in the inner disk of DoAr 33 and consider the two main

mechanisms that have been proposed to explain the origin of

this carbon excess. We also discuss the implications of CH

detection in the context of carbon grain sublimation.

5.1. Depletion of Water

It has been posited that planetesimals and protoplanets

forming beyond the water iceline could trap a signi

fi

cant

fraction of the oxygen content as long as the ef

fi

ciency of

accretion onto larger bodies is higher than the drift rate of icy

pebbles

(

J. R. Najita et al.

2013

)

. Modeling of dust and gas

dynamics also showed that outer disk gaps, if deep enough, can

reduce or stop the inward

fl

ux of icy pebbles and reduce the

oxygen enrichment of inner disks by reducing water ice

sublimation inside the snowline

(

A. Kalyaan et al.

2021

2023

)

These scenarios would raise the C

O as the inner disk accretes

the remaining oxygen-poor gas. Stellar accretion proceeds over

several million years, encompassing the initial phases of rocky

planet assembly and the main phase of giant planet formation

Figure 8.

The top row shows the column density predictions from the thermochemical models for C

, and C

. The bottom row shows slab models that

simulate the conditions predicted by the thermochemical calculations. The parameters used for each slab model are listed in the

fi

gure. Two sets of slab models were

tested, one cospatial with C

(

orange

)

and one cospatial with C

(

green

)

. These models

(

green and orange

)

are shown on top of the data

(

black

)

in the bottom left.

The pink-shaded region in the bottom right panel corresponds to the same composite model as in Figure

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.

(

L. Hartmann et al.

2016

)

. Thus, if a signi

fi

cant population of

planetesimals forms, this is a possible scenario for raising C

However, in this instance, the volatile reservoir in the inner

disk would be mostly dominated by CO and CO

gas, meaning

that, at most, the inner disk could reach a C

O ratio of

∼

−

1. In the case of DoAr 33, we have shown that

1 is needed to reproduce the C

column. Therefore,

even though planetesimal formation in the outer disk could be

invoked to increase C

O, an additional release of carbon is

needed to go beyond C

5.2. Total Carbon Content

To understand the source of the carbon enrichment, we

fi

rst

examine the amount of excess carbon that is consistent with the

observations of DoAr 33. In our thermochemical models, we

assume varying degrees of extra carbon being released into the

system initially in the form of CH

inside the soot line

(

Table

)

. At the midplane in our model the soot line is located

∼

0.08 au

(

Figure

)

, which results in a maximum of

∼

0.5

−

⊕

of carbon released as CH

gas for model

water

dep

co

(

8.5

−

)

, calculated by integrat-

ing the total CH

mass in the disk inside 500 K. Conversely, in

the dust distribution modeling process

(

Appendix

)

,we

considered that the amount of refractory organics carried by the

large dust population accounted for 40% of the total dust mass

(

T. Birnstiel et al.

2018

)

. This results in

∼

1.5

−

⊕

carbon being present as solids inside 0.08 au. Therefore, our

thermochemical models require that, at most, only a third of the

carbon present in refractories is released into the gas.

The presence or absence of carbon grains can strongly

fl

uence the carbon content of forming planetary systems

(

E. A. Bergin et al.

2023

)

. In the context of the solar system,

the sublimation of carbon grains has been posited to explain the

carbon depletion found in meteorites and the Earth

(

E. A. Bergin et al.

2015

; J. Li et al.

2021

)

. It is possible that

the carbon grains that are not sublimated to go into forming

planetesimals. However, we stress our models are not self-

consistent. The SED modeling, which includes predetermina-

tion of the dust composition as part of the dust optical

constants, is a separate step from the thermochemical

calculations, each with independent assumptions of the

elemental carbon content in gas and solids.

Thus, it is entirely possible, given the similarity between the

estimate of the excess carbon in the gas and the total carbon

available in the dust

(

within a factor of 3

)

, that all of the solid-

state carbon inside 0.08 au has been released to the gas. We

note that this solution assumes the presence of a dust trap, as

potentially indicated by the Atacama Large Millimeter

submillimeter Array

(

ALMA

)

continuum and the SED, and

water depletion to lower the water column. Alternately, as

discussed in Section

4.1

, the overall gas

dust ratio may be

lower than assumed by at least an order of magnitude. This

would reduce the water column and allow for solutions with

“

normal or solar

”

water content

(

i.e., the dotted and the solid

lines in Figure

)

. In this instance, the C

column would also

be reduced to near the lower end of the allowed range in the

column density. Here, higher C

O ratios and carbon content

are required to bring the column closer to the middle of our

retrieved estimate

(

∼

few

−

)

. Future work that

combines the effects of carbon destruction on the emission of

both solids and gaseous species is needed to ensure consistency

and derive more accurate estimates. Finally, while carbon

grains represent one potential source term for excess carbon,

other scenarios have been proposed, which are discussed

below.

5.3. Origin of Carbon Enrichment

Currently, there are two potential models for the origin of the

excess carbon. One posits that the carbon supply can originate

in the outer disk via gas accretion of material with C

(

J. Mah et al.

2023

)

, and the other suggests that carbon grain

destruction near the soot line can increase the carbon content of

the inner disk

(

J. Li et al.

2021

)

We will discuss each in turn but note two things:

fi

rst, in our

analysis, the difference in emitting temperature and area found

for C

and C

suggests that the carbon-rich chemistry in

this source is not necessarily restricted to the innermost disk.

On the same note, ALMA has observed many disks that appear

to have C

1 beyond 10 au

(

e.g., A. D. Bosman et al.

2021b

)

, and temperatures are hot enough to sublimate carbon

grains only in the inner 1 au for T Tauri disks

(

J. Li et al.

2021

)

;

thus, these solutions are likely not exclusive. However, they do

have different implications for the composition of planets.

5.3.1. Carbon-rich Gas Advection

The model presented by J. Mah et al.

(

2023

)

suggests that

the accretion of CH

-rich gas into the inner disk can explain the

excess carbon observed in low-mass sources. This scenario

works if pebbles in the outer disk are large enough to be close

to the drift limit, with fragmentation velocities

(

frag

)

higher

than 5 m s

−

. In this case, the system gets an early in

fl

ux of

water-rich pebbles crossing the water iceline, enhancing the

gas. The in

fl

ux of water-rich pebbles eventually runs out, and

the oxygen-enriched gas gets accreted by the star within

∼

2 Myr. Only then, the system is left with carbon-rich gas that

can be advected into the inner disk.

In this concept, CH

is released into the gas with an

abundance that exceeds both CO and CO

providing

(

)

(

2CO

)

1. This

leads to an important point. In protostellar envelopes, the CO

and CO

ice abundances have been determined via Spitzer

spectroscopic surveys and are

∼

60% relative to water

(

K. I. Öberg et al.

2011a

)

. Assuming a water ice abundance

of order

∼

−

(

relative to H

)(

K. M. Pontoppidan et al.

2004

)

, then the typical abundance of CO and CO

ice are of

order

∼

−

relative to H

. This does not account for the gas-

phase CO and CO

content, and thus the gas & ices carry

25% of the cosmically available carbon, which is 4.28

−

relative to H

(

M.-F. Nieva & N. Przybilla

2012

)

. In sum, if CO

and CO

are present at levels expected in protostellar gas

(

inside their respective ice lines

)

, then to achieve C

1, one

must destroy nearly

all

carbon-rich grains, which hold an

equivalent amount of elemental carbon as the gas

ices. This

seems implausible beyond 10 au as the destruction mechanisms

for carbon grains in the outer disk are highly ineffective

(

D. E. Anderson et al.

2017

; L. Klarmann et al.

2018

)

However, the ALMA observations of C

1 also coincide

with many systems where the CO gas-phase abundance has

been reduced

(

A. Miotello et al.

2023

)

. Simulations account for

the presence of strong C

H emission in these disks via the

release of CH

from grains at the level of a few percent of the

solid-state carbon inventory

(

A. D. Bosman et al.

2021a

)

These results empirically infer that

(

)

(

The Astrophysical Journal,

977:173

(

24pp

)

, 2024 December 20

Colmenares et al.