Rapid destruction of protoplanetary discs due to   externalphotoevaporation in star-forming regions

Rhana B. Nicholson; Richard J. Parker; Ross P. Church; Melvyn B.; Davies; Niamh M. Fearon; Sam R. J. Walton

arXiv:1902.11094·astro-ph.GA·April 10, 2019

Rapid destruction of protoplanetary discs due to externalphotoevaporation in star-forming regions

Rhana B. Nicholson, Richard J. Parker, Ross P. Church, Melvyn B., Davies, Niamh M. Fearon, Sam R. J. Walton

PDF

TL;DR

This study uses N-body simulations to show that external photoevaporation by massive stars rapidly disperses protoplanetary discs, especially in dense, substructured, and collapsing star-forming regions, affecting planet formation timelines.

Contribution

It provides new insights into how initial conditions of star-forming regions influence disc dispersal rates due to external photoevaporation, highlighting faster dispersal than previously observed.

Findings

01

Discs in dense, substructured regions are dispersed within 1 Myr.

02

Half of large discs are photoevaporated within 1 Myr in certain regions.

03

Lower-density regions allow longer disc survival, but still experience rapid dispersal.

Abstract

We analyse N-body simulations of star-forming regions to investigate the effects of external far and extreme ultra-violet photoevaporation from massive stars on protoplanetary discs. By varying the initial conditions of simulated star-forming regions, such as the spatial distribution, net bulk motion (virial ratio), and density, we investigate which parameters most affect the rate at which discs are dispersed due to external photoevaporation. We find that disc dispersal due to external photoevaporation is faster in highly substructured star-forming regions than in smooth and centrally concentrated regions. Sub-virial star-forming regions undergoing collapse also show higher rates of disc dispersal than regions that are in virial equilibrium or are expanding. In moderately dense ( $\sim$ 100 M $_{⊙}$ pc $^{- 3}$ ) regions, half of all protoplanetary discs with radii $\geq$ 100 AU are…

Tables6

Table 1. Table 1: The time taken for half of stars within the star-forming region to lose the gas from their 100 AU protoplanetary discs in a 1000 M ⊙ sub-virial ( α vir = 0.3 subscript 𝛼 vir 0.3 \alpha_{\rm vir}=0.3 ) region with two different initial densities; 10 and 100 M ⊙ pc -3 . Four different spatial distributions are analysed; D = 1.6 𝐷 1.6 D=1.6 (highly sub-structured), D = 2.0 𝐷 2.0 D=2.0 (moderately sub-structured), D = 3.0 𝐷 3.0 D=3.0 (smooth), and a Plummer sphere distribution. The highest and lowest values from the 20 different runs are included.

	Half life of protoplanetary discs
Fractal dimension ( $D$ )	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3
1.6	2.12 $\pm$ ${}_{1.11}^{0.51}$ Myr	0.87 $\pm$ ${}_{0.49}^{0.50}$ Myr
2.0	2.60 $\pm$ ${}_{0.62}^{1.36}$ Myr	0.67 $\pm$ ${}_{0.22}^{0.21}$ Myr
3.0	3.62 $\pm$ ${}_{0.89}^{1.68}$ Myr	0.65 $\pm$ ${}_{0.16}^{0.10}$ Myr
Plummer Sphere	3.85 $\pm$ ${}_{1.34}^{3.70}$ Myr	0.84 $\pm$ ${}_{0.29}^{0.90}$ Myr

Table 2. Table 2: The average percentage of 100 AU discs remaining after 10 Myr within a sub-virial ( α vir = 0.3 subscript 𝛼 vir 0.3 \alpha_{\rm vir}=0.3 ) star forming region from 20 realisations of each simulation. The amount of substructure is varied from highly substructured ( D = 1.6 𝐷 1.6 D=1.6 ) to a smooth and centrally concentrated Plummer sphere. The highest and lowest values from the 20 different runs are included. Two different initial densities (10 and 100 M ⊙ pc -3 ) are considered.

	Percentage of discs remaining after 10 Myr
Fractal dimension ( $D$ )	$ρ =$ 10 M_⊙ pc^-3	$ρ$ = 100 M_⊙ pc^-3
1.6	16.40 $\pm$ ${}_{9.8}^{4.58}$ %	5.99 $\pm$ ${}_{2.67}^{2.88}$ %
2.0	17.75 $\pm$ ${}_{4.04}^{17.03}$ %	2.27 $\pm$ ${}_{0.92}^{4.04}$ %
3.0	21.60 $\pm$ ${}_{10.88}^{16.63}$ %	1.35 $\pm$ ${}_{0.42}^{0.63}$ %
Plummer Sphere	29.70 $\pm$ ${}_{14.98}^{18.95}$ %	2.81 $\pm$ ${}_{1.63}^{2.36}$ %

Table 3. Table 3: The time taken for half of stars within the cluster to lose the gas within their 100 AU protoplanetary discs in a 1000 M ⊙ , moderately sub-structured ( D = 2.0 𝐷 2.0 D=2.0 ) star forming region for two different initial densities; 10 and 100 M ⊙ pc -3 . Three different virial ratios are analysed: α vir = 0.3 subscript 𝛼 vir 0.3 \alpha_{\rm vir}=0.3 (sub-virial, or collapsing), α vir = 0.5 subscript 𝛼 vir 0.5 \alpha_{\rm vir}=0.5 (virial equilibrium), and α vir = 0.7 subscript 𝛼 vir 0.7 \alpha_{\rm vir}=0.7 (super virial, or expanding).

	Half life of protoplanetary discs
Virial Ratio ( $α_{vir}$ )	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3
0.3	2.60 $\pm$ ${}_{0.62}^{1.36}$ Myr	0.67 $\pm$ ${}_{0.22}^{0.21}$ Myr
0.5	3.10 $\pm$ ${}_{1.15}^{2.73}$ Myr	0.68 $\pm$ ${}_{0.35}^{0.66}$ Myr
0.7	3.53 $\pm$ ${}_{1.40}^{1.72}$ Myr	0.63 $\pm$ ${}_{0.23}^{0.57}$ Myr

Table 4. Table 4: The average percentage from 20 runs of simulations of 100 AU discs remaining after 10 Myr within a moderately substrutured ( D = 2.0 𝐷 2.0 D=2.0 ) cluster. The bulk motion (virial ratio) of the star-forming region is varied, from collapsing (sub-virial, α vir = 0.3 subscript 𝛼 vir 0.3 \alpha_{\rm vir}=0.3 ) to expanding (super virial, α vir = 0.7 subscript 𝛼 vir 0.7 \alpha_{\rm vir}=0.7 ). The highest and lowest values from the 20 different runs are included.

	Percentage of discs remaining after 10 Myr
Virial Ratio ( $α_{vir}$ )	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3
0.3	17.75 $\pm$ ${}_{4.04}^{17.03}$ %	2.27 $\pm$ ${}_{0.92}^{4.04}$ %
0.5	29.77 $\pm$ ${}_{17.95}^{15.73}$ %	4.16 $\pm$ ${}_{2.06}^{15.14}$ %
0.7	32.67 $\pm$ ${}_{16.14}^{10.9}$ %	5.00 $\pm$ ${}_{3.23}^{8.04}$ %

Table 5. Table 5: The time taken for half of stars in a star forming region to lose the gas within their 100 AU protoplanetary discs in a 1000 M ⊙ , moderately sub-structured ( D = 2.0 𝐷 2.0 D=2.0 ) region for two different initial densities; 10 and 100 M ⊙ pc -3 and two different masses of disc, 10 per cent and 1 per cent. Three different virial ratios are analysed: α vir = 0.3 subscript 𝛼 vir 0.3 \alpha_{\rm vir}=0.3 (sub-virial, or collapsing), α vir = 0.5 subscript 𝛼 vir 0.5 \alpha_{\rm vir}=0.5 (virial equilibrium), and α vir = 0.7 subscript 𝛼 vir 0.7 \alpha_{\rm vir}=0.7 (super virial, or expanding). The highest and lowest values from the 20 different runs are included.

	Half life of cluster protoplanetary discs
	Disc mass = 0.1 M_⋆		Disc mass = 0.01 M_⋆
