Diagnosing 0.1-10 au scale morphology of the FU\,Ori disk using ALMA and VLTI/GRAVITY

TL;DR

This study combines high-resolution ALMA and VLTI/GRAVITY observations to analyze the inner 10 au of the FU Ori disk, revealing a nearly face-on, viscously heated accretion structure with possible dust settling, advancing understanding of young stellar object disks.

Contribution

First combined millimeter and near-infrared high-resolution observations of FU Ori's inner disk, revealing detailed structure and accretion processes at 0.1-10 au scales.

Findings

Inner 0.4 au disk has a triangular spectral shape with water and CO features.

Viscous heating explains the spectral energy distribution from 0.4 to 10 au.

Evidence suggests vertical dust settling within 10 au, pending further confirmation.

Abstract

We report new Atacama Large Millimeter/submillimeter Array Band 3 (86-100 GHz; $\sim$80 mas angular resolution) and Band 4 (146-160 GHz; $\sim$50 mas angular resolution) observations of the dust continuum emission towards the archetypal and ongoing accretion burst young stellar object FU Ori, which simultaneously covered its companion, FU Ori S. In addition, we present near-infrared (2-2.45 $\mu$m) observations of FU Ori taken with the General Relativity Analysis via VLT InTerferometrY (GRAVITY; $\sim$1 mas angular resolution) instrument on the Very Large Telescope Interferometer (VLTI). We find that the emission in both FU Ori and FU Ori S at (sub)millimeter and near infrared bands is dominated by structures inward of $\sim$10 au radii. We detected closure phases close to zero from FU Ori with VLTI/GRAVITY, which indicate the source is approximately centrally symmetric and therefore is likely viewed nearly face-on. Our simple model to fit the GRAVITY data shows that the inner 0.4 au radii of the FU Ori disk has a triangular spectral shape at 2-2.45 $\mu$m, which is consistent with the H$_{2}$O and CO absorption features in a $\dot{M}\sim$10$^{-4}$ $M_{\odot}\,yr^{-1}$, viscously heated accretion disk. At larger ($\sim$0.4-10 au) radii, our analysis shows that viscous heating may also explain the observed (sub)millimeter and centimeter spectral energy distribution when we assume a constant, $\sim$10$^{-4}$ $M_{\odot}\,yr^{-1}$ mass inflow rate in this region. This explains how the inner 0.4 au disk is replenished with mass at a modest rate, such that it neither depletes nor accumulates significant masses over its short dynamic timescale. Finally, we tentatively detect evidence of vertical dust settling in the inner 10 au of the FU Ori disk, but confirmation requires more complete spectral sampling in the centimeter bands.

Tables

Table 1. Table 1: Fluxes measurements from ALMA

	FU Ori		FU Ori S
Frequency	Flux	Flux error	Flux	Flux error
(GHz)	(mJy)	(mJy)	(mJy)	(mJy)
( Band 3 )
86.001	2.1	0.088	1.2	0.13
87.863	2.3	0.13	1.4	0.12
98.196	2.8	0.11	1.5	0.11
100.001	3.0	0.11	1.6	0.14
( Band 4 )
146.002	6.2	0.18	3.5	0.15
147.863	6.2	0.17	3.5	0.18
158.196	7.0	0.17	4.0	0.19
160.002	7.4	0.20	4.0	0.19

Table 2. Table 2: (Sub)millimeter spectral indices

	FU Ori	FU Ori S
Frequency range (GHz)	Spectral index ($\alpha$)
29-37 GHz	1.6$\pm$0.4	1.4$\pm$0.4
29-100 GHz	2.5$\pm$0.05	2.7$\pm$0.05
86-160 GHz	2.0$\pm$0.07	1.9$\pm$0.07
146-232 GHz	2.0$\pm$0.05	2.0$\pm$0.05
218-346 GHz	2.9$\pm$0.2	2.2$\pm$0.2

Table 3. Table 3: Models for the VLTI/GRAVITY data

#	VLTI-compact	VLTI-extended	$\Delta x_e$	$\Delta y_e$	${\chi }_r^2$
	uniform disk	2D Gaussian
	diameter ($mas$)	FWHM ($mas$)	($mas$) to East	($mas$) to North
0	0.0	$\infty$	0	0	3.1
1	0.0	$4.76\pm 0.04$	0	0	1.5
2	$1.14\pm 0.01$	$\infty$	0	0	1.25
3	$1.06\pm 0.01$	$8.4\pm 0.2$	0	0	1.18
4	$1.06\pm 0.01$	$7.9\pm 0.1$	$0.6\pm 0.1$	$1.0\pm 0.1$	1.03

Table 4. Table 4: Fluxes of the best model for the VLT/GRAVITY data

$\lambda$	VLTI-compact	VLTI-extended	Ratio
($\mu$m)	(Jy)	(Jy)	$\pm 0.001$
2.00	4.40	0.13	0.030
2.05	5.34	0.15	0.029
2.10	5.92	0.18	0.031
2.15	6.15	0.21	0.034
2.20	6.00	0.24	0.040
2.25	5.40	0.27	0.051
2.30	4.79	0.31	0.065
2.35	4.20	0.35	0.083
2.40	3.73	0.39	0.104
2.45	3.44	0.43	0.125

Table 5. Table 5: Parameters for SED fittings

FU Ori
	Free-free emission
	$T_{\text{e}}$	EM	${\Omega }_{\text{ff}}$		Obscured by
	(103 K)	(cm-6pc)	(10-14 sr)
	7.0${}_{-1.3}{}^{}{}_{}{}^{+1.4}$	1.8${}_{-4.7}{}^{}{}_{}{}^{+4.7}$$\times$107	1.4${}_{-0.21}{}^{}{}_{}{}^{+0.20}$		none
	Dust components
Comp.	$T_{\text{dust}}$	${\Sigma }_{\text{dust}}$	${\Omega }_{\text{dust}}$	$a_{\text{max}}$	Obscured by	$M_{\text{dust}}$
	(K)	(g cm-2)	(sr)	(mm)		($M_{\oplus }$)
