The Northern Arc of $\epsilon$ Eridani's Debris Ring as Seen by ALMA

TL;DR

This paper reports ALMA observations of the nearby debris disc around $$ Eridani, revealing a narrow, ring-like structure with potential planetary resonances and tentative clump detections, advancing understanding of exoplanetary debris systems.

Contribution

First ALMA imaging of $$ Eridani's debris disc, revealing its narrow structure and possible planetary influences, providing new insights into its morphology and composition.

Findings

Disc is narrow with a width of 11-13AU.

Detected potential clumps, some likely background galaxies.

Confirmed stellar chromospheric emission at the star.

Abstract

We present the first ALMA observations of the closest known extrasolar debris disc. This disc orbits the star $\epsilon$ Eridani, a K-type star just 3.2pc away. Due to the proximity of the star, the entire disc cannot fit within the ALMA field of view. Therefore, the observations have been centred 18" North of the star, providing us with a clear detection of the northern arc of the ring, at a wavelength of 1.3mm. The observed disc emission is found to be narrow with a width of just 11-13AU. The fractional disc width we find is comparable to that of the Solar System's Kuiper Belt and makes this one of the narrowest debris discs known. If the inner and outer edges are due to resonances with a planet then this planet likely has a semi-major axis of 48AU. We find tentative evidence for clumps in the ring, although there is a strong chance that at least one is a background galaxy. We confirm, at much higher significance, the previous detection of an unresolved emission at the star that is above the level of the photosphere and attribute this excess to stellar chromospheric emission.

Tables

Table 1. Table 1: Details of the observations. PWV stands for precipitable water vapour.

ID	Observation Start	Antennas	PWV (mm)
1	2015-01-17 00:37:39	34	4.0
2	2015-01-17 20:09:10	35	4.2
3	2015-01-17 21:40:16	35	4.2
4	2015-01-18 00:12:15	35	4.2
5	2015-01-18 02:32:14	34	4.2
6	2015-01-18 03:50:01	34	4.2

Table 2. Table 2: Free parameters, their priors and the results of the fitting. Note that the uncertainties on the flux densities given here do not take into account the 5% calibration uncertainties. a The maximum likelihood occurs for 48 AU, but inner edges as far in as the star still fit the data.

	Run A		Run B		Run C
Parameter	Prior	Fit	Prior	Fit	Prior	Fit
$R_{in}$ (AU)	Uniform	62.6${}_{-1.5}{}^{}{}_{}{}^{+0.9}$	Uniform		-	-
	5780		057
$R_{out}$ (AU)	Uniform	75.9${}_{-0.9}{}^{}{}_{}{}^{+1.0}$	Uniform	75.0${}_{-0.7}{}^{}{}_{}{}^{+0.7}$	-	-
	100		100
$\gamma$	Uniform	2${}_{-3}{}^{}{}_{}{}^{+3}$	Uniform	7${}_{-1}{}^{}{}_{}{}^{+1}$	-	-
	-810		320
$R_{mid}$ (AU)	-	-	-	-	Uniform	69.4${}_{-0.4}{}^{}{}_{}{}^{+0.5}$
					090
$\Delta R$ (AU)	-	-	-	-	Uniform	11.3${}_{-1.2}{}^{}{}_{}{}^{+1.4}$
					0
$F_{\nu }$ (mJy)	Uniform in ln	8.3${}_{-0.6}{}^{}{}_{}{}^{+0.6}$	Uniform in ln	10${}_{-1}{}^{}{}_{}{}^{+1}$	Uniform in ln	9.2${}_{-0.8}{}^{}{}_{}{}^{+0.8}$
	14		14		14
$I$ (°)	Uniform in cosine	34${}_{-2}{}^{}{}_{}{}^{+2}$	Uniform in cosine	33${}_{-2}{}^{}{}_{}{}^{+2}$	Uniform in cosine	33${}_{-2}{}^{}{}_{}{}^{+2}$
	01		01		01
$F_{cen}$ (mJy)	Uniform in ln	0.86${}_{-0.07}{}^{}{}_{}{}^{+0.07}$	Uniform in ln	0.82${}_{-0.07}{}^{}{}_{}{}^{+0.07}$	Uniform in ln	0.85${}_{-0.07}{}^{}{}_{}{}^{+0.07}$
	-16		-16		-16
$\Omega$ (°)	Uniform	-4${}_{-3}{}^{}{}_{}{}^{+3}$	Uniform	-3${}_{-3}{}^{}{}_{}{}^{+3}$	Uniform	-4${}_{-3}{}^{}{}_{}{}^{+3}$
	-2020		-2020		-2020
${\chi }_{red}^2$		1.100		1.097		1.100

Table 3. Table 3: Fitting results for run A when fitting to the visibilities rather than the image.

Parameter	Fit
$R_{in}$ (AU)	63.0${}_{-0.7}{}^{}{}_{}{}^{+0.6}$
$R_{out}$ (AU)	76.3${}_{-0.7}{}^{}{}_{}{}^{+0.8}$
$\gamma$	1${}_{-2}{}^{}{}_{}{}^{+2}$
$F_{\nu }$ (mJy)	8.3${}_{-0.6}{}^{}{}_{}{}^{+0.6}$
$I$ (°)	33${}_{-3}{}^{}{}_{}{}^{+2}$
$F_{cen}$ (mJy)	0.81${}_{-0.04}{}^{}{}_{}{}^{+0.04}$
$\Omega$ (°)	-4${}_{-3}{}^{}{}_{}{}^{+3}$

Equations

$$\ln \mathcal{L}=-{\chi }^2/2$$

$${\chi }^2=\sum _{i=1}^N{\left(\frac{O_i-M_i}{S_{ncr}\sigma }\right)}^2$$

$${\chi }_{red}^2=\frac{1}{N}\sum _{i=1}^N{\left(\frac{O_i-M_i}{\sigma }\right)}^2$$

$$\delta a=a_{pl}C{\mu }^{2/7}$$

$$n(F_{\nu })\mathrm{d}\log F_{\nu }=A{\varphi }_{\star }{\left(\frac{F_{\nu }}{S_{\star }}\right)}^{\alpha +1}\exp \left(\frac{-F_{\nu }}{S_{\star }}\right)\ln 10\mathrm{d}\log F_{\nu }$$

Full text

††pubyear: 2017

1 Introduction

At a distance of 3.2 pc, the K2V star $\epsilon$ Eridani is the closest extrasolar system known to host a debris disc – a circumstellar disc of dust and planetesimals. This proximity has resulted in it being one of the most intensely studied debris discs since the spatial resolution achievable allows for a detailed analysis of the radial and azimuthal structure of the disc. This is particularly interesting as such analysis can be used to predict the presence of planets in a system (see e.g. Moro-Martín, 2013).