Disc Radius (AU)	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3
10	> 50% remaining	3.92 $\pm$ ${}_{1.87}^{5.91}$ Myr	3.31 $\pm$ ${}_{0.9}^{5.16}$ Myr	0.84 $\pm$ ${}_{0.27}^{0.22}$ Myr
50	3.94 $\pm$ ${}_{0.82}^{5.96}$ Myr	1.04 $\pm$ ${}_{0.24}^{0.40}$ Myr	1.22 $\pm$ ${}_{0.35}^{1.18}$ Myr	0.28 $\pm$ ${}_{0.09}^{0.07}$ Myr
100	2.60 $\pm$ ${}_{0.62}^{1.36}$ Myr	0.67 $\pm$ ${}_{0.22}^{0.21}$ Myr	0.71 $\pm$ ${}_{0.21}^{0.96}$ Myr	0.14 $\pm$ ${}_{0.05}^{0.05}$ Myr
200	1.55 $\pm$ ${}_{0.44}^{1.34}$ Myr	0.36 $\pm$ ${}_{0.11}^{0.09}$ Myr	0.39 $\pm$ ${}_{0.12}^{0.58}$ Myr	0.06 $\pm$ ${}_{0.02}^{0.03}$ Myr
1000	0.37 $\pm$ ${}_{0.11}^{0.55}$ Myr	0.06 $\pm$ ${}_{0.02}^{0.02}$ Myr	0.15 $\pm$ ${}_{0.05}^{0.20}$ Myr	0.02 $\pm$ ${}_{0.01}^{0.02}$ Myr

Table 6. Table 6: The average percentage of 100 AU discs remaining after 10 Myr within a moderately substructured ( D = 2.0 𝐷 2.0 D=2.0 ) star forming region for two different initial densities, ( 10 10 10 and 100 100 100 M ⊙ pc -3 ), with two different initial disc masses, 10 per cent and 1 per cent the mass of the host star. The bulk motion (virial ratio) of the star forming region is varied, from collapsing (sub-virial, α vir = 0.3 subscript 𝛼 vir 0.3 \alpha_{\rm vir}=0.3 ) to expanding (super-virial, α vir = 0.7 subscript 𝛼 vir 0.7 \alpha_{\rm vir}=0.7 ). The highest and lowest values from the 20 different runs are included.

	Percentage of discs remaining after 10 Myr
	Disc mass = 0.1 M_⋆		Disc mass = 0.01 M_⋆
Disc Radius (AU)	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3	$ρ =$ 10 M_⊙ pc^-3	$ρ =$ 100 M_⊙ pc^-3
10	77.29 $\pm$ ${}_{4.71}^{12.79}$ %	39.49 $\pm$ ${}_{14.17}^{14.34}$ %	28.81 $\pm$ ${}_{9.09}^{19.13}$ %	3.99 $\pm$ ${}_{2.18}^{6.65}$ %
50	34.44 $\pm$ ${}_{7.57}^{18.08}$ %	6.73 $\pm$ ${}_{3.62}^{5.97}$ %	4.79 $\pm$ ${}_{1.43}^{9.34}$ %	0.40 $\pm$ ${}_{0.23}^{1.53}$ %
100	17.75 $\pm$ ${}_{4.04}^{17.03}$ %	2.27 $\pm$ ${}_{0.92}^{4.04}$ %	2.27 $\pm$ ${}_{1.13}^{3.91}$ %	0.08 $\pm$ ${}_{0.08}^{0.59}$ %
200	7.12 $\pm$ ${}_{1.82}^{10.79}$ %	0.74 $\pm$ ${}_{0.36}^{1.95}$ %	1.09 $\pm$ ${}_{14.17}^{1.81}$ %	0.00 $\pm$ ${}_{0.00}^{0.38}$ %
1000	1.09 $\pm$ ${}_{0.9}^{1.73}$ %	0.0 $\pm$ ${}_{0.00}^{0.38}$ %	0.29 $\pm$ ${}_{0.29}^{0.51}$ %	0.00 $\pm$ ${}_{0.00}^{0.04}$ %

Equations8

p (m) \propto (\frac{m}{μ})^{- α} (1 + (\frac{m}{μ})^{1 - α})^{- β} .

p (m) \propto (\frac{m}{μ})^{- α} (1 + (\frac{m}{μ})^{1 - α})^{- β} .

N \propto M_{cl}^{- 2},

N \propto M_{cl}^{- 2},

\dot{M}_{FUV} \approx 2 \times 1 0^{- 9} r_{d} M_{⊙} y r^{- 1},

\dot{M}_{FUV} \approx 2 \times 1 0^{- 9} r_{d} M_{⊙} y r^{- 1},

\dot{M}_{EUV} \approx 8 \times 1 0^{- 12} r_{d}^{3/2} \frac{Φ _{i}}{d ^{2}} M_{⊙} y r^{- 1},

\dot{M}_{EUV} \approx 8 \times 1 0^{- 12} r_{d}^{3/2} \frac{Φ _{i}}{d ^{2}} M_{⊙} y r^{- 1},

Peer Reviews

No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.

Videos

No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.

Full text

Rapid destruction of protoplanetary discs due to external photoevaporation in star-forming regions

Rhana B. Nicholson,1,2 Richard J. Parker,2 Ross P. Church,3 Melvyn B. Davies,3 Niamh M. Fearon4 and Sam R. J. Walton1

1Astrophysics Research Institute, Liverpool John Moores University, 146 Brownlow Hill, Liverpool, L3 5RF, UK

2Department of Physics and Astronomy, The University of Sheffield, Hicks Building, Hounsfield Road, Sheffield, S3 7RH, UK

3Lund Observatory, Department of Astronomy and Theoretical Physics, Box 43, SE-221 00, Lund, Sweden

4School of Physics and Astronomy, The University of Manchester, Oxford Road, Manchester, M13 9PL, UK E-mail: [email protected] (RBN)Royal Society Dorothy Hodgkin Fellow

(Accepted XXX. Received YYY; in original form ZZZ)

Abstract

We analyse $N$ -body simulations of star-forming regions to investigate the effects of external far and extreme ultra-violet photoevaporation from massive stars on protoplanetary discs. By varying the initial conditions of simulated star-forming regions, such as the spatial distribution, net bulk motion (virial ratio), and density, we investigate which parameters most affect the rate at which discs are dispersed due to external photoevaporation. We find that disc dispersal due to external photoevaporation is faster in highly substructured star-forming regions than in smooth and centrally concentrated regions. Sub-virial star-forming regions undergoing collapse also show higher rates of disc dispersal than regions that are in virial equilibrium or are expanding. In moderately dense ( $\sim$ 100 M*⊙* pc*-3*) regions, half of all protoplanetary discs with radii $\geq$ 100 AU are photoevaporated within 1 Myr, three times faster than is currently suggested by observational studies. Discs in lower-density star-forming regions ( $\sim$ 10 M*⊙* pc*-3*) survive for longer, but half are still dispersed on short timescales ( $\sim$ 2 Myr). This demonstrates that the initial conditions of the star forming regions will greatly impact the evolution and lifetime of protoplanetary discs. These results also imply that either gas giant planet formation is extremely rapid and occurs before the gas component of discs is evaporated, or gas giants only form in low-density star-forming regions where no massive stars are present to photoevaporate gas from protoplanetary discs.

keywords:

protoplanetary discs – ISM: photodissociation region (PDR) – open clusters and associations: general – methods: numerical

††pubyear: 2018††pagerange: Rapid destruction of protoplanetary discs due to external photoevaporation in star-forming regions–Rapid destruction of protoplanetary discs due to external photoevaporation in star-forming regions

1 Introduction

Protoplanetary discs (‘Proplyds’) are thin Keplerian discs around pre-main-sequence stars (Shu et al., 1994) and are the birth places of planets. Proplyds form as a result of angular momentum conservation during the gravitational collapse of clouds when stars are forming. Within the discs, dust can coagulate to form a range of objects from pebbles to planets.

The evolution and dispersal of protoplanetary discs controls the planet formation process, and observations suggest that disc lifetimes are between $\approx 3-5$ Myr (e.g. Zuckerman et al., 1995; Haisch et al., 2001; Pascucci et al., 2006; Richert et al., 2018). Internal processes remove mass from the protoplanetary disc and, after several Myr, disc accretion slows significantly to the point that these processes begin removing more mass than can be replaced, leading to very rapid disc dispersal (Clarke et al., 2001; Owen et al., 2011).