1	370${}_{-49}{}^{}{}_{}{}^{+63}$	45${}_{-12}{}^{}{}_{}{}^{+10}$	4.1${}_{-0.56}{}^{}{}_{}{}^{+0.64}$$\times$10-14	2.4${}_{-0.32}{}^{}{}_{}{}^{+0.40}$	Comp. 2	510${}_{-190}{}^{}{}_{}{}^{+210}$
2	140${}_{-26}{}^{}{}_{}{}^{+20}$	0.63${}_{-0.099}{}^{}{}_{}{}^{+0.12}$	7.1${}_{-1.6}{}^{}{}_{}{}^{+1.5}$$\times$10-14	0.21${}_{-0.019}{}^{}{}_{}{}^{+0.020}$	none	12${}_{-3.9}{}^{}{}_{}{}^{+6.0}$
3	55${}_{-8.6}{}^{}{}_{}{}^{+6.8}$	0.13${}_{-0.021}{}^{}{}_{}{}^{+0.024}$	2.1${}_{-0.28}{}^{}{}_{}{}^{+0.28}$$\times$10-12	0.0017${}_{-0.00048}{}^{}{}_{}{}^{+0.00047}$	envelope	75${}_{-20}{}^{}{}_{}{}^{+25}$
4	690${}_{-12}{}^{}{}_{}{}^{+18}$	0.0095${}_{-0.0018}{}^{}{}_{}{}^{+0.0019}$	1.8${}_{-0.32}{}^{}{}_{}{}^{+0.34}$$\times$10-15	5.2${}_{-1.3}{}^{}{}_{}{}^{+1.5}$	none	4.7${}_{-1.6}{}^{+2.0}\times$10-3

Equations

$$\begin{split}\tau_{\nu}^{\mbox{\scriptsize ff}}=&\\ &8.235\times 10^{-2}\left(\frac{T_{\mbox{\scriptsize e}}}{\mbox{K}}\right)^{-1.35}\left(\frac{\nu}{\mbox{GHz}}\right)^{-2.1}\left(\frac{\mbox{EM}}{\mbox{pc\,cm${}^{-6}$}}\right),\\ \end{split}$$

$$F_{\nu }^{\text{mbox}ff}={\Omega }_{\text{mbox}ff}(1-e^{-{\tau }_{\nu }^{\text{mbox}ff}})B_{\nu }(T_{\text{mbox}e}),$$

$$\begin{array}{rl}{F}_{\nu }=& \sum _i{F}_{\nu }^i{e}^{-\sum _j{\tau }_{\nu }^{i,j}},\end{array}$$

$$T_{\text{mbox}viscous}(r)={\left[\frac{3G{M}_{\ast }\dot{M}}{8\pi \sigma r^3}\left(1-\sqrt{\frac{R_{\ast }}{r}}\right)\right]}^{\frac{1}{4}},$$

Full text

stars: individual (FU Ori) — protoplanetary disks

††software: CASA (McMullin et al., 2007), Numpy (Van Der Walt et al., 2011), emcee (Foreman-Mackey et al., 2013)

1 Introduction

Understanding the physical mechanisms of protostellar accretion is fundamentally important in studies of star formation. Optical and near infrared surveys have shown that young stellar objects (YSOs) are 10-100 times underluminous with respect to the expected luminosity from steady accretion (Kenyon & Hartmann, 1995; Evans et al., 2009), which indicates that YSOs may accrete episodically (Dunham & Vorobyov, 2012). If episodic accretion is a widespread phenomenon during the YSO phases, it should manifest observationally. This is consistent with the discoveries of two types of YSOs in outburst: the FU Orionis (FUor) and the EX Lupi (EXor) objects, which are characterized by a rapid large increase in their optical and infrared (OIR) brightness (for reviews, see Hartmann & Kenyon, 1996; Herbig, 2007; Audard et al., 2014).

FUors have outburst durations of decades to centuries (Hartmann & Kenyon, 1996). During the outburst state, their optical brightness can increase by $\sim$4 magnitudes or more. Models suggest that the accretion rates of these YSOs vary from 10*-7* M*⊙* yr*-1* in the low (T Tauri) accretion state to 10*-4* M*⊙* yr*-1* in the high (FUors) accretion state (Hartmann & Kenyon, 1996). While accretion processes in quiescent T Tauri stars are generally understood as magnetospheric streams from the inner disk (e.g., Koenigl, 1991; Calvet et al., 2000), how the gas and dust reservoirs immediately around the FUors are different (or altered) compared with quiescent T Tauri stars is not yet well-understood. This limits our understanding of the outburst triggering mechanisms and the consequences of them.

To shed light on this issue, we have performed high angular resolution observations towards the archetypal FU Orionis object, FU Ori, using the Atacama Large Millimeter Array (ALMA) and the General Relativity Analysis via VLT InTerferometrY (GRAVITY) instrument of the Very Large Telescope Interferometer (VLTI). Throughout this manuscript, we assume the distance of FU Ori to be $d\sim$416 pc, according to the parallax measurement published in the second data release of the Gaia space telescope (Gaia Collaboration et al., 2018). According to the prior-assisted parallax distances measurements of Bailer-Jones et al. (2018), we quote a nominal $\pm$2% distance uncertainty, which will not qualitatively affect our analysis. The observations are introduced in Section 2, and the results are presented in Section 3. By jointly analyzing these new observations with the previous (sub)millimeter observations of the ALMA, the Submillimeter Aray (SMA)111The Submillimeter Array is a joint project between the Smithsonian Astrophysical Observatory and the Academia Sinica Institute of Astronomy and Astrophysics, and is funded by the Smithsonian Institution and the Academia Sinica (Ho et al., 2004) and the NRAO222The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. Karl G. Jansky Very Large Array (JVLA), and the infrared spectra taken with the Spitzer and Herschel333Herschel is an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA. space telescopes, our interpretation and the further discussion about the physical implications are provided in Section 4. Our conclusion is given in Section 5. We refer to Berger et al. (2012) for a review of the convention and terminology for the optical and infrared interferometry technique.

2 Observations

We introduce the archival Spitzer and Herschel spectra and the VLTI/GRAVITY observations in Sections 2.1 and 2.2. We provide detailts of our ALMA observations in Section 2.3.

2.1 Spitzer and Herschel spectra