$\epsilon$ Eri was one of the first four debris discs to be discovered (Gillett, 1986). Early submillimetre and millimetre observations suggested that it was slightly extended compared to the beam (Chini et al., 1991). With the introduction of the Submillimetre Common-User Bolometer Array (SCUBA) to the James Clerk Maxwell Telescope (JCMT), Greaves et al. (1998) were able to image the disc finding it to be a ring at $\sim 60$ AU with an orientation close to face-on and resolving the inner edge of a debris disc for the first time. They also noticed ‘clumpy’ structures in the disc, that they proposed could be due to resonant trapping of dust by planets. This inspired a number of authors to numerically determine properties of the planetary system. Ozernoy et al. (2000) modelled the interaction between a planet at 60 AU and the disc and predicted that if the clumps are created by resonances with a planet they should orbit the star at a rate of 0.6-0.8°yr*-1*. Quillen & Thorndike (2002) used the features of the disc to predict that the planet would be on an eccentric orbit with a lower semi-major axis of $\sim$40 AU and this would result in the resonant features orbiting at a faster rate and Deller & Maddison (2005) note that at least one extra planet is needed to explain the lack of sub-mm emission closer to the star. Follow-up observations with SCUBA (Greaves et al., 2005) suggest that the clumps appear to be orbiting at a rate of 2.75°yr*-1* anti-clockwise (Poulton et al., 2006), but leave open the possibility that many of the clumps are actually background galaxies. Observations with the Institut de Radioastronomie Millimétrique (IRAM) 30m telescope (Lestrade & Thilliez, 2015), the Submillimeter Array (SMA) (MacGregor et al., 2015) and SCUBA-2 on the JCMT (Holland et al., 2017) show some tentative signs of asymmetries, but observations with the Caltech Submillimeter Observatory (CSO) (Backman et al., 2009), the Herschel Space Observatory (Greaves et al., 2014b) and the Large Millimeter Telescope (LMT) (Chavez-Dagostino et al., 2016) seem to show smoother structure. Greaves et al. (2014b) and Chavez-Dagostino et al. (2016) also find that there are a large number of background galaxies to the East of $\epsilon$ Eri, some of which would have been behind the disc at the time of the earlier observations due to $\epsilon$ Eri’s high proper motion of almost 1″yr*-1* in a westward direction.

In addition to the main belt seen at long wavelengths, the 24-160 $\mu$m spectral energy distribution (SED) from Spitzer suggests a warmer ring at 20 AU and possibly another at 3 AU (Backman et al., 2009). Reidemeister et al. (2011) propose an alternative model with parent planetesimals just located at the main belt, which then produce dust that fills the inner part of the system via drag forces due to the stellar radiation and stellar wind. Observations with Herschel marginally resolve a ring at $\sim$14 AU rather than 20 AU with the same orientation as the outer ring (Greaves et al., 2014b). Recent observations from the LMT at 1.1mm (Chavez-Dagostino et al., 2016) go deep enough to clearly show dust emission in between the star and the outer ring. The detection of dust interior to the main ring in the mm calls into question the Reidemeister et al. (2011) model as they show that the drag forces only have a strong effect on grains with sizes $\lesssim 10\,\mu$m that are expected to be very faint at long wavelengths. Su et al. (2017) show that the (Reidemeister et al., 2011) model is inconsistent with observations from the Stratospheric Observatory for Infrared Astronomy (SOFIA), whilst the models of Backman et al. (2009) and Greaves et al. (2014b) are consistent with them.

In this paper we present the first Atacama Large Millimeter/submillimeter Array (ALMA) observations of $\epsilon$ Eri. This provides us with the highest resolution image of the disc. We use a forward modelling technique to determine the geometric properties of the disc. We compare these results with prior observations and discuss the implications for planets in the system and the properties of the star.

2 Observations

$\epsilon$ Eri was observed by ALMA in cycle 2 as part of the project 2013.1.00645.S (PI: A. Jordán). The data were taken in band 6 (1.34 mm) on the 17th and 18th January 2015 (see table 1). The total time on target was 4.4 hours. The proximity of this star to the Earth results in the disc having an angular diameter far larger than the primary beam of ALMA and so we decided to focus our pointing at a position of $\alpha=03:32:54.8535$; $\delta=-09.27.11.438$ (J2000), 18″ to the North of the star, coincident with the expected position of the northern ansa of the disc.

The data were reduced using the standard observatory calibration in CASA version 4.3.1 (McMullin et al., 2007), which includes water vapour radiometry correction, system temperature, complex gain calibration and flagging. The correlator was configured such that one spectral window, centred at the CO(2-1) line at 230.538 GHz, had 3840 channels of width 0.5 MHz (0.6 kms*-1*) and the other three spectral windows, centred at 232, 219 and 217 GHz, had 128 channels of width 16 MHz (21 km s*-1*). For the continuum imaging we can, therefore, use the full bandwidth of all four channels adding up to 7.875 GHz. The observations were taken whilst the array was in a compact configuration with baselines between 15 and 350 m.

The dirty image (the inverse Fourier transform of the observed visibilities), created using natural weighting and multi-frequency synthesis, is shown on the left in figure 1. The middle plot of figure 1 shows the dirty beam, which is the Fourier transform of the sampling function. The extensive (positive and negative) substructure is due to the lack of short baselines and resulting lack of sensitivity to extended emission. The resulting synthesised beam has a size 1.6″$\times$1.1″ with a position angle of 92° E of N. As the apparent size of the disc is very large and the dirty beam has extensive substructure, we measure the RMS, $\sigma$, in a region 90″ offset from the primary beam to avoid any influence from the emission. We find this to be $\sigma=14$ $\mu$Jy beam*-1*. Emission is clearly seen where the star is expected to be, despite it being far from the phase centre. The peak flux density is 0.82$\pm$0.07 mJy111Including a 5% calibration uncertainty (see section C.4.1 ALMA Partnership et al., 2016). after correcting for the primary beam. A Gaussian fit around the location of star shows a peak at $\alpha=03:32:54.861\pm 0.004$; $\delta=-09.27.29.415\pm 0.003$. This is offset by 0.11″ from the expected location, which is consistent with the expected astrometric accuracy of 0.10″ (section 10.6.3 of ALMA Partnership et al., 2016). The Gaussian fit also shows no significant deviations from the beam shape, i.e. there is no sign of resolved emission at the location of the star.

