[a]Cédric Lorcé

# 3D structure of hadrons and energy-momentum tensor

###### Abstract

The three-dimensional spatial structure of hadrons is encoded in their form
factors. Via appropriate Fourier transform, the latter describe how charge,
energy, linear and angular momentum, but also pressure are distributed inside
these systems. Electromagnetic form factors of the nucleon have been measured
for a long time, but it is only recently that some gravitational form factors
have been extracted from experimental data, creating a great deal of
enthusiasm in the hadronic community. Here we summarize some of the recent
developments on the interpretation of these form factors, and provide a quick
overview of what we can learn about the nucleon mass, spin and internal
pressure.

##  1 Introduction

Understanding the inner workings of hadrons requires multidimensional mapping
of the distribution of their constituents [1, 2, 3]. Feynman’s parton picture
has now been promoted to a 5D representation based on the notion of light-
front Wigner (or phase-space) distribution [4], which gives the quasi-
probability to find a quark or a gluon with a given light-front three-momentum
$(\vec{k}_{\perp},xP^{+})$ at the transverse position $\vec{b}_{\perp}$.
Integrating out some of these variables, one obtains various lower-dimensional
distributions which can be probed via different (semi-)inclusive and exclusive
processes [5].

In this contribution, we will focus on the spatial structure as probed by
exclusive scattering experiments. The response of hadrons to electromagnetic
probes has been measured with extreme precision over the past decades,
allowing us to gain insight not only into the electromagnetic structure [6,
7], but also into the mechanical structure [8] of hadrons. More precisely,
elastic lepton-hadron scattering gives direct access to the electromagnetic
current $\langle p^{\prime},s^{\prime}|j^{\mu}|p,s\rangle$, where $p$ is the
four-momentum and $s$ is the polarization of the hadron, whereas deeply
virtual Compton scattering provides indirect information on the energy-
momentum tensor $\langle p^{\prime},s^{\prime}|T^{\mu\nu}|p,s\rangle$.

Thanks to Poincaré symmetry these matrix elements can be described by a
limited number of Lorentz-invariant functions of the four-momentum transfer
$\Delta=p^{\prime}-p$ called form factors (FFs). Their Fourier transform in
the Breit frame (BF) are usually interpreted as static 3D spatial
distributions [9, 10, 11]. Since relativistic recoil corrections hinder their
interpretation as probabilistic densities [12, 13, 14], alternative 2D
distributions on the light front (LF) have been proposed [15, 16, 17]. While
the latter admit a strict probabilistic interpretation, they show features
that are at odds with the expected picture of hadrons at rest [18, 19].

To clarify the situation, the concept of relativistic spatial distribution has
recently been revisited by a number of authors, see e.g. [20, 21, 22, 23, 24,
25, 26, 27]. The phase-space approach [28, 29, 30] differs from these attempts
by relaxing the requirement of strict probabilistic interpretation to a quasi-
probabilistic one. In doing so, it is possible to define relativistic spatial
distributions which depend on the hadron average momentum. In this picture, BF
and LF distributions appear as two particular limits, and the LF distortions
are understood as purely kinematical artifacts associated with spin.

##  2 Phase-space formalism

The LF formalism is closely related to the standard description in the
infinite-momentum frame (IMF) [31], i.e. the frame where the target travels at
nearly the speed of light. From a phase-space perspective, elastic frame (EF)
distributions [32, 28, 29] provide a natural interpolation between slow- and
fast-moving targets.

The expectation value of an operator $O$ in a physical state $|\psi\rangle$
can be written as

$\langle\psi|O|\psi\rangle=\int\frac{\mathrm{d}^{3}P}{(2\pi)^{3}}\,\mathrm{d}^{3}R\,\rho_{\psi}(\vec{R},\vec{P})\,\langle
O\rangle_{\vec{R},\vec{P}},$ (1)

where

$\rho_{\psi}(\vec{R},\vec{P})=\int\frac{\mathrm{d}^{3}q}{(2\pi)^{3}}\,e^{-i\vec{q}\cdot\vec{R}}\,\tilde{\psi}^{*}(\vec{P}+\tfrac{\vec{q}}{2})\tilde{\psi}(\vec{P}-\tfrac{\vec{q}}{2})$
(2)

is the Wigner distribution [33] and $\tilde{\psi}(\vec{p})=\langle
p|\psi\rangle/2\sqrt{p^{0}}$ is the wave packet of the system. The amplitude

$\langle
O\rangle_{\vec{R},\vec{P}}=\int\frac{\mathrm{d}^{3}\Delta}{(2\pi)^{3}}\,e^{i\vec{\Delta}\cdot\vec{R}}\,\frac{\langle
p^{\prime}|O|p\rangle}{2\sqrt{p^{0}p^{\prime 0}}}$ (3)

can then be interpreted as the expectation value of $O$ for a system localized
in phase space with average momentum $\vec{P}=(\vec{p}^{\prime}+\vec{p})/2$
and average position $\vec{R}$.

To obtain a static distribution in the case of a spacetime-dependent operator
$O(x)$, we need to pick up a frame where the energy transfer to the system
$\Delta^{0}=\vec{P}\cdot\vec{\Delta}/P^{0}$ vanishes. This is satisfied, e.g.,
in the BF defined by the condition $|\vec{P}|=0$. In the case of a moving
target, we can choose for convenience the $z$-axis along $\vec{P}$. The
elastic condition $\Delta^{0}=0$ is then satisfied if we project the
distribution onto the transverse plane

$O_{\text{EF}}(\vec{b}_{\perp};P_{z})=\int\mathrm{d}z\,\langle
O\rangle_{\vec{R},P_{z}\vec{e}_{z}}(x)=\int\frac{\mathrm{d}^{2}\Delta_{\perp}}{(2\pi)^{2}}\,e^{-i\vec{\Delta}_{\perp}\cdot\vec{b}_{\perp}}\,\frac{\langle
p^{\prime}|O(0)|p\rangle}{2P^{0}}\bigg{|}_{\Delta_{z}=|\vec{P}_{\perp}|=0}.$
(4)

Because of translation symmetry, this EF distribution depends only on the
relative transverse position $\vec{b}_{\perp}=\vec{x}_{\perp}-\vec{R}_{\perp}$
besides $P_{z}$.

For a relativistic system, the $P_{z}$-dependence cannot in general be
factored out. We however expect that [34, 35]

$\displaystyle\langle
p^{\prime},s^{\prime}|O^{\mu_{1}\cdots\mu_{n}}(0)|p,s\rangle=$ (5)
$\displaystyle\qquad\sum_{s_{B}^{\prime},s_{B}}{D}^{*(j)}_{s_{B}^{\prime}s^{\prime}}(p_{B}^{\prime},\Lambda){D}^{(j)}_{s_{B}s}(p_{B},\Lambda)\,\Lambda^{\mu_{1}}_{\phantom{\mu_{1}}\nu_{1}}\cdots\Lambda^{\mu_{n}}_{\phantom{\mu_{n}}\nu_{n}}\,\langle
p_{B}^{\prime},s_{B}^{\prime}|O^{\nu_{1}\cdots\nu_{n}}(0)|p_{B},s_{B}\rangle,$

where $\langle
p_{B}^{\prime},s_{B}^{\prime}|O^{\nu_{1}\cdots\nu_{n}}(0)|p_{B},s_{B}\rangle$
is the BF matrix element and $\Lambda^{\mu}_{\phantom{\mu}\nu}$ is the matrix
that implements the Lorentz boost from the BF to a generic frame. For a target
with spin $j$, the Wigner spin rotation matrix
${D}^{(j)}_{s_{B}s}(p_{B},\Lambda)$ plays a key role for understanding the
aforementioned strange features seen in the LF spatial distributions in terms
of distortions induced by the boost [30, 36, 37, 38]. We will restrict the
discussion to the case of spin-$1/2$ targets in the following.

##  3 Electromagnetic current

The 3D BF distributions of electric charge and current are defined in the
phase-space formalism as

$J^{\mu}_{B}(\vec{r})=\langle
j^{\mu}\rangle_{\vec{R},\vec{0}}(\vec{x})=\int\frac{\mathrm{d}^{3}\Delta}{(2\pi)^{3}}\,e^{-i\vec{\Delta}\cdot\vec{r}}\,\frac{\langle
p^{\prime}_{B},s^{\prime}_{B}|j^{\mu}(0)|p_{B},s_{B}\rangle}{2P^{0}_{B}}$ (6)

with $p^{\mu}_{B}=(P^{0}_{B},-\vec{\Delta}/2)$ and
$p^{\prime\mu}_{B}=(P^{0}_{B},\vec{\Delta}/2)$ the initial and final BF four-
momenta, and $\vec{r}=\vec{x}-\vec{R}$ the position relative to the center of
the system. They differ from the conventional Sachs distributions [9, 10] by
the normalization factor $2P^{0}_{B}$ instead of $2M$ in the denominator. In a
similar way, the corresponding 2D EF distributions are given by

$J^{\mu}_{\text{EF}}(\vec{b}_{\perp};P_{z})=\int\mathrm{d}z\,\langle
j^{\mu}\rangle_{\vec{R},P_{z}\vec{e}_{z}}(\vec{x})=\int\frac{\mathrm{d}^{2}\Delta_{\perp}}{(2\pi)^{2}}\,e^{-i\vec{\Delta}_{\perp}\cdot\vec{b}_{\perp}}\,\frac{\langle
p^{\prime},s^{\prime}|j^{\mu}(0)|p,s\rangle}{2P^{0}}\bigg{|}_{\Delta_{z}=|\vec{P}_{\perp}|=0}.$
(7)

They reduce in the limit $P_{z}\to 0$ to the projection of the BF
distributions onto the transverse plane
$J^{\mu}_{\text{EF}}(\vec{b}_{\perp};0)=\int\mathrm{d}z\,J^{\mu}_{B}(\vec{r})$.
As expected, the total electric charge
$\mathcal{Q}=\int\mathrm{d}^{2}b_{\perp}\,J^{0}_{\text{EF}}(\vec{b}_{\perp};P_{z})\big{|}_{s^{\prime}=s}=\langle
p,s|j^{0}(0)|p,s\rangle/(2p^{0})$ is independent of $P_{z}$, which confirms
that the momentum-space amplitudes are correctly normalized.

The matrix elements of the electromagnetic current are traditionally
parametrized in terms of the so-called Dirac and Pauli FFs

$\langle
p^{\prime},s^{\prime}|j^{\mu}(0)|p,s\rangle=e\,\overline{u}(p^{\prime},s^{\prime})\left[\gamma^{\mu}\,F_{1}(Q^{2})+\frac{i\sigma^{\mu\nu}\Delta_{\nu}}{2M}\,F_{2}(Q^{2})\right]u(p,s),$
(8)

where $Q^{2}=-\Delta^{2}$. One of the appealing features of the BF is that the
spin structure of the matrix elements is the same as in the non-relativistic
theory [12, 9, 10]. One finds that

$\displaystyle\langle
p^{\prime}_{B},s^{\prime}_{B}|j^{0}(0)|p_{B},s_{B}\rangle$
$\displaystyle=e\,2M\,\delta_{s^{\prime}_{B}s_{B}}\,G_{E}(Q^{2}),$ (9)
$\displaystyle\langle
p^{\prime}_{B},s^{\prime}_{B}|\vec{j}(0)|p_{B},s_{B}\rangle$
$\displaystyle=e\,(\vec{\sigma}_{s^{\prime}_{B}s_{B}}\times
i\vec{\Delta})\,G_{M}(Q^{2}),$

where $\vec{\sigma}$ are the three Pauli matrices. The combination
$G_{E}=F_{1}-\tau F_{2}$ with $\tau=Q^{2}/(4M^{2})$ is called the electric
Sachs FF, and $G_{M}=F_{1}+F_{2}$ is called the magnetic Sachs FF.

The spin structure is more complicated in the EF [30, 21, 37]