The Herschel/PACS and SPIRE spectra were taken from the COPS-DIGIT-FOOSH (CDF) archive, a high-level data product provided to the Herschel Science Archive (see Green et al. 2016a for details). Because the source size at these wavelengths is comparable to the SPIRE beam size, and because of the lack of background subtraction, the spectra of the two modules of SPIRE instruments (SLW and SSW) are often mismatched. To resolve this discrepancy, Green et al. (2016a) apply the Semi-Extended Correction Tool in hipe (Ott, 2010; Wu et al., 2013) to calibrate the SPIRE spectra by modeling the source size. The best-fit source size, 23$\farcs$5 for FU Ori, is then convolved with the beam profile of SPIRE, which is a function of wavelength (Makiwa et al., 2013). Therefore, the resulting SPIRE spectrum represents the emission from different apertures at given wavelengths, which correspond to the convolved sizes of the beam and the source size (Yang et al., 2018). For example, the aperture sizes are 29$\farcs$9, 34$\farcs$5, and 43$\farcs$6 at 250 $\mu$m, 350 $\mu$m, and 500 $\mu$m, respectively. We refer to Section 2.2 of Green et al. (2016a) for a complete description of the data reduction.

2.2 VLTI/GRAVITY observations

FU Ori was observed by VLTI/GRAVITY (Gravity Collaboration et al., 2017) on 2016 November 25 and 26. These observations were part of the consortium Guaranteed Time Observations (Program ID 098.C-0765). The observations were carried out at both medium and high spectral resolution, although only the medium resolution data achieved sufficient signal to noise ratio (S/N) for analysis. The medium spectral resolution setting covered the whole near-infrared K-band with a spectral resolution of $\sim$500.

The telescopes chosen for these observations were the medium Auxiliary Telescopes (AT) configurations. This configuration includes the stations K0-G2-D0-J3, which led to baselines ranging from $\sim$40 to $\sim$100 meters (Figure 1).

The calibrator star observed concurrently was HD 38494, which is a K2 star of unknown luminosity class. Its photometric angular diameter was estimated to be $\theta_{\mathrm{UD}}=0.71\pm 0.06$ mas according to the Jean-Marie Mariotti Center Stellar Diameters Catalogue (JSDC) (Bourges et al., 2017). It is nearly unresolved for our observations. The visibility of the fringes is expected to range from 0.971 to 0.996. The uncertainty of the photometric angular diameter of HD 38494 leads to a bias in the reduced data of at most 0.005 in visibility. We reduced the data using the GRAVITY pipeline (Lapeyrere et al., 2014) version 1.0.11. We note that the direct observables from VLTI/GRAVITY are normalized to the total flux. The reduced spectrum does not clearly present emission lines (Figure 1).

2.3 ALMA observations

We have performed ALMA Bands 3 and 4 observations towards FU Ori, which also covered its $\sim$0$\farcs$5 separation companion, FU Ori S (Project code: 2017.1.00388.S, PI: Hauyu Baobab Liu). The pointing and phase referencing center was R.A. (J2000) = 05${}^{\mbox{\scriptsize h}}$45${}^{\mbox{\scriptsize m}}$22$\fs$375, and Decl. (J2000) = +09*∘04′*12$\farcs$400. The uv distance ranges covered by these observations are $\sim$100 m – 13 km. The correlators were configured to cover four 1.875 GHz wide spectral windows with a 976.562 kHz channel spacing.

The Band 3 observations were carried out on 2017 November 08. The four spectral windows were centered on the sky frequencies of 86.000, 87.863, 98.196, and 100.001 GHz. We observed the quasar J0510+1800 for absolute flux and passband calibrations, and J0547+1223 for complex gain calibrations.

The Band 4 observations were carried out on 2017 November 07. The four spectral windows were centered on the sky frequencies of 146.001, 147.863, 158.196, and 160.001 GHz. We observed the quasar J0510+1800 for absolute flux and passband calibrations, and J0536+0944 for complex gain calibrations.

We manually calibrated and phase self-calibrated these data using the CASA software package (McMullin et al., 2007) version 5.4.0. When performing absolute flux scaling, we assumed that J0510+1800 has a 2.0 Jy absolute flux and a $-$0.30 spectral index at the reference frequency 93.015 GHz; and has a 1.6 Jy absolute flux and a $-$0.4 spectral index at the reference frequency 153.016 GHz. These assumptions were based on interpolating the calibrator grid survey measurements. We produced the Briggs Robust = 0 weighted continuum images from line-free spectral channels using the CASA task clean. For each of the two observed bands, we created images for each of the four spectral windows separately using the multi-frequency synthesis (MFS) method, setting the parameter nterm=1. The four spectral windows in each band achieved comparable root-mean-square (RMS) noise levels and angular resolutions. At Band 3, the spectral window centered at 100 GHz achieved a $\theta_{\mbox{\scriptsize maj}}\times\theta_{\mbox{\scriptsize min}}$=0*′′.082$\times$0′′.075 (P.A.=-79∘) synthesized beam and a 69 $\mu$Jy beam-1* RMS noise level; at Band 4, the spectral window centered at 160 GHz achieved a $\theta_{\mbox{\scriptsize maj}}\times\theta_{\mbox{\scriptsize min}}$=0*′′.047$\times$0′′.043 (P.A.=57∘) synthesized beam and a 72 $\mu$Jy beam-1* RMS noise level. In each band, the synthesized beam sizes at other spectral windows are inversely proportional to their central frequencies.

3 Results

3.1 VLTI/GRAVITY data

Figure 2 shows the reduced VLTI/GRAVITY data on FU Ori. The closure phases (CP) are smaller than $\pm$2.5 degrees. In addition, the overall scatter of the CP is less that 1 degree. This indicates that, on the spatial scales resolved by our VLTI/GRAVITY observations, FU Ori appears approximately centro-symmetric. Overall, the squared visibilities have a fairly high level, ranging from 0.7 to 0.9. The variations of the squared visibilities with wavelength are similar for all baselines, irrespective of baseline lengths and orientations: The squared visibilities are approximately constant from wavelength 2.0 to 2.2 $\mu$m, and then drop by about 0.1 from 2.2 to 2.45$\mu$m.

