11institutetext:  (1) European Southern Observatory, Casilla 19001, Santiago
19, Chile
(2) University of Exeter, School of Physics and Astronomy, Astrophysics Group,
Stocker Road, Exeter, EX4 4QL, UK
(3) University of Michigan, Department of Astronomy, S University Ave, Ann
Arbor, MI 48109, USA
(4) Institut de Planetologie et d’Astrophysique de Grenoble, Grenoble 38058,
France
(5) The CHARA Array of Georgia State University, Mount Wilson Observatory,
Mount Wilson, CA 91023, USA

# Imaging the warped dusty disk wind environment of SU Aurigae with MIRC-X

Aaron Labdon 11 Stefan Kraus 22 Claire L. Davies 22 Alexander Kreplin 22
Sebastian Zarrilli 22 John D. Monnier 33 Jean-Baptiste le Bouquin 44
Narsireddy Anugu Benjamin Setterholm 5533 Tyler Gardner 33 Jacob Ennis 33
Cyprien Lanthermann 55 Theo ten Brummelaar 55 Gail Schaefer 55 Tim J. Harries
22

(accepted May 18, 2019)

###### Abstract

Context. T Tauri stars are low-mass young stars whose disks provide the
setting for planet formation, one of the most fundamental processes in
astronomy. Yet the mechanisms of this are still poorly understood. SU Aurigae
is a widely studied T Tauri star and here we present original state-of-the-art
interferometric observations with better uv and baseline coverage than
previous studies.

Aims. We aim to investigate the characteristics of the circumstellar material
around SU Aur, constrain the disk geometry, composition and inner dust rim
structure.

Methods. The MIRC-X instrument at CHARA is a 6 telescope optical beam combiner
offering baselines up to 331 m. We undertook image reconstruction for model-
independent analysis, and fitted geometric models such as Gaussian and ring
distributions. Additionally, the fitting of radiative transfer models
constrains the physical parameters of the disk.

Results. Image reconstruction reveals a highly inclined disk with a slight
asymmetry consistent with inclination effects obscuring the inner disk rim
through absorption of incident star light on the near-side and thermal re-
emission/scattering of the far-side. Geometric models find that the underlying
brightness distribution is best modelled as a Gaussian with a FWHM of $1.53\pm
0.01\,\mathrm{mas}$ at an inclination of $56.9\pm 0.4^{\circ}$ and minor axis
position angle of $55.9\pm 0.5^{\circ}$. Radiative transfer modelling shows a
flared disk with an inner radius at 0.16 au which implies a grain size of
$0.14\,\mathrm{\mu m}$ assuming astronomical silicates and a scale height of
9.0 au at 100 au. In agreement with literature, only the dusty disk wind
successfully accounts for the NIR excess by introducing dust above the mid-
plane.

Conclusions. Our results confirm and provide better constraints than previous
inner disk studies of SU Aurigae. We confirm the presence of a dusty disk wind
in the cicumstellar environment, the strength of which is enhanced by a late
infall event which also causes very strong misalignments between the inner and
outer disks.

###### Key Words.:

Stars: individual: SU Aurigae – Stars: variables: T Tauri, Herbig Ae/Be –
Techniques: interferometric – Protoplanetary disks

##  1 Introduction

Outflows from protoplanetary systems are one of the key mass loss mechanisms
during the planet formation process. They remove both excess mass and angular
momentum from the system, a crucial process as the final masses and rotation
rates of stars are known to be significantly less than the inital mass of
protostellar cores.

Within young stellar objects (YSOs) there are several different mechanisms of
outflow. Firstly, accretion and magnetospherically driven jets can emerge from
the poles of the star. Secondly, photoevaporative winds caused by the ultra-
violet (UV) disassociation of molecules in the upper layers of the outer disk
can cause significant mass loss. Such photoevaporative winds are usually
associated with higher-mass, hotter objects such as Herbig Ae/Be stars.
Finally, magnetospherically driven dusty disk winds can originate from the
inner disk whereby material is lifted from the disk plane along inclined
magnetic field lines. Magnetospheric winds require the presence of a strong
magnetic field, usually associated with T Tauri stars with convective
envelopes rather than fully radiative interiors. This allows for optically
thick material to exist close enough to the central star to contribute to the
Near-Infrared (NIR) emission exterior to the main disk structure. This model
has been shown to successfully account for the NIR excess of the spectral
energy distribution (SED) and the basic visibility features of AB Aur, MWC 275
and RY Tau (Königl & Salmeron 2011; Petrov et al. 2019). While all these
mechanisms have been observed, it is not fully understood why some YSOs only
appear to exhibit a sub-selection of outflow mechanisms. One of the first
stars observed to have an inner dusty disk wind was the T Tauri star SU
Aurigae (Petrov et al. 2019).

For a full description of the literature surrounding SU Aurigae and the basic
stellar properties, see the previous paper by these authors; Labdon et al.
(2019) (hence forth LA19). In this previous work we studied the circumstellar
environment of SU Aur using interferometric observations from the CHARA/CLIMB
and PTI (Palomar Testbed Interferometer) instruments. The disk was found to be
inclined at $51.2\pm 1.2^{\circ}$ with a position angle of $61.0\pm
1.0^{\circ}$ and was best modelled with a ring-like geometry with a radius of
$0.17\pm 0.02\,\mathrm{au}$. Additionally, radiative transfer modelling of
visibilities and the SED found that the NIR excess could only be reproduced in
the presence of a dusty disk wind, where material is lifted from the disk
along magnetic field lines allowing the reprocessing of additional stellar
radiation.

Table 1: Observing log from 2018 from the CHARA/MIRC-X interferometer.

Date | Beam Combiner | Stations | Pointings | Calibrator (UD [mas])
---|---|---|---|---
2018-09-13 | CHARA/MIRC-X | S1-S2-E1-E2-W1-W2 | 2 |  HD 34499 ($0.256\pm 0.007$)
2018-09-16 | CHARA/MIRC-X | S1-S2-E1-E2-W1-W2 | 1 |  HD 28855 ($0.303\pm 0.008$)
2018-09-17 | CHARA/MIRC-X | S1-S2-E1-W1-W2 | 2 |  HD 40280 ($0.599\pm 0.051$)
2018-10-26 | CHARA/MIRC-X | S1-S2-E1-E2-W1-W2 | 6 |  BD+31 600 ($0.391\pm 0.011$),
|  |  |  |  BD+44 1267 ($0.317\pm 0.008$),
|  |  |  |  BD+43 1350 ($0.318\pm 0.008$),
|  |  |  |  HD 28855 ($0.303\pm 0.008$)

111 All uniform disk (UD) diameters quoted obtained from Bourgés et al.
(2014).

Figure 1: Coverage of the uv plane of the interferometric MIRC-X observations
obtained with the CHARA array Figure 2: Visibilities and closure phases of the
image reconstruction. Black triangles with error bars are the original
calibrated observables (squared visibilities on the left and closure phases on
right), over plotted as blue markers are the model observables of the
reconstructed image. Below each plot is the fit residuals normalised by the
standard deviation as black circles.

Since the publication of LA19, additional relevant pieces of literature have
come to light. Spectroscopic and photometric monitoring of SU Aur by Petrov et
al. (2019) has revealed that a dusty disk wind is the potential source of the
photometric variability in both SU Aur and RY Tau at visible wavelengths. The
characteristic time of change in the disk wind outflow velocity and the
stellar brightness indicate that the obscuring dust is located close to the
sublimation rim of the disk, in agreement with previous theoretical disk wind
models (Bans & Königl 2012; Königl & Salmeron 2011). Recent ALMA and SPHERE
observations by Ginski et al. (2021) reveal a significant disk warp between
the inner and out disks of $\sim 70^{\circ}$. This misalignment is shown to
cause large shadows on the outer disk as it blocks light from the central
star. Their observations also suggest that SU Aur is currently undergoing a
late infall event with significant amounts of material falling inwards from
the outermost regions of the disk. Such events have the opportunity to
significantly impact the evolution of the disk.

This paper presents one of the first 6-telescope optical interferometric
studies of a YSO to date utilising state of the art observations covering a
wider range of baseline position angles and lengths (up to 331 m) (other
firsts include Kraus et al. (2020); Davies et al. (2022)). Three different
modelling methodologies were used to interpret our data and to provide direct
comparisons to LA19. (i) Image reconstruction was used to obtain a model-
independent representation of the data and to derive the basic object
morphology. (ii) Following this geometric model fitting allowed us to gain an
appreciation for the viewing geometry of the disk by fitting Gaussian and ring
models to the data. In addition, more complex geometric modelling was used to
explore the chromaticity of the data. (iii) Finally, we combine interferometry
and photometry to derive physical parameters with radiative transfer analysis,
where our focus is on confirming the presence of a dusty disk wind.

##  2 Observations

The CHARA array is a Y-shaped interferometric facility that comprises six
$1\,$m telescopes. It is located at the Mount Wilson Observatory, California,
and offers operational baselines between $34$ and $331\,$m (ten Brummelaar et
al. 2005). The MIRC-X instrument (Anugu et al. 2020; Kraus et al. 2018), a
six-telescope beam combiner, was used to obtain observations in the near-
infrared H-band ($\lambda=1.63\,\mu m,\Delta\lambda=0.35\,\mu m$) between
September and October 2018. We obtained 11 independent pointings of SU Aur,
using a mixture of 5 and 6-telescope configurations, due to the short delay
line limitations of CHARA. We obtained a maximum physical baseline of $331\,$m
corresponding to a resolution of $\lambda/(2B)=0.70\,\mathrm{mas}$
[milliarcseconds], where $\lambda$ is the observing wavelength and $B$ is the
projected baseline. Details of our observations, and the calibrator(s)
observed for the target during each observing session, are summarised in Table
1. The uv plane coverage that we achieved for the target is displayed in
Figure 1. Our data covers an exceptionally wide range of baseline lengths and
position angles, making the data ideally suited for image reconstruction.

The MIRC-X data were reduced using the standard python pipeline developed at
the University of Michigan by (J.B. le Bouquin, N. Anugu, T. Gardner). The
measured visibilities and closure phases were calibrated using interferometric
calibrator stars observed alongside the target. Their adopted uniform
diameters (UDs) were obtained from JMMC SearchCal (Bonneau et al. 2006, 2011),
and are listed in Table 1.

Considering the short timescale over which the observations were taken the
effect of time dependencies/variability of the object is thought to be
minimal. However, care was taken to check for time dependencies in the
visibilities of baselines of similar length and position angle. Variability in
the H band is known to be minimal, so any time dependencies in the visibility
amplitudes is likely geometric. However, no significant time dependencies were
discovered.

Figure 3: TOP LEFT: Image reconstruction resultant bootstrapped image,
including beam size and orientation. The coloured contours represent
significance flux levels of 1$\sigma$ (white), 3$\sigma$ (green) and 5$\sigma$
(blue)