$\displaystyle\langle
p^{\prime},s^{\prime}|j^{0}(0)|p,s\rangle\big{|}_{\text{EF}}$
$\displaystyle=e\,2M\,\gamma\,\Bigg{[}\left(\delta_{s^{\prime}s}\,\cos\theta+\frac{(\vec{\sigma}_{s^{\prime}s}\times
i\vec{\Delta}_{\perp})_{z}}{|\vec{\Delta}_{\perp}|}\,\sin\theta\right)G_{E}(Q^{2})$
(10)
$\displaystyle\qquad\qquad\quad+\beta\left(-\delta_{s^{\prime}s}\,\sin\theta+\frac{(\vec{\sigma}_{s^{\prime}s}\times
i\vec{\Delta}_{\perp})_{z}}{|\vec{\Delta}_{\perp}|}\,\cos\theta\right)\sqrt{\tau}\,G_{M}(Q^{2})\Bigg{]},$
$\displaystyle\langle
p^{\prime},s^{\prime}|j^{3}(0)|p,s\rangle\big{|}_{\text{EF}}$
$\displaystyle=e\,2M\,\gamma\,\Bigg{[}\beta\left(\delta_{s^{\prime}s}\,\cos\theta+\frac{(\vec{\sigma}_{s^{\prime}s}\times
i\vec{\Delta}_{\perp})_{z}}{|\vec{\Delta}_{\perp}|}\,\sin\theta\right)G_{E}(Q^{2})$
$\displaystyle\qquad\qquad\quad+\left(-\delta_{s^{\prime}s}\,\sin\theta+\frac{(\vec{\sigma}_{s^{\prime}s}\times
i\vec{\Delta}_{\perp})_{z}}{|\vec{\Delta}_{\perp}|}\,\cos\theta\right)\sqrt{\tau}\,G_{M}(Q^{2})\Bigg{]},$
$\displaystyle\langle
p^{\prime},s^{\prime}|\vec{j}_{\perp}(0)|p,s\rangle\big{|}_{\text{EF}}$
$\displaystyle=e\,(\sigma_{z})_{s^{\prime}s}\,(\vec{e}_{z}\times
i\vec{\Delta}_{\perp})\,G_{M}(Q^{2}),$

as a result of the Wigner rotation in (LABEL:EFLorentzTrans-Spinj). The
Lorentz boost parameters are given by $\gamma=P^{0}/P^{0}_{B}$ and
$\beta=P_{z}/P^{0}$, and the Wigner rotation angle $\theta$ satisfies

$\cos\theta=\frac{P^{0}+M(1+\tau)}{(P^{0}+M)\sqrt{1+\tau}},\qquad\sin\theta=-\frac{\sqrt{\tau}P_{z}}{(P^{0}+M)\sqrt{1+\tau}}.$
(11)

In the limit $P_{z}\to 0$ one recovers (9) with the restriction
$\Delta_{z}=0$, while in the IMF limit $P_{z}\to\infty$ one finds using the LF
components $a^{\pm}=(a^{0}\pm a^{3})/\sqrt{2}$

$\displaystyle\langle
p^{\prime},s^{\prime}|j^{+}(0)|p,s\rangle\big{|}_{\text{IMF}}$
$\displaystyle=e\,2P^{+}\,\bigg{[}\delta_{s^{\prime}s}\,F_{1}(Q^{2})+\frac{(\vec{\sigma}_{s^{\prime}s}\times
i\vec{\Delta})_{z}}{2M}\,F_{2}(Q^{2})\bigg{]},$ (12) $\displaystyle\langle
p^{\prime},s^{\prime}|j^{-}(0)|p,s\rangle\big{|}_{\text{IMF}}$
$\displaystyle=e\,2P^{-}\,\bigg{[}\delta_{s^{\prime}s}\,G_{1}(Q^{2})+\frac{(\vec{\sigma}_{s^{\prime}s}\times
i\vec{\Delta})_{z}}{2M}\,G_{2}(Q^{2})\bigg{]},$ $\displaystyle\langle
p^{\prime},s^{\prime}|\vec{j}_{\perp}(0)|p,s\rangle\big{|}_{\text{IMF}}$
$\displaystyle=e\,(\sigma_{z})_{s^{\prime}s}\,(\vec{e}_{z}\times
i\vec{\Delta}_{\perp})\,G_{M}(Q^{2}),$

where $G_{1}=(G_{E}-\tau G_{M})/(1+\tau)$ and $G_{2}=-(G_{E}+G_{M})/(1+\tau)$.
After Fourier transform, these coincide with the 2D LF distributions [15, 18,
19].

As a result of the combination of a four-vector transformation with a Wigner
spin rotation, the unpolarized charge distribution receives a magnetic
contribution when $\vec{P}\neq\vec{0}$. Since the induced magnetic
contribution in the neutron is large and opposite to the electric
contribution, a negatively charged region appears at the center of the LF
charge distribution of the neutron, see Fig. 1. This is purely an effect of
perspective which shows that, although probabilistic, the LF distributions
cannot be interpreted as intrinsic densities.

Figure 1: Spin-averaged EF charge distributions $\rho^{EF}_{E}\equiv
J^{0}_{\text{EF}}(\vec{b}_{\perp};P_{z})$ for selected values of the nucleon
average momentum. The small panels on the right correspond to the electric (or
convective) and magnetic (or magnetization) contributions. Figure adapted from
[30].

##  4 Energy-momentum tensor

The energy-momentum tensor (EMT) is a fundamental operator in quantum field
theory, for it measures key properties of physical systems like, e.g., mass,
spin and pressure. A direct access to its matrix elements requires
gravitational scattering, which is far too small to be measured in practice.
Some exclusive processes like deeply virtual Compton scattering are described
in terms of non-local correlators, known as generalized parton distributions,
which can be related to some of the gravitational FFs that parametrize the EMT
matrix elements [39]. Gravitational FFs can therefore in principle be
constrained from experimental data. Applying the same formalism as for the
matrix elements of the electromagnetic current, one can define EMT
distributions in the BF [11] and generalize them to both the EF and the LF
[29].

The nucleon matrix elements of the quark EMT can conveniently be parametrized
as

$\displaystyle\langle p^{\prime},s^{\prime}|T^{\mu\nu}_{q}(0)|p,s\rangle$
$\displaystyle=\overline{u}(p^{\prime},s^{\prime})\left[\frac{P^{\mu}P^{\nu}}{M}\,A_{q}(Q^{2})+\frac{\Delta^{\mu}\Delta^{\nu}-g^{\mu\nu}\Delta^{2}}{4M}\,D_{q}(Q^{2})+Mg^{\mu\nu}\bar{C}_{q}(Q^{2})\right.$
(13)
$\displaystyle\qquad\qquad\quad+\left.\frac{P^{\\{\mu}i\sigma^{\nu\\}\lambda}\Delta_{\lambda}}{M}\,J_{q}(Q^{2})-\frac{P^{[\mu}i\sigma^{\nu]\lambda}\Delta_{\lambda}}{M}\,S_{q}(Q^{2})\right]u(p,s),$

where the first four terms describe the symmetric part
$T^{\\{\mu\nu\\}}=(T^{\mu\nu}+T^{\nu\mu})/2$ while the fifth term accounts for
the antisymmetric part $T^{[\mu\nu]}=(T^{\mu\nu}-T^{\nu\mu})/2$. A similar
parametrization holds for the gluon EMT. Poincaré symmetry implies
conservation of linear and angular momentum, as well as mechanical
equilibrium, which amount to the following constraints on the gravitational
FFs

$A_{q}(0)+A_{g}(0)=1,\qquad
J_{q}(0)+J_{g}(0)=\tfrac{1}{2},\qquad\bar{C}_{q}(Q^{2})+\bar{C}_{g}(Q^{2})=0.$
(14)

In the BF, one finds that

$\displaystyle\langle
p^{\prime}_{B},s^{\prime}_{B}|T^{00}_{q}(0)|p_{B},s_{B}\rangle$
$\displaystyle=2MP^{0}_{B}\,\delta_{s^{\prime}_{B}s_{B}}\left\\{A_{q}(Q^{2})+\bar{C}_{q}(Q^{2})+\tau[D_{q}(Q^{2})-B_{q}(Q^{2})]\right\\},$
(15) $\displaystyle\langle
p^{\prime}_{B},s^{\prime}_{B}|T^{\\{0k\\}}_{q}(0)|p_{B},s_{B}\rangle$
$\displaystyle=2P^{0}_{B}\,(\vec{\sigma}_{s^{\prime}_{B}s_{B}}\times
i\vec{\Delta})^{k}\,J_{q}(Q^{2}),$ $\displaystyle\langle
p^{\prime}_{B},s^{\prime}_{B}|T^{[0k]}_{q}(0)|p_{B},s_{B}\rangle$
$\displaystyle=-2P^{0}_{B}\,(\vec{\sigma}_{s^{\prime}_{B}s_{B}}\times
i\vec{\Delta})^{k}\,S_{q}(Q^{2}),$ $\displaystyle\langle
p^{\prime}_{B},s^{\prime}_{B}|T^{ij}_{q}(0)|p_{B},s_{B}\rangle$
$\displaystyle=2MP^{0}_{B}\,\delta_{s^{\prime}_{B}s_{B}}\left\\{\frac{\Delta^{i}\Delta^{j}}{4M^{2}}\,D_{q}(Q^{2})-\delta^{ij}[\bar{C}_{q}(Q^{2})+\tau
D_{q}(Q^{2})]\right\\},$

where $B_{q}=2J_{q}-A_{q}$. Determining these matrix elements constitute one
of the main motivations of the physics program of the future Electron-Ion
Collider [5]. Here are some examples of the fundamental physical properties
they encode:

  * •

The amplitude $\langle
p^{\prime}_{B},s^{\prime}_{B}|T^{00}_{q}(0)|p_{B},s_{B}\rangle/(2P^{0}_{B})$
is related via Fourier transform to the spatial distribution of the quark
energy $T^{00}_{q,B}(\vec{r})$. In the forward limit $\Delta\to 0$, it reduces
to the quark contribution to the nucleon mass [40, 41].

  * •

The spatial distribution of quark Belinfante angular momentum is obtained from

$\mathcal{J}^{i}_{q,B}(\vec{r})=\epsilon^{ijk}r^{j}T^{\\{0k\\}}_{q,B}(\vec{r}),$
(16)

while the spatial distribution of quark kinetic angular momentum is given by
[42, 32]

$J^{i}_{q,B}(\vec{r})=L^{i}_{q,B}(\vec{r})+S^{i}_{q,B}(\vec{r})$ (17)

with the orbital and spin contributions defined as

$\displaystyle
L^{i}_{q,B}(\vec{r})=\epsilon^{ijk}r^{j}T^{0k}_{q,B}(\vec{r}),\qquad
S^{i}_{q,B}(\vec{r})=\frac{1}{2}\langle\overline{\psi}\gamma^{i}\gamma_{5}\psi\rangle_{\vec{0},\vec{0}}(\vec{r}).$
(18)

  * •

Because of spherical symmetry about the center of the nucleon, the quark
stress tensor is described by two functions [11, 43, 29]

$T^{ij}_{q,B}(\vec{r})=\delta^{ij}\,p_{q}(r)+\left(\frac{r^{i}r^{j}}{r^{2}}-\frac{1}{3}\,\delta^{ij}\right)s_{q}(r).$
(19)

Since the BF distributions are static, EMT conservation amounts to
$\nabla^{i}T^{ij}_{B}(\vec{r})=0$ and leads to a relation between the total
isotropic pressure $p(r)$ and pressure anisotropy $s(r)$

$\frac{\mathrm{d}}{\mathrm{d}r}\bigg{(}p+\frac{2}{3}\,s\bigg{)}+\frac{2}{r}\,s=0.$
(20)

Mechanical equilibrium of the nucleon implies the von Laue condition
$\int\mathrm{d}^{3}r\,p(r)=0$. Stability arguments further suggest that the
quantity $D_{q}(0)+D_{g}(0)=M\int\mathrm{d}^{3}r\,r^{2}\,p(r)$ should be
negative [44].

The first extraction of the pressure distribution inside a proton from
experimental data on deeply virtual Compton scattering has been reported in
[45], and followed by more conservative analyses [46, 47]. Data on
$J/\psi$-photoproduction in the threshold region have also been used to
extract gluon gravitational FFs [48]. Recent calculations of the gravitational
FFs from Lattice QCD can be found for example in [49, 50, 51]. More details
and references can be found in the review [8].

##  5 Conclusions

Electromagnetic and gravitational form factors encode in a Lorentz-invariant
key information about the 3D structure and physical properties of hadrons.
When relativistic recoil effects cannot be neglected, spatial distributions
become frame-dependent. The relativistic notion of spatial distribution is
then naturally formulated in a phase-space approach, where the strict
probabilistic interpretation is lost. From this perspective, it is possible to
interpolate between the Breit frame ($P_{z}=0$) and light-front
($P_{z}\to\infty$) pictures, and to show that the appearance of spatial
distributions in a moving frame is strongly impacted by the spin of the
target.