The fact that the differential visibility variations do not seem to depend on baseline lengths indicates that the intensity distributions may be approximated by a compact component (hereafter VLTI-compact) at the center and a more extended centro-symmetric component (hereafter VLTI-extended)444Note that the two infrared emission components resolved by the VLTI/GRAVITY observations are both more compact than what were detected by JVLA and ALMA. Our terminology is to distinguish them from the spatially more extended (sub)millimeter and centimeter sources.. Including a structure which is nearly resolved out by all baselines leads to the observed $<$1.0 squared visibilities. In this case, the detected values of the squared visibilities depend mostly on the flux ratios of the VLTI-compact and the VLTI-extended components.

To give a qualitative sense, if we define VLTI-unresolved as having a visibility higher than 0.99 and VLTI-resolved as having a visibility amplitude less than 0.01, then for our longest, $\sim$100 meter baseline, a uniform disk with $\lesssim$0.3 mas diameter is VLTI-unresolved. For our shortest, $\gtrsim$30 meter baseline, a two dimensional Gaussian with $\sim$17 mas full width at half maximum is VLTI-resolved. The visibility amplitude of a 1 mas compact uniform disk ranges from 0.93 to 0.99 for the baselines ranging from 60 m to 100 m. The visibility amplitude of a FWHM$=$8 mas Gaussian ranges from $\sim$0.0001 to $\sim$0.15.

3.2 ALMA data

Figure 3 shows the ALMA 100 GHz (Band 3) and 160 GHz (Band 4) images. These ALMA observations detected FU Ori and FU Ori S (Reipurth & Aspin, 2004; Wang et al., 2004) at high significances. However, they only marginally spatially resolved the structures. Their fluxes determined by fitting two-dimensional Gaussians are summarized in Table 1. By quoting the previous JVLA observations at 29-37 GHz (Liu et al., 2017) and the ALMA observations at $\sim$225 GHz (Pérez et al. submitted) and at $\sim$346 GHz (Hales et al., 2015), the derived (sub)millimeter spectral indices ($\alpha$) at various frequency ranges are summarized in Table 2. Figures 4 and 5 summarize the spectral energy distributions (SEDs) of these two protostars at wavelengths from 2 $\mu$m to 33 mm (9 - 1.5$\times$105 GHz).

Both FU Ori and FU Ori S show spectral index values lower than 2.0555 The spectral index $\alpha$ was measured assuming that the flux $F_{\nu}$ around a reference frequency $\nu_{0}$ can be expressed as $F_{\nu}=F_{0}(\nu/\nu_{0})^{\alpha}$ at 29-37 GHz (c.f., Liu et al., 2017); the averaged $\alpha$ are $\sim$2.5 over the frequency range of 29-100 GHz; the averaged $\alpha$ is approximately 2.0 from 100 to 232 GHz, and are higher than 2.0 at higher frequencies. We require multiple emission components with distinct physical properties to fit the complex submillimeter spectral slopes in the observed SEDs. Our detailed SED models for all data presented in Figures 4 and 5 are described in Section 4.

4 Discussion

In Section 4.1 we introduce a simple geometric model to interpret the VLTI/GRAVITY observations. In addition, we have generated simple radiative transfer models to interpret the SEDs of FU Ori and FU Ori S. In Section 4.2 we introduce how we produced the spectra for individual dust or free-free (i.e., from ionized gas) emission components in our radiative transfer model. In Section 4.3, we introduce how we integrate each of the emission components to the abstracted geometric models to reproduce the integrated SEDs, and how we optimized the model free parameters using the Markov chain Monte Carlo (MCMC) method. We discuss the physical implications of our models in Section 4.4.

4.1 Interpreting VLTI/GRAVITY data

The fact that the observed squared visibilities in the VLTI/GRAVITY data vary with wavelength (Figure 2, right panel) implies that the flux ratio of the VLTI-compact and the VLTI-extended components has a wavelength dependence.

We can quantify this dependence by fitting the data. For simplicity, we assumed that the VLTI-compact and the VLTI-extended components have constant sizes over the wavelength range covered by the VLTI/GRAVITY observations. In addition, we assumed that the VLTI-compact component is a uniform disk, while the VLTI-extended component is a two dimensional Gaussian of which the aspect ratio is $\sim$1. We then performed chi-squared fits to determine the sizes of the two components, and to determine the flux ratios as a linear interpolation between eight equally spaced wavelengths between 2.0 and 2.45 $\mu$m (R$\sim$40). We tried various combinations of sizes, allowing the two components to be unresolved, partially resolved, or fully resolved. Given that the observed closure phases are less than 2 degrees, to avoid over-fitting, the two components were concentric in most trials. However, in one of the trials, we also explored how much their centers can deviate.

Our best-fit geometric models are summarized in Table 3. The model in best agreement with data is a uniform disk of $\sim$1 mas in diameter, and a two dimensional Gaussian with FWHM$\sim$8 mas (i.e., solid angle $\sim$1.7$\times$10*-15* sr). The fit is further improved if the VLTI-extended component is slightly offset to the Northeast by $\sim$1.2 mas (Figure 6, right panel).

Assuming that the VLTI/GRAVITY detections arose predominantly from the circumstellar disk, this spatial offset can be interpreted either as a disk that is geometrically thick (e.g., flared) and is slightly inclined (e.g., Figure 5 of Zhu et al., 2008), or as a disk that includes substructures or companions (e.g., Malbet et al., 2005). Using the flux ratios from Table 4, we were able to reproduce the observed slight closure phase signal which increases with wavelength, and the wavelength-dependent variations of the squared visibilities (Figure 2). By implementing an absolute flux scaling, the spectral shape of the VLTI-compact component in our best fit model (Figure 6) appears fully consistent with the viscous accretion disk model of Calvet et al. (1991), which assumed a $\sim$10*-4* $M_{\odot}$ yr*-1* mass accretion rate. The triangular shape of the spectrum presented in Figure 6, following the framework of Calvet et al. (1991), is due to the absorption of the water band and the first-overtone vibration-rotation CO band against the bright continuum emission from the viscously heated mid-plane. Calvet et al. (1991) suggested that these absorption features are predominantly produced at radii of 0.1-0.3 au (i.e., 0.96$\pm$0.48 mas angular diameter assuming $d\sim$416 pc) around the host protostar, which is consistent with the angular sizes in our model fits.

4.2 Individual emission components