. TOP RIGHT: Radiative transfer image produced using TORUS including a dusty
disk wind. BOTTOM MIDDLE: Reconstruction of simulated data created using the
best fit TORUS image above, reconstruction parameters are equivalent to those
of top left image. Colours are normalised intensity.

Figure 4: Flattened reconstructed image, created from angular slices though
the image from the central star. In the full image of Figure 3, north is up
and east is left. R is the radial distance from central star, $\Theta$ is the
polar angular direction.

##  3 Image Reconstruction

Image reconstruction techniques require broad and roughly circular uv coverage
along as many baseline lengths as possible. Fortunately, the data from these
observations lends itself to this process as the uv plane has been well
sampled, though some small gaps remain in the position angle coverage. This
technique is useful for interpretation of non-zero closure phases, indicative
of asymmetric distributions, in a model-independent way. Our closure phase
values are shown in Figure 2. There are many different algorithms with which
to reconstruct images from interferometric data, but the process described
here involved the use of the $SQUEEZE$ algorithm (Baron et al. 2010).
$SQUEEZE$ employs an MCMC approach to image reconstruction and was chosen due
to the wide range of available regularisation options and its ability to
implement $SPARCO$ a semi-parametric approach for image reconstruction of
chromatic objects (Kluska et al. 2014).

Table 2: . Best fit parameters for the simple geometric models investigated.

Parameter | Explored Parameter Space | Gaussian | Ring | Skewed Ring
---|---|---|---|---
$R$ [mas] | $0.0-10.0$ | – | $0.83\pm 0.01$ | $0.17\pm 0.16$
$FWHM$ [mas] | $0.0-15.0$ | $1.52\pm 0.01$ | – | $0.75\pm 0.04$
$INC$ [∘] | $0.0-90.0$ | $56.9\pm 0.4$ | $57.4\pm 0.4$ | $56.9\pm 0.5$
$PA$ [∘] | $0.0-360.0$ | $55.9\pm 0.5$ | $56.8\pm 0.4$ | $55.8\pm 0.5$
$f_{disk}$ | $0.0-1.0$ | $0.43\pm 0.01$ | $0.32\pm 0.01$ | $0.43\pm 0.01$
$\chi^{2}_{vis}$ |  | $11.63$ | $13.87$ | $11.62$
$\chi^{2}_{cp}$ |  | $6.05$ | $6.05$ | $6.01*$
222 (*) The closure phase quoted is the achieved when allowing the skewed ring
to become asymmetric. While the software did detect an asymmetry, it failed to
constrain its location in the disk. PA is the minor-axis position angle of the
disk and is measured from north ($PA=0^{\circ}$) towards east.

In the $SQUEEZE$/SPARCO routine, the object is modelled as an unresolved
central star with an extended, model-independent, environment (Kluska et al.
2014). Both components have different spectral behaviours and so differing
spectral indices. Additionally, the type and weight of the regularisation was
explored, $SQUEEZE$ allows for a very wide range of regularisation algorithms
to be implemented. The regularisation plays the role of the missing
information by promoting a certain type of morphology in the image. Total
variation (TV) was found to most reliably reproduce the best image, TV aims to
minimise the total flux gradient of the image and is useful to describe
uniform areas with steep but localised changes. These regularisations are
considered to be the best ones for optical interferometric image
reconstruction (Renard et al. 2011).

The size and number of pixels also plays an important role in image
reconstruction. One cannot simply use the maximum number of pixels of the
smallest size to obtain better resolution, they have to be chosen to match uv
plane sampling. It was found that a quadratic smoothing regularisation with a
weight of $1\times 10^{5}$ and $252\times 252$ pixels of $0.1$ mas in size
provides the best-fit image reconstruction when utilising exact Fourier
transform methods. The optimal regularisation parameters were determined using
the L-curve method.

The final image is shown in Figure 3 (top left panel) with the 1,3 and
5-$\sigma$ significance levels shown is white, green and blue contours
respectively. The inclination of the disk appears to be greater than that
found by LA19 with a similar minor-axis position angle. There also appears to
be a central bulge along the minor disk axis likely caused by the over
brightness of the star along this axis. The brightness distribution shows a
brighter structure along the north-west of the outer disk, parallel to the
major axis of the disk. This is consistent with the asymmetry found by LA19
and is indicative of a highly inclined disk where the far side of the inner
rim in directly exposed to the observer, while the nearside is obscured by
flaring in the outer disk. There are smaller significant structures to the
south-east of the disk also, we interpret these to be the shadowed near-side
of the rim due to their smaller extent than the northern features. We are not
confident in the exact shape of these objects given their irregularity. In
order to highlight the radial brightness distribution across the disk, Figure
4 shows a flattened profile, the elongation of the bulge can be seen in the NW
and SE directions with the extended rim material at radii out to 1 mas.

In order to quantitatively measure the size and orientation of the emitting
region, a simple ellipse was fitted to 3, 4, and 5-$\sigma$ flux significance
contours. The averaged results find an inclination of $46^{\circ}\pm 6$ and a
position angle of $53^{\circ}\pm 4$. Fitting an ellipse to lower significance
levels is not possible due to their irregular shape. the position angle of the
ellipses are in good agreement with the values derived in Section 4. The
inclination shows slight deviations from these values.

The chromaticity of the object is measured using two variables in the $SPARCO$
implementation. $f_{*}^{0}$, the stellar-to-total-flux ratio at the central
wavelength; and $d_{env}$, the spectral index of the extended environment.
Only the spectral index of the extended environment is needed as
interferometric data is only sensitive to the relative difference in spectral
index between the star en circumstellar material. $d_{env}$ is found be
$1.7\pm 0.8$ which corresponds to a temperature of $1257^{+234}_{-231}$
assuming the objects NIR emission is in a Rayleigh-Jean regime. This is within
the range of sublimation temperatures of typical disk astronomical silicates,
as expected at such small disk radii. $f_{*}^{0}$ is found to be $0.55\pm
0.6$, consistent with values measure from the SED and those found from
geometric modelling (see Section 4). The large errors associated with
$f_{*}^{0}$ and $d_{env}$ are a result of some degeneracy between the
parameters, see Kluska et al. (2014) for a full description of the procedure.

The visibility and closure phase fits of the image reconstruction are shown in
Figure 2 (top panels) along with the residuals of the fit (bottom panels). The
combined visibility and closure phase reduced chi-squared $\chi^{2}_{red}$ of
the image reconstruction was found to be $4.38$.

##  4 Geometric Modelling

In order to understand the geometry of the system one must consider the
application of simple geometric models. In this section we explore several
different approaches to modelling our data with both a non-chromatic ’grey’
models and techniques which explore the chromaticity.

Figure 5: Geometric model images, corresponding to best fit parameters
described in Table 2. TOP LEFT: Gaussian model brightness distribution. TOP
RIGHT: Ring model brightness distribution. BOTTOM LEFT: Skewed Ring model
brightness distribution. BOTTOM RIGHT: Full reconstructed image, repeated from
Figure 3.

###  4.1 Basic Geometric Models

The fitting of Gaussian and ring like distributions to the interferometric
variables allows highly accurate estimations of the characteristic size,
inclination and position angle of the object. In all models the central star
is modelled as a point source, which is an acceptable assumption given the
expected angular diameter of the star. The disk parameters are then fitted in
the RAPIDO (Radiative transfer and Analytic modelling Pipeline for
Interferometric Disk Observations) framework (Kreplin et al. 2018), available
in-house at University of Exeter. RAPIDO utilises the Markov chain Monte Carlo
(MCMC) sampler _emcee_ to produce a fit and error estimate (Foreman-Mackey
2016). Three disk models were employed, a standard Gaussian brightness
distribution which is characterised by it’s full-width-half-maximum (FWHM).
Along with two ring models, a sharp ring with a width fixed to 20% of the disk
radius ($R$) and a ’skewed’ ring with a more diffuse radial profile produced
by convolving with a Gaussian with a $FWHM$. The skewed ring is also capable
of modelling azimuthal modulation or disk asymmetries, a detailed description
of this model can be found in Lazareff et al. (2017). In addition to the model
specific parameters, we also fitted the inclination ($INC$), minor-axis
position angle ($PA$) and disk-to-total flux ratio ($f_{disk}$). As we see no
evidence of time variability in the data we are able to fit all data
simultaneously. The results from the simple geometric model fitting are shown
in Table 2.

Out of the geometric models tested, the Gaussian model is considered to be the
best fit. Even though the skewed ring produced a slightly small $\chi^{2}$
value for the closure phase and visibility measurements, we do not consider
this significant given the additional degrees of freedom in the model. In
addition, as shown in Figure 5, the best fit skewed ring is tending towards a
Gaussian distribution, which a very small inner radius $0.17\,\mathrm{mas}$
with a wide FWHM ring width $0.75\,\mathrm{mas}$. The skewed ring fails to
reproduce the small but complex asymmetries seen in the image reconstruction.

Figure 6: Spectral energy distribution of SU Aurigae. Green points are
photometric data from a variety of instruments. Red line is Spitzer IR data.
Black dashed line is direct radiation from the stellar photosphere. Blue line
is the best TORUS computed radiative transfer model inclined at $56^{\circ}$.
Orange line is the SED computed from the simple temperature gradient models
described in Section 4.2. Figure 7: Temperature gradient profile across the
inner au of SU Aur. The dusty disk extends down an inner radius, consistent
with the expected sublimation radius and temperature. Interior to this
expected a hot dust free inner disk which is not detected in the continuum
observations.

The best fit Gaussian model finds a disk of FWHM $1.52\pm 0.01\,\mathrm{mas}$
which is inclined at $56.9\pm 0.4^{\circ}$ and a minor-axis position angle of
$55.9\pm 0.5^{\circ}$. In addition, we find that $43\pm 1\%$ of the total flux
originates from the disk in the H band. This is consistent with measurements
based on the infrared excess of the spectral energy distribution (SED) in LA19
and also with the flux ratio found in the image reconstruction algorithm.

The primary limitation of the simple geometric models described above is that
they are intrinsically ’grey’ in nature. Meaning they contain no spectral
information, hence the large $\chi^{2}$ values obtained in the fitting
process. In order to better model the spectral dependency of the visibility,
more complex temperature gradient models that are able to account for
observing wavelength must be employed.

###  4.2 Temperature Gradient Models

A physically correct model can be applied by considering the temperature
gradient of the disk. A temperature gradient model (TGM) allows for the
simultaneous fitting of interferometric and photometric observables. The
origin of the photometric data used is describe in Table 4 in Appendix A. It
is built up by several rings extending from an inner radius $R_{\mathrm{in}}$
to an outer radius $R_{\mathrm{out}}$. Each ring is associated with
temperature and hence flux. Therefore, a model SED can be computed by
integrating over the resulting blackbody distributions for each of the
concentric rings. Such a model allows us to not only to build up a picture of
the temperature profile, but also approximate the position of the inner
radius. The TGM is based upon a $T_{R}=T_{0}(R/R_{0})^{-Q}$ profile where
$T_{0}$ is the temperature at the inner radius of the disk $R_{0}$, and $Q$ is
the exponent of the temperature gradient (Kreplin et al. 2020; Eisner &
Hillenbrand 2011). A TGM represents an intrinsically geometrically thin disk.
A point source is used at the centre of each model to represent an unresolved
star, which is a reasonable approximation given the expected angular diameter
of $0.05\,\mathrm{mas}$ (Pérez et al. 2020). Also included in the fitting of
photometric parameters is a treatment of interstellar extinction based on
(Fitzpatrick 1999) with an $E_{(B-V)}=0.5$ (Bertout et al. 2007).