Proplyd host stars do not form in isolation, but rather in clusters and associations with stellar densities that exceed that of the Galactic field by a few orders of magnitude (Lada & Lada, 2003; Bressert et al., 2010). Tens to thousands of stars can form in these regions that are a fraction of a parsec in size (Clarke et al., 2000). Observations of young star-forming regions have revealed that stars form in filamentary structures (André et al., 2010; Wang et al., 2014), resulting in hierarchical spatial distributions. The net motion of stars within these regions indicate that the structures are often collapsing (i.e. sub-virial) (Peretto et al., 2006; Foster et al., 2015; Kuhn et al., 2018).

The initial densities of star-forming regions are difficult to determine, and span a wide range (Bressert et al., 2010; King et al., 2012), but many are thought to be at least $\sim$ 100 M*⊙* pc*-3* at the epoch of star formation (Parker, 2014). Parker (2014) shows that two regions with similar present-day densities at present times may have originally had very different initial densities because initially dense regions expand much faster than lower-density regions due to two-body relaxation. They compare the present-day stellar densities and amount of spatial substructure in seven star-forming regions, including the Orion Nebula Cluster and Upper Scorpius – which both contain massive stars that could act as external photoionising sources – to infer the likely range of initial stellar densities in each of these star-forming regions and all are consistent with having an initial density in the range 10 – 1000 M*⊙* pc3.

External processes, such as close stellar interactions, can also cause proplyds to be truncated or destroyed, as well as disrupting the orbits of fledgling planets (Armitage, 2000; Bonnell et al., 2001; Scally & Clarke, 2001; Adams et al., 2006; Olczak et al., 2008; Parker & Quanz, 2012; Rosotti et al., 2014; Vincke et al., 2015; Portegies Zwart, 2016; Winter et al., 2018a). The density of the star forming region will affect the rate of stellar interactions, with stars in low-density environments experiencing fewer dynamical interactions than higher density environments (Wright et al., 2014; Bressert et al., 2010). Furthermore, star-forming regions can contain massive stars (> 15 M*⊙*), whose intense far ultra-violet (FUV) and extreme ultra-violet (EUV) radiation fields are significantly higher than those in the interstellar medium (Armitage, 2000; Adams et al., 2004; Fatuzzo & Adams, 2008). This high-energy radiation heats the gaseous material of the upper layers of the disc until the thermal energy of the heated layer exceeds the gravitational potential of the disc, causing it to escape as a photoevaporative wind (Hollenbach et al., 1994; Johnstone et al., 1998). This mass loss will affect the evolution of protoplanetary discs, and reduce the reservoir of material available to form gas giant planets (Haworth et al., 2018a).

The effects of external photoevaporation appear to be observed in nearby star-forming regions, such as the Orion Nebula Cluster (ONC) (McCaughrean & O’dell, 1996; Eisner et al., 2016; Eisner et al., 2018) and $\sigma$ Orionis (Ansdell et al., 2017). The ONC has been preferentially observed due to its proximity to Earth ( $\sim$ 415 pc) and because the discs can be viewed in silhouette due to the bright nebulous background. Studies have also shown that 80 – 85% of stars within the ONC host protoplanetary discs (Bally et al., 2000; Lada et al., 2000), making the ONC a favorable target for studying disc evolution and dispersal (though see Clarke, 2007, for an alternative interpretation which posits that the ONC proplyds are merely ionisation fronts of material left over from discs that are almost destroyed).

Due to two-body and violent relaxation, initially highly substructured star-forming regions can evolve to smooth and centrally concentrated clusters after only a few Myr (Allison et al., 2010; Parker et al., 2014). Furthermore, two clusters that presently have similar densities may have had very different initial densities because initially very dense clusters expand faster than lower density counterparts. As mentioned before, the initial density will affect the rate at which protoplanetary discs are disrupted and destroyed due to stellar interactions. However, how much the initial density and substructure of a star-forming region affects the rate of protoplanetary disc dispersal due to external photoevaporation has yet to be studied.

Previous studies into the effects of external photoevaporation on protoplanetary discs in star-forming regions have tended to calculate the background UV radiation without directly calculating the disc mass-loss (Armitage, 2000; Adams et al., 2004). Scally & Clarke (2001) did calculate mass-loss rates in simulations specifically tailored to match the ONC, but assumed rather low stellar densities ( $\sim$ 40 M⊙ pc*-3*), whereas Parker (2014) suggests that the initial density of the ONC may have been much higher (>100 M⊙ pc*-3*).

These previous studies of external photoevaporation used spherically smooth spatial distributions with primordial mass segregation to model the environment of the ONC as observed today. However, observations of star forming regions show that stars form in highly substructured filamentary environments, where the stars are moving with subvirial velocities. The initial net motion and spatial structure of a star-forming region will govern its future evolution, and by extension the degree to which planet formation is hindered. These initial conditions lead to dynamical mass segregation on timescales of the age of the ONC (Allison et al., 2010; Parker et al., 2014), negating the requirement for primordial mass segregation (Bonnell & Davies, 1998).

Here, we focus on initial conditions that more closely reflect observations of young star forming regions (Cartwright & Whitworth, 2004), and determine how much external photoevaporation affects the evolution of protoplanetary discs.. Therefore, we do not centrally concentrate our massive stars, but randomly distribute them in our simulated star-forming regions. We run suites of simulations that cover a range of initial conditions, with varying initial density, spatial distribution and net bulk motion (virial ratio). We then calculate and compare the mass-loss rates due to external photoevaporation for each set of initial conditions.

The paper is organised as follows. In Section 2 we describe our $N$ -body simulations, protoplanetary disc assumptions and our external photoevaporation prescription; in Section 3 we present our results; we provide a discussion in Section 4 and we conclude in Section 5.

2 Method

In this section we describe our method to select low-mass star-forming regions containing massive stars, before describing the subsequent $N$ -body and stellar evolution of these regions.

2.1 Creating low-mass star-forming regions

We adopt two different masses ( $100$ or $1000$ M*⊙*) for our star-forming regions and populate these regions with stars drawn randomly from the initial mass function (IMF) parameterised in Maschberger (2013), which has a probability density function of the form:

[TABLE]

Here, $\mu=0.2$ M⊙ is the average stellar mass, $\alpha=2.3$ is the Salpeter (1955) power-law exponent for higher mass stars, and $\beta=1.4$ is used to describe the slope of the IMF for low-mass objects (which also deviates from the log-normal form; Bastian, Covey & Meyer, 2010). Finally, we sample from this IMF within the mass range $m_{\rm low}=0.1$ M⊙ to $m_{\rm up}=50$ M⊙.

We use a “soft-sampling” technique for sampling the IMF (Elmegreen, 2006), which implies that the only formal limit on the most massive star that can form is that of the upper limit of the IMF (Parker & Goodwin, 2007). From an ensemble of Monte Carlo simulations, we find that typically, a $1000$ M*⊙* star-forming region will contain five massive stars (M ${}_{\star}>15$ M*⊙*).

Low-mass star-forming regions rarely contain massive stars; however, if the only limit on the mass of the star that can form is the total mass of the star-forming region itself then occasionally we would expect a low-mass star-forming region to contain one or more massive stars and such regions are observed (e.g. $\gamma^{2}$ Vel, a low-mass region containing at least two massive stars, Jeffries et al., 2014). Note that we are not explicitly attempting to model the $\gamma^{2}$ Vel star-forming region, which harbours a dense (sub)cluster within a more diffuse region. Instead, we are pointing out that observational examples of low-mass star-forming regions such as this and $\sigma$ Orionis occasionally contain massive stars and their photoevaporative effects on discs in these low-mass star-forming regions could be important.

In order to demonstrate the importance of low-mass star-forming regions, let us consider the demographics of star-forming regions from randomly sampling the underlying probability distributions. The observed mass function of star-forming regions follows a power-law of the form

[TABLE]

where $N$ is the number of star-forming regions with mass $M_{\rm cl}$ (Lada & Lada, 2003), which implies that there many more low-mass star-forming regions compared to high-mass regions.

We follow the procedure in Parker & Goodwin (2007) and Nicholson & Parker (2017) and sample $1\times 10^{6}$ star-forming regions in the mass range 50 – $10^{5}$ M⊙. Of these, $\sim$ 1200 lie in the mass range 1000 $\pm$ 10 M⊙, and $\sim$ 15 000 lie in the mass range 100 $\pm$ 10 M⊙.

We then randomly populate our star-forming regions with stars drawn from Eqn. 1, until the total mass of stars equals or exceeds the chosen star-forming region mass from Eqn. 2.