The debris disc is seen as an arc of mostly low signal-to-noise emission at $\sim$20″ from the star, with a peak directly North of the star of 5$\sigma$ significance. We also show a CLEAN image222CLEAN (Högbom, 1974) is an algorithm that fits point sources to the image. The CLEAN image is a restored image created by convolving the CLEAN model with the synthesised beam and adding the residuals. with a 2.5″ taper on the right in figure 1 to emphasise the detection at the expense of resolution. This also shows that the peak at the northern ansa is not simply due to the higher sensitivity near to the phase centre. No spectral lines were detected down to an RMS sensitivity per beam of 0.98 mJy beam*-1* in a 1.27 km s*-1* wide channel.

3 Modelling

The modelling procedure used here is based on that of Booth et al. (2016). To find the best fit parameters of the disc and their uncertainties, a Markov Chain Monte Carlo (MCMC) routine is run making use of emcee (Foreman-Mackey et al., 2013). For each likelihood calculation, a model disc distribution is created using a specified radial distribution (see the following sections) that is inclined from face on by an angle $I$ and given a position angle, $\Omega$, measured anti-clockwise from North. The disc grains are assumed to act like black bodies i.e., they have a temperature that is $\propto R^{-0.5}$. This assumption is made as we are working at a single wavelength, and so the temperature profile is degenerate with the optical depth profile. From the distribution of particles, the sky image is then determined using a line of sight integrator method (Wyatt et al., 1999). The total flux density of the disc at this wavelength is scaled to $F_{\nu}$ – whilst we only observe part of the disc, we model the whole disc and all disc fluxes given in this paper refer to the full ring not just the part visible in the image. A central component, $F_{cen}$, is added to the sky image at the fitted position of the star (see section 2). Since we only observe one side of the disc, it is difficult for us to place any constraints on the centre of the ring and so we assume that it is coincident with the location of the star. Whilst Greaves et al. (2014b) showed that an offset between the ring centre and the star could explain the brightness asymmetry that they see, no observations have detected any offset yet, with MacGregor et al. (2015) placing a 3$\sigma$ upper limit on the offset of 2.7″. The resulting sky image is then multiplied by the primary beam, convolved with the dirty beam and compared with the dirty image to produce a likelihood given by

[TABLE]

[TABLE]

$O_{i}$ and $M_{i}$ represent pixels in the observed image and model image respectively. $N$ is the total number of pixels used in the calculation. We note that the primary beam image we use is a Gaussian model of the beam. The actual primary beam may deviate slightly from this, which adds extra uncertainty to pixels far from the phase centre. We, therefore, only include pixels where the primary beam power is $>$20% of the peak in the primary beam. $S_{ncr}$ is the noise correlation ratio equivalent to the square root of the number of pixels per beam. This is required because the high resolution of the image compared to the beam size has the side effect of introducing correlated noise. The MCMC is run with 120 walkers and 2000 timesteps for each model tested.

3.1 Single component

Observations at millimetre wavelengths are typically dominated by large grains that are not strongly affected by transport processes and so should be coincident with their ‘birth-ring’. We, therefore, start by considering a single component disc to fit the main belt. This has sharp edges at $R_{in}$ and $R_{out}$ and an optical depth profile, $\tau$, that varies with radius as $R^{\gamma}$. Preliminary tests showed good fits in two regions of parameter space, which we, therefore, split into two MCMC runs: run A, which is constrained to be a narrow belt (with $R_{in}>57$ AU), and run B, which is a wide belt ($R_{in}<57$ AU) with a steeply rising profile ($\gamma>3$). The rising slope of the latter run could indicate a low-level effect of drag forces on the grains as in the model of Reidemeister et al. (2011). As in MacGregor et al. (2015), we also consider a Gaussian function (run C) centred at $R_{mid}$ with a full-width-half-maximum, $\Delta R$, to test whether a smoother edge makes for a better fit. Uniform priors have been assumed in all cases. With the exception of the limits noted above for runs A and B, the limits on the priors have been set to be wide enough so as to not influence the results. The free parameters, their priors and the best fit results are shown in table 2. The reduced $\chi^{2}$ is given by

[TABLE]

The best fit model for each case is shown on the left in figure 2. The middle column of this figure shows the convolved image and the last column shows the residuals for each run. They each show a significant residual roughly coincident with the northern ansa (labelled S1), although for run B it is slightly more significant than runs A and C. There is also another significant residual just to the South of this (S2) and a 4.6$\sigma$ residual in the North-West portion of the disc (S3). These will be discussed further in section 5.3. Other than this, there is very little difference between the best fits for each of the runs and the $\chi_{red}^{2}$ values are almost identical meaning that, from the images alone it is not possible to distinguish between these models.

Given the low signal-to-noise ratio of our data, it is also useful to consider the radially averaged profiles. Comparing directly with the observations (figure 1 left) is not very useful as this is a poor representation of the real radial profile due to the low level positive and negative structures in the PSF (figure 1 middle) that result from the interferometric nature of the observations. We can reduce these interferometric artefacts by creating a restored image for which the PSF is the synthesised beam. This is created using our best fit results (see left-hand plots in figure 3) in a similar manner to how the CLEAN image shown in figure 1 is a restored image based on the CLEAN model. Since the interferometric artefacts cannot be removed completely without knowing exactly the real flux distribution, each restored image will be slightly different. We then deproject these images based on the best fitting $I$ and $\Omega$ parameters for each model and calculate the weighted average in each annulus (where annuli of 1″ have been used here to provide the best compromise between signal to noise ratio and resolution). As the pixels of our image are smaller than the beam size, the noise in the final image is correlated and so we account for this by adjusting the uncertainty for each weighted average by a factor of $S_{ncr}$. As with the likelihood calculation, we only use pixels at a primary beam level of 20% or higher in this calculation.

The middle column of figure 3 shows the deprojected radial profile for each restored image compared to the model and the right column shows the deprojected radial profile of the residuals. From these it is clear that the main ring has a steep outer edge, although both our sharp outer edge and Gaussian models reproduce it well. The inner edge is also reproduced well by both a sharp or Gaussian inner edge up to $\sim$19″. Interior to this the flux then rises and falls again, possibly suggesting an extra component just interior to the main ring, which we shall investigate further in the next section. However, the significance of this is only 2$\sigma$ and this interior emission and the inner edge are also consistent with a $\tau\propto R^{7}$ slope as shown by model B. This also shows why model B has a slightly lower $\chi_{red}^{2}$ than models A and C.

3.2 Multiple components