The inclination and position angle of the disk are maintained at fixed values
of $56.9\pm 0.4^{\circ}$ and $55.9\pm 0.5^{\circ}$ respectively, from the
fitting of the Gaussian distribution. This was done to reduce the number of
free parameters in the model. The fitting was undertaken using all of
visibility data shown in Figure 3 and all the SED points simultaneously. The
fitting and error computation was once again done using _emcee_ (Foreman-
Mackey 2016).

The results of the temperature gradient modelling are shown in Figure 7. We
find an inner disk radius of $0.15\pm 0.04\,\mathrm{au}$ where the temperature
is equivalent to $2100\pm 200\,\mathrm{K}$ and decreases with an exponent of
$Q=0.62\pm 0.02$. The inner disk radius is the point at which the dusty disk
is truncated due to the sublimation of material, in contrast to more extreme
objects such as FU Ori where the inner disk radius is equivalent to that of
the stellar radius indicating boundary layer accretion (Labdon et al. 2021).
Interior to the sublimation radius is expected a hot dust free inner disk from
which material can be magnetospherically accreted. However, our low spectral
resolution continuum observations are not sensitive to these regions.

An inner disk radii temperature of $2100\pm 200\,\mathrm{K}$ may be considered
high, but is broadly consistent with laboratory sublimation temperatures for
silicate grains. An exponent of $Q=0.62\pm 0.02$ is slightly larger than that
expected from a disk heated by stellar radiation alone, and may indicate the
presence of additional disk heating mechanisms, such as viscous heating
(Pringle 1981; Kenyon & Hartmann 1987; Dullemond & Dominik 2004).

The resultant SED of the inner disk is shown in Figure 6 as the orange curve.
Beyond 5-6 $\mu m$ the TGM fails to fit the shape of the disk accurately. This
is likely due to the flared nature of the disk in contrast to the ’flat’ TGM
model which will have the strongest effect and longer wavelengths. In
addition, the strong disk warp reported by Ginski et al. (2021) would also
introduce some temperature discontinuity, the modelling of which is beyond the
scope of this work.

##  5 Radiative Transfer

In order to provide a more physical model and directly compare to LA19, we
used the TORUS Monte-Carlo radiative transfer code (Harries et al. 2019) to
simultaneously fit the visibility, closure phase and photometric data of the
SU Aurigae system.

The models adopted here are based on the disk models used by LA19, adapted to
account for the higher inclination and different observing wavelength. In
these TORUS simulations, the dust was allowed to vertically settle to the
scale height of the gas component and the dust sublimation radius was left as
a free parameter, allowing the inner rim radius to define itself based on
well-defined rules of the Lucy (1999) iterative method to determine the
location and the temperature structure of the whole disk. This is implemented
whereby the temperature is initially calculated for grid cells in an optically
thin disk structure, with dust added iteratively to each cell with a
temperature lower than that of sublimation, until the appropriate dust to gas
ratio is reached ($0.01$). We confirmed that stellar photosphere models of
Castelli & Kurucz (2004) using these stellar parameters can reproduce the
photometry measurements of SU Aur reasonably well across the visible
continuum. We adopt a silicate grain species with dust properties and
opacities adopted from Draine (2003). For a more detailed description of TORUS
and the algorithms used, see Davies et al. (2018); Labdon et al. (2019).

Table 3: Best-fit parameters resulting from TORUS (Harries 2000) radiative
transfer SED and visibility modelling.

Parameter | Literature value | Reference | Range explored | Best fit value | Notes
---|---|---|---|---|---
$R_{\mathrm{{inner}}}$ | $0.18\leavevmode\nobreak\ \mathrm{au}$ | (1)(2) | $0.1-0.6\leavevmode\nobreak\ \mathrm{au}$ | $0.16\leavevmode\nobreak\ \mathrm{au}$ | Determined based on the gas density and grain size
$R_{\mathrm{{outer}}}$ | $100\leavevmode\nobreak\ \mathrm{au}$ | (1)(3) | - | $100.0\leavevmode\nobreak\ \mathrm{au}$ | Fixed to literature values
$h_{\mathrm{0}}$ | $15.0\leavevmode\nobreak\ \mathrm{au}$ | (3) | $7.0-20.0\leavevmode\nobreak\ \mathrm{au}$ | $9.0\leavevmode\nobreak\ \mathrm{au}$ | Scale height at 100 au
$\alpha_{\mathrm{{disk}}}$ | $2.4$ | (3) | $1.0-3.0$ | $2.3$ | $\mathrm{\alpha_{disk}\,\,is\,\,fixed\,\,at\,\,(\beta_{disk}+1)}$
$\beta_{\mathrm{{disk}}}$ | $1.4$ | (3) | $0.0-2.0$ | $1.3$ |
$\mathrm{Dust:Gas}$ | $0.01$ | (1)(3) | - | $0.01$ | Fixed to literature values
$a$ | $0.14\,\mathrm{\mu m}$ | (1) | $0.1-1.4\leavevmode\nobreak\ \mathrm{\mu m}$ | $0.14\leavevmode\nobreak\ \mathrm{\mu m}$ |
$T_{\mathrm{{sub}}}$ | $2000\leavevmode\nobreak\ \mathrm{K}$ | (4) | - | $2000\leavevmode\nobreak\ \mathrm{K}$ | Fixed to literature values
Dusty disk wind parameter | Literature value | Reference | Range explored | Best fit value |
$R_{\mathrm{{0min}}}$ | $4.5$ | (3) | 2.0-6.0 | $4.5\leavevmode\nobreak\ \mathrm{R_{\sun}}$ | New grid of models
$T_{\mathrm{{wind}\leavevmode\nobreak\ (near\leavevmode\nobreak\ surface)}}$ | $1600\leavevmode\nobreak\ \mathrm{K}$ | (3) | 1200-2400 | $1600\leavevmode\nobreak\ \mathrm{K}$ | New grid of models
$\mathrm{Opening\leavevmode\nobreak\ Angle}$ | $50^{\circ}$ | (3) | $25-65^{\circ}$ | $45^{\circ}$ | New grid of models
$\dot{M}$ | $10^{-7}\leavevmode\nobreak\ M_{\odot}yr^{-1}$ | (3) |  $10^{-10}$ \- $10^{-6}$ | $10^{-7}\leavevmode\nobreak\ M_{\odot}yr^{-1}$ | New grid of models

(1) Akeson et al. (2005); (2) Jeffers et al. (2014); (3) Labdon et al. (2019);
(4) Pollack et al. (1994)

The dusty disk wind model is adapted from Bans & Königl (2012). This mechanism
is based on the presence of a large-scale, ordered magnetic field which
threads the disk along which disk material is flung out. The high magnetic
pressure gradient above the disk surface accelerates the material which is
then collimated through the azimuthal and poloidal field components (Bans &
Königl 2012). These centrifugally driven winds are highly efficient at
distributing density above and below the plane of the disk, carrying angular
momentum away from the disk surface. A full description of the implementation
within the TORUS radiative transfer code can be found in LA19.

The disk model adopted follows the curved inner rim prescription of Isella &
Natta (2005) with a density dependent sublimation radius whereby grains
located in the disk midplane are better shielded due to higher densities and
so can exist closer to the central star than grains in the less dense upper
layers. A full summary of the disk parameters can be found in Table 3. The
only parameters we explored with respect to LA19 were the disk scale height
($h_{0}$), the grain size ($a$) and flaring index $\beta_{disk}$. The key
difference in the models described here compared to LA19 is the grain size
adopted. Here we adopt a smaller grain size of $0.14\,\mu m$, which in turn
leads to a slightly smaller inner radius of $0.16\,\mathrm{au}$. Additionally,
in order to improve the SED fit at longer wavelengths we also adopt a more
modest scale height of $9.0\,\mathrm{au}$. The resultant SED from the
radiative transfer model is shown in figure 6, in addition to the model
stellar photosphere, also calculated within TORUS. The SED is fit well across
the optical and IR and longer mm wavelengths, however is a relatively poor fit
across the 8-40 $\mu m$ range. We attribute this to the disk warp reported in
Ginski et al. (2021) which would result in a physical disk break and
temperature discontinuity. The geometry governing disk warps is little
understood, in particular how it would effect the inner and outer edge of the
warp. The shapes of these rims, and how is the vertical structure of the disk
effected at this point, are not known. Such a model are well beyond the scope
of this paper for 2 reasons: Firstly we lack high spatial resolution data at
M/F-IR wavelengths which might cover the location of the disk warp. The
complex geometry of the warp would require extensive modelling, which is
difficult given the limitations of the radiative transfer code used.

Figure 8: Coloured circles are the MIRC-X squared visibilitites, where the
colour represents the wavelength across the H-band. Black crosses are the
corresponding squared visibilities extracted from TORUS radiative transfer
model images. LEFT: Best fit TORUS model including a strong dusty disk wind
component. RIGHT: The same base disk model, excluding any dusty disk wind
component. Lower panels show the residuals in the fits between data and
models.

The presence of a disk wind is once again required to fit both the
visibilities and the SED. The absence of a disk wind fails to reproduce IR
excess across both the H and K bands, with insufficient NIR disk flux. This
significantly impacts both the SED and interferometric fit. A disk wind is
required to eject more hot dust above the midplane of the disk where it is
directly exposed to stellar radiation which is reprocessed as an IR excess.

This is shown in Figure 8, where the squared visibilities are shown for both a
disk model with and without the dusty wind environment. The disk wind model
provides a far superior fit to the observations, being able to successfully
reproduce the NIR excess. Beginning from the wind parameters found in LA19, we
run a new grid of models exploring the disk wind parameter space. The
parameters explored and their ranges are described in Table 3, along with the
resulting best fit values. The results are broadly similar to those found in
LA19, with exception of a slightly higher temperature of material in the disk
wind, closer to dust temperatures found in the inner disk. Of particular note
is the into-wind accretion rate of $10^{-7}\,\odot{M}$. Considering the
historically accepted on-to-star accretion rate to into-wind accretion ratio
of $0.1$, this level of transport is perhaps unphysically high given the age
of the star, but a discussion of possible mechanisms is included in Section 6.

The final computed image is shown in Figure 3 (middle left panel) and shows
the clear asymmetry originating from the inclination in the asymmetry map in
the same figure (middle left panel). In order to approximate what this
computed image would look like if observed at the same resolution as the
original observations, we computed synthetic visibility and closure phases
based on the radiative transfer images (as described in (Davies et al. 2018)).
Artificial noise and error bars were computed to be representative of the
original data and to ensure an accurate representation. These synthetic
observables were then reconstructed in the same manner as the original data,
as described in Section 3. Care was taken to ensure the consistency in the
reconstruction parameters for both the real and synthetic observables. The
reconstructed TORUS image is also shown in Figure 3 (bottom right panel), and
shows clear similarities with both the original TORUS image and the image
reconstructed from the original data.