Nicholson & Parker (2017) found that $\sim$ 10 per cent of low-mass star-forming regions contain at least one massive star ( $>$ 15 M⊙) when using the “soft-sampling” technique described above, and 1 per cent of low-mass regions contain two massive stars. Furthermore, when taking into account the decreasing probability of forming a high-mass star-forming region (Eqn. 2), the number of low-mass ( $M_{\rm cl}=100\pm 1$ M⊙) regions containing at least one massive star is $\sim$ 3100, which is actually greater than the total number of high-mass ( $M_{\rm cl}=1000\pm 10$ M⊙) star-forming regions (1200). Of these 1200 high-mass regions, $\sim$ 1000 contain at least one massive ( $>$ 15 M⊙) star.

If we translate these numbers into the total number of stars that may be affected by photoevaporation, the average number of stars in our $M_{\rm cl}=100\pm 1$ M⊙ star-forming regions containing at least one massive star is $\sim$ 110, so in the 3100 low-mass regions that contain at least one massive star there are $\sim$ 341 000 stars in total. The average number of stars in our high mass star-forming regions ( $M_{\rm cl}=1000\pm 10$ M⊙) is $\sim$ 1710, and so the 1000 regions that contain at least one massive star host a total of $\sim$ 1 700 000 stars that could be affected by photoevaporation. In short, the fraction of stars originating in low-mass star-forming regions containing at least one massive star is $\sim$ 20 per cent of the total number of stars originating from high-mass regions containing at least one massive star. If we stipulate that the high-mass regions must contain three or more massive stars, only $\sim$ 580 regions out of 1200 fulfil this criteria and host a total of 986 000 stars. The fraction of stars originating in low-mass star-forming regions containing at least one massive star is $\sim$ 35 per cent of the total number of stars originating from high-mass regions containing at least three massive stars.

Crucially, this makes no assumption about the disruption and dissolution of these star-forming regions, and how many stars from each type of region eventually enter the Galactic field. The Galactic potential will influence the destruction of low-mass star-forming regions much more than high-mass regions (Binney & Tremaine, 2008), which take longer to dissolve into the Galactic field (and some remain as long-lived open clusters). Therefore, the majority of planet-hosting Field stars may come from lower-mass regions.

Given their significant contribution to the integrated stellar mass function, we therefore also investigate low-mass star-forming regions ( $100$ M*⊙*) that contain either one or two massive stars – these represent an unusual sampling of the IMF but allow us to investigate the effects of photoevaporation in less populous star-forming regions.

Hence we have three different star-forming region set-ups; a 100 M*⊙* region with one massive star (38 M⊙), a 100 M*⊙* region with two massive stars (42 and 23 M⊙) and one 1000 M*⊙* region with 5 massive stars (43, 33, 26, 17 and 17 M*⊙*). These regions were selected as the median outcomes of Monte Carlo sampling of $1\times 10^{6}$ star-forming regions (Nicholson & Parker, 2017), and then filled with stellar masses drawn from the IMF (Maschberger, 2013). We then selected the median regions in terms of the total number of stars from within the mass ranges of 100 $\pm 1$ M⊙ and 1000 $\pm 10$ M⊙, with the stipulation that they had to contain massive stars. For the 100 M⊙ mass regions, we specifically selected the median region containing one and two massive stars. For the 1000 M⊙ we selected the average cluster that contained three or more massive stars.

The external photoevaporation prescriptions we will adopt in this work are those from Scally & Clarke (2001) which only weakly depend on the adopted stellar IMF, but in a future paper we will assign a FUV flux and an EUV flux to each intermediate/high-mass star based on its mass and then determine the respective fluxes incident on every low-mass star and use the recent FRIED grid of models (Haworth et al., 2018b) to determine mass-loss for individual discs.

2.2 $N$ -body simulations

Our simulations are created with initial substructure by following the box-fractal method in Goodwin & Whitworth (2004). We use a range of fractal dimensions for varying amounts of substructure: $D=1.6$ (highly sub-structured), $D=2.0$ (moderately sub-structured), and $D=3.0$ (smooth). The method also correlates stellar velocities on local scales so that nearby stars have similar velocities, but more distant stars can have a wide range of different velocities.

The velocities are then scaled to the virial ratio $\alpha_{\rm vir}$ , where $\alpha_{\rm vir}=T/|\Omega|$ ; $T$ is the total kinetic energy and $\Omega$ is the total potential energy of the stars. A range of virial ratios are used: $\alpha_{\rm vir}=0.3$ (sub-virial, or collapsing), $\alpha_{\rm vir}=0.5$ (virial equilibrium), and $\alpha_{\rm vir}=0.7$ (supervirial, or expanding). Note that this virial ratio determines the net bulk motion, i.e. whether the star-forming region will collapse or expand. The correlated velocities on local scales mean that the local velocity dispersion can be subvirial, facilitating a violent relaxation process as the star-forming region evolves. Such local subvirial velocity dispersions are observed in the earliest stages of star formation (Larson, 1981; Peretto et al., 2006; Foster et al., 2015).

We create star-forming regions with stellar densities of 100 M*⊙* pc*-3* or 10 M*⊙* pc*-3* for the 1000 M⊙ regions; for the 100 M⊙ regions we set an initial density of 100 M*⊙* pc*-3*. Such densities bracket the range observed in present-day star-forming regions (Bressert et al., 2010) as well as allowing for potentially higher primordial densities (Parker, 2014).

Finally, we created a set of simulations with a Plummer sphere distribution (Plummer, 1911) to facilitate comparisons with previous studies. We use the same IMF from our 1000 M*⊙* simulations to create two clusters with Plummer sphere distributions that have initial densities of 10 and 100 M*⊙* pc*-3*.

We evolve each of our star-forming regions for 10 Myr using the kira integrator within the Starlab environment (Portegies Zwart et al., 2001). Stellar evolution is implemented using the SeBa look-up tables (Portegies Zwart & Verbunt, 1996). No binary or multiple stellar systems are included in these simulations. To gauge the amount of stochasticity in the disc photoevaporation, we run 20 realisations of the same initial conditions, identical apart from the random number seed used to assign the positions and velocities.

2.3 Protoplanetary discs and external photoevaporation

The mass loss rate of discs at a certain distance from a neighboring massive star is determined by the strength of the star’s FUV ( $h\nu$ < 13.6 eV) and EUV ( $h\nu$ > 13.6 eV) fluxes at that distance. Mass loss due to FUV photons is caused by heating the circumstellar disc, which creates an unbound neutral layer that can drift towards the ionization front, where it meets the EUV field. FUV is independent of the distance from the massive star because the only requirement is that the FUV flux is strong enough to heat the disc above its escape velocity. With EUV, the mass loss rate depends on the EUV flux and so is directly dependent on the distance from the massive star(s).

We use the same prescriptions for FUV and EUV photoevaporation as Scally & Clarke (2001):

[TABLE]

where $r_{d}$ is the radius of the disc in astronomical units, $\Phi_{i}$ is the ionizing EUV photon luminosity from each massive star in units of 1049 s*-1* and $d$ is the distance from the massive star in parsecs. The UV photon rate ( $\Phi_{i}$ ) for the massive stars (> 15M*⊙*) is dependent on stellar mass and we use the values from Vacca et al. (1996) and Sternberg et al. (2003).

These photoevaporation rates were derived assuming a disc density profile $\Sigma\propto r_{d}^{-2}$ (Hollenbach et al., 2000; Hartmann, 2009); however, our analysis does not take into account the evolution of the surface density profile if the disc radius were to change significantly.

Observations of star forming regions show that disc radii can extend to several 100s AU (e.g. Ansdell et al., 2018). However, the typical initial radius of proplyds is still debated in the literature and therefore we sample a wide range of initial radii: 10, 50, 100, 200 and 1000 AU. We adopt a single value for each disc radius, focusing primarily on 100 AU discs, and then repeat the analysis for the five values. In reality, the radius of the disc will change due to internal processes such as viscous evolution, and due to internal and external photoevaporation, but we are unable to account for this (and the changing disc density profile) in our $N$ -body simulations.

The initial disc masses are also debated, with theoretical constraints from the Minimum Mass Solar Nebula (MMSN, Weidenschilling, 1977; Hayashi, 1981) and observations (Andrews et al., 2013; Ansdell et al., 2016) suggesting an upper limit of $M_{\rm disc}=0.01$ M*⋆. Following the example in Scally & Clarke (2001), the initial disc masses in our simulation are 10 per cent of the host star mass ( $M_{\rm disc}$ = 0.1 M⋆*). Current observations suggest that disc masses are $\sim$ 1 per cent of the host star mass. We select 10 per cent so that we are sampling the upper range of the disc masses. While we do not account for accretion onto the protoplanetary discs, our discs are large enough in mass that we can neglect the accretion onto the disc as it will be minimal in comparison. For completeness, we ran a set of simulations where the disc masses were 1 per cent of the stellar host mass, which is more consistent with the MMSN estimates.