To evaluate the emission properties of dust, we utilized the DSHARP dust optical constants published in Birnstiel et al. (2018). For simplicity, we assumed a constant 170 K water ice sublimation temperature (Pollack et al., 1994). Therefore, we adopted the default DSHARP optical constants for dust emission sources which are cooler than 170 K, and adopted the ice-free optical constants for those which are warmer. Given that the physical conditions of the observed sources (in particular, FU Ori) may be out of equilibrium in various ways, it is not possible for us to evaluate the detailed form of the grain size distribution function from first principles. Therefore, when evaluating the size-averaged dust absorption ($\kappa_{\nu}^{\mbox{\scriptsize abs}}$) and effective scattering ($\kappa_{\nu}^{\mbox{\scriptsize sca,eff}}$) opacities, we simply assumed the typical power-law grain size distribution with a power-law index $q=$3.5, the minimum grain size $a_{\mbox{\scriptsize min}}=$10*-4* mm, and the maximum grain size $a_{\mbox{\scriptsize max}}$ . Before considering mutual obscuration, the SEDs of individual dust emission components ($F^{\mbox{\scriptsize dust}}_{\nu}$) were evaluated based on the analytic radiative transfer solutions published in Birnstiel et al. (2018). Motivated by the small (or negligible) angular offsets of the unresolved and resolved components in the VLTI/GRAVITY models (Table 3; see discussion in Section 4.1), we considered all dust slabs to be approximately face-on. We note that introducing inclinations of the dust slabs will not change the conclusion from our radiative transfer models qualitatively. However, this would increase the total number of free parameters.

Figure 7 shows examples of the SEDs produced for the dust slabs with dust column density of 50 g cm*-2*, temperature of 100 K, solid angle 1 square arcsecond, and maximum grain sizes of 0.002 mm (top panel), 0.2 mm (middle panel), and 2 mm (bottom panel). In the low frequency, low optical depth regime, the SEDs deviate from the blackbody emission model (i.e., Planck function) as dust grains cannot emit/absorb efficiently at wavelengths which are much longer than $a_{\mbox{\scriptsize max}}$ .

In addition, for $a_{\mbox{\scriptsize max}}$ $=$0.2 mm or 2 mm, the SED deviates from a blackbody curve at higher frequencies. As frequency increases, the spectral indices fall below a blackbody curve, and then become steeper; thus the flux in this frequency regime is below that of a blackbody curve. The dust slabs are optically thick in this frequency regime, and the effects of dust (self-)scattering are not necessarily negligible. We attribute this deviation from blackbody emission in the high frequency regime to the frequency variations of albedo, which was addressed in detail in Liu (2019) and Zhu et al. (2019). For example, in Figure 7, the SED of the $a_{\mbox{\scriptsize max}}$ $=$2 mm dust slab shows a rather flat spectral index at $\sim$20-50 GHz, which is because the albedo increases with frequency; the spectral index is steepened at $\sim$50-1000 GHz because the albedo decreases with frequency.

Following Mezger & Henderson (1967) and Keto (2003), we approximated the optical depth of the free-free emission components $\tau^{\mbox{\scriptsize ff}}_{\nu}$ by

[TABLE]

where EM is the emission measure defined as EM$=$ $\int n_{\mbox{\scriptsize e}}^{2}d\ell$, with $n_{\mbox{\scriptsize e}}$ being the electron number density, and $\ell$ is the linear size scale of the free-free emission component along the line of sight. Fluxes of individual free-free emission components were evaluated based on

[TABLE]

where $\Omega_{\mbox{\scriptsize ff}}$ is the solid angle of the free-free emission component, and $B_{\nu}(T)$ is the Planck blackbody function.

4.3 Abstracted geometric model and integrated SEDs

The overall fluxes ($F_{\nu}$) of FU Ori and FU Ori S were determined from the following formulation:

[TABLE]

where $F_{\nu}^{i}$ is the flux of the dust or free-free emission component $i$, and $\tau^{i,j}_{\nu}$ is the optical depth of the emission component $j$ to obscure the emission component $i$. The abstracted geometric information is provided by $\tau^{i,j}_{\nu}$. We chose this approach instead of fitting analytical solutions of (gaseous) disks because dusty protoplanetary disks are commonly composed of sub-structures (e.g., rings, crescent, etc). Our SED fitting procedure for the spatially unresolved target sources effectively decomposed them into sub-structures of certain projected areas but without explicitly constraining the shapes. In this work we considered a simple implementation, such that $\tau^{i,j}_{\nu}=0$ if the emission component $i$ is not obscured by the emission component $j$; otherwise $\tau^{i,j}_{\nu}=\tau^{j}_{\nu}$.

We tried fitting the observed SEDs with the least number of emission components to minimize the total number of free parameters, for various configurations of $\tau^{i,j}$. the free parameters and the configurations of $\tau^{i,j}$ were varied interactively, informed by the results of previous trials. Our interactive fits focused on matching the interferometric data. Nevertheless, we found that once a good fit for the interferometric data was achieved, the infrared spectrum predicted from the model is also very close to the Spitzer and Herschel observations.

After we obtained an approximate fit, we used MCMC to simultaneously optimize all free parameters (i.e., all the parameters in Table 5 except the column of overall dust masses). We assumed flat priors, which permitted each parameter to vary from half of its initial value to two times of the initial value. To prevent the MCMC routine from sampling large unlikely portions of parameter space, we provided an additional constraint from the FU Ori S 9 GHz non-detection (Figure 5). We forced the logarithmic likelihood to be negative infinity when the integrated flux of FU Ori S at 9 GHz is higher than three times the RMS noise of the observations, a condition in the likelihood function to force rejecting such MCMC samples.

The MCMC fittings were initialized with 84 walkers with 1500 iterative steps each; in the end, the results from the first 500 steps were discarded. The Herschel and Spitzer data have very good signal-to-noise ratios, such that their contribution to the likelihood largely outweighted the contribution from interferometric data. To avoid overfitting the Herschel and Spitzer data without achieving a good fit for the interferometric data, we needed to reduce the weight of Herschel and Spitzer data. This is implemented by artificially assigning the flux errors of the Herschel/SPIRE, Herschel/PACS, and Spitzer/IRS to be 1000, 10, and 1 times the detected fluxes. We have monitored how the likelihood evolved over the MCMC iterations to make sure that the contribution of the Herschel and Spitzer data are on the same order with the rest of the data. We note that during the steps of MCMC, some dust emission components may switch from being based on the default DSHARP optical constants to being based on the ice-free optical constants (i.e., the walkers "walked" from below to above the 170 K dust temperature). Because of this mid-routine shift, it is very difficult to implement fitting methods other than MCMC.