##  6 Discussion

Our extensive observations and analysis of the circumstellar environment of SU
Aurigae have revealed the details of the inner disk in unprecedented detail.
The wide variety of techniques used to analyse our interferometric data allows
us to precisely define the disk characteristics.

Image reconstruction is a crucial, model independent, method of analysis which
is ideally suited to our dataset with extensive uv and baseline coverage. Our
analysis reveals an elliptical shape, indicative of an object with a high
inclination as shown in Figure 3. There appears to be a central bulge to the
disk, however this feature is not thought to be physical but rather a
manifestation of the brightness of the central star combined with the width of
the disk at this point. The more extended material in the image are thought to
be a depiction of the far-side of the disk rim which is un-obscured by the
outer disk, this re-enforced by the fact that the extended material is bright
in the north-east of the image. There is a significant asymmetric feature in
the form of a thin brightness on the north-eastern edge of the disk. The
unique shape of this feature indicates that this is again caused by the high
inclination obscuring the nearside disk rim. The effect of an inclined disk on
the observed brightness distribution is described extensively by Jang-Condell
& Turner (2013) and accurately describes the observations here. By fitting an
ellipse to the significant features in the image, we gain a quantitative
measure of the disk geometry. The results find an inclination of
$46^{\circ}\pm 6$ and a position angle of $53^{\circ}\pm 4$, very similar to
those derived from more reliable geometric model fits described below. The
scale and shape of the image is similar to that of LA19, with a slightly
higher inclined viewing angle. We consider the images of this work to be the
more accurate depiction of SU Aur given the higher quality observations, taken
over a much shorter timescale, with the added detail and resolution this
entails.

Although geometric modelling is much more constrained in the geometries it can
explore, in does provide a more quantitative view of the disk. It was found
that the model which best fit our data was a simple Gaussian distribution with
a point source representing the star. A Gaussian model is consistent with
other work, both on this object by LA19, but also in other YSO studies such as
the survey by Lazareff et al. (2017) who find that little under half of their
51 objects can be modelled by a Gaussian structure. The Gaussian fitted in
this work has a FWHM of $1.52\pm 0.01\,\mathrm{mas}$ ($0.239\pm
0.002\,\mathrm{au}$) at an inclination of $56.9^{\circ}\pm 0.4$ at a minor
axis position angle of $55.9^{\circ}\pm 0.5$ with a stellar-to-total flux
ratio of $0.57\pm 0.01$. The reduced $\chi^{2}$ value for the visibilities is
$11.63$ and $6.05$ for the closure phases, which are equivalent to $0^{\circ}$
for this centro-symmetric model. These values are in agreement with the
literature values of Akeson et al. (2005) who find a K band radius of $0.18\pm
0.04$ au and an inclination of $62^{\circ+4}_{-8}$. Similar values for the
inclination in literature are ~$60^{\circ}$ and ~$50^{\circ}$ found by Unruh
et al. (2004); Jeffers et al. (2014) respectively. The minor axis position
angle derived here is significantly greater than the literature values of
$24^{\circ}\pm 23$ and $15^{\circ}\pm 5$ found by Akeson et al. (2005);
Jeffers et al. (2014). This difference is likely due to either: The poor uv
coverage and lack of longer baselines in previous interferometric studies,
both of which make estimating the position angle and inclination particularly
unreliable. Other non-interferometric studies focus on the outer disk which is
shown by Ginski et al. (2021) to be misaligned compared to the inner disk.

The geometric modelling results are broadly similar to those presented in our
previous work LA19 where an inclination of $50.9\pm 1.0^{\circ}$ and minor
axis position angle of $60.8\pm 1.2^{\circ}$ were found and the data were
marginally better described by a ring-like brightness distribution. The values
and models presented here are considered to be more accurate due high
precision observations and significantly smaller potential for temporal
variations, as the data of LA19 was coalesced over 14 years. These values are
also consistent with observations of the outer disk by Ginski et al. (2021)
where dark shadows are observed in scattered light origination from a
significant disk warp between the inner and outer regions. On larger scales
the near-side of disk is seen to the north-east, in our observations it is
seen to the south west.

Figure 9: Black crosses are the closure phase data obtained with MIRC-X.
Overlaid as coloured points are the TORUS model closure phases, the colours
represent the wavelength of the spectral channels and follows the same
convention as other plots in this work. Below in black points are the
normalised residual errors of the fit.

In modelling the temperature gradient of SU Aur we can gain an appreciation
for spectral dependence of our interferometric variables across the 6 spectral
channels of MIRC-X. Our modelling finds a disk which extends down to $0.15\pm
0.04\mathrm{au}$ where the temperature is equivalent to $2100\pm 200\,K$ and
decreases with an exponent of $Q=0.62\pm 0.02$. The outer edge of this
temperature regime was found to be $0.20\pm 0.03\mathrm{au}$, showing this
prescription only covers the innermost regions of the disk. Modelling outer
regions of the disk is beyond the scope of this paper, as our NIR
interferometric data does not cover emission from these regions and larger
radii are likely effected by the significant disk warp causing a temperature
gradient discontinuity. A temperature gradient exponent of $Q=0.62\pm 0.02$
found here lies between two established models from literature. That of
Pringle (1981) who find that a steady state, optically-thick accretion disk
heated by viscous processes will exhibit an exponent of $0.75$ and of Kenyon &
Hartmann (1987); Dullemond & Dominik (2004) who show that a flared disk heated
by reprocessed stellar radiation alone will exhibit an exponent of $\leq 0.5$.
As such, it is difficult to comment on the heating of the circumstellar
environment of SU Aurigae. But it is possible that the disk is not heated by
stellar radiation alone, but additional heating processes such as viscous
heating may also be present.

The radiative transfer modelling presented in this paper is heavily based on
LA19. For a detailed discussion of the motivation behind certain choices,
particularly in relation to the shape of the inner rim, we recommend the
reader see Section 5 of that paper. In this work we are able to achieve a
similar SED fit to LA19, including with the adoption of a smaller
$0.14\,\mathrm{\mu m}$ grain size, it is noted that the smaller grain size is
in line with older radiative transfer work of SU Aur by Akeson et al. (2005).
The smaller grain size results in a smaller inner rim which now extends down
to $0.13\,\mathrm{au}$ at a sublimation temperature of $2000\,\mathrm{K}$.
This is within the uncertainties of values predicted by the temperature
gradient modelling and is roughly consistent with older literature values of
$0.18\pm 0.04\,\mathrm{au}$ and $0.17\pm 0.08\,\mathrm{au}$ by Akeson et al.
(2005) and Jeffers et al. (2014) respectively. The flaring parameters
$\alpha_{\mathrm{disk}}$ and $\beta_{\mathrm{disk}}$ were fixed such that
$\alpha_{\mathrm{disk}}=\beta_{\mathrm{disk}}+1$ and found to be $2.3$ and
$1.3$ respectively, a more physical representation, than the values depicted
in LA19. Similarly to LA19, a dusty disk wind is required in order to fit both
the SED across the NIR and visibilities as shown in Figure 8.

The TORUS implementation of s dusty wind does not depend on the underlying
launching mechanism, but just prescribes a geometry above and below the disk,
which is populated by dust grains, where it can reprocess stellar radiation to
contribute to the NIR excess. This was first put forward by Bans & Königl
(2012) in the context of a magnetospherically driven disk wind. Figure 8
highlights how a standard (no-wind) model cannot sufficiently fir the
interferometric data, and how the addition of a dusty disk wind can provide a
significantly better fit to the data. The maximum temperature of the dust in
the wind is similar to the temperature of the dust at the sublimation
temperature of the disk $1900\,\mathrm{K}$ and $2000\,\mathrm{K}$
respectively, as expected given the dust launches from close to to the
sublimation rim. In Appendix B each baseline is plotted in a separate panel in
order to explore the chromatic (temperature) gradient of the data. The
discrepancy in the gradient and level of some baselines can be explained by
the simplified heating description within the radiative transfer model. TORUS
produces a disk heated by reprocessed stellar radiation alone, with no
internal disk heating such as viscous heating. We see from the specific
temperature gradient modelling that we might expect some viscous heating
within the disk.

However, the implementation of the dusty disk wind in this scenario is not
completely physical, owing to the high into-wind outflow rate of $1\times
10^{-7}\mathrm{M_{\odot}yr^{-1}}$ required. If one assumes the historically
accepted outflow to accretion ratio of $0.1$, the resulting onto-star-
accretion rate is greater than those typically found in T Tauri stars. In
particular, this contrasts with the measured accretion rate of SU Aur (Pérez
et al. 2020). However, the inflow to outflow ratio is subject of some
discussion with recent works suggesting the ratio may be closer to unity in
some situations (Pascucci et al. 2022), particularly by invoking a
magnetically driven winds originating from a dead-zone close to the
sublimation radius. Based on our data we cannot differentiate such specific
wind launching mechanisms. In addition, the suggestion that SU Aur is
undergoing a late stage infall event Ginski et al. (2021) could be a potential
explanation for such a high level of mass transport through the system.

Late infall events typically occur after the depletion of the protostellar
envelope, when a large amount of material exterior to the disk falls inwards.
The triggering mechanisms behind such events are not well understood as there
are very few observational examples. From simulations it is thought that a
large infall event can not only directly increase the accretion rate(Dullemond
et al. 2019), but also dramatically re-sculpt the disk, creating gaps, spiral
arms and strong misalignments (Kuffmeier et al. 2021). In the case of SU Aur a
large infall of material is thought to have caused a strong disk warp between
the inner and outer disk (Ginski et al. 2021). It is also not unreasonable to
assume this lead to an enhanced accretion rate, which would allow for the
strong disks winds needed in the model presented here.

The simplest way to directly compare the analytical methods employed in this
work, is through the images produced. Figure 3 shows a collage of the images
produced by the geometric modelling, image reconstruction and radiative
transfer methods. It is clearly shown the similarities between the images all
with very similar inclinations and position angles with asymmetries oriented
in the same direction. The asymmetries all appear to be the result of the high
inclination of the disk causing the near-side inner rim to be shadowed by the
flared outer disk. The lowest panel of Figure 3 shows an image reconstructed
from synthetic visibilties and closure phases obtained from the radiative
transfer image. Care was taken to ensure the observables were matched in
baseline length and position angle and detector mode, while the image
reconstruction process was exactly that described in Section 3. The remarkable
similarity between the different images, in particular the real image
reconstruction and the synthetic reconstruction, provides strong reinforcement
that the models adopted are truly representative of the disk structure.

###  6.1 Future Observations