We subtract mass from our discs according to Eqns. 2 and 3. The rate of photoevaporation due to EUV radiation is dependent on distance from the ionising source, $d$ , whereas the photoevaporation rate due to FUV is largely independent of distance from the source (Störzer & Hollenbach, 1999). Following Störzer & Hollenbach (1999) and Scally & Clarke (2001) we apply Eqn. 2 if the disc-hosting stars are within 0.3 pc of the ionising source, noting that this distance is calibrated to models where $\theta^{1}$ C Ori is the most massive star (40 M⊙), which is commensurate with the most massive star in our simulations. However, we note that this is likely an underestimate of the amount of photoevaporation due to FUV fields in star-forming regions (Facchini et al., 2016; Haworth et al., 2018a).

3 Results

We focus on 1000 M*⊙* star-forming regions, which typically contain $\sim$ 5 massive stars (M ${}_{\star}>15$ M*⊙) that act as photoionising sources (c.f. Scally & Clarke, 2001). Our 1000 M⊙* cluster contains 5 massive stars; 43.2, 32.7, 25.7 and two 17 M*⊙* with $\Phi_{i}$ values of 1.1, 0.47, 0.19 and $\sim$ 0.013 respectively. We focus on the results for two different initial stellar densities, $\sim$ 10 M*⊙* pc*-3* and $\sim$ 100 M*⊙* pc*-3* and, apart from the final section, the assumed initial mass for every disc is $M_{\rm disc}=0.1$ M*⋆*.

We present the results from varying different initial properties within the star forming regions, focusing on the spatial distribution (fractal dimension, $D$ ) and net bulk motion (virial ratio, $\alpha_{\rm vir}$ ). We focus on protoplanetary discs that have a radius of 100 AU, however the effects of external photoevaporation on discs with different radii and mass are discussed later in Sections 3.3 and 3.4.

We compare our results to the observed disc fractions in both Haisch et al. (2001) and Richert et al. (2018) using ages from the models in Siess et al. (2000). We discuss the caveats associated with these models in Section 4.

We later present two low mass clusters (100 M*⊙) with an initial density of $\sim$ 100 M⊙* pc*-3* that are subvirial ( $\alpha_{\rm vir}=0.3$ ) and highly substructured ( $D=$ 1.6). Our clusters contain either one (38 M⊙) or two (42 M⊙ and 23 M⊙) massive stars. The corresponding $\Phi$ values are $\Phi=0.76$ and $\Phi=1.01,0.11$ respectively.

3.1 Substructure in star-forming regions

We first present the results from four simulations of star forming regions, where in each simulation the star forming region has a different initial spatial distribution; $D=$ 1.6 (highly substructured), $D=2.0$ (moderately substructured), $D=3.0$ (smooth) and a Plummer sphere spatial distribution. In all simulations the star-forming region is subvirial ( $\alpha_{\rm vir}=0.3$ ).

In Fig. 1 we show the average fraction of stars that have retained their (100 AU) discs from 20 runs of each simulation, where the initial substructure of the star-forming region is varied. We present the results for regions with two different initial stellar densities; 10 and 100 M*⊙* pc*-3* respectively.

Fig. 1(a) shows the results from a star forming region with an initial density of 10 M*⊙* pc*-3*. Within highly sub-structured regions ( $D=1.6$ ), the time taken for half of the stars to lose their discs due to external photoevaporation is 2.12 Myr. In moderately sub-structured regions ( $D=2.0$ ), this time increases to 2.60 Myr. However, the average percentage of remaining discs with time in both remain relatively similar throughout the 10 Myr. For regions with an initially smooth and spherical distribution ( $D=3.0$ ), the time taken for half of the discs to disperse is 3.62 Myr. Discs within Plummer spheres have the longest lifetimes (3.85 Myr), with an average of $\sim$ 29.7 per cent of discs surviving for longer than 10 Myr.

Plummer (1911) models (and other models that describe smooth star clusters such as a King (1966) profile or an Elson et al. (1987) profile) are intended the model dynamically relaxed systems, whereas young star-forming regions have yet to relax. Therefore, even a smooth box fractal ( $D=3.0$ ) contains kinematic substructure, which causes the dynamical evolution of such a region to be more violent than a smooth Plummer sphere. It is therefore unsurprising that fewer discs survive in kinematically substructured fractal regions than in Plummer spheres with a similar density.

Fig. 1(b) shows the results for star forming regions with an initial density of 100 M*⊙* pc*-3*. For discs in the highly sub-structured regions ( $D=1.6$ ), the time taken for half of the stars within the cluster to lose their protoplanetary discs is 0.87 Myr. The majority of discs within the highly substructured region ( $D=1.6$ ) are dispersed after 10 Myr, with $\sim$ 6 per cent surviving for the length of the simulation. The majority of discs within smooth, spherical regions are also dispersed within a short time frame, with only $\sim$ 3 per cent remaining after 10 Myr.

In Table 1 we summarise the average time taken for half of the stars in a star forming region to lose their discs for each spatial distribution and Table 2 summarises the percentage of discs remaining at the end of the 10 Myr simulation.

In the low density simulations, regions with more spatial substructure photoevaporate discs faster than smoother regions of a comparable density (Fig. 1(a)) . The reason for this is that the more substructured regions are initially further from dynamical equilibrium than the smooth regions, and low-mass stars undergo more close interactions with the high mass stars as the regions relax.

Interestingly, in the high-density simulations (Fig. 1(b)), whilst the fraction of discs remaining after 10 Myr is lower than in the low density simulations, the initially more substructured star-forming regions contain more discs than the smooth regions after 10 Myr (and their disc-destruction half-life is longer, see Table 1). We attribute this to the higher ejection rate of massive stars in dense, substructured star-forming regions (Parker et al., 2014), which means that some of the ionising sources are no longer near the majority of the proplyds as early as 1 Myr after the start of dynamical evolution.

3.2 Virial ratio

We explore how changing the net bulk motion of the star-forming region affects the rate of disc dispersal due to external photoevaporation. We run simulations of our star-forming region with three different virial ratios; 0.3 (subvirial, or collapsing), 0.5 (virial equilibrium), and 0.7 (supervirial, or expanding). We keep the fractal dimension constant, adopting $D$ = 2.0 and as before we analyse simulations with two different initial densities; 10 M*⊙* pc*-3* and 100 M*⊙* pc*-3*.

In Table 3 we summarise the average time taken for half of the stars in each region to lose their (100 AU) discs for a given bulk virial ratio. Table 4 shows the percentage of discs remaining at the end of the 10 Myr simulation. Fig. 2 shows the the average fraction of stars that have retained their discs from the 20 runs of each simulation for a star forming region where we vary the initial bulk motion (virial ratio).

Fig. 2 (a) shows the average mass loss rate in a star forming region with an initial density of 10 M*⊙* pc*-3*. The time taken for half of the stars within a collapsing (sub-virial) star-forming region to lose their discs is 2.60 Myr. In regions that are expanding (supervirial), this time increases to 3.53 Myr. The percentage of discs within the subvirial region after 10 Myr is 17.75%, in comparison to discs within an expanding region where the percentage rises to 32.67%.

The initial net bulk motion of low density star-forming regions affects the amount of discs that are photoevaporated due to external radiation, with subvirial regions evaporating more discs at a faster rate than either virialised or supervirial regions.

Fig. 2(b) shows the results for a star-forming region with an initial density of 100 M*⊙* pc*-3*. The time taken for half the stars within a collapsing region to lose their discs is 0.67 Myr. This time is similar for regions in virial equilibrium and expanding regions (0.68 and 0.63 Myr respectively). The lower disc half-life for the supervirial regions could again be due to massive stars being ejected in the (sub)virial regions. The percentage of discs remaining after 10 Myr in sub-virial star forming regions is 2.27% whereas in regions where the net motion is expansive, this is increased to 5.00%.

3.3 Disc radii

Here we present the rates of disc dispersal for different initial disc radii in a star-forming region with two different initial densities (10 and 100 $M_{\odot}$ pc*-3*). The region has a fractal dimension of $D=2.0$ (moderately substructured) and a viral ratio of $\alpha_{\rm vir}=0.3$ (sub-virial).