We found at least four dust emission components are required to fit the JVLA, ALMA, and VLTI/GRAVITY data points for FU Ori (for more discussion see Section 4.4). Therefore, we also adopted a four emission component fit for FU Ori S.

In addition, we included an extended common envelope component which is required to fit the far-infrared fluxes detected by the Herschel space telescope. The common envelope component has an extended angular scale, such that it is filtered out by all interferometric observations presented in this work. We note that the envelope component must be included since the previous Herschel photometric imaging observations have spatially resolved complicated structures on sub-parsec scales, which connect to FU Ori and FU Ori S (Green et al., 2013). We used a simplified parametric model for the envelope (Table 5); detailed modeling of the envelope is beyond the scope of our present study.

During the iterations, we found that we can obtain a reasonably good fit to the JVLA and ALMA data of FU Ori S by including only two dust emission components and a free-free emission component. We tentatively assign one additional, $\sim$140 K dust emission component to FU Ori S to better explain the Herschel or Spitzer spectra at (5-10)$\times$103 GHz. Qualitatively, the fact that we need this extra component to explain the mid-far infrared spectra indicates that the dust components in our models are not isothermal. We considered whether each dust component in our models should be allowed to have a small (e.g., 10%-20%) temperature range, which could yield better fits to the infrared spectra. However, our ability to measure any temperature variation is fundamentally limited by the wavelength-dependent aperture used to extract the Herschel/SPIRE spectra. This is because we applied semi-extended source correction to align the two SPIRE modules (Wu et al., 2013; Green et al., 2016a). Residual artifacts from this process can bias our SED fits, although we have mitigated this by artificially lowering the weighting of the Herschel data. Nevertheless, we do not consider it to be meaningful to use a further detailed parameterization for dust temperature profiles to improve the fittings to our infrared spectra.

Parameters of our best fit model are summarized in Table 5. The SEDs of the individual components after incorporating the effect of obscuration, and the integrated SEDs from all emission components, are presented in Figures 4 and 5.

4.4 Model parameters and their physical implications

In Sections 4.4.1 and 4.4.2 we discuss qualitatively the fitting parameters for FU Ori and FU Ori S. The overall geometric picture and the physical implications are discussed in Section 4.4.3.

4.4.1 FU Ori model

Qualitatively, the fact that the spectral index $\alpha$ of FU Ori is $\sim$2 at 86-232 GHz and is $\sim$3 at frequencies higher than 232 GHz (Table 2) indicates that fluxes at intermediate frequency (e.g., $\sim$150 GHz) are a mix of one emission component with $\alpha>$3 and the other emission source with $\alpha<$2: The $\alpha>$3 component (hereafter FUOri_dust3) becomes more prominent at higher frequencies, while the $\alpha<$2 source becomes more prominent at lower frequencies; the observed spectral indices at the intermediate frequencies are weighted averages from these two components.

In order to fit the JVLA data at 29-33 GHz and the ALMA data at 86-160 GHz (Table 1; Figure 5), we need to realize the $\alpha<$2 source by combining at least two dust components: a $\sim$400 K component with high dust column density and $a_{\mbox{\scriptsize max}}$ $\sim$2 mm (hereafter FUOri_dust1), obscured by a $\sim$130 K component with modest dust column density and $a_{\mbox{\scriptsize max}}$ $\sim$0.2 mm (hereafter FUOri_dust2). The dust temperature of FUOri_dust1 is consistent with the high dust brightness temperature observed at $\sim$33 GHz (Liu et al., 2017). The FUOri_dust1 component, whether or not it is mixed with some free-free emission, naturally explains the $<$2.0 spectral index at 29-37 GHz, and the $\sim$2.5 spectral index in between 29-100 GHz (Table 2), due to the albedo effect introduced in Section 4.2 (see also Figure 7). Being obscured by FUOri_dust2 makes the spectral index of FUOri_dust1 much lower than 2.0 at $\sim$100-150 GHz. FUOri_dust2 has $a_{\mbox{\scriptsize max}}$ $\sim$0.2 mm because this $a_{\mbox{\scriptsize max}}$ value yields a high albedo at $\sim$200 GHz. In this case, FUOri_dust2, which is optically thick at $\sim$200 GHz and has a rather flat spectral distribution at this frequency, can scatter off the emission from FUOri_dust1 without contributing much of the emission. This is critical to fit the steeper spectral indices at 232-345 GHz including the optically thin dust component FUOri_dust3. If one or both of FUOri_dust1 and FUOri_dust2 contributes more emission at $\sim$200 GHz, it becomes impossible to reproduce the steep spectral index observed at 232-345 GHz.

The $a_{\mbox{\scriptsize max}}$ of FUOri_dust3 is not well-constrained by the data presented in this paper. Consistent with previous reports of near-infrared scattered light (e.g., Liu et al., 2016; Takami et al., 2018), we presume that the $a_{\mbox{\scriptsize max}}$ of FUOri_dust3 is on the order of $\sim$2 $\mu$m. We cannot accurately determine the dust masses of these two components due to the uncertainties of the dust mass opacities.

Finally, by introducing another FUOri_dust4 component we can simultaneously fit the resolved VLTI-extended component (Table 3) in the VLTI/GRAVITY data and the higher frequency part of the Spitzer spectrum (Figure 4). The dust temperature of FUOri_dust4 ($\sim$700 K) is higher than that of FUOri_dust1 ($\sim$400 K), indicating that FUOri_dust4 is likely the closest component to the host protostar. The $a_{\mbox{\scriptsize max}}$ value of FUOri_dust4 is not constrained by the observations presented in this paper. The resolved VLTI-compact component (Table 3) is hotter than the dust sublimation temperature and therefore is not considered in our dust models. The thermal radiation from the VLTI-compact component may heat the VLTI-extended component (c.f., Zhu et al. 2007). We note that an excellent fit to the 9 GHz observations with the free-free emission component was not necessary, because that particular measurement was impacted by poorly-characterized delay errors, and is rather uncertain (c.f., Liu et al. 2017; Pérez et al. submitted). Emission at 9 GHz may also include a non-thermal emission contribution, which we do not have sufficient data to constrain.

