The Evolution of Inner Disk Gas in Transition Disks

Investigating the molecular gas in the inner regions of protoplanetary disks provides insight into how the molecular disk environment changes during the transition from primordial to debris disk systems. We conduct a small survey of molecular hydrogen (H$_2$) fluorescent emission, using 14 well-studied Classical T Tauri stars at two distinct dust disk evolutionary stages, to explore how the structure of the inner molecular disk changes as the optically thick warm dust dissipates. We simulate the observed HI-Lyman $\alpha$-pumped H$_2$ disk fluorescence by creating a 2D radiative transfer model that describes the radial distributions of H$_{2}$ emission in the disk atmosphere and compare these to observations from the Hubble Space Telescope. We find the radial distributions that best describe the observed H$_2$ FUV emission arising in primordial disk targets (full dust disk) are demonstrably different than those of transition disks (little-to-no warm dust observed). For each best-fit model, we estimate inner and outer disk emission boundaries (r$_{in}$ and r$_{out}$), describing where the bulk of the observed H$_2$ emission arises in each disk, and we examine correlations between these and several observational disk evolution indicators, such as n$_{13-31}$, r$_{in,CO}$, and the mass accretion rate. We find strong, positive correlations between the H$_2$ radial distributions and the slope of the dust SED, implying the behavior of the molecular disk atmosphere changes as the inner dust clears in evolving protoplanetary disks. Overall, we find that H$_2$ inner radii are $\sim$4 times larger in transition systems, while the bulk of the H$_2$ emission originates inside the dust gap radius for all transitional sources.