Although we cannot say for certain from these observations whether there are any further components between the main ring and the star, given that other observations by Spitzer (Backman et al., 2009), Herschel (Greaves et al., 2014b), LMT (Chavez-Dagostino et al., 2016) and SOFIA (Su et al., 2017) have detected such emission, it is worthwhile to use our data to determine any upper limits on such emission. Specifically, we wish to check two of the models suggested to explain the LMT observations. Chavez-Dagostino et al. (2016) suggest that the excess they see interior to the main ring could be fit by either a narrow ($\sim$10 AU) ring at around 20 AU or a wide distribution with a $R^{-3.5}$ radial profile and an inner edge around 14 AU. To test these we assume that the main belt is described as in the best fit for run A and that the interior dust has the same inclination and position angle as the main belt. The geometric parameters of the interior component are allowed to vary slightly from those in Chavez-Dagostino et al. (2016) with the first model (run D) having prior limits of $10<R_{in,int}<30$ and $R_{in,int}<R_{out,int}<30$ and the second model (run E) having prior limits of $0<R_{in,int}<20$, $40<R_{out,int}<75$ and $-5<\gamma_{int}<-2$.

Chavez-Dagostino et al. (2016) find the total flux density of an interior component (at 1.1 mm) to be 1.7-3.3 mJy, depending on the contribution from the stellar chromosphere (see section 5.2). Extrapolating this flux density to 1.3 mm using a spectral slope of $\lambda^{-3}$ (based on their spectral slope fit) means we expect $F_{1.3mm}=$1-2 mJy. For our data, in both cases these extra components have a total flux density that is consistent with 0 mJy (see figure 4). Run D seems especially unlikely as a narrow belt should be relatively easy to detect in our observations and models with a total flux density of this component $>$0.8 mJy are ruled out at the 99.7% level. If a narrow belt does exist at this distance, then the spectral slope must be steeper than the $\lambda^{-3}$ fit of Chavez-Dagostino et al. (2016) and they must have underestimated the contribution of the stellar chromosphere. A wider disc as in run E can more easily be missed by interferometer observations and so the limits on this are less constraining with only total flux densities $>$8.5 mJy being ruled out at the 99.7% level and a peak in the posterior probability distribution around 0.5 mJy.

From the plots in figure 3, there is a tentative sign of emission just interior to the main belt with a peak of 2$\sigma$ at 17″. Given the uncertainties, this could be explained by a rising power law as in run B. Alternatively, this could be a sign of a faint belt just interior to the main belt. We, therefore, conduct one further run (F) that is similar to run D, except that the priors force it to be further from the star: $35<R_{in,int}<70$, $R_{in,int}<R_{out,int}<70$ and $-10<\ln F_{\nu}<1$. This results in a final fit where $R_{in,int}=51^{+3}_{-4}$ AU, $R_{out,int}=59^{+5}_{-5}$ AU and $F_{\nu,int}=0.7^{+0.4}_{-0.4}$ mJy. The plots for this best fit of this run are shown in figure 5. Whilst this model does a better job of fitting the radial profile than a single component alone, it offers just a very slight improvement on the single component models with a $\chi^{2}_{red}=1.086$ and so further observations are necessary to determine the dust distribution at this inner edge.

As noted in section 3.1, there are three prominent residuals seen after finding a best fit model. Given their prominence, it is possible that they may be affecting the conclusions drawn about the radial profile of the ring. S3 is found to have a negligible influence on the radial profile, likely because it is a point source. S1 and S2, however, have a stronger influence. If we consider just a 10° sector north of the star then there are clear peaks at 16″ and 21″. If we consider the azimuthal average not including these central 10°, then the flux just interior to the main ring is diminished. Of course, removing the part of the image with the highest sensitivity, also means the uncertainties then increase and the best fit for run F is still consistent within the 1$\sigma$ uncertainties. Observations of the rest of the ring at a constant sensitivity are necessary to determine whether the emission just interior to the main ring is only present to the north of the star or extends all around the star.

4 CO mass limit

CO has not previously been detected in the $\epsilon$ Eri system (Yamashita et al., 1993; Coulson et al., 2004) and, as noted in section 2, there are no obvious signs of CO in the data cube produced during the reduction process for our data. This is unsurprising as very few debris disc systems have shown signs of gas and those that have are typically young ($\lesssim$40 Myr), A star systems (e.g. Lieman-Sifry et al., 2016). Nonetheless, CO gas has been tentatively detected around a couple of F star systems observed with ALMA by averaging around the disc (Marino et al., 2016, 2017) and so we follow a similar method here.

First, we need to make some assumptions about how the gas is distributed. Given the age of the system, we will assume that any gas is second generation, released from exocometary ice through the collisional cascade of planetesimals in the disc (Zuckerman & Song, 2012; Matrà et al., 2015; Kral et al., 2016). Therefore we assume the CO to be co-located with the dust of the best-fit model for run A, i.e. in a uniform ring between 63 and 76 AU. The channels containing the CO (2-1) line can be found by considering the barycentric radial velocity that the gas would have if it is moving at Keplerian velocity. The largest Keplerian velocity will be for gas at the inner edge and at the ansa. Here the Keplerian velocity is 1.9 km s*-1*, although since we do not know in which direction the disc is rotating, this velocity could either be positive or negative. $\epsilon$ Eri itself has a radial velocity of 16.43$\pm$0.09 km s*-1* (Nidever et al., 2002) and so we consider the channels between 14.5 and 18.3 km s*-1* as potentially including the CO emission and spectrally integrate over these (three) channels. We then correct the image for the primary beam and spatially integrate over the ring, but to avoid including too much noise only consider pixels that are within the 50% power level of the primary beam. Finally this sum is divided by the number of pixels per beam to account for correlated noise. The uncertainty on the individual pixels is given by the primary beam corrected RMS. Then the uncertainties of the pixels are added in quadrature and divided by the square root of the number of pixels per beam to give the uncertainty in the region considered. This result is then extrapolated to the full ring by multiplying by the ratio of the ring area to the total area of the pixels considered in the calculation giving a final result of 0.13$\pm$0.11 Jy km s*-1*. We have, therefore, not detected any CO gas in the disc and can place a 3$\sigma$ upper limit on the CO (2-1) line flux density within the main belt of 0.33 Jy km s*-1*.