## References

  * [1] A. V. Belitsky, X. d. Ji and F. Yuan, _Quark imaging in the proton via quantum phase space distributions_ , _Phys. Rev. D_ 69 (2004), 074014 [arXiv:hep-ph/0307383 [hep-ph]].
  * [2] S. Meissner, A. Metz and M. Schlegel, _Generalized parton correlation functions for a spin-1/2 hadron_ , _JHEP_ 08 (2009), 056 [arXiv:0906.5323 [hep-ph]].
  * [3] C. Lorcé, B. Pasquini and M. Vanderhaeghen, _Unified framework for generalized and transverse-momentum dependent parton distributions within a 3Q light-cone picture of the nucleon_ , _JHEP_ 05 (2011), 041 [arXiv:1102.4704 [hep-ph]].
  * [4] C. Lorcé and B. Pasquini, _Quark Wigner Distributions and Orbital Angular Momentum_ , _Phys. Rev. D_ 84 (2011), 014015 [arXiv:1106.0139 [hep-ph]].
  * [5] R. Abdul Khalek, A. Accardi, J. Adam, D. Adamiak, W. Akers, M. Albaladejo, A. Al-bataineh, M. G. Alexeev, F. Ameli and P. Antonioli, et al. _Science Requirements and Detector Concepts for the Electron-Ion Collider: EIC Yellow Report_ , _Nucl. Phys. A_ 1026 (2022), 122447 [arXiv:2103.05419 [physics.ins-det]].
  * [6] V. Punjabi, C. F. Perdrisat, M. K. Jones, E. J. Brash and C. E. Carlson, _The Structure of the Nucleon: Elastic Electromagnetic Form Factors_ , _Eur. Phys. J. A_ 51 (2015), 79 [arXiv:1503.01452 [nucl-ex]].
  * [7] H. Gao and M. Vanderhaeghen, _The proton charge radius_ , _Rev. Mod. Phys._ 94 (2022) no.1, 015002 [arXiv:2105.00571 [hep-ph]].
  * [8] V. D. Burkert, L. Elouadrhiri, F. X. Girod, C. Lorcé, P. Schweitzer and P. E. Shanahan, _Colloquium: Gravitational Form Factors of the Proton_ , _Rev. Mod. Phys._ 95 (2023) no.4, 041002 [arXiv:2303.08347 [hep-ph]].
  * [9] F. J. Ernst, R. G. Sachs and K. C. Wali, _Electromagnetic form factors of the nucleon_ , _Phys. Rev._ 119 (1960), 1105.
  * [10] R. G. Sachs, _High-Energy Behavior of Nucleon Electromagnetic Form Factors_ , _Phys. Rev._ 126 (1962), 2256.
  * [11] M. V. Polyakov, _Generalized parton distributions and strong forces inside nucleons and nuclei_ , _Phys. Lett. B_ 555 (2003), 57 [arXiv:hep-ph/0210165].
  * [12] D. R. Yennie, M. M. Lévy and D. G. Ravenhall, _Electromagnetic Structure of Nucleons_ , _Rev. Mod. Phys._ 29 (1957) no.1, 144.
  * [13] M. Burkardt, _Impact parameter dependent parton distributions and off forward parton distributions for $\zeta\to 0$_, _Phys. Rev. D_ 62 (2000), 071503 [arXiv:hep-ph/0005108].
  * [14] R. L. Jaffe, _Ambiguities in the definition of local spatial densities in light hadrons_ , _Phys. Rev. D_ 103 (2021) no.1, 016017 [arXiv:2010.15887 [hep-ph]].
  * [15] M. Burkardt, _Impact parameter space interpretation for generalized parton distributions_ , _Int. J. Mod. Phys. A_ 18 (2003), 173 [arXiv:hep-ph/0207047].
  * [16] G. A. Miller, _Transverse Charge Densities_ , _Ann. Rev. Nucl. Part. Sci._ 60 (2010), 1 [arXiv:1002.0355 [nucl-th]].
  * [17] G. A. Miller, _Defining the proton radius: A unified treatment_ , _Phys. Rev. C_ 99 (2019) no.3, 035202 [arXiv:1812.02714 [nucl-th]].
  * [18] G. A. Miller, _Charge Density of the Neutron_ , _Phys. Rev. Lett._ 99 (2007), 112001 [arXiv:0705.2409 [nucl-th]].
  * [19] C. E. Carlson and M. Vanderhaeghen, _Empirical transverse charge densities in the nucleon and the nucleon-to-Delta transition_ , _Phys. Rev. Lett._ 100 (2008), 032004 [arXiv:0710.0835 [hep-ph]].
  * [20] J. Y. Panteleeva and M. V. Polyakov, _Forces inside the nucleon on the light front from 3D Breit frame force distributions: Abel tomography case_ , _Phys. Rev. D_ 104 (2021) no.1, 014008 [arXiv:2102.10902 [hep-ph]].
  * [21] J. Y. Kim and H. C. Kim, _Transverse charge distributions of the nucleon and their Abel images_ , _Phys. Rev. D_ 104 (2021) no.7, 074003 [arXiv:2106.10986 [hep-ph]].
  * [22] A. Freese and G. A. Miller, _Unified formalism for electromagnetic and gravitational probes: Densities_ , _Phys. Rev. D_ 105 (2022) no.1, 014003 [arXiv:2108.03301 [hep-ph]].
  * [23] E. Epelbaum, J. Gegelia, N. Lange, U. G. Meißner and M. V. Polyakov, _Definition of Local Spatial Densities in Hadrons_ , _Phys. Rev. Lett._ 129 (2022) no.1, 012001 [arXiv:2201.02565 [hep-ph]].
  * [24] J. Y. Panteleeva, E. Epelbaum, J. Gegelia and U. G. Meißner, _Definition of electromagnetic local spatial densities for composite spin-1/2 systems_ , _Phys. Rev. D_ 106 (2022) no.5, 056019 [arXiv:2205.15061 [hep-ph]].
  * [25] Y. Li, W. b. Dong, Y. l. Yin, Q. Wang and J. P. Vary, _Minkowski’s lost legacy and hadron electromagnetism_ , _Phys. Lett. B_ 838 (2023), 137676 [arXiv:2206.12903 [hep-ph]].
  * [26] A. Freese and G. A. Miller, _Convolution formalism for defining densities of hadrons_ , _Phys. Rev. D_ 108 (2023) no.3, 034008 [arXiv:2210.03807 [hep-ph]].
  * [27] A. Freese and G. A. Miller, _Synchronization effects on rest frame energy and momentum densities in the proton_ , _Phys. Rev. D_ 108 (2023) no.9, 094026 [arXiv:2307.11165 [hep-ph]].
  * [28] C. Lorcé, _The relativistic center of mass in field theory with spin_ , _Eur. Phys. J. C_ 78 (2018) no.9, 785 [arXiv:1805.05284 [hep-ph]].
  * [29] C. Lorcé, H. Moutarde and A. P. Trawiński, _Revisiting the mechanical properties of the nucleon_ , _Eur. Phys. J. C_ 79 (2019) no.1, 89 [arXiv:1810.09837 [hep-ph]].
  * [30] C. Lorcé, _Charge Distributions of Moving Nucleons_ , _Phys. Rev. Lett._ 125 (2020) no.23, 232002 [arXiv:2007.05318 [hep-ph]].
  * [31] J. B. Kogut and D. E. Soper, _Quantum Electrodynamics in the Infinite Momentum Frame_ , _Phys. Rev. D_ 1 (1970), 2901.
  * [32] C. Lorcé, L. Mantovani and B. Pasquini, _Spatial distribution of angular momentum inside the nucleon_ , _Phys. Lett. B_ 776 (2018), 38 [arXiv:1704.08557 [hep-ph]].
  * [33] E. P. Wigner, _On the quantum correction for thermodynamic equilibrium_ , _Phys. Rev._ 40 (1932), 749.
  * [34] M. Jacob and G. C. Wick, _On the General Theory of Collisions for Particles with Spin_ , _Annals Phys._ 7 (1959), 404.
  * [35] L. Durand, P. C. DeCelles and R. B. Marr, _Lorentz Invariance and the Kinematic Structure of Vertex Functions_ , _Phys. Rev._ 126 (1962), 1882.
  * [36] C. Lorcé and P. Wang, _Deuteron relativistic charge distributions_ , _Phys. Rev. D_ 105 (2022) no.9, 096032 [arXiv:2204.01465 [hep-ph]].
  * [37] Y. Chen and C. Lorcé, _Pion and nucleon relativistic electromagnetic four-current distributions_ , _Phys. Rev. D_ 106 (2022) no.11, 116024 [arXiv:2210.02908 [hep-ph]].
  * [38] Y. Chen and C. Lorcé, _Nucleon relativistic polarization and magnetization distributions_ , _Phys. Rev. D_ 107 (2023) no.9, 096003 [arXiv:2302.04672 [hep-ph]].
  * [39] X. D. Ji, _Gauge-Invariant Decomposition of Nucleon Spin_ , _Phys. Rev. Lett._ 78 (1997), 610-613 [arXiv:hep-ph/9603249 [hep-ph]].
  * [40] C. Lorcé, _On the hadron mass decomposition_ , _Eur. Phys. J. C_ 78 (2018) no.2, 120 [arXiv:1706.05853 [hep-ph]].
  * [41] C. Lorcé, A. Metz, B. Pasquini and S. Rodini, _Energy-momentum tensor in QCD: nucleon mass decomposition and mechanical equilibrium_ , _JHEP_ 11 (2021), 121 [arXiv:2109.11785 [hep-ph]].
  * [42] E. Leader and C. Lorcé, _The angular momentum controversy: What’s it all about and does it matter?_ , _Phys. Rept._ 541 (2014) no.3, 163 [arXiv:1309.4235 [hep-ph]].
  * [43] M. V. Polyakov and P. Schweitzer, _Forces inside hadrons: pressure, surface tension, mechanical radius, and all that_ , _Int. J. Mod. Phys. A_ 33 (2018) no.26, 1830025 [arXiv:1805.06596 [hep-ph]].
  * [44] I. A. Perevalova, M. V. Polyakov and P. Schweitzer, _On LHCb pentaquarks as a baryon- $\psi$(2S) bound state: prediction of isospin-$\frac{3}{2}$ pentaquarks with hidden charm_, _Phys. Rev. D_ 94 (2016) no.5, 054024 [arXiv:1607.07008 [hep-ph]].
  * [45] V. D. Burkert, L. Elouadrhiri and F. X. Girod, _The pressure distribution inside the proton_ , _Nature_ 557 (2018) no.7705, 396.
  * [46] K. Kumerički, _Measurability of pressure inside the proton_ , _Nature_ 557 (2018) no.7759, E1.
  * [47] H. Dutrieux, C. Lorcé, H. Moutarde, P. Sznajder, A. Trawiński and J. Wagner, _Phenomenological assessment of proton mechanical properties from deeply virtual Compton scattering_ , _Eur. Phys. J. C_ 81 (2021) no.4, 300 [arXiv:2101.03855 [hep-ph]].
  * [48] B. Duran, Z. E. Meziani, S. Joosten, M. K. Jones, S. Prasad, C. Peng, W. Armstrong, H. Atac, E. Chudakov and H. Bhatt, et al. _Determining the gluonic gravitational form factors of the proton_ , _Nature_ 615 (2023) no.7954, 813 [arXiv:2207.05212 [nucl-ex]].
  * [49] C. Alexandrou, S. Bacchio, M. Constantinou, J. Finkenrath, K. Hadjiyiannakou, K. Jansen, G. Koutsou, H. Panagopoulos and G. Spanoudes, _Complete flavor decomposition of the spin and momentum fraction of the proton using lattice QCD simulations at physical pion mass_ , _Phys. Rev. D_ 101 (2020) no.9, 094513 [arXiv:2003.08486 [hep-lat]].
  * [50] D. A. Pefkou, D. C. Hackett and P. E. Shanahan, _Gluon gravitational structure of hadrons of different spin_ , _Phys. Rev. D_ 105 (2022) no.5, 054509 [arXiv:2107.10368 [hep-lat]].
  * [51] G. Wang et al. [$\chi$QCD], _Proton momentum and angular momentum decompositions with overlap fermions_ , _Phys. Rev. D_ 106 (2022) no.1, 014512 [arXiv:2111.09329 [hep-lat]].


# Centimeter-sized Grains in the Compact Dust Ring around Very Low Mass Star
CIDA 1

Jun Hashimoto Astrobiology Center, National Institutes of Natural Sciences,
2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan  Subaru Telescope, National
Astronomical Observatory of Japan, Mitaka, Tokyo 181-8588, Japan  Department
of Astronomy, School of Science, Graduate University for Advanced Studies
(SOKENDAI), Mitaka, Tokyo 181-8588, Japan  Hauyu Baobab Liu Department of
Physics, National Sun Yat-Sen University, No. 70, Lien-Hai Road, Kaohsiung
City 80424, Taiwan, R.O.C.  Ruobing Dong Department of Physics & Astronomy,
University of Victoria, Victoria, BC, V8P 5C2, Canada  Beibei Liu Institute
for Astronomy, School of Physics, Zhejiang University, 38 Zheda Road, Hangzhou
310027, China  Takayuki Muto  Division of Liberal Arts, Kogakuin University,
1-24-2, Nishi-Shinjuku, Shinjuku-ku, Tokyo 163-8677, Japan  Leiden
Observatory, Leiden University, P.O. Box 9513, NL-2300 RA Leiden, The
Netherlands  Department of Earth and Planetary Sciences, Tokyo Institute of
Technology, 2-12-1 Oh-okayama, Meguro-ku, Tokyo 152-8551, Japan  Yuka Terada
Institute of Astronomy and Astrophysics, Academia Sinica, 11F of Astronomy-
Mathematics Building, AS/NTU No.1, Sec. 4, Roosevelt Rd, Taipei 10617, Taiwan,
ROC  Department of Astrophysics, National Taiwan University, Taipei 10617,
Taiwan, R.O.C.