INTRODUCTION
Protoplanetary disks (PPDs) provide the raw materials for the formation of stellar systems (Brown et al. 2009;Woitke et al. 2009;Dullemond & Monnier 2010). Planet formation occurs near the midplane of a PPD, where column densities and optical depths are high (Trilling et al. 2002;Armitage et al. 2003), making it difficult to directly observe the material involved in the formation process (Kominami & Ida 2002). Current understanding of the formation and evolution of planetary systems in gaseous disks comes from studies of molecular content above or near disk midplanes, which place limits on the composition and density distribution of the gas and dust content in the inner (r ≤ 10 AU) planetforming regions (Agúndez et al. 2008;Carr & Najita 2008, 2011Salyk et al. 2008Salyk et al. , 2011aWoitke et al. 2009;Willacy & Woods 2009;Heinzeller et al. 2011;Najita et al. 2011).
"Transition" disks refer to a class of PPDs with an optically thick outer zone but an inner region significantly depleted of dust grains (Sato & Nakagawa 1999;Calvet et al. 2002;Salyk et al. 2009) and are traditionally identified by the deficiency in near-infrared (IR) flux and steep rise of far-IR flux in the observed SED (Strom et al. 1989;Calvet et al. 2002Calvet et al. , 2005Espaillat et al. 2007a). Several theories ex-ist for how dust gaps are opened in transition disks, including photoevaporation (Hollenbach et al. 1994;Alexander et al. 2006Alexander et al. , 2013Alexander & Armitage 2007;Gorti et al. 2009), dynamical clearing by protoplanetary systems (Calvet et al. 2002;Rice et al. 2003;Dodson-Robinson & Salyk 2011), and dust grain growth (Tanaka et al. 2005).
Discoveries of significant quantities of gas left within the dust gaps of transition disks (see Najita et al. 2003;Rettig et al. 2004;Salyk et al. 2007) and sharp "walls" between the thin and thick dust disk regions ) support the possibility of transition disks being carved out by giant planet formation and evolution (Salyk et al. 2009;Dodson-Robinson & Salyk 2011;Dong et al. 2014). The remnant gas disks provide constraints on the processes that create the final structure of planetary systems, such as the transfer of gas from the PPD to circumplanetary disks, potentially leading to growth of protoplanets (Lubow et al. 1999;Lubow & D'Angelo 2006;Ayliffe & Bate 2010;Beck et al. 2012). Additionally, the molecular atmosphere of transition disks may respond to the dynamical perturbations caused by the presence of giant protoplanets and can lead to potentially observable effects, such as line asymmetries and distortions in near-IR CO emission profiles (Regály et al. 2010). The strength of molecular emission originating from the inner radii of PPDs is dependent on the gas temperature, density, and degree of grain growth (Salyk et al. 2011b). Molecular line surveys therefore provide the opportunity for a broad examination of the gas distributions in circumstellar environments .
Molecular hydrogen (H 2 ) has been measured to be ∼ 10 4 times more abundant than any other molecule in the inner disks of PPDs (France et al. 2014a). Depending on the density, H 2 can survive at temperatures up to 5000 K (Williams 2000) and selfshields against UV radiation, making it robust to both collisional-and photo-dissociation (Beckwith et al. 1978;Beckwith & Zuckerman 1982;Beckwith et al. 1983). Molecular hydrogen provides a diagnostic for the spatial and structural extent of the warm molecular surface of PPDs (Ardila et al. 2002;Herczeg et al. 2004;Yang et al. 2011). While photo-excited H 2 does not interact strongly with evolving protoplanets, it traces the underlying distribution of gas at planet-forming radii (Ardila et al. 2002;Herczeg et al. 2004Herczeg et al. , 2006France et al. 2012a). However, H 2 has proven difficult to observe in PPDs: cold H 2 (T ∼ 10 K) does not radiate efficiently because it has no permanent dipole (Sternberg 1989), so IR ro-vibrational transitions are weak, making them difficult to observe from the ground. Therefore, studies of molecular material in disks typically rely on other tracers available in the near-and mid-IR, such as CO and H 2 O, to estimate the molecular disk environment and mass of the underlying H 2 reservoir in disks.
The strongest transitions of H 2 are found in the FUV (912 -1700Å), where dipole-allowed electronic transitions are primarily photo-excited ("pumped") by Lyα photons generated near the protostellar surface Schindhelm et al. 2012b). Warm H 2 (T 1500K) has a significant population in excited vibration (v = 1, 2) and rotation quantum states of the ground electronic band (X 1 Σ + g ) (Shull 1978). When a Lyα photon interacts with a warm H 2 molecule in the correct ground-state population [v,J], the H 2 molecule absorbs the photon, exciting it to vibration levels (v ′ → 0-4) of the first electronic band (B 1 Σ + u ). Since molecular hydrogen has strong (A ul ∼ 10 8 s −1 ; see Abgrall et al. 1993) electronic transitions in the FUV, the excited H 2 "immediately" decays back to the ground state, emitting a fluorescent photon, observed as an FUV emission line. The probability for an H 2 excitation-to-ground state transition to emit a photon with wavelength λ depends on the branching ratio of the allowed transitions to the ground electronic state. The brightest H 2 emission lines arise from excited states [v ′ , J ′ ] = [1,4], [1,7], [0,1], and [0,2], which have absorption coincidences with Lyα within 0 and 600 km s −1 of the Lyα line center, large oscillator strengths, and relatively low energy ground-state levels (Herczeg et al. 2002(Herczeg et al. , 2005. The set of emission lines produced in the [v ′ , J ′ ] → [v ′′ , J ′′ ] decay is refered to as a progression. Previous work on FUV fluorescent H 2 emission utilized basic profile fitting or small-sample parametric sets to estimate inner disk diagnostics, such as column density and temperature of the radiating molecular populations (see Herczeg et al. 2004;France et al. 2012a,b). In this study, we create 2D radiative transfer models of PPD atmospheres to reproduce observed FUV H 2 emission lines. The models simulate a disk with radial temperature and density distributions, which depend on physical parameters of the stellar system, such as the disk inclination angle and stellar Lyα radiation profile (taken from Schindhelm et al. 2012b). Using the four strongest H 2 progressions, we compare radiative transfer emission line models to the spectra of 14 CTTSs (8 primordial, 6 transition disks) observed with the Hubble Space Telescope (HST )/Cosmic Origins Spectrograph (COS) and Space Telescope Imaging Spectrograph (STIS). The goal of this modeling work is to examine the relationship between the evolution of warm dust in PPDs and the radial distribution of H 2 in the disk atmosphere. We aim to understand how the spatial distribution of warm H 2 relate to the structure of the dust disk and other well-studied molecular disk tracers, such as carbon monoxide (CO) and water (H 2 O).
In §2, we present the targets, observations, and selection criteria of H 2 emission features explored in this work. In §3, we describe the forward modeling process for estimating the warm H 2 disk radiation fields, and in §4 we analyze how the best-fit models are determined and define metrics used to quantify the evolution of H 2 radiation for each PPD. In §5, we discuss how the modeled radiation distributions of fluorescing H 2 evolve in PPDs, comparing our results with observable warm dust disk evolution, mass accretion rates, and additional inner disk molecular tracers. Finally, we summarize how the gas disk structure correlates with the dissipation of warm dust grains as PPDs evolve to debris disks in §6.

OBSERVATIONS AND H 2 EMISSION LINE SELECTION
We sample a large collection of HST -COS and HST -STIS (for TW Hya) FUV H 2 data to understand the relative changes in the radiation distributions of H 2 arising from the inner regions of primordial and transition disks. The observations were obtained through the DAO of Tau guest observing program (PID 11616; PI -G. Herczeg), the COS Guaranteed Time Observing program (PIDs 11533 and 12036; PI -J. Green), and HST Program GTO-8041 (PI -J. Linsky). The observations have been presented in previous literature (for examples, see Herczeg et al. 2006;Ingleby et al. 2011;Yang et al. 2011;France et al. 2012bFrance et al. , 2014bSchindhelm et al. 2012a;Ardila et al. 2013).
The medium-resolution G130M and G160M FUV modes of COS (Green et al. 2012) were utilized for all targets except TW Hya, which was observed with the E140M mode (1170 -1710Å) with the 0.5 ′′ × 0.5 ′′ aperture of STIS at a resolving power of 25,000 (see Herczeg et al. 2006). The point-source resolution for each mode on COS is ∆v ≈ 17 km s −1 with 7 pixels per resolution element (Osterman et al. 2011) and ∆v ≈ 12 km s −1 for the STIS E140M observing mode of TW Hya (Leitherer 2001). The COS data were smoothed by 3 pixels for analysis. The one-dimensional spectra of COS were produced using the CALCOS COS calibration pipeline, which were aligned and coadded using a custom software procedure (Danforth et al. 2010). The STIS data were reduced using the CALSTIS STScI reduction pipeline (Lindler 1999), with calibration lamp spectra obtained during observations to assign wavelength solutions. An example of the continuous far-UV spectrum of V4046 Sgr is shown in Figure 1.
Stellar properties, such as mass, accretion rate, and inclination angle are used to constrain the underlying model framework.
The HST -COS FUV spectra of all CTTSs from 1300 -1600Å reveal a suite of H 2 fluorescence features linked to Lyα-pumping. We chose to use the strongest transitions from the electronically-excited progressions [v ′ , J ′ ] = [0,1], [0,2], [1,4], and [1,7] for the purposes of studying the underlying general characteristics of the bulk gas disk. We sample 3 emission features from each progression. This gives us access to strong, non-blended emission lines that are well-defined from the FUV continuum, while balancing the CPU time required for detailed line profile analysis. We selected H 2 emission features by locating the strongest transitions for each progression, outlined by Herczeg et al. (2002). See Table 1 for the full outline of transitions chosen.

MODELING ANALYSIS
We create models of warm H 2 in PPDs to constrain the radial distribution of fluorescent H 2 emission in disk atmospheres. Our aim is to understand the relative changes in the H 2 distributions as we observe various stages of dust disk evolution. The fluorescent emission line shape and intensity depend on the physical conditions of the gas, while the observed line width depends predominantly on the disk inclination. We construct a physical model of the disk structure, motivated by the disk modeling analysis done by Rosenfeld et al. (2012).
The models make several basic assumptions on the disk properties: (a) the disk material orbits in Keplerian rotation around a central point mass, representing the stellar mass; (b) the H 2 fluorescence occurs in a shallow, warm layer on the disk surface; and (c) the level populations of warm H 2 that absorb the incident stel-lar Lyα radiation field are in local thermodynamic equilibrium (LTE). (a) implies that the gas disk mass is a small fraction of the stellar mass (M d /M ⋆ ≪ 1). Several studies have shown that the disk mass to stellar mass ratio (M d /M ⋆ ) < 1%, making this assumption plausible ). In the case of a binary system (i.e., V4046 Sgr), both stellar masses are represented as one central mass point. For (b), Herczeg et al. (2004) find that the warm H 2 disk layer interacting with the stellar Lyα to produced the observed fluorescence lines corresponds to mass column density of ∼ 10 −5 g cm −2 , which is a much smaller mass column density predicted to be within 1 AU by D' Alessio et al. (1999). This suggests that the Lyα-pumped fluorescent emission originates from a tenuous layer of warm H 2 on the disk surface and supports a purely radial thermal distribution T (r). For (c), combination of collisional excitation and radiative de-excitation is assumed to be in equilibrium to keep the H 2 gas near the disk surface at warm temperatures (T > 1000 K; Nomura et al. 2005Nomura et al. , 2007. Previous studies of FUV H 2 emission have argued both for and against this assumption (Ardila et al. 2002;Herczeg et al. 2006). LTE conditions keep the assumed parameters straightforward and allow us to model the H 2 ground-state populations as a "snapshot" of the disk atmosphere as it was observed.
The warm H 2 atmosphere is described by the surface density and temperature distribution of gas, which characterizes how much of the warm H 2 is populating excited ground-states [v,J]. We reference these physical quantities in cylindrical coordinate positions in the disk (r,φ,z ). If we consider that a parcel of warm H 2 gas on the disk surface is characterized by its radial position, vertical height from the disk midplane, and velocity distribution (r, z, v φ (r)), the velocity of the gas parcel, v φ (r), is de- where G is the gravitational constant and M ⋆ is the central stellar mass. The mass density at the warm H 2 disk surface is a function of the radial and vertical height in the disk, where Σ(r) is the radial surface density distribution of H 2 , and H p is the pressure scale height as a function of radius, defined as: where c s is the sound speed, Ω is the angular velocity of the gas, k is the Boltzmann constant, T (r) is the radial temperature profile of the warm H 2 disk atmosphere, µ is the "mean molecular weight" of the gas, and m H is the mass of a hydrogen atom. The temperature distribution of the disk atmosphere is approximated as a power-law function: where T 1AU is the temperature of the warm H 2 at r = 1 AU, and q is the temperature gradient. We assume a radial surface density for a static accretion disk, represented by a power-law viscosity profile (see Lynden-Bell & Pringle 1974), where γ is the density gradient, r c is the characteristic radius of the gas in the disk, and Σ c is a normalization factor for the surface density distribution, dependent on the total H 2 mass contributing to the emission lines simulated by these models. The characteristic radius describes the transition from a power-law dominated density distribution to an exponentially-dominated density fall-off in the disk (Lynden-Bell & Pringle 1974;Hartmann et al. 1998). It is important to note that Σ(r) contains a normalization factor (Σ c ), which normalizes to the disk midplace density. Our models only attempt 1442.87 [1,7] (1 − 6)R(6) 1215.726 14 0.9 34.8 1467.08 (1 − 6)P (8) 1.3 1500. 45 (1 − 7)R(6) 1.7 1524.65 (1 − 7)P (8) 1.9 1556.87 (1 − 8)R(6) 1.3 1580.67 (1 − 8)P (8) 1.1 1431.01 [1,4] (1 − 6)R(3) 1216.070 99 1.0 28.9 1446.12 (1 − 6)P (5) 1.4 1489.57 ( a Transitions are from the Lyman-excited to ground electronic states of the H2 band system, Velocity from line center of the pumping transition of Lyα. c Einstein coefficient, describing the spontaneous decay rate from the electronicallyexcited Lyman band, taken from Abgrall et al. (1993). d Oscillator strengths from Abgrall et al. (1993).
to describe the behavior of the disk atmosphere, where the warm, tenuous H 2 resides. As a consequence, the functionality of Σ(r) serves as a structural layout of the radial H 2 disk atmosphere. Since we normalize Σ(r) with a factor describing the disk midplane density, the solutions of Σ(r) describe the radial distributions of warm H 2 , but the resulting H 2 mass estimates are not meaningful.
The level populations of warm, ground state H 2 contributing to the emission line are assumed to be in LTE and are determined using the Boltzmann equation, where X H2 is the fraction of the total H 2 gas mass contributing to the fluorescence observed in the FUV, g [v,J] is the statistical weight of the level population, Z [v,J] (T) is the partition function describing the likelihood that the warm H 2 is in state [v,J], and E [v,J] is the energy of warm H 2 in ground state [v, J].
The radial distribution of molecular hydrogen has two normalization factors (X H2 and Σ c ) that are not independent of disk conditions and are defined by their product in n [v,J] (r, z). The product of these factors describes the total mass of warm H 2 available for photoexcitation to state [v ′ ,J ′ ] (M H2 ), which is obtained by integrating the distribution over (r,φ,z): M H2 = X H2 Σ c 2πr 2 c / (2 − γ). The radiative transfer calculation required to reproduce the observed fluorescent H 2 emission happens in two steps: 1) the warm H 2 in ground state population . We provide a graphical representation of the H 2 disk atmosphere model. The disk contours represent the warm H 2 optical depth (τ λ (r, z)) to stellar Lyα radiation being pumped to state [v ′ ,J ′ ] = [1,4]. The dashed line marks off the approximate location of τ ′ λ ≈ 1, which is where the H 2 disk atmosphere becomes optically thick to the penetrating Lyα photons. The stellar Lyα radiation (purple arrow) is absorbed by the by the warm H 2 , which is excited to state [v ′ ,J ′ ] and emits a photon (λ H 2 ; red arrow) to decay back to ground state [v ′′ ,J ′′ ].
the excited electronic (Lyman band) state by the absorption of an incident stellar Lyα with wavelength λ Lyα , and 2) the excited H 2 molecule decays back to some ground electronic state [X : v ′′ ,J ′′ ], emitting a FUV photon with wavelength λ H2 . Molecular hydrogen has an absorption cross section (σ H2 ) defined by the area around the molecule that can intersect an incoming photon with  Values were chosen to reproduce the desired H2 features (Herczeg et al. 2004;France et al. 2012b). The only parameters without aforementioned constraints were z/r, γ, q and rc because literature estimates of these values were not known. γ and q were constrained by the power-law functionality role they play in the models, and rc was estimated around <r [1,7] > calculated by France et al. (2012b). the appropriate energy for photo-excitation: where λ Lyα is the rest frame wavelength of the stellar Lyα line profile needed to excite the warm H 2 in ground state [X : v, J] up to energy level [B : v ′ , J ′ ], and A lu is the probability that H 2 in population [X : v, J] will be "pumped" to electronic state [B : v ′ , J ′ ]. Note that, for the remainder of this paper, we will omit the ground state H 2 (X) and excited state (B) level branch denominations from the vibration and rotation state discussion.
Assuming an absorption coefficient κ λ (r, z) = σ H2 n [v,J] (r, z), the optical depth of H 2 in ground state [v,J] is described as: For every vertical and radial position in the disk atmosphere that we sample τ λ (r, z), we calculate the amount of the Lyα radiation that will be available for absorption by the warm H 2 . To correct for line absorption overlap of shared Lyα photons, we adopt an effective optical depth τ ′ λ (r, z) (Liu & Dalgarno 1996;Wolven et al. 1997 which corrects for the absorption, scattering, and shielding of Lyα photons. Figure 2 shows a schematic of 1,4] and outlines the radiative transfer process in the disk. We model the emission line flux of each λ H2 produced from the cascade of transitions from energy level [v ′ ,J ′ ] as: All model-to-data reduced-χ 2 statistics for simultaneous emission line fitting, transitioning from excited state [v ′ ,J ′ ]. All χ 2 statistics are calculated between v obs = [-250, 250] km/s. The largest source of errors in the χ 2 statistics come from the linear estimation of the FUV background continuum beneath the emission line. Because the models do not attempt to find the background continuum levels beneath each emission line, extraction of the FUV continuum had to be done manually. Targets with lower signal-to-noise have more uncertainty in the continuum flux, so the χ 2 statistics become large as the errors in the continuum dominate the fitting. Only the [1,4] progression show decent fits for all targets (with χ 2 < 15), so we focus on the relative results of the [1,4] progression emission lines for the remainder of the Discussion section.   [1,4]. The black line represents the observed H 2 fluorescent emission feature, which includes error bars every 5 bins. Each emission line observed has an intrinsic background continuum from the stellar source (see France et al. 2014b), so this continuum was subtracted from the line before model comparisons were made to the observations. The red line is the modeled emission of 1489.57Å from the DM Tau disk model. The blue line is the convolution of the modeled emission line with the COS LSF. This procedure was applied to all modeled emission lines for all targets when comparing the modeled data with FUV observations. The reduced-χ 2 was calculated after the model emission lines were convolved with the COS LSF.
where η represents the coverage fraction of H 2 in the Lyα radiation field (Herczeg et al. 2004 . All reduced-χ 2 values for each progression were calculated by simultaneously fitting each observed emission line profile to those estimated in a given model parameter set. The minimum reduced-χ 2 for each progression is assumed to best represent the H 2 fluorescence distribution. The χ 2 shown in the top left of each emission line box represents the fitting of all emission lines from a given progression with one set of model parameters. Most of the observed emission lines for all targets have single-peaked line profiles (see Brown et al. 2013), but all the best-fit modeled emission lines show a "double horned" profile signature to Keplerian motions. Pontoppidan et al. (2011) points out that the single-peaked profile behavior is suggestive of a low velocity (< 3 km/s) molecular wind located within a few AU of the central star and are typically modeled with an azimuthal velocity vector that is slow relative to Keplerian motion (also see Bast et al. 2011). Since our models leave out the azimuthal velocity component of the H 2 disk gas motions, it is expected that our resulting emission line profiles do not reproduce the line cores of the fluorescent features.
whole suite of transitions available from the progression, and the source function (S λ (r, z)) is defined as the Lyα emission line flux with wavelength λ Lyα , F Lyα (r, z). We calculate how F Lyα (r, z) changes as a function of radial position in the disk. Assuming that the accretiongenerated Lyα flux originates at the stellar surface, we express the ratio of the original F Lyα,⋆ to the flux the warm H 2 disk atmosphere receives at r, To correctly incorporate the Lyα radiation field, we use reconstructued stellar Lyα profiles created by Schindhelm et al. (2012b) and France et al. (2014b), which describe the stellar-Lyα flux seen by the disk surface of each target. After calculating the FUV H 2 fluorescence flux at each disk grid point in our model, we radiate the H 2 emission isotropically, some fraction of which is intercepted by the observer. We calculate the distance of each gas parcel radiating in the disk from the observer s(r, z), based on radial and angular positions of the disk gas parcel, distance to the target, and disk inclination angle. The final modeled emission line flux produced for a fluorescence transition of H 2 is expressed as: Using a total of 6 parameters to represent the physical conditions of the warm, ground-state H 2 populations in the disk atmosphere (z/r, γ, q, T 1AU , r char , M H2 ), Equation 12 characterizes the resulting emission line profiles from H 2 radiating from the disk. All free parameters were allowed to vary over a rough grid of controlled values to create a data cube representing the density distributions, temperature profiles, and radial radiation fields of inner disk H 2 around a given stellar target; see Table 1 for the full list of parameters explored in this study. The resulting models simulate the emission profiles produced for a given fluorescence transition λ H2 , with emission flux as a function of orbital velocity. The radial velocity component of the emission line is determined by v φ (r) of the emitting gas at a given radius in the disk, projected into the sight line of the observer. This model framework was used to describe the observed velocity field of single and binary systems, both close-in and extended. We caution the reader regarding the results of the close-in binary systems (e.g. V4046 Sgr), as the binary potential affects the inner disk velocity-radial relationship differently than a point mass. Therefore, the innermost H 2 modeled for these close-in binary systems may not be accurate, but the outer disk emission distributions will remain unaffected.
Synthesized spectra of each H 2 emission line are compared to HST observations. Each model is convolved with either the HST -COS line spread function (LSF) (Kriss 2011) or a normalized Gaussian distribution with FWHM characterized by the STIS E140M mode spectral resolving power (R∼25,000 for TW Hya; see Herczeg et al. 2006 for more information) prior to comparison with the observed emission line profiles. The FUV continuum level is estimated around each emission feature with a linear fit to the HST -COS data, which is subtracted from the observations before model-to-data comparisons are made. An example of an H 2 emission line, with native and convolved models laid over the HST -COS observed emission line, is shown in Figure 3.

ANALYSIS
The goal of the model-to-data comparison is to find the combination of model parameters that best reproduce the observed fluorescent emission line profiles that cascade from the same excited state [v ′ ,J ′ ]. A reduced-χ 2 statistic is computed when comparing the observed FUV H 2 emission features to the entire data cube of models created for a target. We analyze the reduced-χ 2 statistic data cube for three cases when comparing the modeled emission lines to the observations: (1) fitting individual emission lines; (2) simultaneously fitting all H 2 emission lines fluorescing from the same excited energy level [v ′ ,J ′ ]; (3) fitting only the red wings of the emission lines. The first point was used to set the initial range of temperature and density model parameters of warm H 2 in each disk surface. The third was explored to mitigate the potential influence of a warm molecular wind component that was unresolved at the spectral resolving power of HST -COS. The results of (3) proved inconclusive, which found no significant differences between the red and blue wing line shapes, suggesting that the models are not sensitive to an unresolved warm H 2 disk wind. We focus on the results of (2), which best describe the generalized behavior of the warm H 2 disk atmosphere populations. We simultaneously fit 3 observed fluorescent H 2 transitions for each progression as the most representative of the H 2 radiation distributions in each PPD. Table 1 shows the minimum reduced-χ 2 statistics for  Figure 5. Using the best-fit progression model for GM Aur, we use Equation 12 (integrated over Bmn, which represents the total H 2 flux produced from each progression) to reproduce the observed spectrum. Each progression peaks at different radii, but the overall shape and radial extent of the distributions indicate that the bulk of the radiation for all progressions originates within the same disk annuli.
all targets when simultaneously fitting the 3 progression emission lines from excited state [v ′ ,J ′ ]. Not all minimum reduced-χ 2 simultaneous progression fits for [0,1] and [0,2] were "good", however (i.e., some sources displayed reduced-χ 2 > 25). Many of the strongest lines from [0,1] and [0,2] share similar λ H2 , which makes complex line profiles that depend on the shape of the stellar-Lyα profile illuminating the warm H 2 disk populations to these excited states. The [1,7] and [1,4] progressions are more reliable tracers of the warm H 2 disk atmosphere, and the brightest emission lines in our survey cascade from the [1,4] progression. For the vast majority of the targets, the largest stellar Lyα fluxes pump the warm H 2 disk populations to the [1,4] energy level. This makes the line profile flux fitting more accurate for the [1,4] progression, providing the overall best model fits to the observe FUV emission.
We will focus our discussion around the inner disk diagnostics of the best-fit [1,4] progression for all targets. This progression has good reduced-χ 2 fits (≤ 15) and by-eye model-to-data comparisons for every target in our survey. Figure 4 shows an example of minimum reduced-χ 2 modeled progression lines to those observed with HST -COS for GM Aur. Figure 5 presents the resulting radial radiation distribution for each best-fit progression for GM Aur. While each progression peaks at somewhat different radii, the majority of the radiation distributions originate within similar annuli of the disk. This behavior is typical for all PPD targets that have good minimum reduced-χ 2 fits for all or most progressions.

Uncertainty Estimation and Parameter
Degeneracies Errors in each best-fit parameter per progression are determined after marginalizing the minimum reducedχ 2 parameter fits over all free parameters. Uncertainties are measured as the range of values that encompass 68% of the distribution area, representing the 1-σ uncertainties for a Gaussian distribution. The modeled parameter space was crudely varied over a large range of values for each free variable, so a Gaussian distribution was fit over each marginalized best-fit parameter uncertainty space, and the FWHM of each Gaussian fit was calculated as the uncertainty in each model parameter. Figure 6 displays the reduced-χ 2 marginalized parametric space for each variable in our modeling framework, with filled contours representing the 2-σ uncertainty in the parameter space. Since each parametrized uncertainty is taken within the 1-σ error contours of each marginalized distribution, the uncertainties outlined in red represent the 1-σ errors in the model parameters.
There are noticeable degeneracies amongst several of the parameters; for example, the total mass of emitting H 2 and vertical position of the disk atmosphere (M H2 , z/r) show a trend that requires more mass contributing to the emission lines as the disk height above the disk mid plane increases. This trend makes sense -to produce the same amount of flux in the modeled emission lines, the total mass of H 2 contributing to the emission must increase as the density of H 2 decreases with vertical disk height above the mid plane. The optical depth of the disk atmosphere must remain the same to output the same observed emission line flux, and this relationship between the free parameters maintains the required optical depth. What is important to note is that the models produced are used as a means to describing the H 2 emission flux arising from the inner disk atmosphere. Despite the degeneracies in several parameter pairings relating to the total flux, the radiation distribution of H 2 emission is unaffected by these degeneracies.
We note that our choice in using the reconstructed stellar LyA flux incident on the disk from Schindhelm et al. (2012b) may exacerbate degeneracies in the disk parameters. The Schindhelm et al. (2012b) reconstructed Lyα profiles rely on the same H 2 emission features explored in this study, but we remind the reader that the stellar Lyα flux incident on the H 2 disk scales with the re-emitted H 2 flux (see Eqn 12) and has no effect on the modeled distribution of H 2 flux in each disk. The disk parameters may respond to an inaccurate Lyα flux, but the degeneracies in the disk parameters (for example, the response of M H2 and z/r to the total H 2 flux) would scale to best describe the H 2 radiation that recreate the observed emission profiles. Therefore, the reconstructed Lyα profiles will not change the radial behavior of the best-fit H 2 flux models. Figure 5 presents an example of the radiation distributions of H 2 fluorescence flux produced in the disk for each  Figure 7. The inner and outer radial boundaries which define where 90% of the total radiation is arising from the disk. The black line represents the normalized radial distribution of the [1,4] progression emission to the total amount of flux produced by the [1,4] progression for SU Aur. The green vertical lines show the radial boundaries that encapsulate 90% of the total emission. For r in , we start at the outermost radius and integrate inward to smaller radii the disk until 95% of the total [1,4] progression flux is accounted for. Likewise, rout is defined by starting at the innermost radius defined in our models and integrating the progression flux out until 95% of the total emission flux is accounted for. The resulting annulus between r in and rout represents the ring of disk the majority of the observed FUV H 2 fluorescent emission originates.

The Radial Extent of H 2 Emission
progression explored in this study. We focus our analysis on the [1,4] radiation distributions for all targets in our survey to define inner and outer radial disk boundaries, which describe where the bulk (90%) of the emitting H 2 atmosphere resides. We define the 90% emitting region as follows: ≤ 0.95 for r < r out We use r in and r out to evaluate the evolutionary behavior of the H 2 radiation. Figure 7 presents a schematic of how the inner and outer radial boundaries encapsulate 90% of the total H 2 flux produced in the disk atmosphere. We analyze potential evolutionary characteristics of the molecular disk atmosphere by comparing the FUV H 2 radiation distributions to other dust and molecular disk observables.

Case Study: Model Robustness using AA Tau
We explore how robust our modeling framework is at identifying where the fluorescing H 2 resides in PPDs. We compare two epochs of HST -COS data on AA Tau (2011 and, where the 2013 observations occur during a "dimming" event from X-ray to near-IR wavelengths. Based on the duration of the dimming, Bouvier et al. (2013) suggest an obscuration at r > 8 AU; this hypothesis is strengthened by the gas-to-dust ratio (N H /A v ) of the absorber and the evolution of the FUV H 2 emission (Schneider et al. 2015). We utilize the line profile changes between AA Tau HST -COS FUV observing epochs to determine how those changes relate to radial H 2 radiation distributions in the disk.
There are noticeable differences between the observed FUV H 2 emission line profiles of the 2011 and 2013 AA  (Schneider et al. 2015). On the right: the comparison of the total [1,4] progression flux radiating from the disk of AA Tau in the 2011 and 2013 observations. The 2013 models predict that the observed H 2 fluorescence emission originates from further out in the disk (the peak of the radiation located at r peak = 2.50 AU) than the 2011 radiation distribution (r peak = 0.75 AU), a consequence of the inner disk "shadowing" produced by the extra absorber on the AA Tau sightline (Bouvier et al. 2013).
Tau epochs. The 2013 emission lines are narrower with slightly larger peak fluxes than the same H 2 emission lines observed in 2011 (Schneider et al. 2015). This suggests that less flux is contributed from the innermost disk. The modeling results for the [1,4] progression are shown in Figure 8. Each AA Tau epoch was modeled independently, and the models reproduce the same rest wavelength emission lines. Figure 8 also shows the radiation distributions of [1,4] fluorescence for each epoch in the AA Tau disk. The 2011 emission includes a large contribution from material inside 1 AU (r in,2011 = 0.08±0.01 AU; r peak,2011 = 0.75 AU; r out,2011 = 4.17±2.04 AU), while the 2013 [1,4] emission "appears" to have shifted outward in the disk (r in,2013 = 0.15±0.02 AU; r peak,2013 = 2.50 AU; r out,2013 = 7.59±2.75 AU). Our models indicate that the inner radius of detectable H 2 fluorescence from the [1,4] progression has moved outward radially in the disk as the "extra absorber" moved into our field of view in the AA Tau disk. Schneider et al. (2015), using an independent modeling technique to estimate the radial origins of H 2 fluorescence in the AA Tau disk, come to a similar conclusion: the observed 2013 H 2 emission within ∼ 1 AU is reduced compared to 2011. Additionally, Schneider et al. (2015) find that the outer radial extent of the H 2 fluorescence luminosity doesn't change significantly between epochs, which is a result consistent within the errors on our modeled r out estimates of the AA Tau epochs.
The "extra absorber" obscures the inner disk H 2 fluorescence in the 2013 HST -COS FUV spectrum, making AA Tau appear as a disk with a deficit of inner disk emission -effectively, a pseudo-transition disk. Our modeling framework was capable of identifying the change in emis- Primordial Transitional Figure 9. The normalized modeled radiation field distribution of H 2 fluorescence cascading from the [v ′ ,J ′ ] = [1,4] energy level for all targets. Each radiation distribution was calculated from the minimum reduced-χ 2 model parameters that best reproduce the observed H 2 emission lines. The top plot represents radiation distributions for all primordial disk targets, and the bottom plot shows the distributions for transition disk objects. The two disk evolution types appear to show an evolving H 2 FUV radiation field; primordial disks generally radiating more inward in the disk, with the bulk of the radiation occurring within r 1 AU, and transition disk H 2 radiation starting at larger radii (r ∼ 0.1 AU) extending to larger radii (r ∼ 10 AU). sion line profiles between the 2011 and 2013 AA Tau observations and found that the bulk of the 2013 AA Tau [1,4] radiation in the disk originated at larger radii than the 2011 H 2 fluorescence. We expect our models are therefore capable of distinguishing between H 2 fluorescence evolution in differing disk types.

DISCUSSION
We created 2D radiative transfer models to simulate observed HST -COS and -STIS FUV H 2 emission lines to understand where the majority of the radiation arises in PPDs. We use the best-fit model results to define the inner and outer radii of warm H 2 (r in , r out ) and examine if and how the molecular distributions change as PPDs evolve. We compare r in and r out to other dust and molecular tracers that help describe the evolutionary state of the PPDs. Table 1 provides a detailed list of inner disk observables for each target, including dust cavity radius (r cavity ) and inner disk CO radius (r in,CO ). We also look at where the theoretical snow lines in the disks exist and how these radii relate to the H 2 disk emission.

Radiation Distribution of Modeled H 2 Fluorescent
Emission Figure 9 presents the normalized radial distributions of warm H 2 transitioning from excited state [1,4] for all targets. We modeled 6 primordial disks (AA Tau, BP Tau, DF Tau A, HN Tau A, RECX-11, and RECX-15) and 8 transition disks (CS Cha, DM Tau, GM Aur, LkCa 15, SU Aur, TW Hya, UX Tau A, and V4046 Sgr) to compare the radial distribution of warm H 2 in the disk atmospheres as the dust disk evolves. The H 2 radial distributions of the different dust disk stages appear qualitatively different. The primordial disk population (top plot in Figure 9) generally starts radiating significantly in the very inner disk (r 0.05 AU), and the radiation only extends out to a few AU, consistent with the simple estimates of the average H 2 emitting radius presented by France et al. (2012b). The generalized transition disk radiation behavior (bottom plot) starts further out in the disk (r ∼ 0.1 AU) and extend significantly further out into the disk (r ∼ 10 AU). These different behaviors suggest structural changes in any of all of the following: the spatial distributions of warm H 2 in populations [v,J]; the degree of Lyα penetration into the disk by clearing H 2 from the inner disk atmosphere; or the evolution of the disk surface temperature distribution. This evolving radiation structure is also observable in the line profiles of the [1,4] progression, as seen in Figure 10. As the PPDs in our survey evolve from primordial to transition disks, the majority of the observed H 2 emission migrates to larger radii.
We compare estimates of r in and r out to investigate the idea that the radial distributions of fluorescing H 2 migrate outward in the disks as PPDs evolve. Figure 11 presents a comparison of r in and r out , which shows the annulus of H 2 emission extending further out into the disk as the inner disk radius moves outward. A line can be fit to represent the relationship between the inner and outer radiating disk radii for our survey targets: log 10 (r out (H 2 )) = 0.79 log 10 (r in (H 2 )) + 1.39, where both log 10 (r in (H 2 )) and log 10 (r out (H 2 )) are in units of AU, and the coefficients [1.39 ± 0.22, 0.79 ± 0.21] are computed from a χ 2 minimization (χ 2 min =0.896) of a linear function between log 10 (r in (H 2 )) and log 10 (r out (H 2 )). The Spearman rank correlation coefficient between r in and r out indicates a statistically significant correlation between the variables (ρ = 0.70) with a small probability that the sample is randomized (n = 5.5 × 10 −3 ), providing additional evidence that support the migration of the radial H 2 emission as PPD warm dust dissipates from the inner disk.

Comparison to Dust Evolution
We compare results from our modeled H 2 [1,4] progression radial distributions with dust disk evolution diagnostics to gain insight into how the molecular inner disk environment of PPDs changes as dust grains clear. We identify PPD evolution using observed color-color changes in the near-to mid-IR SED slope of the disk, which provides an estimate of the degree of warm dust clearing (see Espaillat et al. 2014). We interpret the slope of each tar-  Figure 11. We present the relation between the estimated r in and rout quantities, determined from the best-fit modeled radiation distributions for all targets. The inner radial boundary (r in ) is defined as the inner radius of H 2 radiation in the disk that marks where at least 90% of the total radiation is accounted for in the outer disk. Likewise, the outer radial boundary (rout) is defined as the outermost radius of H 2 radiation that encompasses 90% of the total amount of radiation accounted for in the inner disk. The blue diamonds with error bars represent each modeled r in and rout, and the black dashed line represents a linear fit to the data. The Spearman rank correlation coefficient ([ρ, n] = [0.70, 5.5 × 10 −3 ]) between the two radial quantities suggest a strong increasing trend between them, indicating that the whole emitting region is moving outward.  (Furlan et al. 2009). In the top plot: Each blue triangle with error bars represents each target point in our survey. The black dashed line represents the best-fit linear correlation between r in and n 13−31 . In the bottom plot: Each green triangle with error bars represents each target point in our survey. The black dashed line represents the best-fit linear correlation between rout and n 13−31 . In both plots, a clear increasing trend is seen in the radial H 2 emission boundaries as the warm dust disk content evolves.
get SED with the observable n 13−31 (Furlan et al. 2009): which is dominated by longer wavelength continuum emission from the optically-thick dust in the disk and is sensitive to the degree of dust settling towards the disk a AllṀ values taken from Ingleby et al. (2013).

b
(1) Akeson et al. (2005); (2) Andrews & Williams (2007) midplane (D'Alessio et al. 2006). For many targets in this work, n 13−31 were available in Furlan et al. (2009), but for targets not included in the Furlan et al. (2009) survey, we calculate n 13−31 with known or modeled disk SEDs (for example, an intricate model of V4046 Sgr SED was found by Rosenfeld et al. 2013). We interpret the results of n 13−31 as follows: if n 13−31 < 0, the inner dust disk is optically thick, essentially a primordial disk; if n 13−31 ≥ 0, the disk dust is optically thin, indicative of dust clearing or settling and evidence for PPD evolution into the transition state (Lada 1987;Strom et al. 1989;Andre & Montmerle 1994). Table 1 provides a list of n 13−31 values for all targets in this survey.
A comparison of the [1,4] emission boundaries (r in , r out ) to n 13−31 is made in Figure 12. The top figure shows the relationship between r in and n 13−31 , and the bottom figure shows r out versus n 13−31 . The triangles in both plots represent each target in our survey, and the black dashed line in each plot shows the linear correlation between r in versus n 13−31 and r out versus n 13−31 . It is apparent that the molecular inner and outer disk emission radii show a positive correlation with the dust disk evolution: the Spearman rank correlation coefficient for r in versus n 13−31 is 0.72 (n = 4.0 x 10 −3 ), and ρ = 0.69 for r out versus n 13−31 (n = 6.9 x 10 −3 ). Both correlation coefficients suggest a strong increasing trend in the radial outward migration of the FUV H 2 radiation as the warm dust disk evolves in the disk samples. The linear correlation between r in and n 13−31 is expressed as: log 10 (r in (H 2 )) =(0.19 ± 0.07) × n 13−31 − (1.05 ± 0.08), and the linear correlation between r out and n 13−31 is expressed as: log 10 (r out (H 2 )) =(0.25 ± 0.06) × n 13−31 + (0.52 ± 0.07).
We note that, for all transition disks in this study, r out The purple x-points represent all targets with mass accretion rates > 10 −9 M ⊙ yr −1 , while the red diamonds represent the RECX targets (which have low mass accretion rates for primordial PPD targets). The black line is a negative correlation fit through all the purple points, suggesting that the mass accretion rate decreases as r in increases. Since the accretion luminosity, more specifically the stellar Lyα flux produced by the accretion, is directly related to the H 2 emission observed, it is important to note that r in is not necessarily correlated to the flux produced by the mass accretion rate. Instead, r in is sensitive to the observed emission line width, which is independent of the stellar incident flux.
is found to be within the dust gap radius. One interpretation of this result, paired with the correlation between r in and n 13−31 , is that the H 2 FUV radiation observed from the inner PPD atmosphere lags behind the dust disk evolution.
This does not automatically mean that the molecular content of the disk is clearing, and we need further evidence of evolution with other inner disk molecular trac-ers before we can make this distinction. France et al. (2012b) outlined the conditions needed in the H 2 disk atmosphere to produce Lyα-pumped H 2 fluorescence. The opacity of absorbing H 2 in ground-state [v,J] must be large, with excitation temperatures T exc > 1500K, and the mass accretion rate (Ṁ ) onto the proto-star must be large enough to produce enough Lyα photons to stimulate the molecules. The mass accretion rate implies there is a reservoir of material in the inner regions of PPDs that feeds onto the proto-star, and a decrease iṅ M over time (e.g., Muzerolle et al. 2000) strongly suggests that the inner disk material is being depleted. Figure 13 shows the relationship betweenṀ and r in (H 2 ), with purple points representing r in (H 2 ) andṀ for all targets except the RECX targets, which are represented at red diamonds. All mass accretion rates are taken from Ingleby et al. (2013). Figure 13 shows a negative correlation betweenṀ and r in (H 2 ), with Spearman rank correlation [ρ,n] = [-0.80, 1.9 x 10 −3 ] (not including the RECX targets), suggesting that the H 2 atmosphere may be physically thinned or in different ground-state populations not suitable for Lyα-pumping in the very inner disk regions of evolved PPDs. The outlier points in Figure 13, RECX-11 and RECX-15, appear to have abnormally low mass accretion rates given the evolutionary stage of the disks (Ingleby et al. 2011), and more targets of varying evolution may be needed to understand if this result is universal among a large sampling of PPDs. It is important to note that r in is primarily dervied from the observed line widths of H 2 emission profiles, so determination of r in is largely independent of the incident FUV flux.
The link betweenṀ and r in (H 2 ) suggests that the inner disk is clearing of material as the mass accretion rate declines. One explanation for this correlation is that the warm H 2 atmosphere dissipates with the small dust grains. Dust grains present in the disk atmospheres of primordial disks may give warm H 2 a formation site to replenish molecules lost to photo-dissociation and stellar accretion (see Augason 1970;Habart et al. 2004;Fleming et al. 2010). As the dust grains clear out and settle towards the disk midplane or evaporate from the inner disks of evolving PPDs, the warm H 2 atmosphere no longer has a formation site to maintain the molecular reservoir. Via accretion and photo-dissociative processes with FUV continuum photons between 912 -1120Å, the leftover warm H 2 will continue to disperse, even as the accretion flux decreases. This leaves an optically thin (N (H 2 ) 10 18 cm −2 ) path for stellar Lyα to reach the warm H 2 material at larger disk radii (r > 3 AU).
The migration of r out (H 2 ) with increasing n 13−31 also suggests that neutral hydrogen (HI) is being cleared from the inner disks of transitional PPDs. Photo-excitation via stellar Lyα drives the H 2 fluorescence observed in the disk atmospheres, and as the emitting H 2 is observed further out in the disk, there must be new paths open for stellar UV radiation to reach the outer disk material. In primordial disks, HI re-processes and scatters incident stellar Lyα down into the inner disk (Fogel et al. 2011) while H 2 self-shields the radiation from penetrating to the outer disk, preventing the stellar Lyα from reaching the outer disk effectively. If H 2 and HI column densities in the inner disk become optically thin in transitional disks, more stellar Lyα can irradiate molecular material in the outer disk and may explain the observed correlation between r out (H 2 ) and n 13−31 . This suggests that HI clearing from the inner disk may happen over a similar timescale as the characteristic dust dissipation (Wyatt 2008;Ribas et al. 2014) and mass accretion quenching (Fedele et al. 2010). This inner-to-outer disk dissipation is in agreement with the UV switch model, which describes the dispersal of inner disk gas cut off from the gas reservoir of the outer disk, due to selective photoevaporation of material out to r ∼ 5 -10 AU (Clarke et al. 2001;Alexander et al. 2006). Observations of other outer-disk molecules photo-excited by Lyα radiation provide additional evidence for the loss of HI in the inner disks of transitional objects. For example, Schindhelm et al. (2012a) observe FUV-CO fluorescence, also powered by stellar Lyα-pumping, at T exc ∼ 500K, in transitional phase objects with an average emission radius R CO ∼ 1 -5 AU. This indicates that less HI and H 2 column is present in the inner disk to shield the stellar Lyα flux from reaching the cooler CO material at intermediate radii in transition systems. Figure 14 shows a 1D radial comparison of dust and molecular tracers determined in our targets. We present the locations of the outer radiation boundary for H 2 FUV emission, as determined from our models (r out, [1,4] ; green triangles), and the observed dust cavity walls of the transitional disk populations (r cavity ; blue squares). For all transitional disks, we find r out, [1,4] inward of r cavity , meaning that the H 2 population observed in all transition PPDs radiates where the dust is optically thin, suggesting that the H 2 populations remain optically thick even after the dust grains have dissipated. Studies like van der Marel et al. (2015) also find a substantial depletion of the dust-to-gas ratio inside the dust cavities of well-studied transition disks, which is consistent with our findings.

Near-IR CO Emission and Comparison to Snow
Line Radii Figure 14 includes radial estimates of the inner radiation boundary for H 2 FUV emission (r in, [1,4] ; blue x's) and the inner radius of near-IR CO emission, determined from LTE models presented by Salyk et al. (2011a) (r in,CO ; black diamonds). The inner disk emission radii of FUV H 2 and near-IR CO appear to be roughly co-spatial, which is a result also found by France et al. (2012b) when comparing the observed FWHMs of FUV H 2 fluorescence emission and near-IR CO emission. An extensive study by Brown et al. (2013) concluded there is a correlation between the near-IR CO P(8) equivalent width and dust disk dispersal in transitional disks, suggestive of outer radial origins of the CO emission as PPD dust evolves. We have shown that r in, [1,4] increases with n 13−31 and decreases witḣ M , providing further evidence that the inner gas disk environment becomes optically thin as disks evolve towards the transition stage.
We note the disk locations of possible theoretical snow lines in PPDs and these radii coincide with the H 2 fluorescence in Figure 14 Figure 14. A radial comparison of the inner and outer extent of FUV H 2 emission (this work), the innermost radius of near-IR CO emission (Salyk et al. 2011a), and dust cavity locations in transition disk targets (see Table 1 for references). The light blue shaded area from 1 AU ≤ r DISK ≤ 3 AU represents the theoretical water-ice snow line for the presence of water-ice at the midplane of primordial and transitional PPDs (Baillié et al. 2015).
inward in the disk as the protostellar system ages (?). Several independent studies (e.g. Meijerink et al. 2009;Mandell et al. 2012) conclude that the location of the water-ice snow line in PPDs are expected to be found within r ∼ 1 -3 AU for all PPD states. Baillié et al. (2015) shows that the evolution of the water-ice snow line at all stages of PPD evolution (from ages 10 6 -10 7 yr) only varies by ∼ 0.5 AU. Observations of H 2 O and OH (which is thought to be a bi-product of H 2 O photo-dissociation) in the near-and mid-IR are also consistent with these condensation radii (Malfait et al. 1998;Carr et al. 2004;Mandell et al. 2008;Salyk et al. 2008). Figure 14 includes a shaded blue region that represents the assumed generalized H 2 O snow line radii in PPDs, located between r DISK = 1 -3 AU. With the exception of BP Tau and DF Tau A, all targets have outer H 2 emission radii that extend to within or outward of the water-ice snow line.

CONCLUSIONS AND FUTURE WORK
We have created 2D radiative transfer models of FUV H 2 fluorescence emission in PPDs and compared them with observations made with HST -COS and STIS. We analyze the radial distribution of H 2 emission produced by parametrized models, which are determined using a reduced-χ 2 statistic, to understand how the emitting H 2 regions changes as PPD dust disks evolve. We summarize our findings and interpret the evolutionary behavior of the molecular disk atmosphere as the inner dust disks of PPDs disperse: 1. The modeled H 2 radial distributions differ between primordial and transitional disks. Primordial disks have the majority of the total H 2 flux arising from the innermost disk radii and less produced outside ∼ 1 AU. For transitional disks, the total H 2 flux migrates to larger disk radii, producing less flux in the innermost disk and more out to r ∼ 10 AU.
2. We see a positive correlation between the resulting inner and outer emission radii of FUV H 2 (r in and r out ), which supports the result described in conclusion 1. This can mean: a) that the physical structure (i.e., temperature) of the warm molecular disk atmosphere changes as PPDs evolve, b) the warm, ground-state H 2 populations [v,J] change, resulting in evolving regions of the disks where the warm H 2 atmosphere will reprocess the stellar Lyα radiation field, or c) H 2 is being destroyed in the inner disk and not re-formed, owing to the lack of dust grains; the latter point allows stellar Lyα to penetrate to larger r out .
3. We observe positive correlations between r in , r out , and n 13−31 , suggesting that r in corresponds with the loss of warm, small dust grains in the innermost disk. We find a negative correlation between r in andṀ , providing evidence that the warm H 2 inner disk atmosphere may be physically thinned or cleared as the PPDs evolve, possibly by the loss of a molecular formation site as the dust grains dissipate from the atmosphere. Using the observed dust cavity radii of the transitional disk targets, we compare r out to r cavity and find that, for all transition disk targets, r out is found inward of r cavity . This indicates that the warm H 2 disk (for r > r in ) remains optically-thick where the warm dust grains are optically-thin in the disks. This suggests that the physical mechanism that clears or settles the inner disk dust either does not have the same effect on the molecular disk atmosphere, or there is a time lag for the gas disk to respond to the changes observed in the dust distribution.
4. We examine where the emitting H 2 originates in the disk relative to warm CO and the theoretical location of water-ice snow lines. Inner disk CO is roughly co-spatial with r in for all targets, which could point to the dispersal of the warm molecular disk atmospheres of evolving disk systems. With the exception of a few primordial disk targets, all targets have emitting H 2 regions that encapsulate the theoretical water-ice snow line. If disk clearing mechanisms, such as disk photoevaporation via EUV/X-ray photons, are primarily responsible for the final dispersal of the gas disk at the end of the PPD lifetime, it is important to examine late-type PPDs to monitor molecular disk clearing as transitional disks evolve to debris disks.