To determine an upper limit on the CO mass from this we follow the method of Matrà et al. (2015), considering both the case of low and high gas density. If the gas density is low then the excitation of the gas will be dominated by radiation rather than collisions and the CMB will be the dominant radiation source for low-energy rotational transitions at these wavelengths. Under these assumptions, the upper limit on the line flux density given above corresponds to a 3$\sigma$ CO mass upper limit of 9.6$\times 10^{-7}$ M⊕. If, instead, the gas density is high, then the CO excitation is collision dominated and the gas can be treated as in local thermodynamic equilibrium (LTE). The mass is then dependent on the kinetic temperature of the gas; the minimum mass for our flux upper limit is achieved for a temperature of 16.6 K, corresponding to the energy of the upper level (J=2) of the transition. These assumptions then give a 3$\sigma$ CO mass upper limit of 2.6$\times 10^{-8}$ M⊕. In other words, the radiation dominated regime gives a conservative upper limit, but in a secondary gas scenario we expect CO to be released with other species that will contribute to its collisional excitation (as is the case for $\beta$ Pic, Matrà et al., 2017), meaning that this upper limit is likely closer to the LTE value.

By assuming that the bodies that make up the collisional cascade contain a similar fraction of CO as those in the Solar System and comparing the production rate of CO from a collisional cascade to the photodissociation timescale, Kral et al. (2017) show that $\epsilon$ Eri should have roughly 10*-9* M⊕ of CO. Our non-detection is, therefore, consistent with the steady-state production of CO.

5 Discussion

5.1 The cold ring

The high resolution and sensitivity of our observations gives us greater precision of the parameters of the cold ring than previously available. Here we will discuss how these compare to previous observations and what further information we can learn from them.

5.1.1 Radial profile

We find that the main belt is centred at 69 AU with a width of between 11 and 13 AU depending on the model. By considering the radial profiles in figure 3, it is clear that both the inner and outer edges are steep. At the inner edge, however, there are signs of more emission suggesting either a shallower slope or a separate component just inside the main belt with an unresolved gap between the two. Whilst our outer edge seems to be in agreement with previous observations, other observations have found the inner edge and peak to be somewhat closer in than found here with values between 35 and 65 AU for the inner edge and 57 and 64 AU for the peak (Greaves et al., 1998, 2005; Backman et al., 2009; Lestrade & Thilliez, 2015; MacGregor et al., 2015). This discrepancy could be because the lower resolution of previous observations made it impossible to distinguish the main belt, which is only around 4″ wide, from the emission just inside it, strengthening the likelihood that the tentative detection of dust inside of the main belt is real.

Previous observations (Backman et al., 2009; Greaves et al., 2014b; Chavez-Dagostino et al., 2016; Su et al., 2017) have found evidence for emission interior to 15″. We do not find evidence for this in our data meaning that the dust at these distances from the star must be spread out in a wide disc and/or faint at this wavelength.

Such a large clearing interior to a narrow ring strongly implies the presence of planets that have cleared the region interior to the belt. Through scattering, planets clear the chaotic zone around them out to $\delta a$. Wisdom (1980) analytically derived this to be

[TABLE]

where $a_{pl}$ is the semi-major axis of the planet (in AU), $C$ is a coefficient found to be 1.3 by the approximate scaling theory of Wisdom (1980) (although other authors prefer values as high as 2; e.g. Quillen, 2006; Chiang et al., 2009), $\mu$ is the ratio of the planet mass to the stellar mass (0.82 M⊙ for $\epsilon$ Eri; Baines & Armstrong, 2012). Assuming the inner edge is defined solely by the outermost planet in the system (at least two planets are necessary to explain the extent of the clearing; Deller & Maddison, 2005), we use this to calculate the mass of the planet as a function of its semi-major axis as shown in figure 6. The mass required for a given semi-major axis can be reduced by taking into account the eccentricity of the planetesimals (Mustill & Wyatt, 2012), the eccentricity of the planet (Pearce & Wyatt, 2014) and collisions in the disc (Nesvold & Kuchner, 2015). Given the width of the ring, the planetesimals must have an eccentricity less than 0.1 and so we use this value in equation 10 of Mustill & Wyatt (2012) to show the effect of a disc of eccentric planetesimals in figure 6. We know that there cannot be a highly eccentric planet in the system otherwise an offset between the disc centre and the star would be seen but this offset is found to be less than 9 AU (MacGregor et al., 2015) and so the eccentricity must be less than 9 AU/$a_{pl}$. We use this relation along with equations 9 and 10 of Pearce & Wyatt (2014) to show the effect of an eccentric planet in figure 6. These estimates are all purely dynamical. Nesvold & Kuchner (2015) show that collisions are also important in determining the effect of a planet on the disc as collisional destruction of grains is enhanced in mean motion resonances, thus widening the gap. They find that the power law index for $\mu$ in equation 4 is then dependent on the age of the system relative to the collisional timescale. $\epsilon$ Eri is estimated to be between 400 and 800 Myr (Mamajek & Hillenbrand, 2008) and so we assume an age of 600 Myr when demonstrating the effect of collisions in figure 6.

These equations only provide us with a relationship between the mass of the planet and its semi-major axis. By considering the shape of the inner edge it is possible to further constrain these parameters as a low mass planet close to the belt results in a sharper edge than a high mass planet further from the belt (Chiang et al., 2009; Mustill & Wyatt, 2012). The effect of drag forces also needs to be taken into account as these may be strong in the case of $\epsilon$ Eri as shown by Reidemeister et al. (2011). Whilst the effect on large grains is limited, they find that it can be enough to make the slope of the inner edge shallower (A. Krivov, private communication), producing a rising profile as in run B. On the other hand, if there is an extra ring as in run F, then this would be a strong sign of clearing due to a planet either within the gap or through resonances. If the extra emission is not a complete ring but concentrated (see section 3.2) then this could be a sign of Trojans as demonstrated in the simulations of Nesvold & Kuchner (2015). Given this uncertainty in the inner edge, further observations are necessary before any accurate predictions can be yielded from such an analysis.

One way to potentially determine the planet’s semi-major axis is to consider its effect on both the inner and outer edge of the belt. From our results we find a fractional belt width of $\Delta R/R=0.17$, making this one of the narrowest discs known (Lestrade & Thilliez, 2015). For comparison, the classical Kuiper belt in the Solar System extends between semi-major axes of $\sim$40 and 48 AU (e.g. Gladman et al., 2001) giving it a fractional width of 0.18. It is notable that the inner and outer edges of the Kuiper belt correspond to the 3:2 and 2:1 mean motion resonances with Neptune, potentially a result of Neptune’s migration (Hahn & Malhotra, 2005; Levison et al., 2008). If the belt we see around $\epsilon$ Eri is also bounded by the 3:2 and 2:1 resonances of an unseen planet, then the inner and outer edges should correspond to the same planet semi-major axis. Using the parameters for the narrow belt fit, the inner edge would be at the 3:2 resonance of a planet at 47.8$\pm$1.2 AU and the outer edge would be at the 2:1 resonance of a planet at 47.8$\pm$1.0 AU, perfectly matching and suggesting the chance of a planet at this distance. However, it may just be coincidence and there are also a couple of caveats to consider. Firstly, the resonances effect the distribution of the planetesimals as a function of their semi-major axis. The planetesimals are likely to have a non-zero eccentricity and so the radial distribution will be wider than the semi-major axis distribution. Secondly, during the migration of Neptune that created these boundaries for the classical Kuiper belt, Neptune also scattered objects onto eccentric orbits creating the scattered disc extending much further out than 48 AU and resulting in a radial distribution of planetesimals that is broader than when only the classical Kuiper belt is considered. Nonetheless, the radial distribution of the Kuiper belt is found to have a FWHM that is only marginally broader than the locations of these resonances (Vitense et al., 2010). If there is a planet at 48 AU then we can use figure 6 to estimate its mass to be between 0.4 and 5 $M_{J}$.