###### Abstract

We examined the grain size in the dust ring encircling the 0.19 $M_{\sun}$ T
Tauri star CIDA 1 using the Karl G. Jansky Very Large Array (JVLA) at multiple
centimeter wavelengths, with a spatial resolution of 0$\farcs$2–0$\farcs$9\.
We detected distinct partial-ring structures at these wavelengths around CIDA
1. Based on spatial distributions and spectral indexes, we determined that
these centimeter emissions originated from dust, rather than free-free or
synchrotron emissions. To estimate the maximum grain size ($a_{\rm max}$)
within the ring, we compared the observed spectral energy distribution (SED)
with SEDs calculated for different $a_{\rm max}$ values using radiative
transfer calculations. Our findings indicate an $a_{\rm max}$ value of
approximately 2.5 cm in the ring, assuming the dust opacity can be
approximated by the DSHARP models. These results suggest that grain growth
took place within the CIDA 1 ring, potentially facilitating more efficient
planet formation through pebble accretion scenarios involving centimeter-sized
pebbles.

protoplanetary disks — planet–disk interactions — Circumstellar dust — Dust
continuum emission — stars: individual (CIDA 1)

††software:  astropy (Astropy Collaboration et al., 2013), Numpy (van der Walt
et al., 2011), CASA (v6.4.3; McMullin et al., 2007), ††facilities: JVLA, ALMA

##  1 Introduction

The initial stage of planet formation in protoplanetary disks is grain growth.
Understanding this process and determining the maximum grain size are
essential because subsequent planet formation mechanisms such as streaming
instability (e.g., Johansen et al., 2014) and pebble accretion (e.g., Liu &
Ji, 2020) depend on them (Carrera et al., 2015; Morbidelli et al., 2015; Ormel
& Liu, 2018; Drazkowska et al., 2022).

Two methods exist for determining the sizes of dust particles based on
millimeter/centimeter (mm/cm) continuum emission. The first approach involves
observing polarized dust continuum emission in the mm/cm range, resulting from
scattering (Kataoka et al., 2015). The polarization fraction in the dust
emission reaches its maximum when the dust grains grow to a size of
approximately $a_{\rm max}\sim\lambda/2\pi$, where $\lambda$ represents the
observing wavelength. This method requires a signal-to-noise ratio (SNR) of
approximately 100–1,000 in Stokes $I$ and is well-suited for studying bright
objects like DSHARP disks (Andrews et al., 2018). The second approach involves
constructing the spectral energy distribution (SED) (e.g., Testi et al., 2003;
Pérez et al., 2012; Testi et al., 2014; Pérez et al., 2015; Tazzari et al.,
2016; Liu, 2019; Ueda et al., 2020, 2021). This method relies on the fact that
grains emit thermal radiation most efficiently at wavelengths similar to their
sizes (e.g., Draine, 2006). By analyzing the SED, information about the grain
sizes can be inferred.

Polarization observations in the (sub-)millimeter range have provided insights
into the maximum grain size ($a_{\rm max}$) in protoplanetary disks, estimated
to be around 100 $\mu$m (e.g., Kataoka et al., 2016b, a; Ohashi et al., 2018,
2020; Bacciotti et al., 2018; Dent et al., 2019). The efficiency of pebble
accretion, crucial for planet growth, exhibits a size dependency (Liu & Ormel,
2018; Ormel & Liu, 2018). When accreting pebbles of 100 $\mu$m, planet growth
generally occurs at a slower pace compared to mm to cm-sized pebbles (e.g.,
Morbidelli et al., 2015). Consequently, the formation of giant planets with
sub-mm dust grains may require longer timescales, potentially exceeding the
typical disk lifetime of a few million years (Mamajek, 2009). Thus the grain-
size problem remains a prominent topic in the field of planet formation.
However, it is important to note that ALMA polarization observations have
exhibited a bias towards bright disks, leaving uncertainties regarding the
ubiquity of dust grains with sub-mm $a_{\rm max}$ in disks.

Some objects have indeed exhibited the presence of millimeter-/centimeter-
sized grains within their disks (e.g., Carrasco-González et al., 2019; Ueda et
al., 2021; Macías et al., 2021; Hashimoto et al., 2022; Zhang et al., 2023).
In this study, we present another example: CIDA 1 (2MASS J04141760$+$2806096),
which contains dust grains of millimeter to centimeter sizes within its disk.
CIDA 1, located at $d=134.6$ pc in the Taurus star-forming region, has a mass
of 0.19 $M_{\sun}$ and an effective temperature of 3,197 K (Kurtovic et al.,
2021; Gaia Collaboration et al., 2016, 2022). It has been extensively observed
with ALMA in bands 3, 4, 6, and 7 (Ricci et al., 2014; Simon et al., 2017;
Pinilla et al., 2018, 2021; Kurtovic et al., 2021), revealing a ring/gap
structure at a radius of around 20 au in the dust continuum image at bands 4
and 7 (Pinilla et al., 2021). Recent hydrodynamic simulations indicate that a
giant planet with a minimum mass of 1.4 $M_{\rm Jup}$ could account for the
observed ring/gap structure around CIDA 1 (Curone et al., 2022). The spectral
index $\alpha_{\rm band4\mathchar 45\relax 7}$ between ALMA bands 4 and 7 is
measured to be 2.0 $\pm$ 0.2, and the continuum ring at bands 4 and 7 is
optically thin (Pinilla et al., 2021), suggesting that the dust grains within
the CIDA 1 ring have grown to sizes on the order of millimeters to centimeters
(e.g., Draine, 2006).

##  2 Observations and results

The JVLA observations  with B configuration were conducted as part of program
ID 23A-124 (PI: J. Hashimoto), and the details are provided in Table 1. The
data were calibrated using the Common Astronomy Software Applications (CASA)
package (McMullin et al., 2007), following the calibration scripts provided by
JVLA. Due to weak emission from the object, self-calibration of the
visibilities was not performed. The stellar position was corrected by proper
motion (8.285, -23.607) mas/yr (Gaia Collaboration et al., 2016, 2022). The
corrected ICRS coordinate for CIDA 1 is (4h 14m 17.62395s, 28d 6m 9.11174s).
Subsequently,  multi-term multi-frequency synthesis imaging with nterm = 2
(Rau & Cornwell, 2011) was performed using the CASA tclean task, as summarized
in Table 2.  We have attempted phase-only self-calibration on the Ka-band
data, which has the highest signal-to-noise ratio (SNR); see below. However
approximately 62% of solutions on the solution intervals of 100s, 300s, 1,000s
were flagged with a minimum SNR of 2.5. Thus we have decided not to apply
self-calibration.

Table 1: JVLA observations  with B configuration Observations |
---|---
Observing date (UT) | 2023.Jan.13, 14, 15, 20,
| Feb.3, Apr.16
Project code | 23A-124 (PI: J. Hashimoto)
Central frequency (GHz) | 44 (Q), 33 (Ka),
| 15 (Ku), 10 (X)
Continuum band width (GHz) | 8 (Q, Ka), 6 (Ku), 4 (X)
Time on source (min) | 48.5 (Q), 156.0 (Ka),
| 20.7 (Ku), 6.2 (X)
Number of antennas | 27 (Jan.13, 14, 20, Feb.3),
| 26 (Jan.15, Apr.16)
Baseline lengths (km) | 0.243 to 11.1
Bandpass calibrator | J0319$+$4130
Flux calibrator | 3C147
Phase calibrator | J0403$+$2600
Table 2: Imaging parameters | X band (10 GHz) | Ku band (15 GHz) | Ka band (33 GHz) | Q band (44 GHz)
---|---|---|---|---
For imaging in Figure 1 without splitting observing bands
Robust clean parameter | 2.0 | 2.0 | 0.0 | 2.0
Beam shape | 0$\farcs$89$\times$0$\farcs$74 | 0$\farcs$58$\times$0$\farcs$51 | 0$\farcs$17$\times$0$\farcs$16 | 0$\farcs$20$\times$0$\farcs$18
r.m.s. noise ($\mu$Jy/beam) | 11.7 | 4.5 | 5.4 | 18.3
For measuring flux in Table 3 with  or without splitting observing bands
Split frequency (GHz) | 8–12 | 12–15 | 15–18 | 12–18 | 29–33 | 33–37 | 29–37 | 40–48
| (no split) | (split) | (split) | (no split) | (split) | (split) | (no split) | (no split)
Robust clean parameter | 2.0 | 2.0 | 2.0 | 2.0 | 2.0 | 2.0 | 2.0 | 2.0
Beam shape | 0$\farcs$89$\times$0$\farcs$74 | 0$\farcs$63$\times$0$\farcs$55 | 0$\farcs$54$\times$0$\farcs$46 | 0$\farcs$58$\times$0$\farcs$51 | 0$\farcs$29$\times$0$\farcs$25 | 0$\farcs$26$\times$0$\farcs$22 | 0$\farcs$27$\times$0$\farcs$23 | 0$\farcs$20$\times$0$\farcs$18
r.m.s. noise ($\mu$Jy/beam) | 11.7 | 6.1 | 6.7 | 4.5 | 5.4 | 5.2 | 3.7 | 18.3

Figure 1 showcases the JVLA images obtained at Q band (40–48 GHz), Ka band
(29–37 GHz), Ku band (12–18 GHz), and X band (8–12 GHz). To achieve a
comparable spatial resolution of 0$\farcs$18 to the Q-band image, we utilized
Briggs Robust $=$ 0.0 for the Ka-band image. Our observations revealed
significant signals at Q, Ka, and Ku bands, with peak flux densities of  87.6
(4.8 $\sigma$), 36.9 (6.8 $\sigma$), and 19.0 (4.2 $\sigma$) $\mu$Jy/beam,
respectively. The integrated flux measurements, exceeding 3 $\sigma$, are
presented in Table 3. To improve the signal-to-noise ratio (SNR), Briggs
Robust $=$ 2.0 was applied for the Ka-band flux (refer to Table 2). As the
beam size of approximately 0$\farcs$4 is similar to the size of the CIDA 1
ring, the peak flux density from the Ku-band image was utilized as the total
flux.

Figure 1:  The JVLA images (color) show the CIDA 1 ring overlaid with the ALMA
band-7 data (dashed ellipses; $r=0\farcs 16$; $i=37.4$ deg; P.A. $=10.8$ deg;
Pinilla et al. 2021). Panels (a) to (d) depict the Q-band (40–48 GHz), Ka-band
(29–37 GHz), Ku-band (12–18 GHz), and X-band (8-12 GHz) images, respectively.
The stellar position is marked by a black star, corresponding to the ICRS
coordinate of (4h 14m 17.62395s, 28d 6m 9.11174s). The proper motion (pmRA
$=8.285$ mas/yr, pmDEC $=-23.607$ mas/yr; Gaia Collaboration et al., 2016,
2022) has been corrected. The synthesized beams are illustrated as white
ellipses in the lower left corner.

The signals at Q- and Ka-bands are clearly spatially resolved. The locations
of strong signals, depicted in white in Figure 1, well trace the ring
structure detected by ALMA sub-millimeter emission (Pinilla et al., 2021) as
indicated by the dashed black ellipse. Although the ring appears asymmetric at
Q- and Ka-bands, with a brighter region in the northwest, the contrast between
the peak position and the opposing side within the ring are measured to
$1.40\pm 0.43$ and $1.61\pm 0.40$ at Q- and Ka-bands, respectively. The
observed asymmetry at the Ka-band might be genuine. However, due to the low
SNR in the contrast, we refrain from further discussing the asymmetry. It is
imperative to conduct observations with higher S/N ratios to investigate
potential asymmetries within the CIDA 1 ring.

The previously detected central point source in ALMA band-4 and 7 images
(Pinilla et al., 2021) is not robustly recovered in our JVLA observations
(Figure 1), possibly due to the lower spatial resolution. To investigate the
central source in Q- and Ka-band data, we first checked their CLEAN
components. Figures 2(a) and (c) show the CLEAN components of Q- and Ka-band
data. There is a CLEAN component with 4.2 $\mu$Jy/pixel in the central region
of the Q-band model. We then performed tclean again, with a mask on the ring
region indicated by the dashed ellipse in Figure 1. Figures 2(b) and (d) show
the resulting Q- and Ka-band residual images. We find that the peak flux of
the central source in both residual images is less than 1 $\sigma$ (Table 2).
High resolution and high sensitivity observations are needed for further
investigations of the central point source.