4.4.2 FU Ori S model

The spectral index of FU Ori S is $\sim$ 2 over a broad frequency range of 86-346 GHz (Table 2). An optically thick dust component (FUOriS_dust1) with $a_{\mbox{\scriptsize max}}$ $\sim$0.2 mm can explain the slightly smaller than 2.0 spectral index at $\sim$ 90 GHz. Mixing FUOriS_dust1 with an optically thinner dust component (FUOriS_dust3) and a free-free emission component can better fit the observations at 29-37 GHz and at 346 GHz (Figure 5). To reproduce the resolved SEDs for FU Ori S, there is no need of assuming mutual obscurations of the emission components since its spectral index at 232-345 GHz is not as steep as that of FU Ori (Table 5).

4.4.3 Outbursting versus quiescent disks?

By assuming a geometrically thin, axisymmetric, Keplerian, hot inner disk around the center of FU Ori, Calvet et al. (1991) argued that the observed CO linewidths at near infrared bands are consistent with an inclination of $\sim$20*∘-60∘. Based on analyzing the squared visibilities from near and mid-infrared interferometric observations, Malbet et al. (2005), Zhu et al. (2008) and Quanz et al. (2006) suggested that the inclination of the disk is $\sim$50∘*. However, being an accretion outburst object, FU Ori may not be in equilibrium. The assumptions of geometrically thin, axisymmetry, and the Keplerian velocity fields all need to be tested by resolved observations. A great advantage of VLTI/GRAVITY over the previous generations near- or mid-infrared interferometry is that we can anchor the hypothesis of small inclination angle based on the resolved small closure phases (Figure 2; Section 4.1). Thus on the spatial scales of a few au, the morphology and the gas kinematics of the FU Ori disk may be more complicated than previously assumed, which can be further resolved by future observations with better uv coverage.

If we assume an approximately face-on projection of FU Ori (and FU Ori S), then the abstracted geometry we introduced during the SED fits (Section 4.3) follows the picture in Figure 8. Overall, we interpret the observational data for FU Ori as the following: a $>$1000 K hot inner disk at 0.1-0.3 au radii (0.24-0.72 mas) which produces water and CO absorption features at near-infrared bands (unresolved by VLTI/GRAVITY); a $\sim$700 K, not very optically thick dust component with $\sim$3 au radius ($\sim$7 mas; FUOri_dust4, resolved by VLTI/GRAVITY); a very optically thick and a modestly optically thick dust component with up to $\sim$10 au radii (FUOri_dust1,2); an optically thin, cooler dust component on tens of au scales (FUOri_dust3), and some free-free emission. Assuming that the gas-to-dust mass ratio is $\sim$100, the mass surface density of the component FUOri_dust1 (Table 5) is reasonably consistent with the hydrodynamic simulations presented in Zhu et al. (2010) and Bae et al. (2014). However, the vertical thermal profile of FU Ori in its inner 10 au region appears opposite to the typical model of passive disks dominated by radiative heating (e.g., some T Tauri disks; Kama et al., 2009; Tapia & Lizano, 2017). The role of viscous heating in dust thermal dynamics is presently uncertain as it is difficult to observationally constrain gas volume density and viscosity at the disk mid-plane.

FU Ori S can be interpreted as an optically thick dust component with $\sim$10 au radius (FUOriS_dust1) and an optically thin, cooler dust component on tens of au scales (FUOriS_dust3), potentially with contributions from free-free emission.

The qualitative difference between the inner $\sim$10 au region of the FU Ori and the FU Ori S disks, in particular the thermal profile, may be related to the thermal and magnetorotational instabilities triggered during the outburst of FU Ori. For a physical picture, we refer to Figures 1 and 2 of Zhu et al. (2009). The two dimensional hydrodynamic simulations of Zhu et al. (2009) have demonstrated that during outburst, within $\lesssim$10 au of the protostar, viscous energy dissipation is sufficient to heat gas at the disk mid-plane to a considerably higher temperature than the gas at the disk surface. It is not yet very clear to us whether or not this can also explain the vertical dust temperature profile of the FU Ori disk.

For an order-of-magnitude estimate, we compare our fits of dust temperatures (Table 5) with the simplest analytic models of the dust temperature profiles (c.f., Chiang & Goldreich 1997 and references therein) in Figure 9. We caution that many of the underlying assumptions of the simplest analytic models (e.g., axisymmetry, steady or stationary disk, etc.) contradict the observations of the FU Ori disk which is likely asymmetric and may be undergoing instabilities over a broad spatial scale. Such a comparison can serve as a sanity check for whether or not a certain heating mechanism can potentially provide a sufficiently high heating rate to explain the observed dust radiation temperatures. However, the comparison is not yet sufficient for verifying or strictly falsifying a certain scenario.

To assess how the FU Ori disk can be heated due to viscous dissipation, we quoted the effective radiation temperature profile of a steady-state viscous disk $T_{\mbox{\scriptsize viscous}}(r)$ assuming that the disk is very dense and is optically thick such that dust and gas can be thermalized via inelastic collisions, and that there is no radiative heating (Pringle, 1981):

[TABLE]

where $G$ is the gravitational constant, $M_{*}=$0.5 $M_{\odot}$ is the assumed host protostellar mass, $\dot{M}$ is the mass accretion rate which we assumed to be 10*-8*, 10*-6*, and 10*-4* $M_{\odot}$ yr*-1*, $\sigma$ is the Stephen-Boltzmann constant, and $R_{*}$ is the stellar radius which we assumed to be 2 $R_{\odot}$. These profiles, which may be regarded as lower limits to the dust temperature in viscous disks, are presented as the blue lines in Figure 9.