Recent improvements at CHARA include the commissioning of the MYSTIC (Michigan
Young STar Imager at CHARA) beam combiner (Monnier et al. 2018), a 6-telescope
K-band combiner capable of working in parallel with MIRC-X for simultaneous
observations. The combination of MIRC-X and MYSTIC would allow us to
investigate the temperature structure of the inner disk in greater detail,
while also obtaining greater precision due to the higher sensitivity of
MYSTIC. In addition, longer wavelength MIR interferometry with the MATTISE
instrument at VLTI would potentially allow us to resolve the warped region of
the disk, in order to confirm the nature of the geometry in this region for
direct comparison to late infall event simulations.

##  7 Conclusions

This interferometric study of SU Aurigae has revealed the complex geometry and
composition of the disk around SU Aurigae. We summarise our conclusions as
follows:

  * •

We reconstruct an interferometric image that confirms the inclined disk
described in literature. We see evidence for an asymmetry in the brightness
distribution that can be explained by the exposure of the inner-rim on the far
side of the disk and its obscuration on the near side due to inclination
effects.

  * •

Our simple geometric model fits show that the circumstellar environment is
best modelled as a Gaussian distribution with a disk of inclination $56.9\pm
0.4^{\circ}$ along a minor axis position angle of $55.9.0\pm 0.5^{\circ}$ and
an FWHM of $1.52\pm 0.01\,\mathrm{mas}$ ($0.239\pm 0.002\,\mathrm{au}$). Such
geometry is consistent with strong disk shadows observed in the outer disk
originating from a disk warp.

  * •

We model the radial temperature profile of the inner disk and find a disk
which extends down to $0.15\pm 0.04\mathrm{au}$ where the temperature is
equivalent to $2100\pm 200\,K$ and decreases with an exponent of $0.62\pm
0.02$.

  * •

A dusty disk wind scenario is still required to account for both the observed
excess in the SED and the observed visibilities. The dusty disk wind scenario
described here lifts material above the disk photosphere, thus exposing more
dust grains to the higher temperatures close to the star responsible for the
NIR excess. The high accretion rate required to reproduce the stellar-to-total
flux ratio could be explained by a late infall event.

  * •

Our best-fit model (dusty disk wind model) suggests that the dust composition
in the disk is dominated by small grains ($0.14\,\mathrm{\mu m}$) with a
sublimation temperature of $2000$ K. Introducing larger grains results in a
worse fit to the SED shape and NIR excess. The disk is also shown to be highly
flared ($9$ au at $100$ au).

###### Acknowledgements.

We acknowledge support from an STFC studentship (No. 630008203) and an
European Research Council Starting Grant (Grant Agreement No. 639889). MIRC-X
received funding from the European Research Council (ERC) under the European
Union’s Horizon 2020 research and innovation programme (Grant No. 639889). JDM
acknowledges funding for the development of MIRC-X (NASA-XRP NNX16AD43G, NSF-
AST 1909165) and MYSTIC (NSF-ATI 1506540, NSF-AST 1909165). This research has
made use of the VizieR catalogue access tool, CDS, Strasbourg, France. The
original description of the VizieR service was published in Ochsenbein et al.
(2000).

## References

  * Abrahamyan et al. (2015) Abrahamyan, H. V., Mickaelian, A. M., & Knyazyan, A. V. 2015, Astronomy and Computing, 10, 99
  * Akeson et al. (2005) Akeson, R. L., Walker, C. H., Wood, K., et al. 2005, ApJ, 622, 440
  * Ammons et al. (2006) Ammons, S. M., Robinson, S. E., Strader, J., et al. 2006, ApJ, 638, 1004
  * Anderson & Francis (2012) Anderson, E. & Francis, C. 2012, Astronomy Letters, 38, 331
  * Andrews et al. (2013) Andrews, S. M., Rosenfeld, K. A., Kraus, A. L., & Wilner, D. J. 2013, ApJ, 771, 129
  * Anugu et al. (2020) Anugu, N., Le Bouquin, J.-B., Monnier, J. D., et al. 2020, arXiv e-prints, arXiv:2007.12320
  * Bans & Königl (2012) Bans, A. & Königl, A. 2012, ApJ, 758, 100
  * Baron et al. (2010) Baron, F., Monnier, J. D., & Kloppenborg, B. 2010, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 7734, Optical and Infrared Interferometry II, ed. W. C. Danchi, F. Delplancke, & J. K. Rajagopal, 77342I
  * Bertout et al. (2007) Bertout, C., Siess, L., & Cabrit, S. 2007, A&A, 473, L21
  * Bianchi et al. (2011) Bianchi, L., Herald, J., Efremova, B., et al. 2011, Ap&SS, 335, 161
  * Bonneau et al. (2006) Bonneau, D., Clausse, J.-M., Delfosse, X., et al. 2006, A&A, 456, 789
  * Bonneau et al. (2011) Bonneau, D., Delfosse, X., Mourard, D., et al. 2011, A&A, 535, A53
  * Bourgés et al. (2014) Bourgés, L., Lafrasse, S., Mella, G., et al. 2014, in Astronomical Data Analysis Software and Systems XXIII. Proceedings of a meeting held 29 September - 3 October 2013 at Waikoloa Beach Marriott, Hawaii, USA. Edited by N. Manset and P. Forshay ASP conference series, vol. 485, 2014, p.223, Vol. 485, 223
  * Castelli & Kurucz (2004) Castelli, F. & Kurucz, R. L. 2004, ArXiv Astrophysics e-prints
  * Cutri et al. (2014) Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2014, VizieR Online Data Catalog, II/328
  * Davies et al. (2014) Davies, C. L., Gregory, S. G., & Greaves, J. S. 2014, MNRAS, 444, 1157
  * Davies et al. (2018) Davies, C. L., Kraus, S., Harries, T. J., et al. 2018, ApJ, 866, 23
  * Davies et al. (2022) Davies, C. L., Rich, E. A., Harries, T. J., et al. 2022, MNRAS, 511, 2434
  * Draine (2003) Draine, B. T. 2003, ARA&A, 41, 241
  * Dullemond & Dominik (2004) Dullemond, C. P. & Dominik, C. 2004, A&A, 417, 159
  * Dullemond et al. (2019) Dullemond, C. P., Küffmeier, M., Goicovic, F., et al. 2019, A&A, 628, A20
  * Eisner & Hillenbrand (2011) Eisner, J. A. & Hillenbrand, L. A. 2011, ApJ, 738, 9
  * Esplin et al. (2014) Esplin, T. L., Luhman, K. L., & Mamajek, E. E. 2014, ApJ, 784, 126
  * Fitzpatrick (1999) Fitzpatrick, E. L. 1999, PASP, 111, 63
  * Foreman-Mackey (2016) Foreman-Mackey, D. 2016, The Journal of Open Source Software, 1, 24
  * Ginski et al. (2021) Ginski, C., Facchini, S., Huang, J., et al. 2021, ApJ, 908, L25
  * Harries (2000) Harries, T. J. 2000, MNRAS, 315, 722
  * Harries et al. (2019) Harries, T. J., Haworth, T. J., Acreman, D., Ali, A., & Douglas, T. 2019, Astronomy and Computing, 27, 63
  * Isella & Natta (2005) Isella, A. & Natta, A. 2005, A&A, 438, 899
  * Jang-Condell & Turner (2013) Jang-Condell, H. & Turner, N. J. 2013, ApJ, 772, 34
  * Jeffers et al. (2014) Jeffers, S. V., Min, M., Canovas, H., Rodenhuis, M., & Keller, C. U. 2014, A&A, 561, A23
  * Kenyon & Hartmann (1987) Kenyon, S. J. & Hartmann, L. 1987, ApJ, 323, 714
  * Kluska et al. (2014) Kluska, J., Malbet, F., Berger, J.-P., et al. 2014, A&A, 564, A80
  * Königl & Salmeron (2011) Königl, A. & Salmeron, R. 2011, The Effects of Large-Scale Magnetic Fields on Disk Formation and Evolution (University of Chicago Press), 283–352
  * Kraus et al. (2020) Kraus, S., Kreplin, A., Young, A. K., et al. 2020, arXiv e-prints, arXiv:2004.01204
  * Kraus et al. (2018) Kraus, S., Monnier, J. D., Anugu, N., et al. 2018, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 10701, Optical and Infrared Interferometry and Imaging VI, 1070123
  * Kreplin et al. (2020) Kreplin, A., Kraus, S., Tambovtseva, L., Grinin, V., & Hone, E. 2020, MNRAS, 492, 566
  * Kreplin et al. (2018) Kreplin, A., Tambovtseva, L., Grinin, V., et al. 2018, MNRAS, 476, 4520
  * Kuffmeier et al. (2021) Kuffmeier, M., Dullemond, C. P., Reissl, S., & Goicovic, F. G. 2021, A&A, 656, A161
  * Labdon et al. (2019) Labdon, A., Kraus, S., Davies, C. L., et al. 2019, A&A, 627, A36
  * Labdon et al. (2021) Labdon, A., Kraus, S., Davies, C. L., et al. 2021, A&A, 646, A102
  * Lazareff et al. (2017) Lazareff, B., Berger, J.-P., Kluska, J., et al. 2017, A&A, 599, A85
  * Lucy (1999) Lucy, L. B. 1999, A&A, 344, 282
  * Mohanty et al. (2013) Mohanty, S., Greaves, J., Mortlock, D., et al. 2013, ApJ, 773, 168
  * Monnier et al. (2018) Monnier, J. D., Le Bouquin, J.-B., Anugu, N., et al. 2018, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 10701, Optical and Infrared Interferometry and Imaging VI, ed. M. J. Creech-Eakman, P. G. Tuthill, & A. Mérand, 1070122
  * Morel & Magnenat (1978) Morel, M. & Magnenat, P. 1978, Astronomy and Astrophysics Supplement Series, 34, 477
  * Ochsenbein et al. (2000) Ochsenbein, F., Bauer, P., & Marcout, J. 2000, A&AS, 143, 23
  * Ofek (2008) Ofek, E. O. 2008, Publications of the Astronomical Society of the Pacific, 120, 1128
  * Pascucci et al. (2022) Pascucci, I., Cabrit, S., Edwards, S., et al. 2022, arXiv e-prints, arXiv:2203.10068
  * Pérez et al. (2020) Pérez, S., Hales, A., Liu, H. B., et al. 2020, ApJ, 889, 59
  * Petrov et al. (2019) Petrov, P. P., Grankin, K. N., Gameiro, J. F., et al. 2019, MNRAS, 483, 132
  * Pollack et al. (1994) Pollack, J. B., Hollenbach, D., Beckwith, S., et al. 1994, ApJ, 421, 615
  * Pringle (1981) Pringle, J. E. 1981, ARA&A, 19, 137
  * Renard et al. (2011) Renard, S., Thiébaut, E., & Malbet, F. 2011, A&A, 533, A64
  * Röser et al. (2008) Röser, S., Schilbach, E., Schwan, H., et al. 2008, A&A, 488, 401
  * ten Brummelaar et al. (2005) ten Brummelaar, T. A., McAlister, H. A., Ridgway, S. T., et al. 2005, ApJ, 628, 453
  * Tóth et al. (2014) Tóth, L. V., Marton, G., Zahorecz, S., et al. 2014, Publications of the Astronomical Society of Japan, 66, 17
  * Unruh et al. (2004) Unruh, Y. C., Donati, J. F., Oliveira, J. M., et al. 2004, MNRAS, 348, 1301