In order to perform SED modeling of the ring in § 3, it is necessary to
estimate the flux contribution of the central source to the integrated flux.
The visibility fitting analysis of Pinilla et al. (2021) using ALMA band-4 and
7 data indicated that the size (FWHM) of the central source is comparable to
or smaller than the beam size of approximately 30-50 mas, suggesting that the
central source is not well spatially resolved in the ALMA images. Therefore,
we used the peak flux density of the central source, which is summarized in
Table 3.

Figure 2:  The inner disk component is explored in Q- and Ka-band data. Left: the CLEAN components. There is a CLEAN component with 4.2 $\mu$Jy/pixel in the central region of Q-band data in panel (a). Right: residual images with a specific mask on the ring region indicated by the dashed ellipse in Figure 1. The peak flux of the central region in both residual images is less than 1 $\sigma$ (Table 2).  |
---|---
Figure 3:  The observed flux densities from CIDA 1, represented by dots, are
compared to our models for interpretation. The flux densities used are listed
above the separator line in Table 3.  The flux uncertainties are explained in
the table caption in Table 3. Left: The integrated flux densities from CIDA 1
are plotted on a logarithmic scale. Right: The flux density of the CIDA 1
inner disk is shown. In all panels, the blue, green, yellow, and pink lines
correspond to our models for the dust emission from the inner disk, an
embedded dust ring with $a_{\mbox{\scriptsize max}}=25$ mm, a spatially
compact embedded dust component with $a_{\mbox{\scriptsize max}}=0.9$ mm, and
an extended dust ring with $a_{\mbox{\scriptsize max}}=0.55$ mm, respectively
(refer to Table 4 for a summary of our model). The black lines represent the
total flux densities of all model components.  Table 3:  Flux densities in
CIDA 1aaFlux densities listed above the separator line are used in our SED
modeling in § 3.2.