Fig. 3 shows the percentage of protoplanetary discs with initial radii ranging between 10 – 1000 AU that have some remaining mass over 10 Myr in a 1000 M*⊙* star-forming region with different initial stellar densities; 10 M*⊙* pc*-3* in Fig 3(a) and 100 M*⊙* pc*-3* in Fig. 3(b).

In the lower-density star-forming regions (Fig. 3(a)), the time taken for half of the 100 AU discs to be completely photoevaporated is 2.60 Myr. Discs with radii of 10 AU have much greater lifetimes, with an average of $\sim$ 77 per cent of discs surviving the full length of the simulation. The majority of discs with very large radii (1000 AU) are still depleted within very short timescales. Disc depletion rates begin to switch off after $\sim$ 4 Myr due to a combination of a large decrease in density of the star forming region, which peaks at $\sim$ 2 Myr, and the death of the most massive star at 4.33 Myr.

Fig. 3(b) shows that the majority (90 per cent or more) of discs with radii $>$ 10 AU are completely photoevaporated before the end of the 10 Myr simulation in moderately dense star-forming regions. The time taken for half of the stars in the region to lose their 100 AU discs is $0.67$ Myr. The vast majority of the largest discs (1000 AU) are photoevaporated completely within 2 Myr, with half of the stars in the region losing their discs within < 0.1 Myr.

We also ran simulations for low mass star-forming regions (100 M*⊙) with an initial density of $\sim$ 100 M⊙* pc*-3*. These low-mass regions contain two massive stars (42 M⊙ and 23 M⊙), which represents an unusual sampling of the IMF (Parker & Goodwin, 2007), but is observed in nature (e.g. $\gamma^{2}$ Vel, Jeffries et al., 2014). In these low mass regions, half of discs with radii of 100 AU dissipated in 0.51 Myr. This is comparable to the 1000 M*⊙* regions with a similar density. We will further explore the effects of different stellar IMFs on disc dispersal in a future paper.

3.4 Disc Masses

We have assumed that the disc masses are 10 per cent of the host star’s mass, which is likely to be an overestimate and various studies suggest that the disc mass is as low as 1 per cent of the host star’s mass (Weidenschilling, 1977; Hayashi, 1981; Andrews et al., 2013).

In Fig. 4 we show the results for a star-forming region with our two different initial densities (10 and 100 M*⊙* pc*-3* respectively), where the initial disc masses are set to $M_{\rm disc}=0.01$ M*⋆*.

Fig. 4(a) shows that on average the time taken for half of the stars within the low density star forming region to lose their 100 AU disc is $0.71$ Myr, less than half of the time taken for discs with 10 per cent of the mass of their stellar host. For discs with a radii of 10 AU, the half life is 3.31 Myr.

The timescale for half of the 100 AU discs to dissipate in the moderately dense (100 M*⊙* pc*-3*) star forming region (see Fig. 4b) is $\sim$ 0.14 Myr. For discs with a radius of 10 AU, the half life is $\sim$ 0.84 Myr. Less than 5 per cent of 10 AU discs survive for more than 3 Myr.

3.5 Mass of star-forming regions

We also ran simulations for two different low mass star-forming regions (100 M*⊙) with an initial density of $\sim$ 100 M⊙* pc*-3*, which were sub-virial ( $\alpha_{vir}$ = 0.3) and substructured ( $D$ = 1.6). These low-mass regions contain one (38 M⊙) or two (42 M⊙ and 23 M⊙) massive stars, which represents an unusual sampling of the IMF (Parker & Goodwin, 2007), but is observed in nature (e.g. $\gamma^{2}$ Vel, Jeffries et al., 2014). Our expectation from randomly sampling the IMF is that 10 per cent of all star-forming regions can host a massive star, and 1 per cent of regions will host two massive stars. We note that the lack of massive star(s) in any star-forming region would preclude disc destruction from photoevaporation, though as discussed in Section 2 it is unclear which type of star-forming region (in terms of total mass, $M_{\rm cl}$ ) contributes the most (planet hosting) stars to the Galactic field.

In both of these low mass regions, half of discs with radii of 100 AU dissipated before $\sim$ 1 Myr (Fig. 5). The time taken for half of the discs to be destroyed in a region with one massive star is 0.95 Myr. This time is reduced to 0.37 Myr for the cluster with 2 massive stars.

At the end of the 10 Myr simulation, 15.5 per cent of discs within the region with one massive star are surviving. Within the region containing 2 massive stars, less than 5 per cent of discs are remaining, double the number of discs remaining in higher mass regions with the same initial conditions.

4 Discussion

The initial conditions of a star-forming region will affect the rate at which protoplanetary discs are photoevaporated due to the radiation from nearby massive stars. The initial substructure and net bulk motion of a star-forming region impacts the rate of disc dispersal.

4.1 Changing the initial conditions of star-forming regions

In our low-density simulations, highly substructured ( $D=1.6$ ) regions disperse half of the protoplanetary discs within 1.51 Myr, more than twice as fast as smooth ( $D=3.0$ ) regions. In simulations with a Plummer sphere distribution, more than 30 per cent of discs remain at the end of the 10 Myr simulation, almost double that of discs within highly substructured clusters.

At these low densities, the degree of substructure matters because a more substructured star-forming region is further from dynamical equilibrium than a smooth region. When this occurs, a low-mass star is likely to have more close encounters with a massive ionising star than in a smooth region.

In moderately dense initial conditions (100 M*⊙* pc*-3*), the difference in the fraction of discs that are photoevaporated between different initial spatial distributions decreases greatly, although regions with a Plummer sphere distribution retain more of their discs than regions with initial substructure. However, the average of all runs indicates that the amount of initial substructure has little effect on the survival rates of discs at these densities and fewer than 50 per cent of discs remain after 1 Myr.

The effect of changing the net bulk motion of the star-forming region has a similar impact on the rate of disc dispersal as the initial substructure has. For low density regions (10 M*⊙* pc*-3*), the difference between the amount of discs surviving within a collapsing and an expanding star forming region is $\sim$ 15 per cent, with the collapsing regions enabling more photoevaporation than in expanding regions. Again, approximately double the number of discs remain in expanding clusters than in collapsing clusters. In moderately dense clusters it is similar, with the difference being $\sim$ 3 per cent.

For low mass star-forming regions (100 M*⊙), disc dispersal rates are similar to those of in higher mass regions. Whilst the UV field strength can vary due to different realisations of the IMF (Fatuzzo & Adams, 2008), these low mass regions show that the mere presence of a high mass star (> 15 M⊙*) will cause disc lifetimes to be shortened dramatically.

Our simulations are set up to mimic the observations of star formation in filaments, where the pre-stellar cores have subvirial motion (Larson, 1981; Foster et al., 2015). The local velocity dispersion is therefore always subvirial to some degree, and because mass-loss due to photoevaporation is so fast (eqns. 2 and 3), most of the photoevaporation occurs during the substructured phase of a star-forming region.

Previous studies investigating the effects of external photoevaporation on disc dispersal rates assumed smooth and centrally concentrated spatial distributions (Scally & Clarke, 2001; Adams et al., 2004, 2006; Winter et al., 2018b), replicating environments like the present-day conditions of the ONC. However, using the present-day spatial and kinematic distributions to model star clusters may not accurately replicate the dynamical history of the star-forming region from which the cluster formed (Parker et al., 2014).

We cannot provide a direct comparison with the work by Fatuzzo & Adams (2008) and Winter et al. (2018b) because these authors assume initially smooth initial conditions, and we have shown that the severity of photoevaporation depends on the degree of initial substructure, as well as the initial positions of the most massive stars. We have distributed the massive stars randomly within the substructure, and after dynamical evolution these massive stars migrate towards the centre of the star cluster as it forms (Allison et al., 2010; Parker et al., 2014, 2016).

Furthermore, we have used the photoevaporation prescription from Scally & Clarke (2001), rather than determine the photoevaporation rate from the EUV/FUV fluxes as a function of the flux in the interstellar medium (the so-called $G_{0}$ value, $1.6\times 10^{-3}$ erg s*-1* cm*-2* Habing, 1968). In comparison with Scally & Clarke (2001), we find that discs are destroyed earlier in ONC-type regions because the initial densities are higher (in line with current observations, Parker, 2014), and the star-forming regions are substructured (Cartwright & Whitworth, 2004).