##  Appendix A Appendix A.

Table of photometry used in the SED fitting procedure.

Table 4:  Photometric values used to construct the SED of SU Aur Wavelength [$\mathrm{\mu m}$] | Flux [Jy] | Reference
---|---|---
0.15 | 1.31E-04 | Bianchi et al. (2011)
0.23 | 0.00219 | Bianchi et al. (2011)
0.42 | 0.233 | Ammons et al. (2006)
0.44 | 0.402 | Anderson & Francis (2012)
0.53 | 0.6 | Ammons et al. (2006)
0.69 | 1.32 | Morel & Magnenat (1978)
0.79 | 1.47 | Davies et al. (2014)
0.88 | 1.75 | Morel & Magnenat (1978)
1.24 | 2.08 | Röser et al. (2008)
1.25 | 2.12 | Ofek (2008)
1.63 | 2.47 | Ofek (2008)
2.17 | 2.71 | Röser et al. (2008)
2.19 | 2.62 | Ofek (2008)
3.35 | 2.6 | Cutri et al. (2014)
3.40 | 2.44 | Bourgés et al. (2014)
4.50 | 1.75 | Esplin et al. (2014)
4.60 | 2.78 | Cutri et al. (2014)
5.03 | 2.58 | Bourgés et al. (2014)
7.88 | 1.99 | Esplin et al. (2014)
8.62 | 2.36 | Abrahamyan et al. (2015)
11.57 | 2.83 | Cutri et al. (2014)
11.60 | 3.52 | Abrahamyan et al. (2015)
18.40 | 6.47 | Abrahamyan et al. (2015)
22.11 | 9.24 | Cutri et al. (2014)
23.90 | 12.8 | Abrahamyan et al. (2015)
61.89 | 12.2 | Abrahamyan et al. (2015)
65.04 | 9.89 | Tóth et al. (2014)
90.06 | 8.8 | Tóth et al. (2014)
140.10 | 10.2 | Tóth et al. (2014)
160.11 | 8.88 | Tóth et al. (2014)
849.86 | 0.074 | Mohanty et al. (2013)
887.57 | 0.071 | Andrews et al. (2013)
1300.90 | 0.03 | Mohanty et al. (2013)
1333.33 | 0.0274 | Andrews et al. (2013)

##  Appendix B Appendix B.

Figure 10: Squared visibilities split by baseline/pointing against wavelength.
Black crosses are the measured visibilities. Red lines represent synthetic
visibilities from the final TORUS radiate transfer model. The legend of each
plot is the length of the individual baseline. Figure 11: Figure 10 continued…
Figure 12: Figure 10 continued… Figure 13: Figure 10 continued… Figure 14:
Figure 10 continued… Figure 15: Figure 10 continued…


# A note on the renormalization group approach to the Central Limit Theorem

Sébastien Ott  Département de Mathématiques, Université de Fribourg, Chemin du
Musée 23, 1700 Fribourg, Switzerland<EMAIL_ADDRESS>
###### Abstract.

A proof of the Central Limit Theorem using a renormalization group approach is
presented. The proof is conducted under a third moment assumption and shows
that a suitable renormalization group map is a contraction over the space of
probability measures with a third moment. This is by far not the most optimal
proof of the CLT, and the main interest of the proof is its existence, the CLT
being the simplest case in which a renormalization group argument should
apply. None of the tools used in this note are new. Similar proofs are known
amongst expert in limit theorems, but explicit references are not so easy to
come by for non-experts in the field.

## A word from the author

I am neither a renormalization group expert, nor a limit theorem expert. So it
is perfectly possible that the results described this note appeared somewhere
else 111Indeed, after the second version of this note appeared on the arXiv, I
received a reference to a paper, [11], containing the proof presented in
section 1. It is historically interesting to note that they introduce (and
use) the Fourier based metrics: the paper is from 1984 and the Fourier based
metrics where believed to have been introduced in [8], more than 10 years
later., or that I missed relevant references to the “Renormalization group
approach to CLT” story. In both cases, I would be extremely grateful to
receive pointers towards the relevant references. After the appearance of (the
first version of) this note on the arXiv, I received several emails with
references and/or comments. The updated version is motivated by the content of
those mails, I hope that the enhanced bibliography (which is still far from
being exhaustive) will help other people interested in these questions as much
as it helped me. Many thanks to all the people who took the time to send me
information!

## Introduction/bibliographical review

This note is about providing a simple renormalization group style proof of the
classical Central Limit Theorem. The CLT is likely to be the most well known
feat of probability theory, so I do not intend to say anything new about it.
The purpose of this work is more to perform a “sanity check” for the
renormalization group approach: if one can not get the method to work in the
simplest instance it should apply to, one has little chances of success in
much more involved situations. This note is not the first result about the
renormalization road to the CLT, and one can distinguish two general
approaches: a “Lyapunov” method, consisting of solving the fixed point
equation and, more or less explicitly, finding a suitable Lyapunov function
for the discrete dynamical system $\mu_{n+1}=T\mu_{n}$ (where $T$ is the
renormalization map); and a “Banach” method, consisting in finding a subset of
the set of probability measures $\Omega$, and a distance $d$ on $\Omega$ such
that 1) $(\Omega,d)$ is complete 2) the renormalization operation is a
contraction on $(\Omega,d)$. On the one side, the “Banach” method yields more
quantitative results and _one does not need to solve the fixed point equation_
to prove convergence towards something. On the other side, the “Lyapunov”
method is more flexible and allows for weaker conditions. Section 1 contains
an example of what I called a “Banach” approach, while Section 2 presents a
“Lyapunov” approach.

### “Lyapunov” method

On the “perturbative side” (perturbation around the fixed point), one can find
the articles [13, 4, 16] which deal with the transform of equation (4) (and
its natural generalizations) that have as fixed points stable laws. Suitable
neighbourhood of the fixed points are then shown to contract to the fixed
point. Reviews/introductory texts (mostly focused on the Gaussian fixed point)
can be found in [20, pages 131-132], [14], [15, Section 10.3]. A non-rigorous
discussion on that problem can be found in [1]. The same set of ideas has been
used to prove CLT for dependent fields, as this is not the main topic of this
note, I will only mention a few early references: [3] deal with hierarchical
spin models, their method influenced several of the papers mentioned
previously; [9, 7] contain reviews of CLT for the “magnetization” of (weakly)
dependent lattice random fields (the second also provides a general discussion
about renormalization and CLT).

Another route has also been studied: the entropy can be used as a Lyapunov
function for the dynamical system underlying the renormalization process, see
[17] for the CLT using entropy, [5] for an “Lyapunov function” proof, and [12]
for a book on the topic.

### “Banach” method

A Banach fixed point argument can be found in [19] (see the slides [18] for a
discussion of the application to the CLT). The proof there is very similar to
the one of Section 1: use of ideal metric (Zolotarev metric -introduced in
[21, 22]-, in [18] versus Fourier based metric here) to prove contraction
properties of the renormalization map.

### Present note

This note contain an example of each methods. As said above, the approaches
taken here are not new: the renormalization group map is the one of [13] and
is probably older, the contraction principle of Section 1 is included in the
paper [10]: equations (6), (7) are morally the contraction used here (which is
a property of certain ideal metrics, see Section 3). The argument presented in
[19, 18] is the same as the one of Section 1 but with a different metric (see
Section 3).

The argument of Section 2 is a simplified version of an argument that has been
communicated to me by Jiwoon Park. It provides an example of the larger
flexibility of the “Lyapunov” method, and applies with only a second moment
condition.

### General comment

It is worth noting that a contraction principle is underlying most approaches
to the CLT, but the latter is usually formulated as follows: one has a
sequence of operators, $(T_{n})_{n\geq 1}$ (which acts by $n$-fold convolution
and rescaling by $\sqrt{n}$), having the normal distribution as fixed point,
and having “law-dependent contraction constant” going to $0$ with $n$, rather
than iterating a fixed transform with uniform (over probability distribution)
contraction constant.

Moreover, the proof of Section 2 can be seen as a slightly convoluted re-
writing of the usual proof of the CLT through pointwise convergence of
characteristic function, I nevertheless find it a nice example of global
asymptotic stability of the fixed point of the renormalization process. The
argument of Section 1 is more quantitative but more restrictive.

##  1\. Renormalization road to CLT: a “Banach” result

As nothing is really new in the content of this text, I will work under a
third moment assumption. The latter can be relaxed to a $2+\epsilon$ moment
assumption: the distance used (denoted $\mathrm{d}_{3}$ in the text) is a
particular case of the Fourier-based distances introduced in [8]. Replacing
this distance by its $\mathrm{d}_{2+\epsilon}$ version handles the extension
($3$ is chosen for aesthetic reasons and to slightly simplify appendix A). A
review on these metrics and their applications can be found in [6].

###  1.1. Framework

Probability measures on $\mathbb{R}$ will be denoted $\nu,\mu$, the
expectation under $\nu$ will be denoted $E_{\nu}$. Inside expected value, $X$
will denote a random variable of the relevant law, and $(X,Y)$ a random vector
of law $\nu\otimes\mu$ when this case is considered. Let
$\mathcal{P}_{r}=\mathcal{P}_{r}(\mathbb{R})$ be the set of probability
measures on $\mathbb{R}$ with finite $r$th absolute moment. Denote

$\mathcal{Q}_{3}:=\\{\nu\in\mathcal{P}_{3}:\ E_{\nu}(X)=0,\
E_{\nu}(X^{2})=1\\},$

the set of centred, reduced probability measures with a third absolute moment.
Equip $\mathcal{Q}_{3}$ with the Fourier-based distance

$\mathrm{d}_{3}(\nu,\mu)=\sup_{\xi\in\mathbb{R}^{*}}\frac{\big{|}\varphi_{\nu}(\xi)-\varphi_{\mu}(\xi)\big{|}}{|\xi|^{3}}$
(1)

where $\mathbb{R}^{*}=\mathbb{R}\setminus\\{0\\}$, and

$\varphi_{\nu}(\xi)=E_{\nu}\big{(}e^{\mathrm{i}X\xi}\big{)},$ (2)

is the characteristic function of $\nu$ ([6] defines it with a $-$ sign in the
exponential).

######  Lemma 1.1.

$\mathrm{d}_{3}$ is a finite distance on $\mathcal{Q}_{3}$, and convergence in
$\mathrm{d}_{3}$ implies weak convergence.

This result can be imported form [8, 6], but a proof is included in Appendix
A. The goal will be to study convergence towards a normal distribution
$\mathcal{N}(0,1)$. Denote $\gamma$ the normal law:

$d\gamma(x)=\frac{1}{\sqrt{2\pi}}e^{-x^{2}/2}dx.$ (3)