Frequency | Integrated | Inner disk | Refsjj(1) This study. (2) Ricci et al. (2014). (3) Pinilla et al. (2021). (4) Andrews et al. (2013).
---|---|---|---
(GHz) | (mJy) | (mJy) |
8–12 |  $<$0.0351b,cb,cfootnotemark:  |  $\cdots$ | (1)
12–15ddThe Ka and Ku-bands are split into higher and lower bands. | $0.0128\pm 0.0061$ee As the area where the flux is detected at more than 3 $\sigma$ is smaller than the beam size, we regarded the peak flux density as the integrated flux density. The photometric error is the RMS noise in Table 2. |  $\cdots$ | (1)
15–18ddThe Ka and Ku-bands are split into higher and lower bands. | $0.0259\pm 0.0067$ee As the area where the flux is detected at more than 3 $\sigma$ is smaller than the beam size, we regarded the peak flux density as the integrated flux density. The photometric error is the RMS noise in Table 2. |  $\cdots$ | (1)
29–33ddThe Ka and Ku-bands are split into higher and lower bands. | $0.1007\pm 0.0095$ff The integrated flux is measured in the area where the flux is detected at more than 3 $\sigma$. The error is calculated as RMS $\times\sqrt{N}$, where RMS is taken from Table 2 and $N$ is the number of beams. As the absolute accuracy of flux in JVLA is 2% (https://science.nrao.edu/facilities/vla/docs/manuals/oss2017A/performance/fdscale), we only take account of random noise. |  $\cdots$ | (1)
33–37ddThe Ka and Ku-bands are split into higher and lower bands. | $0.1390\pm 0.0097$ff The integrated flux is measured in the area where the flux is detected at more than 3 $\sigma$. The error is calculated as RMS $\times\sqrt{N}$, where RMS is taken from Table 2 and $N$ is the number of beams. As the absolute accuracy of flux in JVLA is 2% (https://science.nrao.edu/facilities/vla/docs/manuals/oss2017A/performance/fdscale), we only take account of random noise. |  $\cdots$ | (1)
40–48 | $0.2050\pm 0.0350$b,fb,ffootnotemark:  |  $\cdots$ | (1)
93 |  $2.1\pm 0.21$ggWe assume nominal 10% (1-$\sigma$) absolute flux uncertainty. |  $\cdots$ | (2)
141 |  $5.1\pm 0.51$ggWe assume nominal 10% (1-$\sigma$) absolute flux uncertainty. |  $0.087\pm 0.0087$g,hg,hfootnotemark:  | (3)
225.5 |  $13.5\pm 2.8$ii The error is taken from Andrews et al. (2013). |  $\cdots$ | (4)
339 |  $31\pm 3.1$ggWe assume nominal 10% (1-$\sigma$) absolute flux uncertainty. |  $0.60\pm 0.06$g,hg,hfootnotemark:  | (3)
12–18 | $0.0190\pm 0.0045$b,eb,efootnotemark:  |  $\cdots$ | (1)
29–37 | $0.1303\pm 0.0074$b,fb,ffootnotemark:  |  $\cdots$ | (1)

hhfootnotetext: The peak flux density of the central point source.

##  3 Constraints on Maximum Grain Sizes

###  3.1 Qualitative interpretation for the SED

The flux densities measured at all observed frequencies (Table 3) are
primarily attributed to dust emission rather than free-free or synchrotron
emissions. If there were strong free-free or synchrotron emission, the
intensity at frequencies below 50 GHz would peak at the location of the star
itself (e.g., ZZ Tau IRS; Hashimoto et al., 2022), which is not detected in
Figure 1. We will discuss later the possibility of spatially extended free-
free or synchrotron emission by examining the observed spectral indices.

Weak free-free emission typically exhibits a spectral index close to 0, while
synchrotron emission near young stellar objects (YSOs) often shows negative
spectral indices (e.g., Anglada et al., 1998). In contrast, the observed
spectral indices at 15–33 GHz and 33–44 GHz are  $2.44\pm 0.31$ and  $1.58\pm
0.63$, respectively (calculated using the flux densities from the non-split
JVLA bands listed in Table 3). If there were significant free-free and/or
synchrotron emission at frequencies below 50 GHz, we would expect lower
spectral indices at lower frequencies. It is also worth noting that free-free
and synchrotron emissions are generally very weak in typical Class II YSOs
(e.g., Liu et al., 2014; Galván-Madrid et al., 2014), so our current findings
are not surprising.

The ALMA observations have successfully resolved two sources of dust emission:
the inner disk and the dusty ring (refer to Figure 4 in the Appendix). The
spectral index of the inner disk at 141–339 GHz is $2.06\pm 0.16$ (Table 3),
indicating optically thick dust emission. On the other hand, the dominant flux
densities at (sub)millimeter and centimeter wavelengths can be attributed to
the dusty ring, which can be further divided into multiple dust emission
components through analysis of the spectral energy distribution (SED). Figure
3 presents the SED of CIDA 1, which exhibits a change in slope around 40–50
GHz, separating it into the (sub)millimeter region (with a spectral index of
$3.11\pm 0.26$ at 44–93 GHz) and the centimeter region (with a spectral index
of  $1.58\pm 0.63$ at 33–44 GHz). In the (sub)millimeter region, the spectral
index remains relatively constant at approximately 2.0 (as shown in Table 3;
Ricci et al., 2014; Pinilla et al., 2018), whereas it varies with frequency in
the centimeter region. These features of the SED require an interpretation
that considers dust emission components with different column densities
($\Sigma_{\rm dust}$), temperatures, and maximum grain sizes ($a_{\rm max}$).

The spectral index of  $3.51\pm 2.70$ at 13.5–16.5 GHz, derived from splitting
the JVLA Ku-band flux densities (Table 3), provides a qualitative
characterization of the upper limit for the maximum grain size ($a_{\rm max}$)
of the dominant dust emission source contributing to the detected centimeter
flux densities (further discussed in § 3.2.2). However, it is still possible
that there are faint dust emission sources with larger $a_{\rm max}$ grains.
Through more detailed quantitative modeling, we have determined that the
(sub)millimeter region of the spectral energy distribution (SED) is dominated
by optically thick dust emission with $a_{\rm max}<1$ mm, while the centimeter
region is dominated by dust emission with $a_{\rm max}\gtrsim$1–25 mm (as
described in § 3.2). It should be noted that although Pinilla et al. (2021)
associated the low spectral index of $\sim$2 in the (sub)millimeter region
with optically thin emissions from grown dust grains, larger dust grains with
$a_{\rm max}>1$ mm are unlikely to explain the steeper spectral index of
$3.11\pm 0.26$ at 44–93 GHz presented in Table 3 (refer to, for example,
Figure 3 in Draine, 2006).

Furthermore, the emission from the dust component with $a_{\rm
max}\gtrsim$1–25 mm is suppressed at frequencies above 50 GHz due to
obscuration caused by the dust component with $a_{\rm max}<$1 mm. This
indicates a vertical segregation of dust grain sizes within the system. A
similar vertical segregation of dust sizes was also observed in the previous
analysis of the SED for the inner $\sim$10 au region of FU Ori (Liu et al.,
2021). Our current SED model for CIDA 1, incorporating multiple dust
components, shows qualitative similarities to another recently studied very
low mass (VLM) star, ZZ Tau IRS (Hashimoto et al., 2022). However, in the case
of ZZ Tau IRS, the SEDs can be explained without assuming mutual obscuration
between dust emission components.

###  3.2 SED modeling

We conducted SED modeling using the flux densities listed above the separator
line in Table 3. As explained in § 3.1, attempting to fit all observed flux
densities with a single dust emission slab model would result in a poor fit
and misleading parameters. Instead, we followed the approach of Liu et al.
(2019b, 2021) and employed a multi-layered dust slab model to fit the
integrated (sub)millimeter and centimeter flux densities of CIDA 1 by

$F_{\nu}=\sum\limits_{i}F_{\nu}^{i}e^{-\sum\limits_{j}\tau^{i,j}_{\nu}},$ (1)

where $F_{\nu}^{i}$ is the flux density of the single dust component $i$ and
$\tau^{i,j}_{\nu}$ is the optical depth of the emission component $j$ that
obscures the emission component $i$. Since our observations show no evidence
of free-free and synchrotron emission (§ 3.1), we did not include them in our
SED models. We assumed that the dusty disk of CIDA 1 is symmetric with respect
to the disk midplane, meaning that component $i$ is closer to the midplane
than component $j$. The dust emission component embedded in the midplane is
only obscured by the front side of the other dust emission components. Hence,
we assumed that the $\tau^{i,j}_{\nu}$ is half of the total extinction optical
depth of dust emission component $j$. In our models, we considered mutual
obscuration only when it was necessary to explain the observations. In cases
where the data could be interpreted both with and without assuming mutual
obscuration, we assumed there was no obscuration. We will discuss both
scenarios qualitatively when relevant.

We determined the flux density and total extinction optical depth of each dust
component by utilizing the analytic radiative transfer solutions presented in
Birnstiel et al. (2018) (Equations 10–20). These solutions are based on the
Eddington approximation. The free parameters for each dust component include
the solid angle ($\Omega_{\rm dust}$), dust temperature ($T_{\rm dust}$), dust
mass surface density ($\Sigma_{\rm dust}$), and maximum grain size ($a_{\rm
max}$). To evaluate the size-averaged, frequency-dependent dust absorption and
scattering opacity, we employed the default (i.e., water-ice-coated) DSHARP
opacity table from Birnstiel et al. (2018), which assumes compact dust grains
composed of water ice (20 %), astronomical sillicates (32.91 %), trollite
(7.43 %), and refractory organics (39.66 %). It is worth noting that some
recent observational studies (e.g., Ueda et al., 2023) have raised concerns
about the presence of highly porous dust, which is not favored by the DSHARP
dust opacity table. For further discussion on this matter, we refer to Tazaki
et al. (2019).

In determining the size-averaged dust opacity, we made the assumption that the
dust grain size ($a$) distribution, denoted as $n(a)$, follows a power law
with an exponent of $-q$ (i.e., $a^{-q}$) between the minimum and maximum
grain sizes ($a_{\mbox{\scriptsize min}}$, $a_{\mbox{\scriptsize max}}$), and
is zero elsewhere. The size-averaged dust opacity exhibits a very weak
dependence on the $a_{\mbox{\scriptsize min}}$, which we set to $10^{-4}$ mm.
For our SED models, we adopted a constant dust temperature of
$T_{\mbox{\scriptsize dust}}=18.4$ K111The ring midplane temperature $T_{\rm
mid}$ at $r=20$ au is calculated by a simplified expression for a passively
heated, flared disk in radiative equilibrium (e.g., Dullemond et al., 2001):
$T_{\rm mid}=\left(\frac{\phi L_{*}}{8\pi r^{2}\sigma_{\rm
SB}}\right)^{0.25}$, where $L_{*}$ is the stellar luminosity ($0.19L_{\odot}$;
Herczeg & Hillenbrand, 2014), $\phi$ is the flaring angle (assumed to 0.02),
and $\sigma_{\rm SB}$ is the Stefan–Boltzmann constant., except for the
unresolved inner disk. The value of 18.4 K represents the median expected dust
temperature in the CIDA 1 dust ring. As for the unresolved inner disk, its
dust temperature is not constrained by observations and was fixed to 50 K in
our SED models to ensure consistency with the Rayleigh-Jeans limit, which is
likely the case.

After conducting several trials, we realized that at least four dust
components are required to reproduce the frequency-dependent variation of
spectral indices: an unresolved inner disk and three components within the
dusty ring. These components exhibit qualitative differences in terms of their
maximum grain size $a_{\rm max}$ and dust surface density $\Sigma_{\rm dust}$.
However, it is important to note that the existing observations do not provide
sufficient constraints for the models, as discussed in § 3.2.1, § 3.2.2, and §
3.2.3). Due to this under-constrained nature of the models, the model
parameters are degenerate, making it impossible to perform automatic fits of
the model parameters that are meaningful. Instead, we can only iteratively
determine plausible model parameters that result in predicted flux densities
consistent with the observations. These parameters serve as working hypotheses
for the properties of the inner disk and dusty ring. They can be tested and
refined through future observations and are summarized in Table 4. While
qualitative discussions based on these parameters are feasible, it is
important to exercise caution, as some parameters may not have precise face
values due to the limited constraints provided by the present observations.

Table 4:  Dust models for CIDA 1

Component | $T_{\mbox{\scriptsize dust}}$ | $\Sigma_{\mbox{\scriptsize dust}}$ |  $\Omega_{\mbox{\scriptsize dust}}$aa1 sr $\sim$4.25$\times$1010 square arcsecond, corresponding to $\sim$1.72$\times$1041 cm2 at $d=134.6$ pc of CIDA 1 (Gaia Collaboration et al., 2016, 2022). | $a_{\mbox{\scriptsize max}}$ | Mass
---|---|---|---|---|---
| (K) | (g cm-2) | ($10^{-13}$ sr) | (mm) | ($10^{-4}$ $M_{\sun}$)
Inner DiskbbFlux densities at  X, Ku, Ka, and Q-bands without splitting. | 50.0 | 1.5 | 0.055 | 25 | 0.0071
$a_{\mbox{\scriptsize max}}$=25 mmccThe 3 $\sigma$ upper limit. | 18.4 | 20 | 4.2 | 25 | 7.3
$a_{\mbox{\scriptsize max}}$=0.9 mmccObscured by the $a_{\mbox{\scriptsize max}}=0.55$ mm component which represents the integral dust ring. | 18.4 | 30 | 0.5 | 0.9 | 1.3
$a_{\mbox{\scriptsize max}}$=0.55 mm | 18.4 | 3.0 | 11.5 | 0.55 | 3.0

bbfootnotetext: Dust temperature and $a_{\mbox{\scriptsize max}}$ in the inner
disk were not constrained by observations (see discussion in § 3.1).

####  3.2.1 Inner disk

We approximate the inner disk as an unobscured dust slab with uniform
properties. In addition to the dust temperature, three additional free
parameters ($\Omega_{\rm dust}$, $\Sigma_{\rm dust}$, $a_{\rm max}$) are
required to describe the dust emission within the inner disk. However, due to
the limited availability of independent measurements, the model parameters for
the inner disk exhibit degeneracy. While it is not possible to determine the
physical properties of the inner disk with certainty, modeling its dust
emission remains valuable as it enables us to evaluate the contribution of the
inner disk emission to the overall integrated flux densities.

We set the maximum grain size $a_{\rm max}$ of the inner disk to the largest
value among the components in the dusty ring (25 mm; discussed further in §
3.2.2). The spectral index of the inner disk at 141–339 GHz is $2.20\pm 0.16$
(Table 3), suggesting that the inner disk is marginally optically thick within
this frequency range. Consequently, we can only establish a lower limit for
its dust surface density $\Sigma_{\rm dust}$, which is approximately 1.5 g
cm2. However, the solid angle ($\Omega_{\rm dust}$) of the inner disk is
predominantly determined by the assumed dust temperature, and its face value
should not be employed in scientific discussions. It happens that when $a_{\rm
max}$ is around 25 mm and the optical depth is very high, the frequency-
dependent variation of the dust albedo can account for the observed spectral
index of approximately 2.2 at 141–339 GHz (refer to Figure 3 in Liu 2019).
Thus, our current assumption regarding $a_{\rm max}$ serves as a valid and
intriguing hypothesis that can be tested through future observations with
higher angular resolution at frequencies below 140 GHz. The model indicates
that the inner disk makes a negligible contribution to the integrated flux
density across all observed frequencies.

####  3.2.2 Dusty ring: centimeter region

The calculated spectral indices for the frequency ranges 13.5–16.5 GHz,
16.5–31 GHz, 31–35 GHz, and 35–44 GHz, based on the flux densities listed
above the separator line in Table 3, are  $3.51\pm 2.70$, $2.15\pm 0.44$,
$2.66\pm 0.97$, and $1.70\pm 0.81$, respectively. The drop in spectral indices
observed in the frequency ranges of 16.5–31 GHz (12.6 mm) and 35–44 GHz (7.6
mm) can only be explained when considering the presence of multiple dust
emission components with distinct values of $a_{\rm max}$ and/or $\Sigma_{\rm
dust}$. In fact, to account for the abrupt decrease in spectral indices, it is
crucial to consider the effects of dust scattering, as discussed by Liu 2019.
The anomalously lowered spectral indices due to dust scattering occur within a
narrow frequency range for a specific $a_{\rm max}$ value. Therefore, although
other model parameters may be degenerate, our SED modeling provides a
qualitative estimate of $a_{\rm max}$.

In our current fiducial model, we incorporate two dust components with $a_{\rm
max}$ values of 25 mm and 0.9 mm (Table 4) to reproduce the flux densities at
frequencies below 50 GHz. We assume that there is no mutual obscuration
between these two dust components. The spectral index of the $a_{\rm max}=0.9$
mm component is expected to exhibit an anomalously low value around
wavelengths of approximately 0.9$\times$2$\pi$ mm ($\sim$53 GHz), which
explains the observed low spectral index ($\sim$1.5) at 35–44 GHz. However,
the model tends to overestimate the emission at around 44 GHz (Figure 3). It
is possible that there are additional dust components with $a_{\rm max}$
values close to 0.9 mm, and they may exhibit complex structures of mutual
obscuration. Due to limited signal-to-noise ratio and independent
measurements, we are unable to faithfully reproduce these structures.
Physically, this suggests the presence of a spatial gradient in $a_{\rm max}$
within a region where $a_{\rm max}$ is around 1 mm. The $a_{\rm max}=0.9$ mm
component in our model roughly represents the average dust properties in this
region.

Due to the rapid decrease in spectral index of the $a_{\rm max}=0.9$ mm
component at frequencies below 30 GHz, the presence of another dust component
with a larger $a_{\rm max}$ is required to account for the remaining
centimeter emission. However, the determination of the $a_{\rm max}$ for this
component is challenging due to the degeneracy between its $\Sigma_{\rm dust}$
and $a_{\rm max}$. Nevertheless, it is evident that the $a_{\rm max}$ of this
component must exceed 0.9 mm in order to produce sufficiently strong thermal
dust emission at 12–18 GHz. In our current fiducial model, the $a_{\rm max}$
is approximately 25 mm for this component. It is worth noting that it may be
possible to resolve the $\Sigma_{\rm dust}$–$a_{\rm max}$ degeneracy and
consequently constrain the $a_{\rm max}$ more accurately with follow-up deep
JVLA observations at K-band (18–26 GHz). It is also worth considering the
possibility that there are additional dust components with larger $a_{\rm
max}$, which may exhibit flux densities below the detection threshold at 12–18
GHz in our current observations.

####  3.2.3 Dusty ring: (sub)millimeter region

The spectral index between 44 and 93 GHz is measured to be  $3.11\pm 0.26$,
which is significantly higher than the spectral index between 33 and 44 GHz (
$1.58\pm 0.63$). This sudden change in spectral index creates a dip-like
feature in the SED around 40–50 GHz (Figure 3). This indicates that the flux
density in the 44-93 GHz range is primarily contributed by a dust emission
component with an optical depth of approximately unity (marginally optically
thick/thin) and a dust opacity spectral index ($\beta$) greater than 1.
Consequently, the expected spectral index $\alpha$ is around 3. This dust
emission component, with an optical depth of approximately unity, must obscure
the $a_{\rm max}=25$ mm and $a_{\rm max}=0.9$ mm components that dominate the
centimeter flux densities. Otherwise, the spectral index between 44 and 93 GHz
would be too low to align with the observations, considering that the $a_{\rm
max}=25$ mm and $a_{\rm max}=0.9$ mm components are bright at frequencies
below 50 GHz, possess high column densities (Table 4), and thus have high
optical depths (resulting in a spectral index of approximately 2 between 44
and 93 GHz).

The spectral indices between 93–141 GHz, 141–225.5 GHz, and 225.5–339 GHz are
measured to be $2.13\pm 0.34$, $2.07\pm 0.49$, and $2.04\pm 0.56$,
respectively (Table 4). If we neglect the dust scattering opacity, these
spectral indices of approximately 2 over such a wide frequency range can be
explained by the modified black body emission, as described by Hildebrand
(1983), of an optically thick dust slab (with an optical depth $\tau\gtrsim
10$) in the Rayleigh-Jeans limit. However, considering the assumed dust
temperature of 18.4 K, the Rayleigh-Jeans limit ($h\nu\ll k_{b}T$) is not an
accurate approximation at $\nu=339$ GHz. The non-Rayleigh-Jeans effect causes
the spectral indices to decrease at frequencies above 300 GHz. Moreover, based
on our current understanding of dust opacity, it is not physically feasible to
simultaneously achieve an optical depth of $\tau\gtrsim 10$ at frequencies
above 93 GHz and an optical depth of $\tau\sim$1 between 44 and 93 GHz
because dust opacity is not expected to change by more than one order of
magnitude in such a narrow frequency range (e.g., Figire 11 of Birnstiel et
al., 2018).

To reproduce the observed spectral indices at frequencies above 93 GHz, we
have two options: (1) assuming that a single dust slab with uniform properties
(such as $\Sigma_{\rm dust}$ and $a_{\rm max}$) dominates the emission, or (2)
considering a mixture of emissions from optically thick and optically thin
dust slabs that may have different $a_{\rm max}$ values. Option (1) requires
fewer parameters. In this case, we carefully adjust $a_{\rm max}$ and
$\Sigma_{\rm dust}$ and rely heavily on the frequency-dependent behavior of
dust albedo to match the observed spectral indices. Although it is not always
possible to achieve a perfect match using this approach, our study of CIDA 1
demonstrates its feasibility. Specifically, by assuming $a_{\rm max}\sim 0.5$
mm, the anomalously lowered spectral indices at $a_{\rm max}\times 2\pi$ (as
shown in Figure 3 in Liu 2019) makes it possible to reproduce the spectral
indices of approximately 2 when the optical depth is around 1. Additionally,
the corresponding anomalously elevated spectral indices at higher frequencies
compensate for the non-Rayleigh-Jeans effect.

The alternative approach (case 2) can also explain the observed spectral
indices at frequencies above 93 GHz. However, it is generally not preferable
to rely on fine-tuning parameters for data interpretation. We find that the
dust emission with an $a_{\rm max}$ value of approximately 0.55 mm (as listed
in Table 4) is necessary to account for the observed spectral energy
distribution (SED) between 44 and 141 GHz. If the $a_{\rm max}\sim 0.55$ mm
component does not dominate the emission in the (sub)millimeter range, the
remaining flux densities in that range may be contributed by dust emission
sources that are marginally optically thin and have smaller $a_{\rm max}$
values. Our current model for interpreting the emission above 93 GHz
represents an extreme scenario where we mix dust emission sources with various
$a_{\rm max}$ values, which in a sense overestimates the average $a_{\rm max}$
of the dust components that dominate the emission above 93 GHz.

##  4 Discussions

###  4.1 Comparison with other sources

Previous analyses of the SED using multi-wavelength data have estimated the
maximum grain size $a_{\rm max}$ in the millimeter/centimeter range (e.g.,
Carrasco-González et al., 2019; Ueda et al., 2021; Macías et al., 2021;
Hashimoto et al., 2022; Zhang et al., 2023). Liu (2019) highlighted the impact
of dust scattering on thermal emission from disks, resulting in anomalously
reddened SEDs. This effect has led to an overestimation of
$a_{\mbox{\scriptsize max}}$ and an underestimation of dust masses. This
finding was further supported by subsequent observational studies conducted by
Ueda et al. (2020). A comprehensive exploration of parameters using radiative
transfer modeling has consistently reached the same conclusion by Zhu et al.
(2019). Carrasco-González et al. (2019) conducted a re-analysis of the SED of
HL Tau and demonstrated the significant influence of dust scattering in this
source as well.

Ueda et al. (2021) proposed an interpretation for the observed ALMA
polarization and SED of HL Tau. They suggested that millimeter-sized grains
settle deeper into the midplane of the disk, while 100 $\mu$m-sized grains
remain in the disk atmosphere. This scenario is similar to the case of CIDA 1,
where we assume a multi-layered distribution of dust grains (§ 3). In CIDA 1,
the emission from larger dust grains with $a_{\rm max}$ values ranging from 1
to 25 mm in the disk midplane is obscured by the presence of small dust grains
with an $a_{\rm max}$ of 0.55 mm in the upper layer of the disk.

More recently, Zhang et al. (2023) explored the porous particles to explain
both the polarization and SED of HL Tau. While polarization observations
strongly constrain the $a_{\rm max}$ for compact particles (Kataoka et al.,
2015), porous particles allow for a wider range of $a_{\rm max}$ values that
can explain the observed polarization (Tazaki et al., 2019). In light of this,
Zhang et al. (2023) proposed the presence of porous grains with $a_{\rm max}$
ranging from 1 mm to 1 m in the HL Tau disk at $r\lesssim 60$ au. To test the
potential existence of porous grains in the CIDA 1 ring, conducting multi-
wavelength polarization observations with ALMA could prove valuable. Overall,
the CIDA 1 ring, with its centimeter-sized grains, represents one of the
objects with the largest $a_{\rm max}$ values.

###  4.2 Dust mass

The dust mass of the CIDA 1 ring is approximately 383 $M_{\earth}$ according
to Table 4. Comparing it with previous ALMA studies (e.g., Andrews, 2020), the
dust mass of the CIDA 1 ring appears to be exceptionally large. One possible
reason for our higher estimate of the CIDA 1 ring’s mass is that we assumed it
to be optically thick. Zhu et al. (2019) discussed that a compact disk ($<30$
au), where most of the disk is optically thick in ALMA (sub-)millimeter
observations, can underestimate the actual dust mass by a factor of 10. While
Pinilla et al. (2021) estimated a dust mass (including the central point
source) of approximately 10 $M_{\earth}$ for the CIDA 1 ring, assuming it to
be optically thin, our estimate using the optically thick model for grains
with $a_{\rm max}=$ 0.55 mm, which ALMA is sensitive to, is around 100
$M_{\earth}$. These differences could explain the disparity. Similar to our
estimation, SED analyses based on multi-wavelength data of TW Hya suggested a
dust mass approximately five times higher than those derived from
(sub-)millimeter disk surveys (Macías et al., 2021).

###  4.3 Gravitational stability

Assuming a gas-to-dust mass ratio of the canonical value of 100 (Bohlin et
al., 1978), the total disk mass (gas $+$ dust) of CIDA 1 is approximately 0.12
$M_{\sun}$ (Table 4), while CIDA 1 itself has a mass of 0.19 $M_{\sun}$
(Kurtovic et al., 2021). A large disk-to-star mass ratio exceeding 0.1 can
potentially lead to gravitational instability (e.g., Kratter & Lodato, 2016).
We will now examine the gravitational stability of the system. Theoretical
studies (e.g., Durisen et al., 2007) suggest that if the Toomre (1964)
Q-parameter, defined as

$\displaystyle Q=\frac{c_{\rm s}\Omega_{\rm K}}{\pi G\Sigma_{\rm disk}},$ (2)

where $c_{\rm s}$, $\Omega_{\rm K}$, and $\Sigma_{\rm disk}$ represent the
sound speed, the Keplerian angular velocity, and the disk’s gas surface
density, respectively, is less than unity, the disk may be prone to
gravitational instability. With a midplane temperature of 18.4 K at the ring
located at $r=20$ au (refer to § 3.1) and a gas surface density of $5.3\times
10^{3}$ g/cm-2 (Table 4), assuming the gas-to-dust mass ratio of 100, the
Toomre’s Q-parameter at the CIDA 1 ring is approximately 0.06.

A Q-parameter well below unity indicates that the disk is strongly unstable,
in which spiral arms are expected to be present (e.g., Dong et al., 2015).
However, CIDA 1 does not exhibit such spiral arms in its disk (Pinilla et al.,
2021; Kurtovic et al., 2021). Estimating the Q-value comes with a major
uncertainty, namely the gas-to-dust mass ratio. If the ring represents a dust
trap, the local gas-to-dust ratio could be significantly lower than 100. To
further constrain on the gas mass in the CIDA 1 system, observations of gas-
mass tracers such as HD (Bergin et al., 2013) would provide valuable insights.

###  4.4 Speculation of planet formation around CIDA 1

CIDA 1, with the mass of 0.19 $M_{\sun}$ (Kurtovic et al., 2021), falls into
the category of very low-mass (VLM) stars, typically defined as having a mass
of $\lesssim 0.2M_{\sun}$. According to the core (pebble) accretion scenario,
these systems are expected to be capable of forming only low-mass planets,
typically ranging up to a few Earth masses (Payne & Lodato, 2007; Ormel et
al., 2017; Schoonenberg et al., 2019; Liu et al., 2019a, 2020; Miguel et al.,
2020; Burn et al., 2021).

According to the pebble accretion scenario (e.g., Liu et al., 2019a), the
growth of planetary cores is impeded when a growing core creates a partial gas
gap in the disk, because the local gas pressure bump efficiently traps
accreting pebbles (e.g., Pinilla et al., 2012). Consequently, the planet
becomes isolated from further pebble accretion, and its mass is referred to as
the ‘pebble isolation mass’ ($M_{\rm iso}$). The $M_{\rm iso}$ is determined
by the central stellar mass and the disk aspect ratio ($h_{\rm g}=h/r$) and is
given by the following relationship: $M_{\rm iso}\approx 25(h_{\rm
g}/0.05)^{3}(M_{\star}/M_{\sun})M_{\earth}$ (e.g., Liu et al., 2019a). In
flared disks, where the $h_{\rm g}$ is generally higher in outer disk regions
(e.g., Hayashi et al., 1985), the $M_{\rm iso}$ can reach masses comparable to
giant planets, even for VLM stars.

For CIDA 1, the $h_{\rm g}$ value at the ring position of $r=20$ au is
estimated to be 0.076222The gas scale height $h$ is the ratio of the gas sound
speed to the angular velocity ($h=c_{s}/\Omega$), then the disk aspect ratio
can be written as $h_{\rm g}=h/r\approx 0.03(\frac{M_{\rm
star}}{M_{\sun}})^{-\frac{1}{2}}(\frac{T_{\rm mid}}{300{\rm
K}})^{\frac{1}{2}}(\frac{r}{1{\rm au}})^{\frac{1}{2}}$, where $T_{\rm mid}$ is
a disk midplane temperature with 18.4 K (§ 3.2)., resulting in an $M_{\rm
iso}$ of 16.7 $M_{\earth}$. The dust mass in the CIDA 1 ring is estimated to
be approximately 383 $M_{\earth}$ (Table 4), which could provide sufficient
material for the formation of giant planets (Jang et al., 2022). The presence
of centimeter-sized pebbles in the CIDA 1 ring (§ 3) may further facilitate
the growth of planetary cores (e.g., Morbidelli et al., 2015). Therefore, it
is plausible that giant planets could form within the CIDA 1 ring during the
typical disk lifetime of a few million years.

Another avenue of planet formation in the CIDA 1 ring is gravitational
instability (e.g., Boss, 1997; Mercer & Stamatellos, 2020). When the
gravitationally-unstable disks cool sufficiently fast ($t_{\rm
cool}\Omega\lesssim$ a few, where $t_{\rm cool}$ and $\Omega$ are the cooling
time and orbital frequency), they may fragment (e.g., Gammie, 2001; Rice et
al., 2003; Dong et al., 2016), resulting in bound, self-gravitating objects
(i.e., companions). As discussed in § 4.3, the CIDA 1 system might be
gravitationally unstable but the disk does not show eminent asymmetric
structures. More studies, especially on gas-mass estimates, are needed for
further investigation.

##  5 Conclusion

We conducted observations of the CIDA 1 ring using the JVLA at centimeter
wavelengths ranging from $\lambda=0.7$ cm  (44 GHz) to 3.0 cm  10 GHz. We
successfully detected signals at $\lambda=0.7$ cm  (44 GHz), 0.9 cm  (33 GHz),
and 2.0 cm  (15 GHz). At $\lambda=0.7$ (44 GHz) and 0.9 cm  (33 GHz), we were
able to spatially resolve a partial ring structure. Although a potential weak
asymmetry was observed in the ring, further deep observations are required to
achieve a more reliable detection. Based on the spatial distributions and
spectral indexes of these emissions, we attributed them to dust rather than
free-free or synchrotron emissions.

The emission in the frequency range of 15–340 GHz was analyzed using a simple
SED model consisting of four dust components. It is important to note that our
modeling results depend on the uncertain properties of the dust. The model
includes an inner dust disk component, as well as dust ring components with
maximum grain sizes ($a_{\mbox{\scriptsize max}}$) of 25, 0.9, and 0.55 mm.
The analysis suggests that grain growth in the CIDA 1 ring proceeds to at
least centimeter-sized grains.

While recent ALMA dust polarization observations have suggested the presence
of sub-millimeter-sized dust grains (approximately 100 $\mu$m) in Class II
disks (e.g., Kataoka et al., 2016b), it is important to note that these
observations were biased towards bright disks. Despite being faint in
comparison (with a total flux of approximately 10 mJy at $\nu=225.5$ GHz;
Andrews et al., 2013) compared to the brighter DSHARP disks ($\gtrsim 100$ mJy
at band 6; Andrews et al., 2018), the CIDA 1 ring, with its centimeter-sized
grains, is placed in a unique parameter space in terms of dust size. The exact
origins of grain growth in the CIDA 1 ring remain uncertain. Expanding the
sample size to include more disks with millimeter to centimeter-sized grains
will provide valuable insights into grain growth and shed light on key
physical processes such as the snowline, the orbital timescale, and
turbulence.

The authors thank the anonymous referee for a timely and constructive report.
We thank Paola Pinilla for sharing fits files of ALMA dust continuum images at
band-4 and 7. The National Radio Astronomy Observatory is a facility of the
National Science Foundation, operated under a cooperative agreement by
Associated Universities, Inc. This work has made use of data from the European
Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed
by the Gaia Data Processing and Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the Gaia Multilateral Agreement. This study used the
following ALMA data: ADS/JAO. ALMA #2011.0.00259.S and 2018.1.00536.S. ALMA is
a partnership of ESO (representing its member states), NSF (USA), and NINS
(Japan), along with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic
of Korea), in cooperation with the Republic of Chile. The Joint ALMA
Observatory is operated by ESO, AUI/NRAO, and NAOJ. This study was supported
by JSPS KAKENHI Grant Numbers 21H00059, 22H01274, 23K03463 and 18H05441. T.M.
is supported by Yamada Science Foundation Overseas Research Support Program.