Attempts to detect a planet in the outer parts of the system have been made using the direct imaging technique with NaCo on the Very Large Telescope, the Spitzer Space Telescope and Clio on the MMT (Marengo et al., 2006, 2009; Janson et al., 2007, 2008, 2015; Heinze et al., 2008). Whilst none of these have been successful, they do provide useful upper limits on the magnitude of any planets in the system, which can be converted to mass limits given a planetary evolution model and an age. Janson et al. (2015) provide the most constraining limits on planets in the system, which are shown by the grey region in figure 6. Some uncertainty is shown in these limits for two reasons. Firstly, if we assume the planet has the same inclination as the disc then the limits are more constraining when it is north or south of the star than east or west. Secondly, the limits are more constraining if a young age (0.4 Gyr) is assumed than if an older age (0.8 Gyr) is assumed as younger planets are brighter. Also note that the planet models used by Janson et al. (2015) to convert from luminosity limits to mass limits only work for masses $>0.5$ MJ and so the limits for the 0.4 Gyr case beyond $\sim$50 AU could actually be more constraining than shown here.

Combining the dynamical and observational constraints we see that if a planet is responsible for the location of the inner edge of the disc then it must be further out than $\sim$45 AU and have a mass $<1.3M_{J}$. If it is located at 48 AU, as suggested above, then these limits strongly rule out planets $>1.2M_{J}$ and potentially place limits down to $M_{pl}>0.7M_{J}$ if the system is only 400 Myr old in addition to the lower limit of 0.4 MJ from the dynamical analyses noted above. It is important to note that all of the above only considers the effect of one planet and it is possible for other planets in the system to contribute to shaping the inner edge of a debris disc, thereby meaning that the outermost planet does not need to be so massive.

In section 3.2 we discussed how there may be dust interior to the belt. If this is produced by planetesimals in that region then this suggests that planets have not completely removed comets on planet crossing orbits and we can use equation 2 of Wyatt et al. (2017) to derive an upper limit for the planet mass as a function of semi-major axis. This is also shown in figure 6 assuming a $t_{\star}=0.6$ Gyr.

5.1.2 Orientation to the line of sight

Despite only observing part of the ring, we are still able to precisely measure the orientation of the disc to the line of sight. Our models result in an inclination of around 34$\pm 2$° and position angle of $-4\pm 3$° East of North. Previous work has found a range of different inclinations. The Herschel (Greaves et al., 2014b) and LMT (Chavez-Dagostino et al., 2016) results most closely agree with ours, both finding an inclination of around 30° (with uncertainties of around 5° and 10° respectively). The SMA result shows the largest inclination difference to ours – 17$\pm$14° (MacGregor et al., 2015) – although the large uncertainties mean it is still consistent with our result. All observations so far have consistently given position angles around 0°, with Chavez-Dagostino et al. (2016) also finding a slight rotation West of North with $\Omega=-7\pm 10\degr{}$.

In general, planetary systems are found to be roughly coplanar (Kennedy et al., 2013; Greaves et al., 2014a) although this is not always the case (e.g. Winn et al., 2010) and so it is worth considering how the orientation of the disc compares to that of the planet and the star. $\epsilon$ Eri b was first discovered by Hatzes et al. (2000) using the radial velocity (RV) technique. RV alone cannot determine the orientation of a planet’s orbit, but combining with astrometry can provide this information. Benedict et al. (2006) combine the RV data with astrometry from the Hubble Space Telescope to calculate an inclination of 30.1$\pm$3.8° and Reffert & Quirrenbach (2011) make use of Hipparcos to calculate an inclination of 23$\pm$20°, both of which agree with our result. However, it is not enough to just compare the inclinations as the disc and planet may be inclined in different directions. For a system to be coplanar, the longitudes of ascending node must also match. Benedict et al. (2006) and Reffert & Quirrenbach (2011) determine the longitude of ascending node to be 254$\pm$7° and 282$\pm$20° respectively, considerably different to our position angle meaning that the system is not coplanar. Although the difference may, alternatively, be because there are issues with the RV analysis, that would result in a different orbital fit (Anglada-Escudé & Butler, 2012).

The inclination of the star has been determined spectroscopically to be 30$\pm$20° (Campbell & Garrison, 1985; Saar & Osten, 1997). More recently, the Microvariability and Oscillations of STars telescope has observed the star but different analyses provide different results with Croll et al. (2006) finding 30$\pm$3°, Fröhlich (2007) finding 45${}^{+11}_{-19}$° – albeit with a possible alternative fit of 72° – and Giguere et al. (2016) finding 69.5${}_{-7.6}^{+5.6}$. The analyses are unable to provide any information on the direction in which the stellar axis is tilted. It is therefore not yet possible to say for sure whether the disc, known planet and star are coplanar.

5.1.3 Disc flux density

From our modelling, we find the total flux density of the disc to be between 8 and 10 mJy with uncertainties around 1 mJy (depending on the assumed model). As the SMA observations are at the same wavelength we can compare directly with these. As the SMA is also an interferometer, MacGregor et al. (2015), like us, provide a flux density found through model fitting rather than aperture photometry, which they find to be 17$\pm$5 mJy. Whilst lower than this, our result is still within 2$\sigma$. However, a full SED analysis of all the infrared and millimetre photometry (Chavez-Dagostino et al., 2016) finds an SED fit that is consistent with the SMA result and not consistent with our result.

It is possible that our modelling approach is responsible for this discrepancy. White et al. (2017) fit ALMA data of Fomalhaut using two different approaches. First fitting in the image plane, as we have done in this paper, and second fitting in the visibility plane. They find that the best fit disc flux density from the first approach is lower than both the value found from the second approach and an extrapolation from single dish measurements. To make sure that a similar issue is not the reason for our discrepancy, we re-run model A, this time fitting to the visibilities, and find only a marginally ($<$1%) higher disc flux density (see appendix A).