One obviously has $\gamma\in\mathcal{Q}_{3}$.

######  Remark 1.1.

It is worth mentioning that $(\mathcal{Q}_{3},\mathrm{d}_{3})$ is not expected
to be a complete metric space, but asking for $3+\epsilon$ moments rather than
$3$ guarantees completeness (see [6, Proposition 2.7] and the follow-up
discussion).

###  1.2. Renormalization transform

Denote $\nu*\mu$ the convolution of $\nu$ and $\mu$ (the law of $X+Y$ where
$X\sim\nu,Y\sim\mu$ are independent). Also denote $[\nu]_{\lambda}$ the law of
$\lambda X$, where $X\sim\nu$. One then consider the following renormalization
transformation on probability measures: $T\nu$ is the law of $2^{-1/2}(X+Y)$
where $X,Y$ are independent random variables of law $\nu$:

$E_{T\nu}(f)=E_{\nu\otimes\nu}\Big{(}f\big{(}(X+Y)/\sqrt{2}\big{)}\Big{)}.$
(4)

In words: $T$ maps $\nu$ to the renormalized convolution of $\nu$ with itself
($T\nu=[\nu^{*2}]_{2^{-1/2}}$). Taking the sum is the “coarse graining” part
of a renormalization step and dividing by $\sqrt{2}$ is the “rescaling” part.

######  Lemma 1.2.

If $\nu\in\mathcal{Q}_{3}$, then $T\nu\in\mathcal{Q}_{3}$.

###### Proof.

From the definition, one has that if $\nu$ has a second moment, and
$E_{\nu}(X)=0$, then $E_{T\nu}(X)=0$, and $E_{T\nu}(X^{2})=E_{\nu}(X^{2})$.
Then,

$E_{T\nu}(|X|^{3})=\frac{1}{2^{3/2}}E_{\nu\otimes\nu}(|X+Y|^{3})\leq\frac{16}{2^{3/2}}E_{\nu}(|X|^{3})<\infty.$

∎

The goal of this note is to study the CLT, so the main interest of this
transformation is

$T\gamma=\gamma,$ (5)

which is a standard consequence of the stability of the Gaussian distribution.

The link with the CLT is as follows: if $X_{1},X_{2},\dots$ form an i.i.d.
sequence of law $\nu\in\mathcal{Q}_{3}$, then $T^{n}\nu$ is the law of

$\frac{1}{2^{n/2}}\sum_{k=1}^{2^{n}}X_{k}\equiv\frac{1}{\sqrt{N}}\sum_{k=1}^{N}X_{k},$

where $N=2^{n}$.

###  1.3. Contraction and CLT

######  Theorem 1.3.

The application $T$ is a contraction on $(\mathcal{Q}_{3},\mathrm{d}_{3})$
with contraction constant $\leq 2^{-1/2}$. In particular (by (5)),

$\mathrm{d}_{3}(T^{n}\nu,\gamma)\leq 2^{-n/2}\mathrm{d}_{3}(\nu,\gamma).$ (6)

###### Proof.

The proof is almost trivial. Let $\mu,\nu\in\mathcal{Q}_{3}$. First, notice
that

$\varphi_{T\nu}(\xi)=E_{\nu\otimes\nu}(e^{\mathrm{i}\xi(X+Y)/\sqrt{2}})=\varphi_{\nu}(\xi/\sqrt{2})^{2}.$
(7)

Then, writing the definition, one has

$\displaystyle\mathrm{d}_{3}(T\nu,T\mu)$
$\displaystyle=\sup_{\xi\in\mathbb{R}^{*}}\frac{\big{|}\varphi_{\nu}(\xi/\sqrt{2})^{2}-\varphi_{\mu}(\xi/\sqrt{2})^{2}\big{|}}{2^{3/2}|\xi/\sqrt{2}|^{3}}$
$\displaystyle=2^{-3/2}\sup_{\xi\in\mathbb{R}^{*}}\frac{\big{|}\varphi_{\nu}(\xi)^{2}-\varphi_{\mu}(\xi)^{2}\big{|}}{|\xi|^{3}}$
$\displaystyle=2^{-3/2}\sup_{\xi\in\mathbb{R}^{*}}\frac{\big{|}\varphi_{\nu}(\xi)-\varphi_{\mu}(\xi)\big{|}}{|\xi|^{3}}\big{|}\varphi_{\nu}(\xi)+\varphi_{\mu}(\xi)\big{|}$
$\displaystyle\leq
2^{-1/2}\sup_{\xi\in\mathbb{R}^{*}}\frac{\big{|}\varphi_{\nu}(\xi)-\varphi_{\mu}(\xi)\big{|}}{|\xi|^{3}}=2^{-1/2}\mathrm{d}_{3}(\nu,\mu),$

as $\left\|\varphi_{a}\right\|_{\scriptscriptstyle\infty}=1$ for any
probability measure $a$. ∎

######  Corollary 1.4 (Central Limit Theorem).

Let $\nu\in\mathcal{Q}_{3}$. Let $X_{1},X_{2},\dots$ be an i.i.d. sequence of
law $\nu$. Then,

$\frac{1}{\sqrt{N}}\sum_{k=1}^{N}X_{k}\xrightarrow{N\to\infty}\mathcal{N}(0,1),$
(8)

where the convergence is $\mathrm{d}_{3}$ distance (and therefore also in
law).

The claim along the sequence $N=1,2,4,8,\dots$ (or any geometric sequence)
follows directly from Theorem 1.3, Lemma 1.2, and the fact that
$\mathrm{d}_{3}$ is a _finite_ distance on $\mathcal{Q}_{3}$ (by Lemma 1.1).
Extending this to arbitrary sequences can be done in several ways, which shall
not be exposed here. It is worth stressing out that extending the result to
arbitrary sequences is never simpler than using stability of Gaussian and
scaling+convolution properties of $\mathrm{d}_{3}$ to obtain the CLT directly.
Indeed,

$\begin{gathered}\mathrm{d}_{3}(\nu_{1}*\nu_{2},\mu_{1}*\mu_{2})\leq\mathrm{d}_{3}(\nu_{1},\mu_{1})+\mathrm{d}_{3}(\nu_{2},\mu_{2}),\\\
\mathrm{d}_{3}([\nu]_{\lambda},[\nu]_{\lambda})\leq\lambda^{3}\mathrm{d}_{3}(\nu,\mu),\end{gathered}$
(9)

where $[\nu]_{\lambda}$ is the law of $\lambda X,X\sim\nu$. The first point
follows from $|ab-cd|\leq|a-c||b|+|b-d||c|$, and the definition. Therefore, by
the properties of the Gaussian,

$\mathrm{d}_{3}([\nu^{*n}]_{n^{-1/2}},\gamma)=\mathrm{d}_{3}([\nu^{*n}]_{n^{-1/2}},[\gamma^{*n}]_{n^{-1/2}})\leq\frac{1}{n^{3/2}}\mathrm{d}_{3}(\nu^{*n},\gamma^{*n})\leq\frac{\mathrm{d}_{3}(\nu,\gamma)}{n^{1/2}}.$

##  2\. A “Lyapunov proof” under a second moment condition

As said in the introduction, the following proof is a simplified version of an
argument that has been sent to me by Jiwoon Park. It is presented here with
his kind permission.

Let

$\displaystyle\mathcal{Q}_{2}=\\{\nu\in\mathcal{P}_{2}:\ E_{\nu}(X)=0,\
E_{\nu}(X^{2})=1\\},$
$\displaystyle\mathrm{d}_{2}(\nu,\mu)=\sup_{\xi\in\mathbb{R}^{*}}\frac{|\varphi_{\nu}(\xi)-\varphi_{\mu}(\xi)|}{\xi^{2}}.$

Then, $\mathrm{d}_{2}$ is a (finite) distance on $\mathcal{Q}_{2}$ (which
metrizes weak convergence). The problem with $\mathrm{d}_{2}$ is that it is
only $2$-ideal (see Section 3), so the argument of the previous Section can
not be used with this metric (but still gives continuity of $T$). Yet one can
still use it as a Lyapunov function.

######  Theorem 2.1.

Let $V:\mathcal{Q}_{2}\to\mathbb{R}_{+}$ be defined by
$V(\nu)=\mathrm{d}_{2}(\nu,\gamma)$. Then, $V$ satisfies

  1. (1)

$V(\gamma)=0$, and $V(\nu)>0$ for any
$\nu\in\mathcal{Q}_{2}\setminus\\{\gamma\\}$.

  2. (2)

$V$ is continuous on $(\mathcal{Q}_{2},\mathrm{d}_{2})$.

  3. (3)

$V$ has bounded level sets.

  4. (4)

$V(T\nu)<V(\nu)$ for every $\nu\in\mathcal{Q}_{2}\setminus\\{\gamma\\}$.

###### Proof.

Items 1, 2, and 3 are obvious as $\mathrm{d}_{2}$ is a distance that metrizes
weak convergence. Only item 4 requires some care. Let $\nu\neq\gamma$. Then,

$V(T\nu)=\sup_{\xi\in\mathbb{R}^{*}}\frac{|\varphi_{\nu}(\xi)-\varphi_{\gamma}(\xi)|}{\xi^{2}}\frac{|\varphi_{\nu}(\xi)+\varphi_{\gamma}(\xi)|}{2}\leq\sup_{\xi\in\mathbb{R}^{*}}\frac{|\varphi_{\nu}(\xi)-\varphi_{\gamma}(\xi)|}{\xi^{2}}\frac{1+e^{-\xi^{2}/2}}{2}.$

Now, by second order Taylor expansion,

$\varphi_{*}(\xi)=1-\frac{\xi^{2}}{2}+\xi^{2}h_{*}(\xi)$

with $h_{*}(\xi)\xrightarrow{\xi\to 0}0$. Let
$a=V(\nu)=\mathrm{d}_{2}(\nu,\gamma)>0$, and $b>0$ be such that
$\sup_{|\xi|<b}\max(|h_{\nu}(\xi)|,|h_{\gamma}(\xi)|)<\frac{a}{4}$. Then, the
last term of the last display is less than or equal to

$\max\Big{(}\sup_{|\xi|\geq
b}\frac{|\varphi_{\nu}(\xi)-\varphi_{\gamma}(\xi)|}{\xi^{2}}\frac{1+e^{-b^{2}/2}}{2},\sup_{\xi\in(-b,b)\setminus\\{0\\}}|h_{\nu}(\xi)-h_{\gamma}(\xi)|\Big{)}\leq\\\
\leq\max\Big{(}\frac{1+e^{-b^{2}/2}}{2}a,\frac{a}{2}\Big{)}<a,$

proving the last item. ∎

The immediate consequence is

######  Corollary 2.2.

Let $\nu\in\mathcal{Q}_{2}$. Let $X_{1},X_{2},\dots$ be an i.i.d. sequence of
law $\nu$. Then,

$\frac{1}{2^{n/2}}\sum_{k=0}^{2^{n}}X_{i}\xrightarrow{n\to\infty}\mathcal{N}(0,1).$
(10)

