{'fraction_non_alphanumeric': 0.02592502259651942, 'fraction_numerical': 0.018421543449878768, 'mean_word_length': 4.289269995446957, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Using Ice and Dust Lines to Constrain the Surface Densities of Protoplanetary Disks

We present a novel method for determining the surface density of protoplanetary disks through consideration of disk 'dust lines' which indicate the observed disk radial scale at different observational wavelengths. This method relies on the assumption that the processes of particle growth and drift control the radial scale of the disk at late stages of disk evolution such that the lifetime of the disk is equal to both the drift timescale and growth timescale of the maximum particle size at a given dust line. We provide an initial proof of concept of our model through an application to the disk TW Hya and are able to estimate the disk dust-to-gas ratio, CO abundance, and accretion rate in addition to the total disk surface density. We find that our derived surface density profile and dust-to-gas ratio are consistent with the lower limits found through measurements of HD gas. The CO ice line also depends on surface density through grain adsorption rates and drift and we find that our theoretical CO ice line estimates have clear observational analogues. We further apply our model to a large parameter space of theoretical disks and find three observational diagnostics that may be used to test its validity. First we predict that the dust lines of disks other than TW Hya will be consistent with the normalized CO surface density profile shape for those disks. Second, surface density profiles that we derive from disk ice lines should match those derived from disk dust lines. Finally, we predict that disk dust and ice lines will scale oppositely, as a function of surface density, across a large sample of disks.


INTRODUCTION
Extrasolar planetary systems display a large diversity in both orbital architecture and the physical characteristics of the planets. This diversity could be the result of late stage planetary collisions (e.g. Inamdar & Schlichting 2016), the properties and evolution of the initial gas disk (e.g. Ginzburg et al. 2016), different initial planetary formation locations in the disk (e.g. Inamdar & Schlichting 2015) or a combination of these factors. The immediate initial conditions of planet formation are encapsulated in the protoplanetary disks that surround young stars. However, many disk characteristics remain largely unconstrained. Recent telescopic advances, particularly the Atacama Large Millimeter Array (ALMA), have enabled exploration of disks with unprecedented spatial resolution. These advances have already given us many insights into the spatial structure of disks (e.g. Nomura et al. 2016;Pérez et al. 2016;Tobin et al. 2016;Andrews et al. 2016;Williams & McPartland 2016). Here, we take advantage of this spatial resolution to propose new observational constraints on disk surface density, a property that is fundamental for understanding both disk evolution and planet formation.
Protoplanetary disk surface density cannot be measured directly because the majority of the disk mass resides in H 2 which is a symmetric particle that does not readily emit. Instead, densities have been inferred using observations of disk dust, CO, and HD. The reliability of each of these tracers has recently been called into question, leaving open the possibility that disk surface densities are entirely unconstrained (Mundy et al. 1996;Andrews et al. 2009Andrews et al. , 2010Isella et al. 2009Isella et al. , 2010Guilloteau et al. 2011;Bergin et al. 2013;Williams & Best 2014;Schwarz et al. 2016).
The first and most commonly used tracer of the surface density is the disk dust mass which is typically derived from resolved continuum observations or spectral energy distribution fitting (Andrews & Williams 2005;Calvet et al. 2002;Guilloteau et al. 2011). Using dust as a tracer of total disk surface density is fallible, however, as it requires an assumed dust-to-gas ratio. This ratio is not well constrained and can be altered from the arXiv:1704.04693v1 [astro-ph.EP] 15 Apr 2017 ISM value of 10 −2 (Williams & Cieza 2011) through processes such as grain growth and particle drift ) and can further have a non-uniform value throughout a disk.
CO gas has also been used as a tracer of the total gas present in disks (e.g. Rosenfeld et al. 2012). However, recent observations have called into question the typically assumed abundance of CO (10 −4 in warm molecular clouds; Lacy et al. 1994), suggesting that the existence of disk processes such as photodissociation may alter this value or that there may be a global depletion of gas phase carbon in disks such as TW Hya (van Zadelhoff et al. 2001;Dutrey et al. 2003;Chapillon et al. 2008;Schwarz et al. 2016). In addition to this unknown, many CO lines are optically thick and are therefore unreliable measures of mass. The use of CO observations thus further requires careful consideration of lower optical depth CO isotopologues to estimate the gas mass to within an order of magnitude (Williams & Best 2014;Ansdell et al. 2016).
More recently, observations of the HD J = 1 − 0 line have been used to probe the gas mass in disks (Bergin et al. 2013;McClure et al. 2016). HD is thought to be a good tracer of the total gas mass as the deuterium to hydrogen fraction is relatively well-known for objects near to the sun (Linsky 1998). However, because HD emits at temperatures above 20 K, HD emission will only trace the warm gas and thus provides a lower limit on the total gas mass present in a disk.
Given the uncertainties that accompany these observational tracers, in this paper we choose to adopt an agnostic point of view in regards to surface density. We develop a novel method to derive this quantity through asserting that the dust line locations are determined by the microphysical process of particle drift. We use physics that has been studied extensively in previous disk models. Our contribution is to suggest a new interpretation of disk observations in light of this physics.
Previous theoretical work indicates that, for the outer regions of evolved disks, drift dominates the processes of growth and collisional fragmentation in determining the maximum particle size at a particular radial scale. Particle growth can be limited by either a lack of total time (disk age), by fragmentation, or by drift. For typically turbulent disks, the disk lifetime allows plenty of time for particles to grow. The particle size at a given radial scale is limited, however, as particles are removed from the outer disk due to particle drift at a smaller maximum size than could be removed by collisional fragmentation. The maximum particle size at a given radius is thus defined as the size for which the growth timescale (t grow ) and the drift timescale (t drift ) are equal (i.e. t drift = t grow ). This is described both numerically and analytically in Birnstiel et al. (2012) and Birnstiel & Andrews (2014) as the "late phase" of disk evolution in which dust growth has progressed such that it is limited by the removal of larger grains via the process of radial drift for roughly sub-centimeter sized particles. Thus, the maximum particle size at a given radial location is limited by particle drift which will remove all particles with a drift timescale less than the age of the system (i.e. t disk ). This indicates that the equilibrium particle size at a particular disk location is limited by the time in which that particle has been able to interact dynamically with the disk gas. The maximum disk radius where we would therefore expect to see emission from a given particle size is defined as the location for which the drift timescale (and thus the growth timescale) is equal to the age of the system (i.e. t drift = t grow = t disk ). These locations can be seen observationally as the point where emission drops off at an observed wavelength (λ obs ) where we assume that the observed particle size is roughly equal to λ obs . Any particle larger than this size would result in a shorter drift and/or growth timescale and would therefore reach an equilibrium location at a shorter radial scale. We note that we expect the equilibrium particle size at a particular radial location to decrease with increasing system age such that we expect larger disk radial scales at the same λ obs for younger disks, limited only by the size of the gas disk itself.
In this paper we thus make the assumption that t drift = t grow = t disk to derive fundamental disk properties. We demonstrate that, if validated, this assertion allows us to recover the surface density distribution for observed disks. To test this central premise of our modeling we use two sets of observations that give empirical information about disk dust and ice lines.
The first set of observations is a collection of recent observations, using ALMA, the Jansky Very Large Array (JVLA) and the sub-millimeter array (SMA), of TW Hya that demonstrate that the disk radial scale is distinctly smaller at longer wavelengths (Menu et al. 2014;Cleeves et al. 2015;Andrews et al. 2012). These observations provide information about the distribution of dust grains throughout the disk. We describe these observations in terms of disk dust lines which refer to the disk radial scale as it corresponds to a particular wavelength. We use these dust lines to derive disk surface densities through equating t disk and t drift as described in Section 3.1.
ALMA also provides direct observations of disk ice lines, either through direct measurements of CO gas line emission (Nomura et al. 2016;Schwarz et al. 2016) or through indirect measurements of N 2 H + , which is only present in large abundance when CO freezes out (Qi et al. 2013. For TW Hya, the observation of N 2 H + yields an ice line location of ∼ 30 AU (Qi et al. 2013) while the emission from the C 18 O J = 3−2 line indicates an ice line of ∼ 10 AU (Nomura et al. 2016). The direct observations of CO emission also give insight into the surface density contributed by CO (Cleeves et al. 2015;Rosenfeld et al. 2012). The CO ice line location depends on the CO surface density through grain adsorption rates and particle drift. These observations therefore provide a constraint on the disk surface density and CO fraction.
After discussing representative parameters for our fiducial disk, TW Hya (Section 2), we explain how we derive disk surface densities from dust lines (Section 3). In Section 4 we provide additional tests of our model through considering the disk ice lines. In Section 5 we present a description of three observational diagnostics of our model and an application to a larger range of disk parameter space. If these diagnostics confirm our interpretation, this will provide a new way to observationally measure disk masses and surface density profiles. We conclude with a paper summary and a discussion of the presented observational diagnostics in Section 6.

PARAMETERS FOR FIDUCIAL DISK TW HYA
We adopt TW Hya as our fiducial protoplanetary disk because it is the nearest observed disk (d = 54 ± 6 pc) that is nearly exactly face-on (i ∼ 7 • ; Qi et al. 2004) and hence boasts a wealth of observational data. TW Hya is a long-lived disk (t disk = 3 − 10 Myr; Barrado Y Navascués 2006;Vacca & Sandell 2011) that is likely an unusually massive representative of a class of evolved protoplanetary disks as the typical disk lifetime is thought to be a few million years (Mamajek 2009). We note that disks are typically assumed to have the same ages as their host stars and estimates of stellar ages are subject to observational uncertainties. For our discussion of TW Hya we use an approximate age of 5 Myr, however, t disk could be treated as a tunable quantity as appropriate.
We assume a temperature structure for TW Hya that is dominated by passive stellar irradiation. We expect that this model holds for the outer disk and note that TW Hya may be irradiation dominated at all but the very inner radii (Dullemond et al. (2007) Figure 3, based on models from D' Alessio et al. (2006)). Our parameterization of the disk midplane temperature follows Chiang & Goldreich (1997) where the canonical temperature profile is: where the coefficient T 0 is a function of stellar luminosity and stellar mass, defined at r 0 = 1 AU, and is determined via: We adopt the following parameters for TW Hya: L = 0.28L , M = 0.8M , and µ = 2.3m H assuming a hydrogen/helium disk composition (Qi et al. 2013;Rhee et al. 2007). Using Equation 2, we derive T 0 ∼ 82 K. We note that our derived midplane temperature profile is in good agreement with the observationally constrained midplane temperature derived in Cleeves et al. (2015) as well as the upper limit on the midplane temperature from Schwarz et al. (2016). We vary the normalization of this temperature profile in Section 5.2 and discuss the effect that this has on our ice line derivations. We also note that a factor of 2 change in temperature normalization (T 0 ) leads to a change of a factor of ∼ 0.7 in our surface density profile derived in Section 3.1.
Spatially resolved CO observations of TW Hya have been well fit by the following surface density profile: which follows from the self-similar solution to the viscous evolution equations as shown in Lynden-Bell & Pringle (1974) and Hartmann et al. (1998). This profile is a shallow power law at small radii and follows an exponential fall off at radii larger than the critical radius, r c . Using an assumed CO abundance of ∼ 10 −4 n H (the standard CO fraction in warm molecular clouds; Lacy et al. 1994), Rosenfeld et al. (2012) derive best fit parameters for TW Hya of r c = 30 AU, γ = 1, and Σ c ∼ 0.5. As we move forward with our discussion of TW Hya we accept the best fit parameters for all values mentioned above, except for Σ c , which relies on an assumed CO abundance. Instead, we treat Σ c as a free parameter. This is motivated by discrepancies between assumed and derived CO abundances in disks. We further note that there is also a potential discrepancy between the measured and derived mass accretion rates for TW Hya.
TW Hya has an average measured accretion rate of ∼ 1.5 × 10 −9 M yr −1 (Brickhouse et al. 2012). As a consistency check, we can use the surface density profile from Equation 3 to derive an approximation of the mass accretion rate using the following expression (Jones et al. 2012 whereṀ is the mass accretion rate, t disk is the age of the disk, and M disk is the disk mass which we take to be the mass of TW Hya interior to the critical radius of 30 AU. This estimate forṀ is a rough approximation for protoplanetary disks under the assumption that the primary mode of disk evolution is accretion. For instance, a rate that is higher than this derivedṀ would quickly evolve the disk past the current state and a lowerṀ would indicate that a process other than accretion drove the disk to evolve into its current state. Thus, if the disk age is a proxy for evolution timescale, Equation 4 gives an accretion rate of ∼ 10 −11 M yr −1 for TW Hya using the best fit parameters from Rosenfeld et al. (2012) -a value that is 2 orders of magnitude smaller than the observational value. This value is inconsistent with observations. However, it is important to note that the accretion rates onto pre-main sequences stars are likely variable or episodic in nature (e.g. Armitage et al. 2001;Salyk et al. 2013;Hein Bertelsen et al. 2016). For episodic accretion to explain this discrepancy, TW Hya would have to be currently undergoing an episode of high accretion -a result that is unlikely given the smooth, axisymmetric nature of the disk and its observed central cavity.
We further note that photoevaporation can also remove mass in the outer disk. Since this process reduces a disk's accretion rate onto its star for a given disk mass, if important, it would make Equation 4 an upper limit foṙ M , making its agreement with observed accretion rates worse. The rate of photoevaporative mass loss for typical fluxes is < 10 −10 M yr −1 (Alexander et al. 2006), which is less than TW Hya's observed accretion rate, so it is likely subdominant. Thus, while the estimate of the mass accretion rate from Equation 4 is not necessarily conclusive, it nevertheless provides a reason to believe that mass in the disk may be higher than indicated by CO observations. TW Hya also has an observational lower limit of total gas mass of 0.05 M from HD measurements of the warm gas in the disk (Bergin et al. 2013). The mass estimate is inconsistent with the mass estimate from the CO observations and gives an accretion rate of ∼ 6 × 10 −9 M yr −1 which is more consistent with the observed rate.
We consider the discrepancy between the measurements of the CO emission, HD gas emission, and observed accretion rate to be additional motivation for treating Σ c as a free parameter.

DUST LINES
The observed extent of TW Hya is wavelength dependent, ranging from a radius of r ≈ 25 AU at a wavelength of λ = 9 mm to r ≈ 130 AU at a wavelength of 1.6 µm as shown in Figure 1. In particular, observations at 0.87 mm show a disk size of approximately 60 AU , at 1.3 mm the disk size is around 50 AU (Cleeves et al. 2015), and at 9 mm the disk size is approximately 25 AU (Menu et al. 2014). We Figure 1. The dominant particle size in the disk as a function of radius which follows a power law relationship. The three black points represent disk sizes derived from observations by Menu et al. (2014), Cleeves et al. (2015), and Andrews et al. (2012Andrews et al. ( , 2016. The upper limit arrow represents measurements of the radial extent of both the CO gas disk and the disk size of the smallest grains which we take to be indicative of the total radial extent of the disk (Debes et al. 2013). This is a theoretical lower limit as we expect that small dust sizes would have a larger radial extent if the gas disk was larger. The error bars shown correspond to our chosen nominal error for the disk radius of ±10 AU, consistent with the observationally sharp cut-off in disk emission (see text). take the observed disk size of approximately 130 AU at 1.6 µm to indicate the total radial extent of the disk as this distance matches the observed radial extent of the CO emission (Qi et al. 2013).
Recent observational work has found that the continuum emission at each wavelength exhibits a markedly sharp decrease over a very narrow radial range such that ∆r/r 0.1 de Gregorio-Monsalvo et al. 2013). Thus, while these are only approximate radial sizes for the TW Hya disk, they are adequate for our preliminary physical interpretation. In our further discussion of TW Hya we allow the radius error bars to vary by ±10 AU and find little change in our theoretical modeling. These error bars are slightly larger than those derived from treating the systematic uncertainty in distance to the TW Hya system alone (d = 54 ± 6 pc or d = 51 ± 6 pc van Leeuwen 2007; Mamajek 2005, which gives uncertainties that roughly range from ± 3 AU at λ obs = 9 mm to ± 7 AU at λ obs = .87mm). We inflate the error bars to account for error in measuring the disk radial scale without modeling the disk visibilities. We will improve on these error estimates in future work.
Models of disk emission from particles with a range of sizes indicate that λ obs roughly corresponds to the primary particle size, s, contributing to the observed emission (Birnstiel & Andrews 2014). We therefore set the particle radius s = λ obs for all that follows and leave a more detailed evaluation of the particle size distributions for future work. We call the cutoff distance for emission at λ obs = s the "dust line" for particles of size s.
In Figure 1 we note that the observed locations of dust lines as a function of particle size (equivalent to the disk size as a function of λ obs ) can be well fit by a power law, with the exclusion of the disk size at the shortest λ obs = 1.6 µm. Because the radial extent of the 1.6 µm grains matches that of the CO gas, they are limited by the total disk size and not by the equilibrium that governs disk sizes at longer λ obs . We conclude that this point acts as a theoretical lower limit for the radius at which you would expect to find these grains (Rosenfeld et al. 2012;Debes et al. 2013). Given their present location, they have not yet had to time to grow or drift due to the long dynamical timescales for micron sized particles in the outer disk.
The presence of dust lines in TW Hya indicate that a physical process is removing particles larger than a set maximum size from the outer disk. Particle growth can be truncated by particle fragmentation, particle drift, or a lack of total dynamical time. If we consider the case of TW Hya, the age of the system (t disk ) is sufficiently large that there is plenty of time for particle growth to occur.
Previous work has shown that for the outer regions of evolved disks, such as TW Hya, particle drift is significant for sub-centimeter sized bodies such that these particles are not able to grow to a size that can be disrupted due to collisional fragmentation (see Birnstiel et al. 2012;Birnstiel & Andrews 2014). This is physically intuitive as small particles in the outer disk drift faster as they grow larger (see Section 3.1).
This drift limited regime reaches an equilibrium such that the radial drift of particles imposes a size limit as large particles are removed faster than they can be replenished due to particle growth. Particles in the outer disk are also removed at a size that is smaller than the limit imposed by collisional fragmentation, for disks of standard turbulence, as the relative turbulent velocities between particles of these sizes is not sufficiently large for efficient destruction (Birnstiel et al. 2012). Thus, the maximum particle size at a given radial location in the disk is given by the size for which the growth timescale is equal to the drift timescale (i.e. t grow = t drift ).
Smaller particles have longer drift timescales (see Section 3.1) which means that the dust lines in a disk will evolve with disk age such that we expect the radial location of a dust line at a specific λ obs to decrease with increasing disk age. This therefore indicates that the dust line of a disk is determined by the time in which Figure 2. Cartoon of our model for disk dust lines (dashed lines). Particle sizes are denoted s1, s2, and s3 where s1 > s2 > s3. Particles of size s1 are present in the disk throughout the yellow region. Particles of size s2 extend throughout the yellow and red regions while particles of size s3 are present throughout all depicted disk regions. At the dust lines for each particle size the growth and drift timescale are equal to the age of the system. When observed at λ obs = s1 only the yellow region of the disk will be observed, while for λ obs = s2 the disk extends radially to the end of the red region and for λ obs = s3 the disk appears to extend to the end of the blue region.
particles are able to drift, which is necessarily the age of the system. We can thus expect the fall-off of emission from a given particle size to correspond to the location where a particle's drift timescale is equal to the lifetime of the disk (t drift = t disk ). Therefore, the system can be described at any given time by an equilibrium state in which t drift = t grow = t disk . This new picture of disks is summarized in Figure 2.
Using the dust line observations of TW Hya in conjunction with our theoretical premises we can now consider the dominant physical processes of growth and drift and use these calculations to derive the total disk surface density for TW Hya as well as the dust-to-gas ratio in the outer disk. We do this through the use of our central modeling premise that t drift = t grow = t disk .

Radial Drift
As discussed above, we assume that the disk radial scale as a function of wavelength is set by the distance that a particle of radius, s = λ obs , can drift in the age of the system. We therefore set the drift timescale t drift equal to the lifetime of the system t disk .
In a protoplanetary disk particle drift occurs because the gas disk orbits at a sub-Keplerian velocity due to an outward pressure gradient (Weidenschilling 1977). The particles in the disk continue to rotate at a Keplerian velocity (v k ≡ Ω k r) and experience a headwind from the gas. This headwind causes the particles to loose angular momentum and drift radially inwards (Weidenschilling 1977;Takeuchi & Lin 2002). The amount of drift that a particle experiences depends on how well-coupled the particle is to the gas, quantified by a dimensionless stopping time: τ s ≡ Ω k t s where Ω is the Keplerian frequency and (2010)). Here ρ is the gas midplane density, ρ s = 2 g cm −3 is the density of a solid particle, s is the particle size, and λ = µ/ρσ coll is the gas mean free path where σ coll = 10 −15 cm 2 .
The radial particle drift velocity is: where η ≈ c 2 s /2v 2 k , c s is the sound speed of the gas, and v k is the Keplerian velocity (see the review by Chiang & Youdin 2010).
We can now derive an equation for drift timescale (t drift = |r/ṙ|) that directly depends on the disk surface density, r, and s. We first note that the particle sizes we consider (s < 1 cm; see Section 3) typically interact with the disk gas via the Epstein law for gas drag and have a dimensionless stopping time less than 1. Once our calculation is complete we verify that the Epstein regime applies for our reconstructed surface density profile and that τ s < 1.
We approximate a dimensionless stopping time in the Epstein regime as τ s = Ωρ s s/ρc s . We can rewrite this quantity by noting that the disk's scale height H = c s /Ω and its surface density in gas Σ g = Hρ to obtain τ s = ρ s s/Σ g . As found in Rosenfeld et al. (2013), increasing temperature with height can inflate the aspect ratio of an observed disk when layers several scaleheights above the midplane are probed. As most of the disk mass is contained within the scale height closest to the midplane, however, the temperature at the midplane sets the volumetric density in the disk to an order of magnitude as described above.
We reformulate the drift timescale through defining the parameter v 0 which physically corresponds to the maximum drift velocity. For a passively irradiated disk, v 0 ≡ ηv k = c 2 s /2v k varies very weakly with radius. We find that c s ∝ r −3/14 where c s = kT /µ using Equation 1 for temperature as a function of radius. For the Keplerian velocity we find that v k ∝ r −1/2 using Kepler's third law. This gives a parameter v 0 that is weakly dependent on radius: v 0 ∝ r 1/14 .
We can now write the drift timescale directly in terms of the surface density, radius, particle size, and the maximum drift velocity v 0 : We can now set the drift timescale equal to the age of the disk (t drift = t disk ) and solve for the gas surface density at a given dust line radius. We thus derive the following equation for disk surface density as a function of radius.
From this equation we see that Σ(r) ∝ s/r with s(r) plotted in Figure 1. This reformulation is thus particularly powerful as it gives a direct scaling between surface density and radius with only a very weak dependence on the temperature profile of the disk via the term v 0 .
We now apply Equation 8 to TW Hya again assuming that the dominant particle size s is given by the wavelength of the observation (λ obs ). The resulting surface density points are well described by a steep power law function of approximately r −4 . Initially, this appears to be an unrealistically steep relation for the disk surface density as typical scalings for the minimum mass solar nebula are Σ ∝ r −1 . However, as these points fall close to the observationally derived exponential fall-off range for the measurements of CO in Rosenfeld et al. (2012), the data match this profile quite well when we allow Σ c to be a floating normalization factor (see Section 2). We find that a surface density normalization of Σ c ≈ 10 2.5 , approximately 3 orders of magnitude larger than Σ c from Rosenfeld et al., adequately matches our derived surface density points. The derived surface density profile with error estimates is shown in Figure 3. We find that the particle sizes we consider in the outer disk are indeed still well within the Epstein regime as was assumed in our derivation of the disk surface density profile.
We note that this surface density estimate depends on the age of the disk which may be uncertain by as much as a factor of 2 (e.g. Pecaut et al. 2012). For TW Hydra the surface density normalization, Σ c , ranges from ∼ 200 at 3 Myr to ∼ 800 at 10 Myr. We note that the upper end of this age range will lead to a derived disk surface density that is Toomre-Q unstable, however, due to the order using the three observed disk sizes r(s) with radius error bars of ± 10 AU and their corresponding errors in surface density. The surface density profile at these radii are well fit by an r −4 power law shown by the solid black line. The normalized surface density profile is shown in blue with the corresponding shading region indicating the roughly normalized surface density profile with the inclusion of the radial error estimates. The normalized surface density profile is an r −1 power law interior to the critical radius (rc) of 30 AU and is then described by an exponential fall off at radii larger than rc. We find that a surface density normalization of Σc ≈ 10 2.5 adequately matches our derived surface density points.
of magnitude nature of this derivation those numbers could be revised through a more detailed calculation. This trend generally holds such that the derived disk surface density increases with an increase in estimated disk age.
We can now convert our surface density profile to a derived gas mass as well as a mass accretion rate for TW Hya. We find that our derived gas mass of approximately 0.05 M is consistent with the observational lower limit of 0.05 M as derived from HD measurements of the warm gas (Bergin et al. 2013). We remark that this agreement implies no low-density cool gas and note that there is room within the errors for there to be a comparable amount of cold and hot mass present within the disk, however, to avoid being Toomre-Q unstable the mass cannot increase by a substantial amount.
Following the discussion in Section 2, with the inclusion of rough radius error estimates, we find a mass accretion rate ofṀ ∼ 4 × 10 −9 − 2 × 10 −8 M yr −1 . Compared to the measured value of ∼ 1.5 × 10 −9 M yr −1 (Brickhouse et al. 2012) our derived accretion rates give a larger value more in line with measured accretion rates for younger systems than the derived disk accretion rate using CO surface density observations as in Rosenfeld et al. (2012). We note that TW Hya has an inner disk gap at 1 AU . It is possible that this gap was formed by the presence of an accreting protoplanet which could be allowing only ∼ 10% of the accretion flux onto the star (Papaloizou & Nelson 2005;Najita et al. 2007;Zhu 2015;Eisner 2015).
We also note that, given our estimated disk age, our derived surface density profile is Toomre-Q unstable for the higher region of our estimated errors. We note, however, that assuming a shorter age of the system (∼ 3 Myr) gives a stable disk surface density profile for all estimated errors.

Particle Growth and the Dust-to-Gas Ratio
Given that TW Hya is an evolved disk in dynamical equilibrium, we can use an approximation of the growth timescale to determine the dust-to-gas ratio in the outer disk through allowing t grow = t disk as discussed in Section 3. We derive our growth timescale by first considering the growth rate in terms of particle size (s): where ρ p is the density of the particles which we take to be ∼ ρ g f d where f d is the dust-to-gas ratio and σ = πs 2 where s is the size of the dominant particle at that radius r. In the outer disk we assume that the relative particle velocities ∆v are due to turbulence which we approximate as: (Ormel & Cuzzi 2007) where t L is the overturn time of the largest eddies which we take to be the orbital period, t η = Re −1/2 t L , τ s,1 is the stopping time of our dominant particle size, τ s,2 is the stopping time of the particle size that contributes to growth, and α is the parameter of ignorance for turbulence in a disk which we take to be a standard value of 10 −3 (e.g. Rafikov 2017). This ∆v is valid for particles that are tightly coupled to the gas. We note that direct observations of disk turbulence (Hughes et al. 2011a;Teague et al. 2016) do not yet probe the disk midplane.
We are now able to solve for the dust-to-gas ratio in the disk through setting the growth timescale (t grow ) equal to the age of the system (t disk ) and solving for f d . We do this while assuming that the particles grow through perfectly efficient sticking with particles that are similar in size. We make this approximation as most of the dust surface density will be in the largest grains such that these particles dominate particle growth. Given this assumption, typical values of ∆v for our TW Hya calculation are ∼ 3 − 10 cm s −1 . For TW Hya this exercise gives an average dust-to-gas ratio of f d ≈ 10 −3 in the outer disk where we observe the disk dust lines.
As a consistency check we can compare this derived dust-to-gas ratio with the observed dust surface density profile described in Andrews et al. (2012). We note that the total dust mass and the dust surface density profile are model dependent quantities with values that vary across the literature (Calvet et al. 2002;Andrews & Williams 2005;Hughes et al. 2008Hughes et al. , 2011bKamp et al. 2013;Hogerheijde et al. 2016). We thus choose the Andrews et al. (2012) model as it is based off of 870 µm emission which is within the range of our dust line observations.
In their model, the dust surface density profile is well described by a shallow power law Σ d ∝ r −0.75 until the emission falls off at roughly 60 AU. This profile is normalized such that Σ d = 0.39 g cm −2 at 10 AU. When we divide this dust surface density profile by our derived total surface density profile we find an average dust-togas ratio of ∼ 10 −2.5 which is roughly consistent with our derived dust-to-gas ratio above.
We further compare our dust surface density points at the three dust line locations using our derived dust-togas ratios with the best fit model from Andrews et al. (2012) in Figure 4. We find a dust surface density that is systematically lower than the best-fit model from the 870 µm emission, however, without error bars on this observationally derived model it is difficult to know if our results are discrepant. We determine an extremely rough estimate of the uncertainty of this model by propagating uncertainties on the mass opacity, flux calibration (∼ 10 % systematic uncertainty), and distance (a few percent uncertainty). The largest source of uncertainty in the Andrews et al. (2012) model is the mass opacity. We adopt an uncertainty on this value of roughly ∼ 99 % as the range in assumed mass opacities for dust at ∼ 870 µm ranges from ∼ 10 −2 − 10 0 in the relevant literature (e.g. Hildebrand 1983;Pollack et al. 1994;Henning & Stognienko 1996;Andrews & Williams 2005;Andrews et al. 2012). The error on the mass opacity thus dominates the uncertainty in the dust measurement. We note that our points represent the dust surface density comprised of the large particles in some set size distribution and that including the contribution of smaller grains to the total surface density should increase our answer by a few tens of percent.
We also compare to dust depletion in a time evolving disk as shown in Birnstiel & Andrews (2014) (see their Figure 4a) where they find a solid depletion of up to two orders of magnitude in their fiducial model for the outer regions of late stage disks. We therefore find that our estimate of dust depletion is consistent with their model results for evolved disks.
In conclusion, by considering particle drift we derive a total disk surface density for TW Hya that is roughly 3 orders of magnitude larger than the surface density de-  Figure 4. The dust surface density of TW Hya (points) determined from our derived dust-to-gas ratio in conjunction with our dust line derived total surface density profile. These points can be compared to the best-fit dust surface density profile (black line) from Andrews et al. (2012) for TW Hya derived using 870 µm emission. We find a dust surface density that is systematically lower than the best-fit model from the 870 µm emission. We provide a rough estimate of the uncertainty surrounding the dust surface density profile from Andrews et al. (2012) (gray shaded region, see text) and find that our results are not discrepant within the assumed error.
rived from CO observations through considering particle drift. We further derive a dust-to-gas ratio of approximately 10 −3 in the outer disk (∼ 25 − 60 AU).

ICE LINES
We can now further our discussion of the disk surface density through a careful treatment of the molecular ice lines. The three ice lines that are most frequently calculated are H 2 O, CO 2 , and CO as these species are considered to be in relatively high abundance (Öberg et al. 2011). To date, the CO ice line is the most readily observed due to its large radial distance from the central star (Qi et al. 2013). We thus primarily focus on the CO ice line in our calculations. However, when we expand our discussion to consider a larger parameter space of disk parameters, we further extend our arguments to the other volatile species as well.
While we consider a passively heated disk (see Section 2), we note that the water ice line may be impacted by accretion heating. We include some plots of the H 2 O ice line for reference, but these should be viewed with caution.
Our discussion focuses on the midplane ice lines as opposed to the surface ice lines frequently discussed in the literature (Blevins et al. 2016;Schwarz et al. 2016). As a result we will compare our derived ice lines to observations that probe the midplane ice line directly (via the C 18 O line; Qi et al. 2015;Dartois et al. 2003) or indirectly (via the N 2 H + ion; Qi et al. 2013Qi et al. , 2015.
There are three pieces of physics that we consider in our ice line calculations: particle adsorption, desorption, and drift.

Volatile Adsorption and Desorption
The classic ice line calculation balances adsorption and desorption flux onto a grain to determine the ice line radius (Hollenbach et al. 2009;Öberg et al. 2011). We refer to this ice line as the 'classical ice line'. Following Hollenbach et al. 2009, these two fluxes are quantified as: where n i is the relevant gas density species, c s is the sound speed, N s,i ≈ 10 15 sites cm −2 is the number of adsorption sites per volatile per cm 2 , ν vib = 1.6 × 10 11 (E i /µ i ) s −1 is the molecules vibrational frequency in the surface potential well, E i is the adsorption binding energy in units of Kelvin, and f s,i is the fraction of the surface adsorption sites that are occupied by species i (which we take to be unity). Finally, we assume that T grain = T , meaning that the dust and gas have the same temperature in the disk mid-plane. Balancing Equations 11 and 12 allows for us to solve for the freezing temperature of a species as a function of radius for a given disk surface density profile. We then locate the classical ice line by finding the disk radius where the molecular freezing temperature is equal to the disk temperature. This self-consistent method allows us to determine how the ice line location changes as a function of both disk surface density and temperature.

The Influence of Particle Drift
Particle drift, as described in Section 3, influences the location of the ice lines (Piso et al. 2015(Piso et al. , 2016. This is because particles that drift faster can cross the ice line before desorbing, thus potentially moving the location of the ice line inwards. The drift ice line location can be calculated by setting the desorption timescale, (where µ i is the molecular weight of the desorbing species) equal to the drift timescale (the radial location of a particle divided by Equation 6) and solving for the desorption distance, r des (as verified by time evolving calculations in Piso et al. 2015). This is done analytically for the small stopping time approximation in the Appendix of Piso et al. (2015). We extend this calculation to the τ s > 1 regime through making no approximations in Equation 6 in regards to τ s . Instead, we use the two distinct stopping time expressions given both Stokes and Epstein regimes (Equation 5). We use the SciPy subroutine fsolve (Jones et al. 2001) to solve for radius through balancing the drift timescale with the desoprtion timescale (t drift = t des ) where these timescales are given by t drift = |r/ṙ|, whereṙ is given by Equation  6, and Equation 13. We consider the radial dependence of each term without approximation, which allows us to derive a completely self-consistent estimate of r des for any arbitrary surface density and/or temperature profile. We compare the difference between the drift ice line and the classical ice line for our new self-consistent solver and the analytic solver from Piso et al. (2015) in Figure  5. We use the disk temperature and surface density profiles from Piso et al. (2015) and find that particles smaller than ∼ 10 cm (τ s ∼ 0.5) do not drift past the CO ice line. Particles larger than ∼ 10 cm do experience drift, with the maximum drift reached at a τ s ≈ 1 that then reaches a near constant value at higher stopping time. We note that, while small particles do not drift to the same radial location as large particles, accretion of disk gas has been shown to pull smaller particles in past the ice line as shown in Piso et al. (2015).
The nearly constant drift distance at particle sizes with a stopping time greater than 1 occurs due to the interplay between two different pieces of physics. The first is that the larger the particle size the longer it takes for that particle to desorb and the further it can drift past the ice line. The second is that the maximum drift velocity occurs for particles with τ s = 1. The decrease in drift velocity past a stopping time of one is offset by the increase in particle size until the regime changes from Epstein to Stokes. The drift timescale for a particle in the Epstein with τ s >> 1 is t drift = |r/ṙ| ∼ τ s /ηΩ ∼ ρ s s/ρc s η. The desorption timescale is simply proportional to the particle size (t des ∝ s) with a desorption constant such that t des ∼ C des s. When we equate these two timescales we find that the dependence on particle size cancels out and that the drift distance is therefore independent of this quantity. We find that this approximation is roughly accurate for particles in the Epstein regime with τ s close to or greater than 1.
Due to this effect, we consider the "drift ice line" to be the location for which a τ s = 1 particle starting at the classical ice line location drifts and then desorbs for situations in which we know the disk surface density. If we do not know the disk surface density we use a large sized particle such that it desorbs at the maximum drift radius for a range of disk parameters. We numerically solve for the drift distance without making any further approximations in all further discussions of drift. The existence of a maximum drift distance found in our self-consistent solution allows us to predict the circumstances in which  Piso et al. (2015) (blue, solid) and our new extended solver (green, dashed). We find great agreement for particles with a stopping time less than unity. Right: The distance that a CO particle is able to drift past the classical CO ice line before desorbing using the analytic method from Piso et al. (2015) (blue, solid) and our new self-consistent solver (green, dashed). Using our extended solver we see that, for this comparison case, particles smaller than ∼ 10 cm do not drift past the ice line and are not shown. Particles larger than ∼ 8 cm do experience drift in this case, with the maximum drift reached at a particle size slightly larger than τs ≈ 1 that reaches a near constant value at larger stopping times.
drift affects the ice line locations.
We find that the radial location of the classical ice line strongly impacts the amount in which drift will play a role. In particular, we find that the process of drift becomes important when the classical ice line is located where larger particles have a stopping time close to 1 as τ s = 1 particles drift the fastest and larger particles take longer to desorb (see Figure 6). The size of the τ s = 1 particle depends on the overall surface density profile of the disk. Drift will therefore matter the most in determining the true ice line location at the disk radial location where the particles with τ s = 1 reach a maximum size. If the classical ice line of a particular molecule is located near this point then drift will affect the location of the ice line radius pushing it inwards. This is shown in Figure 6 through a consideration of our derived surface density profile for TW Hya and the standard profile used for the minimum mass solar nebula (Σ ∝ 1700 g cm −2 × r −3/2 ) (Öberg et al. 2011; Chiang & Youdin 2010).

Application to TW Hya
Given our derived total surface density profile for TW Hya, the last parameter needed to derive the CO ice line location is the CO abundance. We use the measured CO surface density from Schwarz et al. 2016 in conjunction with our derived total surface density profile to uncover an approximate CO abundance of 10 −7 n H from ∼ 10 − 60 AU. This value is consistent with the average upper limit of 10 −6 n H found in Schwarz et al. (2016)  [cm] TW Hydra MMSN Figure 6. More massive disks reach a peak in the size of τs = 1 particles at radii further from their central star as is the case for TW Hya as compared to the minimum mass solar nebula. Shown here is the size of τs = 1 particles in the disk as a function of radius for two different surface density profiles: our TW Hya surface density profile (blue, solid) and the commonly used minimum mass solar nebula surface density profile (green, dashed). The τs = 1 particles are the largest particles that are still well-coupled to the gas. Drift affects the ice line locations the most when the classical ice lines occur close to the peak in these plots.
of the warm gas (Bergin et al. 2013).
This reduction of 3 orders of magnitude from the measured abundance in the ISM (10 −4 n H ) demonstrates that our model requires a global depletion of volatile carbon in TW Hya. Cleeves et al. (2015) infer a similar level of depletion using observations of CO and modeling by Du et al. (2015) supports this conclusion. The blue lines are the dust lines solved by assuming that the drift timescale is equal to the age of the system. The blue lines adequately reproduce the observed disk radial scale of TW Hya at various wavelengths. The solid red line is the classical CO ice line solved by balancing the adsorption and desorption flux onto a grain. This line is in agreement with the observed CO ice line of ∼ 30 AU using N2H + (Qi et al. 2013). The dashed red line is the CO drift ice line for a τs = 1 particle at a radius of ∼ 15 AU which we find to be in close agreement with the ice line derived from C 18 O measurements in Nomura et al. (2016), suggesting that C 18 O is a sensitive probe of the CO drift ice line.
We can now calculate both the classical and maximum drift CO ice lines. For our calculations we use a binding energy of E i /k ∼ 850 K (Aikawa et al. 1996) 1 . For TW Hya, we derive a classical CO ice line of ∼ 30 AU and a drift ice line of ∼ 15 AU.
We find that our derived classical ice line is in agreement with the ice line determined observationally via the N 2 H + ion of ∼ 30 AU (Qi et al. 2013). We further find that our derived drift ice line is in good agreement with the CO ice line as derived through the C 18 O J = 3 − 2 line where a decrement of CO emission is observed until ∼ 10 AU (Nomura et al. 2016). Given the discrepancy between the ice line as determined by N 2 H + and C 18 O and the corresponding agreement we find with our theoretical model, we note that the C 18 O J = 3 − 2 line may be the most sensitive probe of a disk's midplane ice line, particularly when particle drift effectively smears the ice 1 We note that observations discussed in Schwarz et al. (2016) indicate a CO binding energy closer to E i /k ∼ 960 K for TW Hya. Using this binding energy we determine that the classical CO ice line is located at ∼ 23 AU and the drift ice line is located at ∼ 10 AU. These derived values are also reasonably consistent with the observed ice lines of TW Hya (see text).
line out over a relatively large radial scale. This could be because the presence of N 2 H + requires a lack of CO gas in the disk and may thus require that all particle sizes have frozen out (i.e. the classical ice line).
We can now put together our complete model for the dust and ice lines in TW Hya as shown in Figure 7.

OBSERVATIONAL DIAGNOSTICS
We now propose three different observational tests of our model framework beyond the application to TW Hya. The first test is whether or not other disks have dust lines that can be well fit by a normalized form of the physically motivated surface density profile in Equation 3 for disks where r c is constrained by CO observations. This is a particularly powerful test of concept if we find dust lines in other disks that are located in the powerlaw region of the surface density profile interior to the exponential fall-off. The second test is whether we can derive a surface density profile from the disk ice lines that is consistent with the surface density profile derived from disk dust lines. The third test is whether or not the dust and ice lines scale oppositely, as a function of surface density, across a large sample of disks.
Our second and third tests rely on an understanding of both classical (i.e. classical regime) and drift (i.e. drift dominated regime) ice lines as described in Section 4. The particular assumptions and/or prior knowledge of a disk will inform which of the two regimes is relevant as discussed in Section 5.2.

Test 1: Surface Density from Disk Dust Lines
The first test of our model is whether or not other disks have dust lines that can be well fit by a normalized form of the physically motivated surface density profile in Equation 3. In the era of ALMA, multiple disks will have high resolution data at several wavelengths -the ideal data set for this fundamental test. Following the technique described in Section 3 we can convert these dust lines into a surface density measurement and see if these measurements scale with the observationally derived surface density profiles.
In the case of TW Hya, the disk dust lines were near the exponential fall-off point, which unexpectedly (due to the steep slope of the derived power law) matched previous surface density models with the application of a normalization factor. For other disks there is no reason to think that the disk dust lines will necessarily fall in the same exponential fall-off region of the disk surface density profile. If a disk has a critical radius further from the star then it is increasingly likely that we will be able to detect dust lines interior to this point.
This test is thus a powerful test of concept that will allow us to see whether other disk sizes can be well fit through considering drift as the primary driver of disk dust lines and using these derived surface density values to normalize previous observationally fit disk surface density profiles.
• Test 1 Observational Requirements: This test can be best carried out on disks with CO observations, such that there is a derived critical radius, and observations of dust emission at several wavelengths such that the disk dust lines are known.

Ice Line Regimes
It is now helpful for us to distinguish between the drift and classical ice lines from observations of CO or other ice line tracers. Here we provide several rules of thumb to aid in the use of ice lines as useful diagnostics. However, we note that these aids should be used initially and then verified for self-consistent results.
Our proposed tests using these ice lines differ as the classical ice line depends on the CO abundance while the drift ice line does not. This makes the drift ice line particularly useful as it relies on fewer assumptions while the classical ice line retains its usefulness when considered across a sample of disks.
As mentioned in Section 5, there are two regimes that dictate when the drift or classical ice lines are relevant: the drift-dominated regime and the classical regime. In the drift-dominated regime the drift ice line is interior to the classical ice line and thus there should be only be CO freeze out exterior to this point with complete freezeout occurring exterior to the classical ice line, making the drift ice line the most interior detectable ice line in the disk. In the classical regime particle drift does not happen quickly enough for particles to cross the ice line without desorbing and the classical ice line is observed.
To interpret the surface density from the disk ice line (see Section 5.3) we need a priori knowledge of the relevant ice line regime in the disk. Without knowledge of the disk surface density profile we can roughly determine the correct regime through an analysis of the ice line dependence on elemental abundance and disk temperature as well as a comparison of the ice line location with respect to the disk's critical radius.
One simple metric in determining the ice line regime is the knowledge that the ice line location is preferentially dominated by drift for disks that have large particles with a stopping time of unity near the ice line location. As seen in Section 4.2, particles can drift further without desorbing where the τ s ≈ 1 particle size is large. This size increases at a given radius with increased disk critical radii and surface density as demonstrated in Figure  8 where we calculate the τ s = 1 particle size at 30 AU for a disk with the same temperature profile as TW Hya for a range of surface density normalization factors and [cm] Figure 8. The τs = 1 particle size for a disk using our derived temperature profile for TW Hya (see Equation 1 with T0 = 82 K) at a radius of 30 AU (the observed classical ice line radius). We vary the disk critical radius and the surface density normalization and find that the ice line is likely to be dominated by drift for disks with large critical radius and high surface density normalization. disk critical radii. We find that an ice line observed interior to the critical radius is likely dominated by drift for massive disks with large critical radii.
Another parameter that determines the ice line regime is the amount of radiation that the disk receives from its host star (i.e. the stellar luminosity). The importance of this parameter is shown in Figure 9 for generally assumed molecular abundances (n CO = 1.5 × 10 −4 n H , n CO2 = 0.3 × 10 −4 n H , n H2O = 0.9 × 10 −4 n H ; Pontoppidan 2006).
We find that we are in the drift dominated regime for disks with moderate densities and high temperatures across all three molecular species. At very high disk densities we find that the drift ice line is irrelevant and only the classical ice line can be observed, as a large increase in disk density moves the classical ice line inwards such that the disk temperature at the ice line is higher and the desorption timescale is shorter (i.e. the increased temperature near the star overwhelms the effect of the increase in disk density, see Equation 12). As the disk density is not known a priori we posit that the drift ice line is likely relevant for disks with high stellar irradiation.
We further find that, given a molecular depletion as in the case of our derived disk parameters for TW Hya and may be true for the majority of disks, drift sets the true ice line location across a wide range of parameter space as shown in Figure 10. Drift is more important when the CO abundance is smaller because the decrease in CO abundance moves the classical ice line further The fractional difference between the classically derived ice line r and the drift ice line r des as a function of stellar luminosity and disk density (top) as well as disk density and temperature (bottom) for the following molecular abundances: nCO = 1.5 × 10 −4 nH, nCO 2 = 0.3 × 10 −4 nH, nH 2 O = 0.9 × 10 −4 nH (Pontoppidan 2006). Drift is most important for CO2 and generally increases at moderate densities and high temperatures. The y-axis label T0 refers to the temperature normalization for Equation 1 and we can convert from this temperature normalization to stellar luminosity using Equation 2. For stellar luminosities below 10 −2 L drift does not play a role in determining the ice line locations.
from the star while the drift ice line remains unchanged with the change in abundance. Therefore, we find that disks depleted in CO should have ice lines that are determined by drift. This is further demonstrated in Figure  11 where, across a wide range of molecular CO abundances, the drift ice line is constant and interior to the classical ice line given our derived disk parameters for TW Hya. It is therefore a relatively safe assumption to consider the CO ice line location to be determined by drift for disks that are hot, depleted in CO, and/or have an observed ice line interior to a large critical radius and appear to be relatively massive. Disks for which the opposite is true are in the classical regime and their ice lines should be theoretically treated as classical ice lines.

Test 2: Surface Density from Disk Ice Lines
The second test of our model is whether or not the surface density profile derived from a consideration of the disk ice line matches the surface density profile derived from the disk dust lines. We now discuss how empirical knowledge of a disk's molecular ice line can be used to independently determine the total surface density profile of a disk. This can be done as both the drift and classical ice lines depend on the disk surface density.

Drift Ice Line Surface Density Determination
If the disk in question is in the drift dominated ice line regime (see Section 5.2), solving for the disk surface density profile simply involves setting Equation 6 for τ s = 1 (where particle drift reaches a maximum) equal to Equation 13 and solving for surface density given an assumed characteristic particle size. This surface density measurement can then be used to normalize the total disk surface density profile, thus giving us an independent es- timate that can be compared to the profile derived via the disk dust lines.
While this method requires the assumption of a characteristic particle size, we note that a wide range of particle sizes will result in the same drift ice line location (see Figure 5) and thus, the derived surface density profile will not be sensitively affected by this value. For computational purposes we use a particle size of 1 m such that its stopping time exceeds unity for a wide range of disk parameter space. As is the case for all particles close to or larger than a stopping time of unity, these particles should desorb after drifting to the maximum drift ice line.
Consideration of the drift ice line in determining disk surface density has the advantage of not needing a measure of the molecular abundance in the disk. Regardless of CO abundance this location will be constant as it is dependent on disk surface density and not the CO surface density. Thus, a sensitive probe of the CO ice line (i.e. the C 18 O J = 3 − 2 emission) should be able to detect the uptick in CO emission past this point and provide a probe of disk surface density from this measurement alone without further assumption.

Classical Ice Line Surface Density Determination
If the disk in question is in the classical ice line regime (see Section 5.2), the surface density can be derived by balancing Equation 11 and Equation 12. Unlike the drift dominated regime, this calculation requires an assumed molecular abundance. While these uncertainties may diminish the robustness of this test, we can make the simplifying assumption that molecular abundances are roughly constant across a single stellar type. This could be a reasonable assumption if molecular abundances are primarily shifted from ISM values via photochemical processes. If we accept this assumption as true, we can determine the disk surface density profile for disks of a stellar type for which one member has a relatively well determined molecular abundance. We note that, for this method, a factor of 2 change in the ice line location will change the derived surface density by a factor of 3.
The relative behavior of the classical ice line across many disks will also improve this test and may indeed allow the potential for deriving the CO abundance in conjunction with Test 1.
• Test 2 Observational Requirements: This test can be carried out on disks with accurate ice line measurements (either via C 18 O or N 2 H + ). To carry out this test more precisely it is also useful to have observationally derived CO surface density profiles for these disks.

Test 3: Disk Dust and Ice Line Scalings
The third test is whether or not the dust and ice lines scale oppositely across a large sample of disks. Our model predicts that this will happen as disks with dust lines at larger radial scales should have ice lines located at shorter radii across a particular disk temperature profile and molecular abundance. This arises as larger particles become well coupled to the gas as disk surface density increases, thus causing a given particle size to have a dust line located further out in the disk. Conversely, increased surface density increases the adsorption flux onto a grain such that the freezing temperature is hotter and the ice line is moved closer to the star.
To clearly see how the dust and ice lines in a disk scale as a function of surface density alone, we vary the disk surface density for a set disk temperature profile. Using the empirical evidence that the Earth is not formed of water ice, we consider the case of a passively irradiated disk that is normalized such that the H 2 O ice line is outside of 1 AU as shown in Figure 12 as our fiducial temperature profile. We take this temperature profile to be roughly representative of a young sun like star. This profile follows Equation 1 with a derived T 0 ≈ 210 K.
This profile is derived through the use of an assumed canonical H 2 O abundance of 0.9×10 −4 n H (Pontoppidan 2006). While this assumption is sufficient for our illustrative example, we again note that disk molecular abundances may be poorly constrained. However, let us naively assume that molecular abundances vary from ISM values as a function of photochemistry such that they are constant across a particular stellar type and thus a particular temperature profile. Thus, using this temperature profile we can derive the trends that we expect to see as we vary the disk surface density.
We find that for all three of our molecular species: CO, CO 2 , and H 2 O, the radial extent of the disk increases with increasing surface density while the ice line location decreases as expected (see Figure 13). For the CO and H 2 O ice lines, the drift distance past the classical ice line is negligible while for the CO 2 ice line drift plays a small role that decreases in importance with increasing surface density. We note that while the importance of the drift ice line will vary with disk parameters, the general trend should hold.
This trend remains true across a range of molecular abundances as long as they are held constant across a sample of disks. Here we have assumed the following abundance values for our molecular species: n CO = 1.5× 10 −4 n H , n CO2 = 0.3 × 10 −4 n H , n H2O = 0.9 × 10 −4 n H (Pontoppidan 2006). We note that lower molecular abundances result in lower freezing temperatures (see Section 4.1) which, depending on the temperature profile of the disk, pushes the classical ice lines further from the star with the reverse being true for higher molecular abundances. This trend is therefore a key observational diagnostic that will be confirmed if, across a constant molecular abundance, disks with dust lines at larger radii have ice lines closer to their star. We also note that, if drift determines the ice line location (see Section 5.2) we can observe these trends without the need for an assumed molecular abundance.
• Test 3 Observational Requirements: To best perform this test there needs to exist a significant sample of disks with observations of both the dust and ice lines.

SUMMARY & DISCUSSION
We propose a novel method to derive disk surface density through the consideration of disk dust and ice lines. To derive this method we adopt an agnostic point of view in regards to disk surface density which we do in particular response to the uncertainties that accompany typical observational tracers. This method relies on the assumption that, at late stages of evolution, the growth timescale, drift timescale, and the lifetime of the disk are all equal for the dominant particles at a dust line location (s = λ obs ).
While other work finds that these timescales are equal they are often unable to match the disk surface density. We therefore make the assumption that these timescales are equal and use this to determine the surface density profile without evoking other observational tracers. These assumptions allow us to self-consistently derive a disk surface density profile as well as a dust-to-gas ratio in the outer disk without the further assumption of a given molecular abundance. H 2 O Particle locations in disk Classic iceline Figure 13. The disk radial scale for different particle sizes (black, solid), the classically derived ice line locations (blue, dashed) and the drift ice line (red, dotted) for CO (top; with particle sizes labeled), CO2 (middle), and H2O (bottom) as a function of surface density. The dust line locations increase dramatically with increased surface density while the ice line location decreases.
We apply our modeling technique to our fiducial disk TW Hya. We find that our derived surface density profile and dust-to-gas ratio are consistent with the lower limits found through measurements of the HD gas in the disk (Bergin et al. 2013;Schwarz et al. 2016). Using our derived surface density profile we uncover a theoretical estimate of the disk accretion rate that is more closely aligned with the measured accretion onto TW Hya. We further find that our theoretical classical and drift ice lines have clear observational analogues where the classical ice line is predictive of the ice line found via N 2 H + emission (Qi et al. 2013) and the drift ice line is probed by the C 18 O emission (Nomura et al. 2016). We conclude that the ice line derived through observations of C 18 O emission may be more sensitive to the extent of the drift ice line in the disk and thus may be the best observational method to test our model's assumptions. We note here that our method highlights the likelihood of detecting multiple ice lines for a given molecular species as there will be relatively large regions in the disk where some of the material has desorbed and some of the material has not. Furthermore, if the banded structure in the ALMA images of HL Tau and TW Hydra reflect substantial variations in surface density then our model should still hold. However, ice lines may be preferentially found in less dense bands in the disk where drift slows.
We next consider a large range of theoretical disk parameter space and uncover three observational tests of our model. The first test is whether or not the dust lines of other disks, once converted to surface densities, can be matched with a previously derived normalized surface density profile for the disk. For TW Hya, the disk dust lines fall near the exponential cut off region of the disk but there is no reason that this needs to be true across many objects and thus provides a powerful test of the model. In the age of ALMA, dust lines can be determined observationally for a significant sample to protoplanetary disks which will provide the ideal observations necessary for this test.
The second test is whether we can derive a surface density profile from the disk ice lines that matches the surface density profile derived via disk dust lines. The third and final test is whether or not disk dust and ice lines scale oppositely, as a function of surface density, across a large sample of disks. Disk ice lines have been observed for an increasing number of disks such that the last two tests could be carried out in the near future with the aid of facilities such as ALMA.

ACKNOWLEDGEMENTS
This work was done as a part of the Kavli Summer Program in Astrophysics: Exoplanetary Atmospheres 2016. We would like to thank Jonathan Fortney and Pascale Garaud for providing us with the opportunity to work together on this project and for organizing a wonderful workshop. We thank James Owen, Michael Rosenthal, Catherine Espaillat and Laura Pérez for their useful comments and insightful discussion. This material is based upon work supported by the National Science Foundation Graduate Research Fellowship under Grant DGE1339067. RMC acknowledges support from NSF grant number AST-1555385.