The actual reason for our lower total flux density is likely twofold. Firstly, we have only observed part of the ring and assumed there to be no large variations around the ring. Greaves et al. (2014b) and Chavez-Dagostino et al. (2016) both find the South-East of the disc to be brighter than the North. This could be due to an actual asymmetry in the disc that our model does not take into account, perhaps due to apocentre glow (Pan et al., 2016), or due to a background source. Observations of the rest of the ring are necessary to confirm whether this is the case. Secondly, the interferometer set-up used for these observations has a shortest baseline of 15 m meaning that we are not sensitive to emission on scales larger than 18″. Whilst our modelling procedure takes account of this so that the fitted flux density is accurate for the given model, we could easily be missing larger scale components that do not show up in our data and so we are unable to place any strong constraints on them. For instance, we showed in section 3.2 that adding an interior component of up to a few milli-Janskies would still be consistent with the data. Observations with shorter baselines are necessary to determine whether we are missing flux.

5.2 Stellar flux

As mentioned in section 2, the flux density at the star’s location is found to be $F_{cen}=820\pm 68\,\mu$Jy. The flux density of the stellar photosphere at 1.3 mm is expected to be $530\pm 11\mu$Jy based on a Kurucz atmosphere model (MacGregor et al., 2015), showing that there is an excess over the photosphere of $290\pm 55\mu$Jy, confirming at higher significance the excess seen at the star by other millimetre observations (Lestrade & Thilliez, 2015; MacGregor et al., 2015; Chavez-Dagostino et al., 2016). $\epsilon$ Eri is known to be a magnetically active star with a bright chromosphere (see Jeffers et al., 2014, and references therein). Gillett (1986) even tested whether emission from the chromosphere could explain the infrared excesses, but found the expected emission would be far too low at those wavelengths. At longer wavelengths, however, chromospheric emission can be bright enough to be detectable and this offers the best explanation for the detection of an excess above the photosphere at 7 mm with the Australia Telescope Compact Array (MacGregor et al., 2015). On the other hand, we know from shorter wavelength observations in the infrared that there is dust close to the star (Backman et al., 2009; Greaves et al., 2014b), that could be contributing to an unresolved excess.

One way to determine the relative dominance of these is to compare to a star of a similar spectral type that is not known to have a significant debris disc. Liseau et al. (2015) have observed $\alpha$ Cen AB with ALMA, neither of which have shown any previous signs of a debris disc, and found a similar excess above the photosphere at long wavelengths, which they attribute to heated plasma from a chromosphere. $\alpha$ Cen B is a star of spectral type K1V (Torres et al., 2006) and, therefore, we expect it to be similar enough to $\epsilon$ Eri (of spectral type K2V; Gray et al., 2006) to estimate what the chromospheric contribution should be. Liseau et al. (2015) find that $\alpha$ Cen B has excesses above the photosphere of 40% and 220% at 870 $\mu$m and 3 mm respectively. Following Chavez-Dagostino et al. (2016), we interpolate between these, assuming that the excess follows a power law, and find that the likely chromospheric contribution at 1.3 mm is $\sim$70%. This would give a total stellar flux density of around 900 $\mu$Jy, which is slightly above our measured value suggesting that the entirety of excess we see is due to the chromosphere. This is different to the findings of Chavez-Dagostino et al. (2016) who find that, even including the chromosphere, the star makes up only about half of their measured central emission with $\sim$1.3 mJy being due to unresolved dust. Given the difference in beam sizes (11″ with the LMT compared to our 1.4″) this gives extra weight to the hypothesis that there is dust in between the main ring and the star as any dust between about 2 and 20 AU would be unresolved by the LMT but resolved and, potentially, too spread out to be detected in these ALMA observations. Whilst Backman et al. (2009) do show evidence of dust at $\sim$3 AU, the quantity of dust is too low for this to contribute significantly at millimetre wavelengths.

5.3 Azimuthal variations

Previous observations in the sub-mm have noted a possible clumpy structure in the disc (Greaves et al., 1998, 2005; Lestrade & Thilliez, 2015; MacGregor et al., 2015), although the Herschel (Greaves et al., 2014b) and the LMT (Chavez-Dagostino et al., 2016) observations find the disc to be smooth with a slight brightness asymmetry between the North and South ends of the disc (the South end being brighter). Chavez-Dagostino et al. (2016) show that there are a large number of galaxies to the East of the star’s current location, which would have been behind the disc at the time of the previous observations and likely explains at least some of the clumps.

As noted in section 3, our smooth models do leave some significant residuals. There are some residuals at the northern ansa (S1 and S2 in figure 2) and a point source to the northwest of the star (S3). Our further analysis of these shall consider the residuals from our best fit of run F (see figure 5), although the residuals for the best fits of each run are roughly the same. The significant residual at the northern ansa peaks at 50$\pm$14 $\mu$Jy beam*-1* and is extended compared to the beam. A Gaussian fit gives an estimate of the total flux density within $2.76\arcsec\times 2.40\arcsec$ to be 188$\pm$27 $\mu$Jy. The northwest point source is at a position angle of -63° and is located 16″ from the star. If it is associated with the $\epsilon$ Eri planetary system, its deprojected distance would therefore be 18″ or 58 AU, just interior to the main ring. It has a flux density of 210$\pm$50 $\mu$Jy beam*-1* after correcting for the primary beam (this point is at the 28% primary beam level, so we should be careful about over-interpreting this clump).

Given the depth of this map, let us first consider whether any background galaxies are coincident with the ring. Galaxy number counts can be parametrised using a Schechter (1976) function. The number of galaxies, $n$, per logarithmic bin can be found using

[TABLE]

where $A$ is the area of sky in deg2 and $\phi_{\star}$, $S_{\star}$ and $\alpha$ are all constants of the equation. Many authors have investigated galaxy counts at sub-mm wavelengths. Carniani et al. (2015) have recently done this for ALMA observations at 1.3 mm and so their results are most useful to us. They find values of these constants of $\phi_{\star}=1800\pm 400$ deg*-2*, $S_{\star}=1.7\pm 0.2$ mJy and $\alpha=-2.08\pm 0.11$. The significant residuals mentioned above are either within the main ring or just interior to it in the potential faint ring and have signal to noise ratios of at least 3.6$\sigma$. Considering this, we find the possibility of finding a galaxy between 15 and 24″ (deprojected) from the star, $\pm$65° from North of at least 3.6$\sigma$ (here we are considering signal to noise ratio rather than sensitivity since the sensitivity varies across the image) to be 103%. Carniani et al. (2015) show that physical models predict slightly lower galaxy counts than this Schechter function and, as they only go deep enough to detect galaxies brighter than 40 $\mu$Jy, they are not able to distinguish between these models, so the number we predict should probably be taken as an upper limit. Nonetheless, this does show that if any of the significant residuals in our image are real, it is likely that one is a galaxy.