###### Proof.

The proof is the standard proof that global asymptotic stability is implied by
the existence of a suitable Lyapunov function. It is included for the reader
convenience. Let $V(\nu)=\mathrm{d}_{2}(\nu,\gamma)$. The claim is equivalent
to $V(T^{n}\nu)\to 0$ as $n\to\infty$ for any $\nu\in\mathcal{Q}_{2}$. Fix
$\nu\in\mathcal{Q}_{2}$. Define

$D=\\{\mu\in\mathcal{Q}_{2}:\ 0\leq V(\mu)\leq V(\nu)\\}.$

By items 3,2 of Theorem 2.1, $D$ is compact. Moreover, as $V(T\mu)\leq
V(\mu)$, $\nu_{n}\equiv T^{n}\nu\in D$ for all $n\geq 0$. The sequence
$V(\nu_{n})$ is then a bounded decreasing sequence, denote $a\geq 0$ its
limit. By continuity of $V$, every cluster point, $\nu_{*}$, of
$(\nu_{n})_{n}$ has $V(\nu_{*})=a$. Now, if $(\nu_{n_{k}})_{k}$ is a
converging subsequence of $(\nu_{n})_{n}$, so is $(\nu_{n_{k}+1})_{k}$ (by
continuity of $T$). Moreover, the limit of $(\nu_{n_{k}+1})_{k}$ is $T\nu_{*}$
where $\nu_{*}$ is the limit of $(\nu_{n_{k}})_{k}$. So, $a=V(T\nu_{*})\leq
V(\nu_{*})=a$ with equality only if $\nu_{*}=\gamma$ (by Theorem 2.1, item 4).
So either $\nu_{*}=\gamma$ or $a=0$ which imply each other. ∎

##  3\. Additional remarks

### About the metric

The metrics $\mathrm{d}_{3},\mathrm{d}_{2}$ (and their generalizations
$\mathrm{d}_{s}$, see [10]) as well as the Zolotarev metric mentioned in the
introduction and used in [19, 18], belong to a class of metric on probability
measures which all allow the same type of argument. These metric satisfy

$\displaystyle d(\nu*\eta,\mu*\eta)\leq d(\nu,\mu),$ (11) $\displaystyle
d([\nu]_{\lambda},[\mu]_{\lambda})\leq\lambda^{s}d(\nu,\mu),$ (12)

for some $s>0$. Such a metric with equality in (12) is call $s$_-ideal_.
Suppose that the metric is defined as a supremum over a suitable class of test
function of weighted expectation difference:

$d(\nu,\mu)=\sup_{f\in F}w(f)\big{|}E_{\nu}(f)-E_{\mu}(f)\big{|},$ (13)

such metric is apparently said to have a _$\zeta$ -structure_. Suppose that
one can construct a metric $d$ satisfying (11), (12) with $s>2$, and (13).
Further suppose that convergence in the topology induced by $d$ implies weak
convergence. Then one can perform the same argument, the key being

$\displaystyle d(\nu*\nu,\mu*\mu)$ $\displaystyle=\sup_{f\in
F}w(f)\big{|}E_{\nu\otimes\nu}(f(X+Y))-E_{\mu\otimes\mu}(f(X+Y))\big{|}$
$\displaystyle\leq\sup_{f\in
F}w(f)\big{|}E_{\nu\otimes\nu}(f(X+Y))-E_{\nu\otimes\mu}(f(X+Y))\big{|}$
$\displaystyle\quad+\big{|}E_{\nu\otimes\mu}(f(X+Y))-E_{\mu\otimes\mu}(f(X+Y))\big{|}$
$\displaystyle\leq 2d(\nu,\mu).$

It is also clear that the same type of argument as in Section 1 with $s=2$
will not work so easily, so a different idea is required. Limitations of the
method are discussed in [19]. One can also wonder if a “Banach” proof can even
exist. An argument going in this direction is [2]: if
$T:\mathcal{X}\to\mathcal{X}$ is such that $T,T^{2},T^{3},\dots$ have a unique
fixed point, one can construct a metric on $\mathcal{X}$ such that $T$ is a
contraction (with arbitrary contraction constant). But the argument is non-
constructive, uses the axiom of choice, and the metric has no chance of
metrizing weak topology. I would tend to think that the general claim is wrong
but that the following weaker claim could true: one can find a sequence
$\mathcal{Q}_{2}^{(1)}\subset\mathcal{Q}_{2}^{(2)}\subset\dots$ converging to
$\mathcal{Q}_{2}$ (i.e.: $\nu\in\mathcal{Q}_{2}$ implies there is $n\geq 1$
such that $\nu\in\mathcal{Q}_{2}^{(n)}$) and a sequence of metrics
$d_{1},d_{2},\dots$ on those spaces such that
$T\mathcal{Q}_{2}^{(n)}\subset\mathcal{Q}_{2}^{(n)}$ and $T$ is a contraction
on $(\mathcal{Q}_{2}^{(n)},d_{n})$ (and $(\mathcal{Q}_{2}^{(n)},d_{n})$ is
complete, and $d_{n}$ metrizes weak convergence on $\mathcal{Q}_{2}^{(n)}$…).
Any comment on this question would be welcome!

## Acknowledgements

Thanks to Nicolas Curien for suggesting me to contact Ralph Neininger, and to
Ralph Neininger for pointers to [19, 18], and comments on the general picture.
Also thanks to Jiwoon Park for sending me the argument which led to Section 2.
Finally, thanks to all the people who send me information and comments after
the first version went online.

The author is supported by the Swiss NSF grant 200021_182237 and is a member
of the NCCR SwissMAP.

##  Appendix A Fourier based distance

This Section contains a proof of Lemma 1.1. The proof is by no mean new, it is
included for the reader convenience.

First, note that $\mathrm{d}_{3}$ is a distance: symmetry is obvious, and
separation follows from the fact that two probability measures are the same if
and only if they have the same characteristic function. Triangular inequality
follows from triangular inequality for the absolute value, and from $\sup
f+g\leq\sup f+\sup g$.

Then, show that $\mathrm{d}_{3}$ is finite on
$\mathcal{Q}_{3}\times\mathcal{Q}_{3}$. Take $\nu,\mu\in\mathcal{Q}_{3}$. As
$\nu,\mu$ have a third moment, $\varphi_{\nu},\varphi_{\mu}$ are three times
continuously differentiable and they admit a Taylor expansion at $0$:

$\varphi_{*}(\xi)=1-\frac{\xi^{2}}{2}-\mathrm{i}\frac{E_{*}(X^{3})\xi^{3}}{6}+h_{*}(\xi)\xi^{3},$
(14)

with $h_{*}(\xi)\xrightarrow{\xi\to 0}0$ bounded uniformly over $[-1,1]$ (as
it is continuous over $\mathbb{R}$), $*\in\\{\nu,\mu\\}$. So,

$\mathrm{d}_{3}(\nu,\mu)\leq
6^{-1}\sup_{0<|\xi|<1}\big{|}-\mathrm{i}E_{\nu}(X^{3})+6h_{\nu}(\xi)+\mathrm{i}E_{\mu}(X^{3})-6h_{\mu}(\xi)\big{|}+2<\infty.$

Finally, convergence in $\mathrm{d}_{3}$ implies pointwise convergence of
characteristic functions, as well as continuity of the limit at $0$, which is
equivalent to weak convergence by Levy’s continuity theorem.

## References

  * [1] A. Amir.  An elementary renormalization-group approach to the generalized central limit theorem and extreme value distributions.  Journal of Statistical Mechanics: Theory and Experiment, 2020(1):013214, 2020.
  * [2] C. Bessaga.  On the converse of Banach "fixed-point principle".  In Colloquium Mathematicum, volume 7, pages 41–43. Instytut Matematyczny Polskiej Akademii Nauk, 1959.
  * [3] P. M. Bleher and J. G. Sinai.  Investigation of the critical point in models of the type of dyson’s hierarchical models.  Communications in Mathematical Physics, 33:23–42, 1973.
  * [4] I. Calvo, J. C. Cuchí, J. G. Esteve, and F. Falceto.  Generalized central limit theorem and renormalization group.  Journal of Statistical Physics, 141:409–421, 2010.
  * [5] E.A. Carlen and A. Soffer.  Entropy production by block variable summation and central limit theorems.  Communications in mathematical physics, 140:339–371, 1991.
  * [6] J. A. Carrillo and G. Toscani.  Contractive probability metrics and asymptotic behavior of dissipative kinetic equations.  Riv. Mat. Univ. Parma, 6(7):75–198, 2007.
  * [7] M. Cassandro and G. Jona-Lasinio.  Critical point behaviour and probability theory.  Advances in Physics, 27(6):913–941, 1978.
  * [8] G. Gabetta, G. Toscani, and B. Wennberg.  Metrics for probability distributions and the trend to equilibrium for solutions of the boltzmann equation.  Journal of statistical physics, 81:901–934, 1995.
  * [9] G. Gallavotti and G. Jona-Lasinio.  Limit theorems for multidimensional Markov processes.  Commun. math. Phys, 41, 1975.
  * [10] T. Goudon, S. Junca, and G. Toscani.  Fourier-based distances and berry-esseen like inequalities for smooth densities.  Monatshefte für Mathematik, 135:115–136, 2002.
  * [11] G.G. Hamedani and G.G. Walter.  A fixed point theorem and its application to the central limit theorem.  Archiv der Mathematik, 43:258–264, 1984.
  * [12] O. Johnson.  Information theory and the central limit theorem.  World Scientific, 2004.
  * [13] G. Jona-Lasinio.  The renormalization group: A probabilistic view.  Nuovo Cimento B, 26:99–119, 1975.
  * [14] G. Jona-Lasinio.  Renormalization group and probability theory.  Physics Reports, 352(4-6):439–458, 2001.
  * [15] L. Koralov and Y. G. Sinai.  Theory of probability and random processes.  Springer Science & Business Media, 2007.
  * [16] D. Li and Y. G. Sinai.  An application of the renormalization group method to stable limit laws.  Journal of Statistical Physics, 157:915–930, 2014.
  * [17] J. V. Linnik.  An information-theoretic proof of the central limit theorem with Lindeberg conditions.  Theory of Probability & Its Applications, 4(3):288–299, 1959.
  * [18] R. Neininger.  School and Workshop on Probability Theory and Applications: Probabilistic trees, January 2007.
  * [19] R. Neininger and L. Rüschendorf.  A general limit theorem for recursive algorithms and combinatorial structures.  The Annals of Applied Probability, 14(1):378–418, 2004.
  * [20] Y. G. Sinai.  Probability theory: an introductory course.  Springer, 1992.
  * [21] V. M. Zolotarev.  Metric distances in spaces of random variables and their distributions.  Mathematics of the USSR-Sbornik, 30(3):373, 1976.
  * [22] V. M. Zolotarev.  Ideal metrics in the problem of approximating distributions of sums of independent random variables.  Theory of Probability & Its Applications, 22(3):433–449, 1978.