An initially highly substructured star-forming region can become smooth and centrally concentrated within a few Myr due to a combination of violent and two body relaxation. Protoplanetary discs in these highly substructured environments will be photoevaporated at faster rates than discs within initially smooth regions. Even though they will both appear smooth within a few Myr, the percentage of discs remaining, and possibly the population of planets within the regions, will vary greatly.

The initial density of the cluster has the largest effect on the disc dispersal rate due to external photoevaporation. The ‘moderately-dense’ clusters reflect the likely initial densities of many star-forming regions (Parker, 2014). However, we find that these ‘moderately-dense’ environments are very destructive for protoplanetary discs and evaporate discs at rates faster than suggested by observations (compare the black points in Figs. 1 and 2 with our simulated data). Our results suggest that protoplanetary discs (or at least their gas content) would always be significantly depleted in moderately dense (100 M*⊙* pc*-3*) star-forming regions, if those regions contain massive stars.

Haisch et al. (2001) finds that the fraction of disc-hosting stars in young star-forming regions falls to 50 per cent after $\sim$ 3 Myr whereas Richert et al. (2018) find that after only $\sim$ 2 Myr half of the discs remain in their observed regions. In comparison, more than half of the discs in our simulations that are within dense environments are destroyed within $\sim$ 1 Myr. One interesting data point in the Haisch et al. (2001) sample is the ONC. With an age of $\sim$ 1 Myr, the centre of the ONC contains 4 massive stars and a density of $\sim$ 400 M*⊙* pc*-3* (Hillenbrand & Hartmann, 1998). Observations of this part of the ONC suggest that $\sim$ 80 – 85 per cent of stars within the cluster are surrounded by bright ionization fronts, interpreted to be discs, with radii of $\sim$ 1000 AU (Bally et al., 2000; Lada et al., 2000).

The age of our simulated regions, where 80 per cent of stars still possess a 100 AU disc with some mass, is 0.48 Myr – likely to be less than half the age of the ONC. This suggests that the massive stars within the ONC should have destroyed the majority of 100 AU discs. From similar arguments, Clarke (2007) concluded that the possible discs in the ONC with radii >10 AU are likely to be merely ionisation fronts, containing little mass. Our simulations with different initial disc radii show that the radius of the disc will greatly affect the rate at which it is photoevaporated (see also Clarke, 2007) due to the dependence on disc radii within the FUV and EUV photoevaporation prescriptions.

Recent surveys suggest that most stars in the Galactic field host planets, and many of these are gas or ice giants (Mayor et al., 2011). This implies that the majority of planet forming discs were able to survive a significant amount of time in their birth environment. Our simulations suggest that this is only possible in low-density regions that contain no photoionising sources (i.e. massive stars). Therefore, (giant) planet formation must occur on very rapid timescales (<1 Myr), or stars that host giant planets must have formed in very benign environments.

Many observed protoplanetary discs are located in low-mass, low-density star-forming regions (Andrews et al., 2013; Ansdell et al., 2018) and would be unaffected by external photoevaporation. However, many star-forming regions are typically moderate-density ( $\sim$ 100 M*⊙* pc*-3*) environments (Parker, 2014) and our results suggest that the majority of protoplanetary discs in star-forming regions with these densities do not survive for long enough periods of time to form giant planets.

4.2 Caveats

There are several caveats to our results, which we discuss below.

The effects of external EUV radiation on protoplanetary discs can be reduced when thick winds are present, caused by FUV heating of the disc (Alexander et al., 2014). However, the majority of the disc mass loss occurs due to FUV radiation. We repeated our analysis without EUV photoevaporation and find the disc dispersal rates to be similar.

It is possible that we are overestimating the amount of photoevaporation from massive stars. However, recent research suggests that the prescriptions used here are actually underestimating the amount of FUV radiation that discs receive (Facchini et al., 2016; Haworth et al., 2018a). As FUV is the dominant source of external photoevaporation, the protoplanetary discs in our simulations could dissipate on even shorter timescales.

Star formation is an inherently inefficient process, with typically only $\sim$ 30 per cent of the mass of a giant molecular cloud converted into stars. Young star-forming regions are observed to contain a large amount of dust and gas, which could shield the proplyds from significant photoevaporation. At these early stages the stellar density within the substructure is highest, and is therefore when the largest percentage of stars are in closest proximity to the massive stars. However, hydrodynamical simulations of star-forming regions show that massive stars blow large ( $\sim$ pc-scale) cavities within the gas on short time scales (Dale et al., 2013), and so low-mass disc-hosting stars that would be affected by EUV/FUV radiation will likely reside in the cavities blown out by the massive stars. If the gas and dust could shield the disc, this would protect them for a very short period of time (Gorti et al., 2015). Whether this grace period would be long enough to allow gas-rich giant planets to form is uncertain.

Given that most star-forming regions have stellar densities above a few M⊙ pc*-3* (Bressert et al., 2010), external photoevaporation will detrimentally affect protoplanetary discs in any star-forming region that contains massive stars. This implies that star-forming regions that do not contain massive stars are more likely to form giant planets, but we note that massive stars appear necessary in order to deliver short-lived radioisotopes to the young Solar System (Lugaro et al., 2018). The number of massive stars in a star-forming region appears to only be limited by the mass of the star-forming cloud (Parker & Goodwin, 2007), but this also means that low-mass star forming regions ( $<10^{4}$ M⊙) stochastically sample the IMF, meaning that our simulations cannot be described as ‘typical’ star-forming regions.

Quantifying disc dispersal is further complicated by how difficult it is to determine the ages of young stars, especially before 1 Myr (Siess et al., 2000). We use the stellar ages from the Siess et al. (2000) model. However, models of pre-main sequence stellar evolution calculate different ages depending on the physics that is implemented. Of the three models presented in Richert et al. (2018) we use the ages from Siess et al. (2000) so that we are comparing the lower end of cluster ages to our simulations. The average stellar age calculated for the clusters in Richert et al. (2018) is significantly shorter than in more recent models from Feiden (2016). By using these lower age limits, we more than halve the possible average life times of the discs within the observed clusters.

Furthermore, recent work by Bell et al. (2013) suggests that the ages of pre-main sequence stars may be underestimated by a factor of two, meaning that the observed discs (e.g. Haisch et al., 2001) could be a factor of two older. This would make it even more surprising that discs would remain around low-mass stars, if those stars form in regions containing massive stars.

There is also the question of how quickly the photoionising massive stars form. In the competitive accretion models (Bonnell et al., 2001), massive stars gradually gain in mass over $\sim$ 1 Myr (Wang et al., 2010), suggesting high-mass stars form later than low-mass stars (Tan et al., 2014), which would in turn decrease the amount of time low-mass stars spend near the photoionising sources (Dale et al., 2012; Dale et al., 2014). In our simulations all stars form simultaneously, and therefore the disc-hosting low-mass stars do not have this grace period, which would increase disc lifetimes.

The growth of planetesimals into planets can be greatly accelerated by the accretion of cm-scale pebbles. Johansen & Lambrechts (2017) show that once a 10 ${}^{-2}\,M_{\earth}$ planetesimal has formed it can grow to Jupiter mass in 1 Myr when starting as far out as about 15 AU. An initial phase of accreting pebbles forms a 10 $\,M_{\earth}$ core in about 0.8 Myr, which then undergoes runaway gas accretion to reach Jupiter mass. Such processes potentially allow close-in giant planets to be formed even in the relatively hostile conditions that we consider here.

However, photoevaporation by the central star can cause large amounts of mass loss in the inner disc, potentially affecting giant planet formation (Alexander et al., 2014). Grain size also has a significant effect on disc dispersal rates. Mass loss occurs much more quickly when grain growth has occurred because the FUV radiation can penetrate deeper into the disc (Facchini et al., 2016).

Discs that can survive in moderately dense environments have small radii (10 – 50 AU). This is because of the disc radius dependency in the external photoevaporation prescriptions. Internal UV radiation can cause significant mass loss and erosion of the disc within short time scales (1 Myr, Gorti & Hollenbach, 2008). The timescale for internal disc dispersal is very short ( $10^{5}$ yr), with a UV switch being triggered due to the slowing of accretion onto the inner 10 AU of the disc (Clarke et al., 2001), also calling into question the survivability of small discs.