To assess how the FU Ori S disk can be heated due to protostellar irradiation, we scaled the approximate solutions for the surface ($T_{\mbox{\scriptsize s}}(r)$) and interior ($T_{\mbox{\scriptsize i}}(r)$) dust temperature profiles of a radiative equilibrium disk (c.f. Equations 11 and 14a in Chiang & Goldreich 1997) according to the total protostellar luminosity. The upper and lower bounds of the yellow filled area are shown with respect to $T_{\mbox{\scriptsize s}}(r)$ and $T_{\mbox{\scriptsize i}}(r)$ in Figure 9. We assumed an effective stellar temperature $T_{*}=$4000 K and stellar radius $R_{*}=$2 $R_{\odot}$, typical for T Tauri stars. We note that due to the Stefan-Boltzmann law the dust temperatures have a weak dependence on the protostellar luminosity.

In Figure 9, we also overplotted our fits of dust components (c.f., Table 5). The inner and outer radii of these dust components were estimated to be 1% and 100% of their solid angle, assuming a circular geometry in a face-on projection. We found that it is plausible to interpret the observed radiation temperature of FUOri_dust1 based on $T_{\mbox{\scriptsize viscous}}(r)$ given the $\sim$10*-4* $M_{\odot}$ yr*-1* accretion rate of FU Ori. If this is the case, a higher dust temperature at the disk mid-plane than at the surface can be expected, which explains why FUOri_dust1 has a higher temperature than FUOri_dust2 (Figure 8). Moreover, this explains how the 0.1-0.3 au scales hot inner disk with a $10^{-4}$ $M_{\odot}$ yr*-1* mass accretion rate (c.f., Section 4.1) is being replenished by the up to $\sim$10 au scales gas reservoir at a modest rate, such that the hot inner disk neither becomes depleted nor accumulates mass over a short time scales. This may explain the relatively stable mid-infrared and (sub)millimeter fluxes in the previous monitoring observations (Green et al., 2016b; Liu et al., 2018).

The optically thinner components FUOri_dust3, FUOriS_dust3, FUOri_dust4, and FUOriS_dust4 are likely dominated by radiative heating. Radiative heating alone can reasonably explain the observed temperature distributions from FU Ori S.

The comparisons in Figure 9 are uncertain since the accretion rates of FU Ori and FU Ori S are not necessarily constant over all radii. In addition, FU Ori is unlikely to be in equilibrium, and it is not trivial to accurately estimate the disk scale-height and thus the radiative heating. Moreover, these comparisons have ignored other mechanical processes which can potentially be important in asymmetric or unstable systems (e.g., shocks, adiabatic compression, etc; Dong et al. 2016; Sakai et al. 2014). More realistic considerations of dust and gas dynamics, grain growth, and dust heating/cooling would provide better comparison. We additionally hypothesize that, during the outburst, the inner 0.1-10 au disk may expand significantly in the vertical direction, may be partly thermally ionized, and some dust may be sublimated. The morphology of the 0.1-10 au disk may also become porous due to accretion and instabilities, allowing dust to be radiatively heated close to the disk mid-plane at a relatively large range of radii.

Finally, why might we have detected millimeter sized $a_{\mbox{\scriptsize max}}$ from FU Ori (i.e., from component FUOri_dust1) but not from FU Ori S? A tentative hypothesis is that at the quiescent stage, dust grains of millimeter or larger sizes may either be radially trapped in regions too small in projected area to be detected by observations (e.g., Vorobyov et al. 2018; Okuzumi & Tazaki 2019), or areas that are fully obscured due to a combination of very high optical depth and the vertical dust settling. These mechanisms may be particularly efficient if the inner few au regions are effectively dead zones with negligible ionization fraction during the quiescent stage. The instabilities during the outburst may help radially and vertically mix dust grains of various sizes, which make the millimeter-sized grains more easily detectable. That we find tentative evidence of vertical dust settling by comparing the $a_{\mbox{\scriptsize max}}$ values of FUOri_dust1 and FUOri_dust2, may also be because viscous heating is more efficient in heating the vertically settled grown dust from the mid-plane. This may be further tested by a systematic comparison of the (sub)millimeter and radio spectral indices of the inner disks of outbursting and the quiescent T Tauri sources.

5 Conclusion

We have analyzed unpublished archival data from the Guaranteed Time Observations of VLTI/GRAVITY at near infrared K-band (2-2.45 $\mu$m) towards the archetypal accretion outburst young stellar object, FU Ori. In addition, we have performed high angular resolution ALMA observations at 86-100 GHz and 146-160 GHz bands, which simultaneously covered FU Ori and its companion, FU Ori S.

The observed small closure phases by VLTI/GRAVITY indicate that the FU Ori disk may be approximately face-on. In addition, by comparing with the squared visibilities resolved by previous generation near and mid-infrared interferometry, we found that the inner few au region of FU Ori may not be simply an axisymmetric, Keplerian rotating thin disk. Instead, it may have a more complicated morphology, which may be related to the instabilities which occurred during the accretion outbursts. Combined analysis of all existing ALMA, SMA, and JVLA observations along with Spitzer and Herschel infrared spectra also points to an unconventional vertical dust thermal profile in the inner $\sim$10 au region of the FU Ori disks. This consistently suggests a complicated disk morphology in comparison to a quiescent T Tauri disk. The observed thermal profiles in the inner $\sim$10 au region may be explained by a viscously heated disk of which the mass inflow rate is $\sim$10*-4* $M_{\odot}$ yr*-1*, which can explain how the 0.1-0.3 scales hot inner disk detected from infrared observations is being replenished.

This paper is based on data obtained from the ESO Science Archive Facility under request number AMERAND384481. This paper makes use of the following ALMA data: ADS/JAO.ALMA #2011.0.00548.S, #2016.1.01228.S, and #2017.1.00388.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together 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 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 work is based [in part] on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. H.B.L. is supported by the Ministry of Science and Technology (MoST) of Taiwan (Grant Nos. 108-2112-M-001-002-MY3 and 108-2923-M-001-006-MY3). E. Vorobyov acknowledges financial support from the Russian Foundation for Basic Research (RFBR), Russian-Taiwanese project #19-52-52011. R.G.M. acknowledges support from UNAM-PAPIIT Programme IN104319. Y.-L. Yang acknowledges the support of University Continuing Graduate Fellowship from The University of Texas at Austin. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. Y.H. is supported by the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme under grant agreement No 716155 (SACCRED, PI: Á. Kóspál). ALMA, VLTI/GRAVITY