##  Appendix A Zoomed-in Ka-band image

The image with the highest signal-to-noise ratio (SNR) among our observations,
captured in the Ka-band, is overlaid on the contours of the ALMA band-7 image
for a closer view.

Figure 4:  Zoomed-in Ka-band image (0$\farcs$18 $\times$ 0$\farcs$16)
superposed with contours of ALMA band-7 image (54 mas $\times$ 37 mas; Pinilla
et al., 2021). Contour levels are 17.1 $\mu$Jy beam-1 (1 $\sigma$) $\times$[5,
15, 45]. The stellar position is marked by a black star

## References

  * Andrews (2020) Andrews, S. M. 2020, ARA&A, 58, 483, doi: 10.1146/annurev-astro-031220-010302
  * Andrews et al. (2013) Andrews, S. M., Rosenfeld, K. A., Kraus, A. L., & Wilner, D. J. 2013, ApJ, 771, 129, doi: 10.1088/0004-637X/771/2/129
  * Andrews et al. (2018) Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018, ApJ, 869, L41, doi: 10.3847/2041-8213/aaf741
  * Anglada et al. (1998) Anglada, G., Villuendas, E., Estalella, R., et al. 1998, AJ, 116, 2953, doi: 10.1086/300637
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33, doi: 10.1051/0004-6361/201322068
  * Bacciotti et al. (2018) Bacciotti, F., Girart, J. M., Padovani, M., et al. 2018, ApJ, 865, L12, doi: 10.3847/2041-8213/aadf87
  * Bergin et al. (2013) Bergin, E. A., Cleeves, L. I., Gorti, U., et al. 2013, Nature, 493, 644, doi: 10.1038/nature11805
  * Birnstiel et al. (2018) Birnstiel, T., Dullemond, C. P., Zhu, Z., et al. 2018, ApJ, 869, L45, doi: 10.3847/2041-8213/aaf743
  * Bohlin et al. (1978) Bohlin, R. C., Savage, B. D., & Drake, J. F. 1978, ApJ, 224, 132, doi: 10.1086/156357
  * Boss (1997) Boss, A. P. 1997, Science, 276, 1836, doi: 10.1126/science.276.5320.1836
  * Burn et al. (2021) Burn, R., Schlecker, M., Mordasini, C., et al. 2021, A&A, 656, A72, doi: 10.1051/0004-6361/202140390
  * Carrasco-González et al. (2019) Carrasco-González, C., Sierra, A., Flock, M., et al. 2019, ApJ, 883, 71, doi: 10.3847/1538-4357/ab3d33
  * Carrera et al. (2015) Carrera, D., Johansen, A., & Davies, M. B. 2015, A&A, 579, A43, doi: 10.1051/0004-6361/201425120
  * Curone et al. (2022) Curone, P., Izquierdo, A. F., Testi, L., et al. 2022, A&A, 665, A25, doi: 10.1051/0004-6361/202142748
  * Dent et al. (2019) Dent, W. R. F., Pinte, C., Cortes, P. C., et al. 2019, MNRAS, 482, L29, doi: 10.1093/mnrasl/sly181
  * Dong et al. (2015) Dong, R., Hall, C., Rice, K., & Chiang, E. 2015, ApJ, 812, L32, doi: 10.1088/2041-8205/812/2/L32
  * Dong et al. (2016) Dong, R., Vorobyov, E., Pavlyuchenkov, Y., Chiang, E., & Liu, H. B. 2016, ApJ, 823, 141, doi: 10.3847/0004-637X/823/2/141
  * Draine (2006) Draine, B. T. 2006, ApJ, 636, 1114, doi: 10.1086/498130
  * Drazkowska et al. (2022) Drazkowska, J., Bitsch, B., Lambrechts, M., et al. 2022, arXiv e-prints, arXiv:2203.09759.  https://arxiv.org/abs/2203.09759
  * Dullemond et al. (2001) Dullemond, C. P., Dominik, C., & Natta, A. 2001, ApJ, 560, 957, doi: 10.1086/323057
  * Durisen et al. (2007) Durisen, R. H., Boss, A. P., Mayer, L., et al. 2007, in Protostars and Planets V, ed. B. Reipurth, D. Jewitt, & K. Keil, 607.  https://arxiv.org/abs/astro-ph/0603179
  * Gaia Collaboration et al. (2016) Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016, A&A, 595, A1, doi: 10.1051/0004-6361/201629272
  * Gaia Collaboration et al. (2022) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2022, arXiv e-prints, arXiv:2208.00211, doi: 10.48550/arXiv.2208.00211
  * Galván-Madrid et al. (2014) Galván-Madrid, R., Liu, H. B., Manara, C. F., et al. 2014, A&A, 570, L9, doi: 10.1051/0004-6361/201424630
  * Gammie (2001) Gammie, C. F. 2001, ApJ, 553, 174, doi: 10.1086/320631
  * Hashimoto et al. (2022) Hashimoto, J., Liu, H. B., Dong, R., Liu, B., & Muto, T. 2022, ApJ, 941, 66, doi: 10.3847/1538-4357/aca01d
  * Hayashi et al. (1985) Hayashi, C., Nakazawa, K., & Nakagawa, Y. 1985, in Protostars and Planets II, ed. D. C. Black & M. S. Matthews, 1100–1153
  * Herczeg & Hillenbrand (2014) Herczeg, G. J., & Hillenbrand, L. A. 2014, ApJ, 786, 97, doi: 10.1088/0004-637X/786/2/97
  * Hildebrand (1983) Hildebrand, R. H. 1983, QJRAS, 24, 267
  * Jang et al. (2022) Jang, H., Liu, B., & Johansen, A. 2022, A&A, 664, A86, doi: 10.1051/0004-6361/202243368
  * Johansen et al. (2014) Johansen, A., Blum, J., Tanaka, H., et al. 2014, in Protostars and Planets VI, ed. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 547, doi: 10.2458/azu_uapress_9780816531240-ch024
  * Kataoka et al. (2016a) Kataoka, A., Muto, T., Momose, M., Tsukagoshi, T., & Dullemond, C. P. 2016a, ApJ, 820, 54, doi: 10.3847/0004-637X/820/1/54
  * Kataoka et al. (2015) Kataoka, A., Muto, T., Momose, M., et al. 2015, ApJ, 809, 78, doi: 10.1088/0004-637X/809/1/78
  * Kataoka et al. (2016b) Kataoka, A., Tsukagoshi, T., Momose, M., et al. 2016b, ApJ, 831, L12, doi: 10.3847/2041-8205/831/2/L12
  * Kratter & Lodato (2016) Kratter, K., & Lodato, G. 2016, ARA&A, 54, 271, doi: 10.1146/annurev-astro-081915-023307
  * Kurtovic et al. (2021) Kurtovic, N. T., Pinilla, P., Long, F., et al. 2021, A&A, 645, A139, doi: 10.1051/0004-6361/202038983
  * Liu & Ji (2020) Liu, B., & Ji, J. 2020, Research in Astronomy and Astrophysics, 20, 164, doi: 10.1088/1674-4527/20/10/164
  * Liu et al. (2019a) Liu, B., Lambrechts, M., Johansen, A., & Liu, F. 2019a, A&A, 632, A7, doi: 10.1051/0004-6361/201936309
  * Liu et al. (2020) Liu, B., Lambrechts, M., Johansen, A., Pascucci, I., & Henning, T. 2020, A&A, 638, A88, doi: 10.1051/0004-6361/202037720
  * Liu & Ormel (2018) Liu, B., & Ormel, C. W. 2018, A&A, 615, A138, doi: 10.1051/0004-6361/201732307
  * Liu (2019) Liu, H. B. 2019, ApJ, 877, L22, doi: 10.3847/2041-8213/ab1f8e
  * Liu et al. (2014) Liu, H. B., Galván-Madrid, R., Forbrich, J., et al. 2014, ApJ, 780, 155, doi: 10.1088/0004-637X/780/2/155
  * Liu et al. (2019b) Liu, H. B., Mérand, A., Green, J. D., et al. 2019b, ApJ, 884, 97, doi: 10.3847/1538-4357/ab391c
  * Liu et al. (2021) Liu, H. B., Tsai, A.-L., Chen, W. P., et al. 2021, ApJ, 923, 270, doi: 10.3847/1538-4357/ac31b9
  * Macías et al. (2021) Macías, E., Guerra-Alvarado, O., Carrasco-González, C., et al. 2021, A&A, 648, A33, doi: 10.1051/0004-6361/202039812
  * Mamajek (2009) Mamajek, E. E. 2009, in American Institute of Physics Conference Series, Vol. 1158, Exoplanets and Disks: Their Formation and Diversity, ed. T. Usuda, M. Tamura, & M. Ishii, 3–10, doi: 10.1063/1.3215910
  * McMullin et al. (2007) McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in Astronomical Society of the Pacific Conference Series, Vol. 376, Astronomical Data Analysis Software and Systems XVI, ed. R. A. Shaw, F. Hill, & D. J. Bell, 127
  * Mercer & Stamatellos (2020) Mercer, A., & Stamatellos, D. 2020, A&A, 633, A116, doi: 10.1051/0004-6361/201936954
  * Miguel et al. (2020) Miguel, Y., Cridland, A., Ormel, C. W., Fortney, J. J., & Ida, S. 2020, MNRAS, 491, 1998, doi: 10.1093/mnras/stz3007
  * Morbidelli et al. (2015) Morbidelli, A., Lambrechts, M., Jacobson, S., & Bitsch, B. 2015, Icarus, 258, 418, doi: 10.1016/j.icarus.2015.06.003
  * Ohashi et al. (2018) Ohashi, S., Kataoka, A., Nagai, H., et al. 2018, ApJ, 864, 81, doi: 10.3847/1538-4357/aad632
  * Ohashi et al. (2020) Ohashi, S., Kataoka, A., van der Marel, N., et al. 2020, ApJ, 900, 81, doi: 10.3847/1538-4357/abaab4
  * Ormel & Liu (2018) Ormel, C. W., & Liu, B. 2018, A&A, 615, A178, doi: 10.1051/0004-6361/201732562
  * Ormel et al. (2017) Ormel, C. W., Liu, B., & Schoonenberg, D. 2017, A&A, 604, A1, doi: 10.1051/0004-6361/201730826
  * Payne & Lodato (2007) Payne, M. J., & Lodato, G. 2007, MNRAS, 381, 1597, doi: 10.1111/j.1365-2966.2007.12362.x
  * Pérez et al. (2012) Pérez, L. M., Carpenter, J. M., Chandler, C. J., et al. 2012, ApJ, 760, L17, doi: 10.1088/2041-8205/760/1/L17
  * Pérez et al. (2015) Pérez, L. M., Chandler, C. J., Isella, A., et al. 2015, ApJ, 813, 41, doi: 10.1088/0004-637X/813/1/41
  * Pinilla et al. (2012) Pinilla, P., Birnstiel, T., Ricci, L., et al. 2012, A&A, 538, A114, doi: 10.1051/0004-6361/201118204
  * Pinilla et al. (2018) Pinilla, P., Natta, A., Manara, C. F., et al. 2018, A&A, 615, A95, doi: 10.1051/0004-6361/201832690
  * Pinilla et al. (2021) Pinilla, P., Kurtovic, N. T., Benisty, M., et al. 2021, A&A, 649, A122, doi: 10.1051/0004-6361/202140371
  * Rau & Cornwell (2011) Rau, U., & Cornwell, T. J. 2011, A&A, 532, A71, doi: 10.1051/0004-6361/201117104
  * Ricci et al. (2014) Ricci, L., Testi, L., Natta, A., et al. 2014, ApJ, 791, 20, doi: 10.1088/0004-637X/791/1/20
  * Rice et al. (2003) Rice, W. K. M., Armitage, P. J., Bate, M. R., & Bonnell, I. A. 2003, MNRAS, 339, 1025, doi: 10.1046/j.1365-8711.2003.06253.x
  * Schoonenberg et al. (2019) Schoonenberg, D., Liu, B., Ormel, C. W., & Dorn, C. 2019, A&A, 627, A149, doi: 10.1051/0004-6361/201935607
  * Simon et al. (2017) Simon, M., Guilloteau, S., Di Folco, E., et al. 2017, ApJ, 844, 158, doi: 10.3847/1538-4357/aa78f1
  * Tazaki et al. (2019) Tazaki, R., Tanaka, H., Kataoka, A., Okuzumi, S., & Muto, T. 2019, ApJ, 885, 52, doi: 10.3847/1538-4357/ab45f0
  * Tazzari et al. (2016) Tazzari, M., Testi, L., Ercolano, B., et al. 2016, A&A, 588, A53, doi: 10.1051/0004-6361/201527423
  * Testi et al. (2003) Testi, L., Natta, A., Shepherd, D. S., & Wilner, D. J. 2003, A&A, 403, 323, doi: 10.1051/0004-6361:20030362
  * Testi et al. (2014) Testi, L., Birnstiel, T., Ricci, L., et al. 2014, in Protostars and Planets VI, ed. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 339, doi: 10.2458/azu_uapress_9780816531240-ch015
  * Toomre (1964) Toomre, A. 1964, ApJ, 139, 1217, doi: 10.1086/147861
  * Ueda et al. (2020) Ueda, T., Kataoka, A., & Tsukagoshi, T. 2020, ApJ, 893, 125, doi: 10.3847/1538-4357/ab8223
  * Ueda et al. (2021) Ueda, T., Kataoka, A., Zhang, S., et al. 2021, ApJ, 913, 117, doi: 10.3847/1538-4357/abf7b8
  * Ueda et al. (2023) Ueda, T., Okuzumi, S., Kataoka, A., & Flock, M. 2023, arXiv e-prints, arXiv:2305.12598, doi: 10.48550/arXiv.2305.12598
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science and Engineering, 13, 22, doi: 10.1109/MCSE.2011.37
  * Zhang et al. (2023) Zhang, S., Zhu, Z., Ueda, T., et al. 2023, arXiv e-prints, arXiv:2306.00158, doi: 10.48550/arXiv.2306.00158
  * Zhu et al. (2019) Zhu, Z., Zhang, S., Jiang, Y.-F., et al. 2019, ApJ, 877, L18, doi: 10.3847/2041-8213/ab1f8c