The SCUBA (Greaves et al., 2005), IRAM (Lestrade & Thilliez, 2015) and SCUBA-2 (Holland et al., 2017) images all show tentative evidence for a ‘clump’ to the northwest of the star that could be related to the point source that we detect, however, this does not appear in the deeper LMT image (Chavez-Dagostino et al., 2016). For comparison with the previous observations we also plot the azimuthal variation in the residuals for run F in figure 7. Each point is the mean of a 10° annular sector between deprojected radii of 15″ and 24″ (note that this is slightly narrower than the analysis in previous papers since our higher resolution allows us to more accurately define the extent of the ring). This has been limited to sectors where all pixels are above the 20% primary beam level as the noise increases drastically beyond this, which does mean that the northwest point source is not included. The uncertainties are calculated in the same way as for the radial profile (see section 3.1). None of the azimuthal profiles of the previous observations show a rise to the north of the star, although comparing this azimuthal profile to that of MacGregor et al. (2015), the rise we see directly north of the star does look similar to the rise they see $\sim$20° East of North, but the SMA observations were only taken roughly half a year prior to ours and so this would represent an infeasibly fast rotation rate (for comparison the Keplerian rate is 0.6°/yr at 70 AU).

Follow-up observations will be necessary to determine whether any of the features seen here are real and associated with the planetary system.

6 Conclusions

In this paper we present the first observation of $\epsilon$ Eri using ALMA. We focused our observation on the northern ansa of the main ring, 18″ north of the star. Emission at the ansa of the ring is clearly detected at the 5$\sigma$ level, with a faint arc of emission from the ring stretching away from this. Whilst the signal-to-noise ratio along the arc is low, the area it covers compared to the beam size is large, enabling us to determine various parameters of the disc to a high degree of accuracy. We first consider three different single component models to fit the main belt. The best fits for all three have a FWHM between 11 and 13 AU, showing that the ring is narrow. For instance, assuming a sharp edged distribution, the observed disc emission is reproduced well with a model that extends from $63\pm 1$ AU to $76\pm 1$ AU. It is inclined at 34$\pm 2$° from face on with a position angle of -4$\pm 3$° measured east of north. The total flux density of the disc is found to be between 8 and 10 mJy, although caution should be taken here since this is based on a fit to only a small part of the ring and there is evidence for the southern side of the ring to be brighter, which is likely why this is discrepant with previous results.

The edges of the belt are found to be very steep. Considering both dynamical and observational constraints, if a planet is responsible for the location of the inner edge of the disc then it must have a $M_{pl}<1.3$ MJ and $a>45$ AU. It is notable that the inner and outer edges of the Solar System’s Kuiper belt line up with the 3:2 and 2:1 mean motion resonances of Neptune. If the same is true in the $\epsilon$ Eri system, then we would expect a planet in the system to have a semi-major axis of 48 AU and a mass between 0.4 and 1.2 MJ.

We find tentative evidence for extended residual emission at the northern ansa and a point source to the North-West of the star, just interior to the disc. The latter of which could relate to clumps seen in previous observations, although there is a strong probability of it being due to a background galaxy. Further monitoring is necessary to determine whether these features are real and associated with the disc.

We find no strong evidence for emission around 20 AU in our data that has been previously reported in other observations. A narrow component at this distance would need to have $F_{1.34mm}<0.8$ mJy to not show up in our data whereas a flux density an order of magnitude higher would be possible if it were distributed more broadly. We do find tentative evidence for a faint, extra component just interior to the main disc, peaking at around 55 AU. Due to this observation being of such a small part of the ring, it is unclear whether this extra component is entirely due to the residual emission at the northern ansa or extends around the entire ring.

As with previous observations we find the stellar flux density to be higher than would be expected for a simple photosphere. This is consistent with a contribution from the chromosphere.

Acknowledgements

The authors wish to thank Markus Janson for providing the observational limits from Janson et al. (2015), Johan Olofsson and Claudio Caceres for help with fitting visibilities, and G. Fritz Benedict, Aaron Boley and Jacob White for useful discussions. The authors thank the referee for their helpful review. MB acknowledges support from a FONDECYT Postdoctoral Fellowship, project no. 3140479 and the Deutsche Forschungsgemeinschaft (DFG) through project Kr 2164/15-1. AJ and SC acknowledge financial support from the Millennium Nucleus RC130007 (Chilean Ministry of Economy). AJ acknowledges support from FONDECYT project 1130857, BASAL CATA PFB-06 and the Ministry for the Economy, Development, and Tourism’s Programa Iniciativa Científica Milenio through grant IC 120009, awarded to the Millennium Institute of Astrophysics (MAS). MCW acknowledges the support of the European Union through ERC grant number 279973. GMK is supported by the Royal Society as a Royal Society University Research Fellow. JCA acknowledges support from PNP/CNES. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2013.1.00645.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSC 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. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration et al., 2013).

Appendix A Fitting to the visibilities

The modelling in this paper is all done in the image plane. In theory, fitting a model to the dirty image is equivalent to fitting a model to gridded visibilities and should be a good approximation to fitting to the ungridded visibilities. In practice, each method will likely involve some averaging, which could plausibly bias the results. White et al. (2017) model their ALMA data of Fomalhaut using both an image fitting method and visibility fitting method to compare the results. They find that in general the fitted parameters for the two different approaches are very similar accept for the total disc flux, which is 15% lower in the image fitting approach. With this in mind, here we re-run run A but this time fitting to the visibilities.

To do this requires only some minor changes to the set-up. We average the visibilities to one channel for each of the four spectral windows. We create models at each of these four wavelengths, Fourier transform them and align them with the coordinates of the data. We then interpolate this image to provide the expected visibilities at the same $u$, $v$ coordinates as the data, calculating the $\chi^{2}$ values using the weights in the data.

The results are shown in table 3. The differences between the results found when fitting to the visibilities and those found when fitting to the image (table 2) are small and all within the uncertainties. In particular, the total disc flux density is only 0.6% higher. We, therefore, do not find the discrepancy seen by White et al. (2017). Considering the plots (their figures 3 and 5) of the residuals in White et al. (2017), it seems likely that the discrepancy they see is due to their image plane fit not successfully accounting for all of the flux in the ring rather than an inherent issue with fitting in the image plane.