Our disc radii are fixed, but in reality disc radii change with time, often in an inside-out fashion where the initial radius is small (and not as susceptible to photoevaporation) compared to later in the disc’s life. We include several different disc radii to help visualise what happens for different disc initial conditions, but we cannot model the full viscous evolution in our post-processing analysis.

In our simulations, we have elected to keep the stellar IMF constant across different realisations of the spatial and kinematic initial conditions of our star-forming regions. The reasons for this are two-fold. First, we wish to isolate the possible effects of stochastic dynamical evolution (Allison et al., 2010; Parker & Goodwin, 2012; Parker et al., 2014), which could lead to different photoevaporation rates even if the ionising flux from massive stars were kept constant. The uncertainties shown by the shaded regions in Figs. 1, 2, and 5 show this stochasticity for the same initial conditions. Second, the photoevaporation prescriptions we adopt (following Scally & Clarke, 2001) are actually quite insensitive to the mass of the most massive stars (but rather depend on whether the massive stars are present or not).

However, Fatuzzo & Adams (2008) show that the FUV and EUV fluxes can vary if the stellar IMF is extremely top-heavy and contains more massive stars than expected on average. In a forthcoming paper we will calculate the EUV and FUV fluxes in our substructured star-forming regions and use recent the FRIED models of disc photoevaporation from Haworth et al. (2018b) to determine mass-loss based on these models, and whether it depends strongly on stochastic sampling of the stellar IMF.

Similarly, in some of our simulations the massive stars are ejected early on, which is a common occurrence in simulations of dense star-forming regions (Allison et al., 2010; Oh & Kroupa, 2016, Schoettler et al., in prep.). We will also quantify the effects of these ejections on the fraction of surviving protoplanetary discs in an upcoming paper.

The majority of discs observed with ALMA have been located in low-mass, low-density star-forming regions. Current observations suggest that the majority of stars form in moderately dense ( $\sim$ 100 M*⊙* pc*-3*) environments (Parker, 2014). However, the majority of protoplanetary discs in clusters with these densities do not survive for long enough periods of time to form planets, as planet formation is thought to take place over a few million years (Pollack et al., 1996). The fact that the majority of stars have planetary systems around them poses important questions as a result of the discrepancies that seemingly arise. This may indicate that the majority of stars form in low mass clusters where there are few to no high mass stars.

We adopt initial disc masses that are 10 per cent the mass of the host star, which is likely to be a large overestimate. When looking at more realistic values (1 per cent), discs are destroyed on even shorter timescales. However, it should be noted that accretion and internal photoevaporation will have much larger effects on disc mass evolution for these lower mass discs.

5 Conclusions

We have calculated the mass loss due to external photoevaporation of protoplanetary discs in $N$ -body simulations of the evolution of star-forming regions. We ran a suite of simulations where we vary the initial spatial structure, bulk motion and initial density of the regions. We compared our simulations that more closely represent observed star-forming regions (subvirial, substructured) with those of primordially mass segregated, spherical clusters, similar to those used in previous studies of external photoevaporation.

The parameter that most affects rates of disc dispersal is the initial density of the star-forming region. The majority of protoplanetary discs within simulated regions that mimic the conditions in nearby star-forming regions are dispersed due to external photoevaporation within very short time scales. In moderately dense ( $\sim$ 100 M*⊙* pc*-3*) star-forming regions which have moderate levels of substructure ( $D=2.0$ ) and are collapsing ( $\alpha_{\rm vir}=0.3$ ), we find the time taken for half of 100 AU discs to dissipate is significantly shorter (3 times less) than suggested in observational studies (Haisch et al., 2001). Lower density clusters ( $\sim$ 10 M*⊙* pc*-3*) allow discs to survive long enough to match observations of disc lifetimes, although the half-life of 100 AU discs is still less than that found by Haisch et al. (2001).

The initial spatial distribution of the star-forming region also affects the rate of protoplanetary disc dispersal due to external photoevaporation. The degree to which initial substructure affects disc dispersal rates depends on the initial density. In moderately dense ( $\sim$ 100 M*⊙* pc*-3*) regions the effects are washed out, but in lower-density regions ( $\sim$ 10 M*⊙* pc*-3*) we find that the more fractal and clumpy a star-forming region is, the higher the rate of disc dispersal. This is due to violent relaxation and the rapid increase in density (sometimes up to an order of a magnitude) of the star forming region within a short amount of time. As most star forming regions appear to have a high degree of substructure, it is important for future studies of disc dispersal to take the initial conditions into consideration due to external photoevaporation in dense environments.

The virial ratio of the star forming region affects the rate of disc dispersal in a similar way to substructure. Regions that have a low initial density and are collapsing photoevaporate more discs on average than clusters which are expanding. The effects of varying the initial net bulk motion in moderately density clusters is negligible.

The majority of observed stars in the Galactic field host planetary systems, implying their protoplanetary discs survived long enough for formation to take place. There are three possible scenarios to resolve this apparent tension between observations and our simulations:

i) The majority of planets may not form in moderately dense star-forming regions ( $\sim$ 100 M*⊙* pc*-3*); rather, they would form in low density regions with no photoionising massive stars present. Many of the protoplanetary discs have been in observed in these low-density ambient environments (Ansdell et al., 2018), but significant numbers of protoplanetary discs (or at least their remnants) have been observed in dense, hostile regions like the ONC (McCaughrean & O’dell, 1996).

ii) If some planets do form in dense, clustered environments containing massive stars (such as the ONC), then this suggests that giant planet formation must happen on very short time scales (less than 1 – 2 Myr), or be confined to discs with radii significantly smaller than the orbit of Neptune in our Solar System. Johansen & Lambrechts (2017) show that giant planet formation can occur on these timescales once large enough planetesimals have formed. However, internal photoevaporation processes can deplete the inner disc and set limits on the formation time of giant planets (Alexander et al., 2014).

iii) The current calculations of mass-loss in discs due to external photoevaporation are severely overestimating the detrimental effects of EUV and FUV radiation. However, recent research (Facchini et al., 2016; Haworth et al., 2018a) suggests that photoevaporative mass-loss rates caused by FUV radiation may be underestimated, and our calculations also underestimated the effects as we adopt conservatively high initial disc masses.

Acknowledgements

We thank A. R. Williams for useful discussions and S. Habergham-Mawson for comments. RBN is partially supported by an STFC studentship. RJP acknowledges support from the Royal Society in the form of a Dorothy Hodgkin Fellowship. RC and MBD are supported by grant 2014.0017, “IMPACT", from the Knut and Alice Wallenberg Foundation.

Bibliography92

The reference list from the paper itself. Each links out to its DOI / PubMed record.

1Adams et al. (2004) Adams F. C., Hollenbach D., Laughlin G., Gorti U., 2004, Ap J , 611, 360 · doi ↗
2Adams et al. (2006) Adams F. C., Proszkow E. M., Fatuzzo M., Myers P. C., 2006, Ap J , 641, 504 · doi ↗
3Alexander et al. (2014) Alexander R., Pascucci I., Andrews S., Armitage P., Cieza L., 2014, Protostars and Planets VI , pp 475–496 · doi ↗
4Allison et al. (2010) Allison R. J., Goodwin S. P., Parker R. J., Portegies Zwart S. F., de Grijs R., 2010, MNRAS , 407, 1098 · doi ↗
5André et al. (2010) André P., et al., 2010, A&A , 518, L 102 · doi ↗
6Andrews et al. (2013) Andrews S. M., Rosenfeld K. A., Kraus A. L., Wilner D. J., 2013, Ap J , 771, 129 · doi ↗
7Ansdell et al. (2016) Ansdell M., et al., 2016, Ap J , 828, 46 · doi ↗
8Ansdell et al. (2017) Ansdell M., Williams J. P., Manara C. F., Miotello A., Facchini S., van der Marel N., Testi L., van Dishoeck E. F., 2017, AJ , 153, 240 · doi ↗

TL;DR

Contribution

Findings

Abstract

Peer Reviews

Videos

Rapid destruction of protoplanetary discs due to external photoevaporation in star-forming regions

Abstract

keywords:

1 Introduction

2 Method

2.1 Creating low-mass star-forming regions

2.2 NNN-body simulations

2.3 Protoplanetary discs and external photoevaporation

3 Results

3.1 Substructure in star-forming regions

3.2 Virial ratio

3.3 Disc radii

3.4 Disc Masses

3.5 Mass of star-forming regions

4 Discussion

4.1 Changing the initial conditions of star-forming regions

4.2 Caveats

5 Conclusions

Acknowledgements

2.2 $N$ -body simulations