Age spreads and the temperature dependence of age estimates in Upper Sco

TL;DR

This study investigates how the star formation history of Upper Sco causes apparent age differences across spectral types, suggesting that observed age spreads are consistent with prolonged star formation and pre-main sequence evolution physics.

Contribution

It demonstrates that star formation history and pre-main sequence evolution physics can explain the spectral type-dependent age estimates in Upper Sco.

Findings

Star formation history influences age estimates across spectral types.

Constant star formation over several Myr explains age differences.

Prolonged star formation aligns with observed luminosity spreads.

Abstract

Past estimates for the age of the Upper Sco Association are typically 11-13 Myr for intermediate-mass stars and 4-5 Myr for low-mass stars. In this study, we simulate populations of young stars to investigate whether this apparent dependence of estimated age on spectral type may be explained by the star formation history of the association. Solar and intermediate mass stars begin their pre-main sequence evolution on the Hayashi track, with fully convective interiors and cool photospheres. Intermediate mass stars quickly heat up and transition onto the radiative Henyey track. As a consequence, for clusters in which star formation occurs on a similar timescale as the transition from a convective to a radiative interior, discrepancies in ages will arise when ages are calculated as a function of temperature instead of mass. Simple simulations of a cluster with constant star formation over several Myr may explain about half of the difference in inferred ages versus photospheric temperature; speculative constructions that consist of a constant star formation followed by a large supernova-driven burst could fully explain the differences, including those between F and G stars where evolutionary tracks may be more accurate. The age spreads of low-mass stars predicted from these prescriptions for star formation are consistent with the observed luminosity spread of Upper Sco. The conclusion that a lengthy star formation history will yield a temperature dependence in ages is expected from the basic physics of pre-main sequence evolution and is qualitatively robust to the large uncertainties in pre-main sequence evolutionary models.

Tables

Table 1. Table 1 Rejected members included in Pecaut et al. ( 2012 )

2MASS	SpT	${\pi }_{GAIA}\left[mas\right]$
J15545986-2347181	G3V	13.08 $\pm$ 0.68
J15582930-2124039	F0IV	3.45 $\pm$ 0.26
J16123604-2723031	K4	1.03 $\pm$ 0.88

Table 2. Table 2 Median age versus spectral type for Upper Sco

		magnetic
SpT	Dotter	Feiden	Siess
F2	13.3 $\pm$ 0.7	15.8 $\pm$ 0.9	13.9 $\pm$ 1.5
F5	13.6 $\pm$ 1.5	16.3 $\pm$ 1.5	14.9 $\pm$ 1.4
F8	10.8 $\pm$ 2.0	13.3 $\pm$ 3.3	14.0 $\pm$ 1.7
G2	8.6 $\pm$ 2.3	11.4 $\pm$ 4.7	13.1 $\pm$ 2.4
G5	8.7 $\pm$ 3.0	16.0 $\pm$ 6.6	14.8 $\pm$ 4.0
G8	7.9 $\pm$ 1.5	16.1 $\pm$ 3.3	13.5 $\pm$ 2.1
K2	5.5 $\pm$ 1.2	11.4 $\pm$ 2.7	9.5 $\pm$ 1.6
K5	2.7 $\pm$ 0.5	7.3 $\pm$ 1.1	4.1 $\pm$ 0.7
K8	2.2 $\pm$ 0.3	6.4 $\pm$ 0.7	3.4 $\pm$ 0.4
M0	2.2 $\pm$ 0.4	5.9 $\pm$ 0.8	3.0 $\pm$ 0.4
M2	2.9 $\pm$ 0.2	6.2 $\pm$ 0.4	3.7 $\pm$ 0.2

Table 3. Table 3 Parameters

	Smooth-Burst	Smooth-Quench
$t_{\mathrm{start}}$	12, 14, 16, 18, 20	12, 14, 16, 18, 20
$t_{\mathrm{quench}}$ (Myr)	2, 4, 6, 8, 10	2, 4, 6, 8, 10
$A$	20, 40, 60, 80	–

Table 4. Table 4 Ages of Nearby Associations with Feiden magnetic models

Region	$\log L_{4200}/L_{\odot }$	Age	$\log L_{3400}/L_{\odot }$	Age
Upper Sco	$-0.25$	$10.5$	$-0.86$	$6.2$
$ϵ$ Cha	$-0.28$	$11.6$	$-0.89$	$6.7$
$\eta$ Cha	$-0.32$	$13.5$	$-0.93$	$7.6$
TWA	$-0.40$	18	$-1.01$	$9.7$
BPMG	$-0.65$	40	$-1.25$	$21$
Tuc-Hor	$-0.72$	50	$-1.33$	$25$

Equations

$$\log \frac{L}{L_{\odot }}=-0.4(m_J-5\log \frac{d}{10pc}+B{C}_J-M_{bol,\odot })$$

$$\sigma (\log \frac{L}{L_{\odot }})=0.4\sqrt[]{\sigma (m_J{)}^2+[\frac{5\sigma (d)}{dln10}{]}^2+\sigma (B{C}_J{)}^2}$$

$$p(q)=(\frac{1+\beta }{1-q{}_{lo}{}^{1+\beta }})q^{\beta }$$

$$SFR(t)=\{\begin{array}{ll}1& t\le t_{\mathrm{start}}-t_{\mathrm{quench}}\\ e^{-[t-(t_{\mathrm{start}}-t_{\mathrm{quench}}){]}^2}& t_{\mathrm{start}}-t_{\mathrm{quench}}<t\end{array}$$

$$SFR(t)=\{\begin{array}{ll}1(t\le t_{\mathrm{start}}-t_{\mathrm{quench}})& \\ A[t-(t_{\mathrm{start}}-t_{\mathrm{quench}})]e^{-[t-(t_{\mathrm{start}}-t_{\mathrm{quench}}){]}^2}& \\ (t_{\mathrm{start}}-t_{\mathrm{quench}}<t)& \end{array}$$

$$\frac{N_1}{N_0}=\frac{0.886}{t_{\mathrm{start}}-t_{\mathrm{quench}}}$$

$$\frac{N_1}{N_0}=\frac{0.5A}{t_{\mathrm{start}}-t_{\mathrm{quench}}}$$

Full text

1. INTRODUCTION

Nearby clusters of young stars are our laboratories to study the formation and early evolution of stars and planets. Masses and ages of young stars in these clusters are usually measured by isochrone fitting, i.e. comparing their location on Hertzsprung-Russell (HR) diagram and evolution models. In any given cluster, an accurate stellar locus should provide a robust test for star formation theories. However, star formation histories and age spreads are challenging to measure empirically because of observational errors, including unresolved binarity, stellar spots, and membership. These observational errors are exacerbated by uncertainties in the input physics of pre-main sequence evolution, including birthlines, accretion histories, and magnetic activity (see reviews by Soderblom et al., 2014; Hartmann et al., 2016), and by inconsistencies between models of stellar evolution.

One of the nearest young clusters, the Upper Sco Association (Upper Sco), has emerged as a test-bed for pre-main sequence evolution (de Geus et al., 1989; de Zeeuw et al., 1999). The temperatures and luminosities of the low-mass111In this paper, low-mass stars refers to stars from $0.2$ to 1.0 M⊙, with spectral types K to mid-M at the age of Upper Sco. Intermediate-mass stars refers to AFG stars, with masses $1.0$ to 5 M⊙, while B-type and MS-turnoff stars in Upper Sco are 5-15 M⊙ stellar population in Upper Sco indicate an age of $\sim 5$ Myr (e.g. Preibisch et al., 2002; Slesnick et al., 2008; Rizzuto et al., 2015). However, analyses of intermediate and of high mass stars in Upper Sco yield an older age of 10–12 Myr (Pecaut et al., 2012). Subsequent studies of low-mass populations confirm these age discrepancies between intermediate- and high-mass stars and low-mass stars (Herczeg and Hillenbrand, 2015; Rizzuto et al., 2016; Pecaut and Mamajek, 2016). Other regions show similar differences, with intermediate-mass stars and fits of empirical isochrones to empirical isochrones typically appearing two times older than K and M-type stars (e.g. Hillenbrand, 1997; Hillenbrand et al., 2008; Naylor, 2009; Malo et al., 2014; Bell et al., 2015).

Some of this discrepancy may be resolved with improved calculations of pre-main sequence evolution. In the most recent models, prescriptions for convection have been updated from the historical assumption of a single, constant mixing length. Feiden (2016) and MacDonald and Mullan (2017) introduce strong magnetic fields in the stellar interiors to slow the efficiency of convection and thereby also slow the stellar contraction rates, yielding inflated radii and therefore older ages for all low-mass pre-main sequence stars (see, e.g., Somers and Stassun, 2017). These models reduce but do not eliminate temperature-dependent discrepancies in age estimates of Upper Sco, and also introduce challenges when applied to other regions (see Appendix A). Similar reductions in convective efficiency have led to older ages when photospheric spots are considered (Somers and Pinsonneault, 2015) or when the treatment of convection is based on results from 3D MHD simulations (Baraffe et al., 2015). The uncertainty in ages related to accretion history would preferentially make low-mass stars less luminous (e.g. Hartmann et al., 1997; Baraffe et al., 2016), thereby exacerbating the age discrepancy. Uncertainties in stellar atmospheres caused by errors in opacities (e.g. Bell et al., 2013) and boundary conditions (e.g. Chen et al., 2014) also propagate into errors in the evolutionary tracks and in the placement of stars in HR diagrams.

We investigate an alternative explanation, that the conflict in ages between cool and hot stars may be explained with basic stellar evolution and an age spread within the cluster. During their pre-main sequence evolution, low-mass young stars are fully convective and evolve along the Hayashi track, with a roughly constant temperature and a luminosity that decreases as the star contracts. On the other hand, intermediate-mass stars begin their lives fully convective, but the stellar interior heats up, thereby reducing the opacity and allowing for efficient radiative energy transport (Hayashi, 1961). Their evolution then proceeds along the roughly-horizontal Henyey track, during which the stellar photosphere becomes hotter while maintaining a similar luminosity (Henyey et al., 1955; Hayashi and Nakano, 1965). If the star formation timescale in a cluster is similar to the timescale for intermediate mass stars to become radiative (1-10 Myr), then measuring ages as a function of photospheric temperature (e.g. Pecaut et al., 2012; Herczeg and Hillenbrand, 2015) will mean that hotter (spectral types AFG) stars will tend to be the older stars evolving on the Henyey track, while younger stars with the same mass will still be on the Hayashi track and will have cool effective temperatures (K and early M spectral types). For young ($\lesssim 30$ Myr) clusters, the AFG stars will preferentially be among the oldest stars in the cluster, while M stars will sample the full star formation history.

These effects are likely broadly applicable to young regions. Some age spread is required by the presence of embedded objects, objects with disks, and diskless objects within the same cluster (e.g. Evans et al., 2009). A halo of mostly-diskless young stars spread about between dense molecular clouds in the Taurus star forming region suggests either a large age spread or past bursts of star formation (Kraus et al., 2017). While concentrated regions may form stars quickly, the measured age spread should be expected to be longer when the region considered is larger (see, e.g., description of Chamaeleon by Sacco et al. (2017) and age spreads of tens to hundreds of Myr in massive clusters inferred by, e.g., Goudfrooij et al. (2011)).

In this paper, we model young clusters with a range of star formation histories and then apply the results to Upper Sco, which has evidence for a spatial age gradient (Preibisch and Zinnecker, 2007; Pecaut and Mamajek, 2016) and an age spread from HR diagrams, including B-type MS-turnoff stars with ages that range from 5–14 Myr (e.g. Pecaut et al., 2012; Rizzuto et al., 2015). In §2, we first establish the census of probable Upper Sco members from the literature to create a stellar locus. In §3, we use stellar evolution models to simulate HR diagrams to establish how age depends on temperature for a cluster with an age spread. In §4, we then use apply these simulations to Upper Sco to evaluate whether star formation history may lead to a temperature dependence in age estimates. We discuss these results in §5.

2. THE POPULATION OF UPPER SCO

In this section, we compile a list of high-probability members of Upper Sco based on past selection criteria, proper motion, photometry, and distance, and then characterize their luminosities and photospheric temperatures. Our selection criteria excludes members with anomalous proper motions or distances to minimize interlopers. We also exclude stars with disks. These selection criteria should lead to an accurate measurement of the locus of Upper Sco stars in the HR diagram.

The methods described for selecting members and assigning temperatures and luminosities are similar to, and in some cases follow exactly, the methods implemented for Upper Sco by Pecaut and Mamajek (2016). Our sample, adopted from several literature studies, extends to M stars, while Pecaut and Mamajek (2016) concentrates on K stars. We also adopt a single distance for stars without parallax measurements, while Pecaut and Mamajek (2016) use distance estimates from kinematic analyses. Relative to Pecaut and Mamajek (2016), our sample of stars cooler than 4000 K is more complete and less luminous.

2.1. Sample Description

The initial membership list of Upper Scorpius and spectral type measurements is compiled from Rizzuto et al. (2015), Pecaut and Mamajek (2016) and the compilation of Luhman and Mamajek (2012), which included previous surveys of Preibisch and Zinnecker (1999), Ardila et al. (2000), Preibisch et al. (2002), and Slesnick et al. (2006). Photometric data for low mass members of Upper Sco are obtained from the 2MASS $JHKs$ survey (Skrutskie et al., 2006), WISE $W1$, $W2$, $W3$, $W4$ survey (Cutri and et al., 2013), and APASS (Henden et al., 2012). Proper motions in this work are compiled from GAIA DR1 (Gaia Collaboration et al., 2016), UCAC4 (Zacharias et al., 2013), PPMXL (Roeser et al., 2010), and SPM4 surveys (Girard et al., 2011). We focus on members with spectral types that range from F0–M5.

The Gaia Collaboration et al. (2016) DR1 TGAS astrometric catalog includes parallaxes of 126 F, G, and K-type objects that were previously identified as members of Upper Sco and are located between $234$ to $253$ degrees in right ascension and $-34$ to $-14$ degrees in declination. Figure  1 shows that the distance distribution of these members has a mean of 144 pc, with a standard deviation of 17.6 pc, a mean error of 8.6 pc, and a systematic error of $\sim 6$ pc (Arenou et al., 2017). A quadratic subtraction of these two errors indicates that the distance spread is $\sim 15$ pc, consistent with previous estimates (e.g. Slesnick et al., 2008). Three stars with distances smaller than 100 pc (10 mas in parallax) or larger than 200 pc (5 mas in parallax) are rejected as members (Table 1). We adopt distances from the Gaia DR1 parallaxes when available, or the average distance of $145\pm 15$ pc for stars without parallax measurements.

Some members of the older Upper Centaurus-Lupus subgroup, a related association also within the Sco-Cen OB Association, likely contaminate this membership list. The related moving groups have similar space motions and overlap on the sky (de Zeeuw et al., 1999; Pecaut et al., 2012). Some contamination is also expected from the younger $\rho$ Ophiucus star forming region.

Figure 2 shows the proper motions of the candidates in both directions, where $\mu_{\alpha*}\equiv\mu_{\alpha}cos\delta$. The average proper motion of all Upper Scorpius members considered in this paper is $\mu_{\alpha*}=-10.45$ mas yr*-1* and $\mu_{\delta}=-23.21$ mas yr*-1*, with a total standard deviation of 6.41 mas yr*-1*. To reject objects that are moving inconsistently with the common streaming motion, all stars with proper motions that deviate 3$\sigma$ from the mean are rejected as members. The remaining 767 members are kinematic selected members and will be used in the following analysis.

2.2. Stellar properties

In this section, we describe how we combine the literature spectral type information with existing photometry to measure stellar properties and identify stars with disks, which are subsequently excluded from analysis.

2.2.1 Extinction

We adopt the conversion between spectral type, temperature, and color scales from Pecaut and Mamajek (2013). Uncertainties in temperatures are assigned to be 200 K at temperatures $>4100$ K, 125 K at 3750 – 4100 K, and 50 K below 3750 K. Following Pecaut and Mamajek (2016), for stars with $V$-band photometry, extinction is estimated from the color excesses $E(V-J)$, $E(V-H)$, $E(V-K)$ relative to the expected photospheric color. for other sources the extinction is estimated from the $E(J-K)$ color excess. The color excesses are then converted to extinction using $A_{J}/A_{V}$ = 0.27, $A_{H}/A_{V}$ = 0.17, $A_{K}/A_{V}$ = 0.11 (Fiorucci and Munari, 2003). The adopted V band extinction is the mean value of $A_{V}$ estimated from these colors.

The top panel of Figure 3 shows the $J-K$ color for low mass members of Upper Sco. The black solid line is the relation between effective temperature and intrinsic color, which is derived from (Pecaut and Mamajek, 2013) and the dashed black line is the color reddened by $A_{V}$ = 2 mag. The median extinction is $A_{V}\sim 0.60$ mag. All objects reddened by $A_{V}>2$ mag are excluded from further analysis. All extinctions $A_{V}<0$ mag are set to $A_{V}=0$ mag. Use of negative extinctions would have a negligible affect on the empirical isochrone of Upper Sco. The median $A_{V}$ is constant versus spectral type to within $0.25$ mag, which leads to systematic errors of $<7$% in luminosity.

2.2.2 Disk identification

Disks produce excess emission in the infrared from dust and in the optical from accretion onto the star. Since correcting photometry for these processes is challenging, we instead identify and exclude all disk sources in our sample. In this work, disks are identified based on a mid-IR color excess, with $E(K-W3)>0.5$ (see also, e.g., Luhman and Mamajek, 2012), as shown in Figure  4. These results are consistent with the disk identification in Upper Sco by Luhman and Mamajek (2012) and Rizzuto et al. (2015). Uncertainty in the $K-W3$ color leads to a higher rejection rate of diskless stars at $<3000$ K.

2.2.3 Luminosity

Luminosities are calculated from 2MASS $J$-band photometry and the bolometric correction of Pecaut and Mamajek (2013) obtained for the appropriate spectral type, as follows:

[TABLE]

where the bolometric luminosity of the Sun, $M_{bol,\odot}$, is 4.755 mag (Mamajek 2012). The error of luminosity, $\sigma L$, is estimated as

[TABLE]

Estimating the first two terms is straightforward, as $\sigma$(mJ) is the observation error of the $J$ band magnitude and $\sigma$(d) is distance spread (15 pc, see §2). To estimate the errors of bolometric corrections, we perform 1000 Monte Carlo simulations for each star. In each simulation, we assign the same error in temperature as discussed in previous sections, and derived a set of bolometric corrections for each star. The standard deviation is the adopted error in the third term in Equation 2. The error in bolometric magnitude then propagates to the error in luminosity.

2.3. HR diagram of Upper Sco

Figure 5 shows the temperatures and luminosities of the diskless members of Upper Sco in an HR diagram. The empirical isochrone is calculated by dividing the sample into $\log T_{\rm eff}$ ranges from 3.5 (3000 K) to 3.9 (7000 K) into 20 bins of $\sim 0.02$ dex. Bins with fewer than three stars are ignored. All stars that are located either above the 0.1 Myr isochrone or below the main sequence are rejected as non-members. Errors of temperatures and luminosities are assigned as described above. For each star, the adopted value and uncertainty for age is evaluated from ages obtained in $10^{4}$ Monte Carlo simulations in which these 1 $\sigma$ uncertainties are randomly applied to the effective temperature and luminosity. The empirical isochrone is then defined as the median of all luminosities within the bin. As has been established in other studies (e.g. Pecaut et al., 2012; Rizzuto et al., 2015; Herczeg and Hillenbrand, 2015; Pecaut and Mamajek, 2016), the average age of hotter stars is older than that of younger stars.

Figure 6 compares the empirical isochrone with stellar evolution tracks222An overview of many pre-main sequence evolutionary models can be found in Stassun et al. (2014). The Siess and Dotter tracks use standard mixing length theory, while Feiden introduces magnetic fields into the stellar interior. from Siess et al. (2000, hereafter Siess), Dotter et al. (2008, hereafter Dotter) and the magnetic tracks of Feiden (2016, hereafter Feiden). The isochrone of Upper Scorpius is not parallel to any of the theoretical isochrones, as should be expected for co-eval members of a stellar cluster (if pre-main sequence evolutionary models are accurate). The age measurement strongly depends on temperature in all models for pre-main sequence evolution (see median age of each spectral bin in Table 2). The empirical isochrone calculated here for Upper Sco is very similar to that of Pecaut and Mamajek (2016) at temperatures hotter than 4000 K and is several tenths of dex fainter at cooler temperatures.

This sample suffers from several biases. First, the low-mass stellar population is incomplete, with significant spatial biases. The disk fraction versus luminosity within a small range of subclasses indicates that the sample may be preferentially missing faint targets at spectral types later than M3. Early-type members were identified by de Zeeuw et al. (1999) using Hipparcos cover a large area on the sky, while many searches for low-mass stars focus on smaller regions. The sample also suffers from contamination of stars that belong to nearby, related regions. Finally, our exclusion of disk sources may preferentially exclude some of the youngest sources in the region. None of these biases affects our comparisons of inferred ages of low-mass and intermediate-mass stars.

3. POPULATION SYNTHESIS OF YOUNG STAR CLUSTERS

If a population of young stars forms over some duration of time, then stars with hotter photospheres will be preferentially older than stars with cooler photospheres. In the Dotter evolutionary models, at 5 Myr the mass range of F type (here defined as $6000-7500$ K) stars is 1.95 $\sim$ 2.14 M*⊙, while at 10 Myr same temperature range corresponds to masses of 1.50 $\sim$ 1.79 M⊙*. By integrating over the Chabrier initial mass function (Chabrier, 2003), a 10 Myr old population will have five times more F-type stars than a 5 Myr old population. On the other hand, the 10 Myr old population will have 1.5 times fewer K-type stars (here defined as 4000-5000 K) than a 5 Myr old population. In this context, a star cluster that forms over 10 Myr will have a measured age that depends on the spectral type or temperature range being measured. The median age estimated for F type stars will be biased to older stars, while the age estimated for K type stars will be biased to younger stars.

In this section, we quantify the size of this effect and characterize the age spread and star formation history that would be needed to reproduce the ages of the population in Upper Sco. We first simulate HR diagrams in a cluster to derive model isochrones to compare with the observed isochrone of Upper Sco under some simple assumptions. We then introduce some speculative changes that are able to entirely reproduce the empirical isochrones and age spread of Upper Sco.

3.1. Simulation setup

Stellar populations are simulated using the pre-main sequence evolutionary models of Dotter et al. (2008). A sample of $10^{4}$ stars with masses $0.1$ to 2.0 M⊙ are generated with masses randomly selected from the Chabrier initial mass function, which describes Galactic disk stars (Chabrier, 2003). We then introduce to each star (or binary system) an age randomly generated from the star formation history.

The main sources of uncertainty are errors in spectral types, photometry, and distance, and in unresolved multiplicity. In the previous section, we assign to each member in Upper Sco errors in both effective temperature and luminosity. The mean errors in $\log T$/K and $\log L/L_{\odot}$ are assumed to be in a Gaussian distribution with a $1\sigma$ scatter of 0.01 and 0.05 dex, respectively (see also a similar implementation in Jose et al., 2016a).

Unresolved binaries will lead to incorrect effective temperatures and luminosities, resulting in inferred ages that are too young (e.g. Preibisch and Zinnecker, 1999; Slesnick et al., 2008; Reggiani et al., 2011). To estimate the effect of binarity (higher order systems are ignored), the effective temperature of multiple systems are assumed to be that of the primary star, although the combined temperature will be lower than the temperature of the primary alone (e.g., see discussion in Herczeg and Hillenbrand, 2014). The luminosity of binary systems are estimated by summing the primary and companion luminosity, which will lead to small differences with our observed luminosity measurements from the $J$-band photometry.

We follow the quantification of the multiplicity fraction and mass ratio of Lu et al. (2013), which was estimated from the empirical results of Preibisch and Zinnecker (1999); Kouwenhoven et al. (2007); Kobulnicky and Fryer (2007); Lafrenière et al. (2008). The multiplicity fraction adopted here follows power law with the form $MF(m)=0.44m^{0.51}$ given by Lu et al. (2013) for stars with masses less than $5$ M⊙, and where m is the primary mass. All systems in our work have either one or two stars. The mass ratio $q\equiv$ $\frac{M_{2}}{M_{1}}$ between the companion star and primary star also follow a power law,

[TABLE]

where $\beta$ ranges from -0.6 to 1.4. In this work we adopt a flat companion mass ratio distribution ($\beta$ = 0), following Reggiani and Meyer (2013).

The effect of pairing function is carefully discussed in previous work (see, e.g. Kouwenhoven et al., 2009). Here we adopt the primary-constrained random pairing method, where the primary mass is derived from Chabrier IMF and the companion mass is determined by $p(q)$. All companion stars with masses lower than 0.1 M⊙ will not contribute to significant changes in the HR diagram and are excluded from further analysis.

3.2. The Star Formation History of a Synthetic Stellar Population

In this subsection, we describe how the star formation history is implemented into our population synthesis. We first calculate stellar populations for the Smooth-Quench model, which has a constant star formation over some time period. We then calculate populations in the Smooth-Burst model, which adds a burst of star formation at the end of the Smooth-Quench model. The Smooth-Quench model is the simplest possible characterization of a star formation history. While ad hoc additions to the Smooth-Quench model could become hopelessly complicated, adding a burst at the end of the Smooth-Quench model maximizes possible age differences between stars of different temperatures, thereby probing the limits of age differences versus spectral type. The Smooth-Burst model also may represent a reasonable scenario of a past supernova triggering a burst in star formation (see §4). The accelerating star formation model of Palla and Stahler (2000) is not investigated in this work but would yield results that are intermediate between the Smooth-Quench and Smooth-Burst models.

These models are summarized in Figure 7. The Smooth-Quench model is described by two parameters: the time since the association first started forming stars, $t_{\rm start}$, and the time since the star formation was quenched, $t_{\rm quench}$. The decrease in star formation after quenching is described by e${}^{-t^{2}}$. The star formation rate is then analytically described by

[TABLE]

In addition to $t_{\rm start}$ and $t_{\rm quench}$, the Smooth-Burst model adds a burst of strength $A$ relative to the quiescent star formation rate at time $t_{quench}$. The star formation rate is then described as

[TABLE]

For the Smooth-Quench model, by integrating over the star formation history, the ratio of the number of stars form at the quench/burst process N1 and the number of stars form during the constant star formation rate $N_{0}$ is

[TABLE]

and for the Smooth-Burst model

[TABLE]

where tstart and tquench are both in the unit of Myr. Thus, if star formation occurs over 10 Myr ($t_{\rm start}-t_{\rm quench}=10$), then an amplitude of the burst $A=20$ would yield an equal number of stars formed during the burst and during the smooth period of star formation.

In both of these constructions, the IMF, binary formation, and stellar evolution are assumed to be constant and independent of star formation rate. In our models, these parameters are calculated by generating $10^{4}$ stars over $10^{4}$ possible ages during the star formation history and randomly assigning these stars a mass from the Chabrier IMF and binarity as described above.

3.3. Simulation Results

Our simulations are run with several sets of parameters, listed in Table  3. For each simulated population, the median ages of stars are estimated in different temperature ranges, with an empirical isochrone calculated following the same steps in §2. The empirical isochrones of the simulated HR diagrams from these models are then simplified into two measureable parameters: the age difference between F5 and M0 stars ($\sim 6500$ K and $3900$ K, respectively) for the entire cluster and the luminosity spread of M0 stars. Binaries have only a minor effect in these simulations.

3.3.1 Age Differences versus Photospheric Temperature

When star formation occurs over some duration, hotter stars will have older ages than cooler stars. Figure 8 shows the age versus temperature for a cluster that started forming stars 16 Myr ago and ceased forming stars 4 and 8 Myr ago. For this plot, the start time of star formation is fixed at t${}_{\rm start}=16$ Myr, with the quench time $t_{\rm quench}$ fixed at 4 Myr (top panel) and $8$ Myr (bottom panel). Figure 9 shows these same differences, but parameterized by the age difference between F5 and M0 stars to show results from additional models. The time since the last star formation, $t_{\rm quench}$, is fixed in each of the panels.

In the most basic Smooth-Quench model, the age differences between F5 and M0 stars are largest when star formation occurred recently and decrease with time. The recent star formation maximizes the number of intermediate-mass stars that are still fully convective. If no star formation has occurred in 8 Myr, then the differences between ages of F5 and M0 stars would be $\sim 2-3$ Myr. The duration of star formation ($t_{\rm start}-t_{\rm quench}$) increases these differences. In the limit where $t_{\rm start}-t_{\rm quench}=0$, then there would be no difference in age between F5 and M0 stars since all star formation would occur instantaneously.

Moving to the burst model, a modest burst ($A=20$) increases the age difference by $\sim 4$ Myr. Much larger bursts do not increase this effect significantly. Again, as $t_{\rm quench}$ increases, the intermediate mass stars born during the burst will have enough time to evolve into radiative phases, thereby reducing any apparent age differences between F5 and M0 stars. Larger $A$-values will decrease the median ages hot and cold stars. When $t_{\rm quench}$ is large (long duration since the burst), larger $A$ values will result in smaller age differences.

The age differences become larger as the duration of star formation, $t_{\rm start}-t_{\rm quench}$, increases. For fixed $t_{\rm start}$ and $A$, larger $t_{\rm quench}$ will result in smaller age difference. For a single burst, we should not expect a large age difference between stars of different spectral type. A lengthy star formation process ($t_{\rm start}-t_{\rm quench}$) is the most significant factor contributing to the spectral type-dependent age measurement.

3.3.2 Luminosity Spread

Any age spread during star formation will also lead to a spread in observed luminosities. In this subsection, we describe how different parameters in Smooth-Quench and Smooth-Burst models affect the luminosity spread of M0 type stars. Earlier spectral types are affected by stellar evolution and are a poor probe of age spread with additional analysis, while M0 stars are fully convective and contract along the Hayashi track for the duration of the first 20 Myr of evolution simulated here.

Figure 10 show the effects of $A$ and $t_{\rm start}$ on the luminosity spread. The luminosity spreads are always largest when star formation is recent because the luminosity of an M0 star decreases with $\log$ age, while we are measuring age spread in linear space. In the Smooth-Quench model, for example the age spread (standard deviation of all ages) will be $\sim 0.29\times(t_{\rm start}-t_{\rm quench})$. While this age spread is constant, the luminosity spread decreases with age.

All other factors are second-order effects in the luminosity spread. A modest burst ($A=20$) will increase the luminosity spread. In a large burst ($A=80$), the number of cool stars born during the burst overwhelm the number of stars formed during the smooth star formation, thereby reducing the age spread.

4. APPLYING STAR FORMATION HISTORIES TO UPPER SCO

The Upper Sco region covers a projected distance of 30 pc on the sky, with a sound crossing time of $\sim 85$ Myr (see discussion in Slesnick et al., 2008). The lack of any known Class 0 or I stars indicates that all star formation within Upper Sco has ceased. de Geus (1992) and Preibisch and Zinnecker (1999) suggest that star formation in Upper Sco may have been triggered by a supernova-driven shock wave from a massive star in the parent Sco-Cen OB Association (see also Preibisch and Zinnecker, 2007; Slesnick et al., 2008). In this context, smooth star formation could occur over a long interval, followed by a supernova that triggers a burst and subsequent quenching of star formation.

4.1. Comparing Synthetic Populations to the Known Upper Sco Members

Our simulations of the star formation history of a cluster correspond to smooth star formation with no trigger (the Smooth-Quench model) and smooth star formation followed by a supernova trigger (Smooth-Burst model). Figure  11 compares of the results of these simulations to the stellar population in Upper Sco for the Smooth-Burst model with parameters $t_{\rm start}=16$ Myr, $t_{\rm quench}=4$ Myr and $A=60$, and for the Smooth-Quench model with $t_{\rm start}=16$ Myr and $t_{\rm quench}=4$ Myr.

The Smooth-Quench model exhibits a temperature-dependent age measurement up to $\sim 6$ Myr between KM stars and F stars, thereby accounting for half of the $\sim 12$ Myr discrepancy (see Table 2). This simplest model demonstrates that an age spread of more than a few Myr in a young ($<15$ Myr) population will lead to a significant temperature dependence in age estimates. This qualitative result is required by the evolution of intermediate mass stars from the Hayashi track to the Henyey track. Cool stars are young and have a median age that changes slowly. For stars with $\log T_{\rm eff}/K=3.6-3.7$, the median age increases rapidly. At temperatures higher than $\log T_{\rm eff}/K>3.7$, stars tend to be older and have median ages that change slowly. These simulation results are consistent with the age spreads of 4 – 5 Myr over a relatively large temperature range of cool stars (e.g. Herczeg and Hillenbrand, 2015), and with the $\sim 3$ Myr age spread for G type and F type stars (Pecaut et al., 2012).

Figure 11 shows that the Smooth-Burst model may explain entirely the empirical temperature-age relationship of Upper Sco. The burst increases the population of cool stars at the end of the star formation process, thereby maximizing the possible differences in ages between cool and hot stars.

Some of the age difference between K-M type and F type stars may be attributed to the incomplete physics in the calculations of pre-main sequence evolution of low-mass stars (e.g. Feiden, 2016), although these explanations may introduce problems into other clusters. However, the age difference between G type stars ($9\pm 3$ Myr) and F type stars ($13\pm 1$ Myr), as calculated by Pecaut et al. (2012), are not easily explained by the same changes that are invoked to move the K-M type stars to older ages. The Smooth-Burst model recovers this age difference between G and F type stars without invoking any shortcomings in stellar evolution models or in empirical temperature and luminosity measurements.

Figure 12 compares the simulated temperature-dependent luminosity spread to the measured luminosity spread within each temperature bin. The luminosity spread in the simulations and observations are consistent at the level of 1$\sigma$. The luminosity spread in the simulation peaks at slightly cooler temperatures than the observed luminosity spread, and also yields luminosity spreads that are systemically smaller than the measured spreads. Both of these uncertainties may be caused by observational errors, including sample incompleteness and mistaken membership, and also by errors in the stellar evolution models. A full evaluation of the star formation history will require a more complete membership list that is free of any biases with spectral type.

4.2. Revised Star Formation History of Upper Sco

The star formation in Upper Sco may have been triggered $\sim 5$ Myr ago by the shock wave from a supernova explosion in Upper Centaurus Lupus (de Geus, 1992; Preibisch and Zinnecker, 1999). This star formation may have then been halted by a shock wave from another supernova that exploded in Upper Sco, destroying the molecular cloud. About 1 Myr ago, this shock wave would have arrived at the $\rho$ Oph molecular cloud and triggered star formation. The scenario of the shock wave from Upper Centaurus Lupus triggering star formation in Upper Sco is challenged by the $\sim 11$ Myr old age of Upper Sco, as measured from intermediate- and high-mass stars (Pecaut et al., 2012; Mamajek et al., 2013), but is consistent with the younger age estimated for low-mass stars. The conflict in the age measurements of different spectral types in the same cluster is an impediment toward our understanding of star formation history.

In this work, we find that when a lengthy star formation is introduced before the trigger event, stars of different spectral types give different age estimations. In this context, we speculate about a plausible revision to the sequentially trigger star formation proposed by Preibisch and Zinnecker (1999), with a possible scenario outlined as follows (see also Figure 13): the star formation in Upper Sco started about $16$ Myr ago with a constant star formation rate with time. About 12 Myr ago, the most massive star in Upper Centaurus Lupus exploded, creating a shock wave that arrived at Upper Sco 5 Myr ago, triggering a burst of star formation. About 3 Myr later, the most massive star in Upper Sco may have destroyed the parent molecular cloud quenched star formation. This star would need to have formed early during the initial constant phase of star formation, because the formation of a $20$ M⊙ star, which explodes in $\sim 10$ Myr, is far more likely than a $60$ M⊙ star, which explodes in $\sim 3$ Myr (Ekström et al., 2012). The shock wave from this supernova explosion may have then arrived at $\rho$ Oph, triggering star formation there. Both supernova and hot star winds are expected to have this type of positive and negative feedback from simulations (e.g. Dale et al., 2012) and also observations of clusters in massive star forming complexes (e.g. Bik et al., 2012; Rivera-Ingraham et al., 2013; Jose et al., 2016a).

Compared with the scenario proposed by Preibisch and Zinnecker (1999), the similar history includes weak star formation with a long duration before the trigger event. This constant star formation may be related to that of Upper Centaurus Lupus and Lower Centaurus Crux, which have older average ages than Upper Sco (e.g. Pecaut and Mamajek, 2016). The physics behind a lengthy star formation may be related to the cooling down of molecular cloud. The depletion time, i.e. the lag time between the formation of molecular time and the start time of the burst of star formation, is an important timescale for the understanding of molecular cloud evolution. As a molecular cloud cools, star formation can occur in over-dense regions; the star formation is local so the star formation rate would be low. The over-dense regions may be distributed randomly over the molecular cloud, with small bursts that appears to be a constant star formation with time. Although the estimation on depletion time requires larger samples and more detailed simulations, the temperature-dependent age estimate may provide a lower limit on past star formation.

5. Discussion and Conclusions

This paper investigates the effects of stellar evolution in analyses of ages and histories of young stellar clusters. Ages of young associations are often estimated in specific spectral type bins (e.g. Pecaut et al., 2012; Herczeg and Hillenbrand, 2015) or estimated by fitting empirical isochrones to clusters (Bell et al., 2015), with membership lists that are incomplete and may be biased to hotter or cooler stars for different regions. However, if star formation proceeds over millions of years, then the different pre-main sequence tracks of intermediate-mass and low-mass stars will lead directly to a spectral type (temperature) dependence in the star formation history. These differences arise when ages are measured in spectral type (or temperature) bins rather than in mass, which is a pragmatic approach because mass is not a measurable quantity.

In this context, in Upper Sco about half of the age discrepancy between F and KM stars may be attributed to the star formation history. A star formation burst at the end of the star formation epoch will maximize these differences and could entirely explain the locus of stars in the HR diagram of Upper Sco. These simulations also explain the differences in median age between F and G stars in Upper Sco, where any age discrepancy is less likely to be introduced by problems in models of pre-main sequence evolution. The accelerated star formation scenario of Palla and Stahler (2000) would be intermediate between the two scenarios modeled here. Isochrone fitting may also suffer from related problems when the known members of different clusters have a different distribution in spectral type or mass. Although the use of Dotter models maximizes differences between intermediate-mass and low-mass stars, other models yield qualitatively similar results. Very low-mass star and brown dwarfs in young associations have much younger ages than assessed from the Li Depletion Boundary (e.g. Malo et al., 2014), which cannot be explained by age spreads and is instead likely caused by errors in evolutionary tracks.

These effects are only present if star formation occurs over some duration of time, and are most significant if the timescale for star formation in the cluster are similar to the timescales for intermediate mass stars to transition from convective to radiative interiors. These timescales correspond to the expected age and possible age spread of Upper Sco. In some cases, most star formation may occur in a single burst, as may be the case for the Orion Nebula Cluster (Reggiani et al., 2011) and for regions of triggered regions of star formation within more massive clusters (e.g. Rivera-Ingraham et al., 2013; Jose et al., 2016a), although even in these cases some past of star formation would affect age measurements for hotter stars without measurable differences for cooler stars. If star formation occurs over an extended timescale of $\sim 5-10$ Myr, then the early-to-mid M-dwarfs will provide the most robust historical record of star formation rate versus time, including average age of stars in the cluster and the luminosity/age spread. This mode of star formation may be quite common, based on inferences of age spreads in high-mass star forming regions (e.g. Bik et al., 2012; Kiminki et al., 2015; Jose et al., 2016b). The early spectral types will provide a record of the earliest epoch of star formation (e.g. Román-Zúñiga et al., 2015; Kiminki et al., 2015; Wu et al., 2016) but will be biased against the most recent epochs of star formation.

As a significant caveat to the application of the simulations to Upper Sco, observational uncertainties in the effective temperature may be especially large. A dedicated study of a single T Tauri star indicates that spots may radically change effective temperatures of low-mass young stars (Gully-Santiago et al., 2017). Spots are likely common and large on low-mass stars through the $\sim 125$ Myr age of the Pleiades (Fang et al., 2016), which may inflate radii (Somers and Stassun, 2016). Some of the age difference between cool and hot stars may be driven by this observational uncertainty and deficiencies in the physics of the pre-main sequence models. Updated stellar models of Somers and Pinsonneault (2015) include spots but are not yet available as a grid and may be challenging to apply to large datasets with stars that have a wide range of spot properties.

Pre-main sequence evolution is still an unsolved problem, with unsettled differences between different models and between models and observations. The age discrepancies between stars of different temperatures in Upper Sco have helped to motivate a new generation of models. Our results here suggest that, while such re-evaluation is important, it is unclear whether these improvements should be judged by comparing coeval isochrones to the population of Upper Sco. A complete census of members in Upper Sco, without any bias versus spectral type, will be necessary to test whether the age discrepancy versus spectral type is driven primarily by an age spread or is the result of large errors in pre-main sequence evolution.

6. Appendix A: Applying magnetic Feiden models to young star-forming regions

The Feiden (2016) evolutionary models introduce magnetic fields in the stellar interior, which reduce the efficiency of convection and thereby slows contraction rates. The available models include large fluctuations at young ages ($<10$ Myr), before the results settle into a stable configuration. We adapt these models by assuming that the evolution at slightly older ages can be extrapolated to young ages with linear fits. This extrapolation then allows these evolutionary tracks to be applied to young associations.

Table 4 applies these models to several associations to obtain a consistent set of ages for these tracks. The known members of each cluster were fit with an empirical isochrone by Herczeg and Hillenbrand (2015), which was then scaled to each region to obtain the expected median luminosities for stars with temperatures of 4200 K, 3800 K, and 3400 K. In this appendix, the ages of 4200 K and 3400 K stars are calculated here using the magnetic Feiden tracks, which were not yet available in 2015. Further details in this simplistic approach can be found in Herczeg and Hillenbrand (2015).

In the magnetic Feiden models, the inferred ages of a 4200 K low mass star are 40 Myr in the Beta Pic Moving Group and 50 Myr in the Tuc-Hor Association, older than the respective Lithium depletion boundary ages of $\sim 24$ and $40$ Myr (Binks and Jeffries, 2014; Malo et al., 2014; Kraus et al., 2014). Stars with spectral types later than M3 have much younger inferred ages for all regions, including Upper Sco, because of discrepancies between observations and predictions from pre-main sequence evolution (see Discussion in Herczeg and Hillenbrand, 2015).

We thank the anonymous referee for constructive comments, which helped improve the quality and clarity of the paper. We also thank Adam Kraus, Michael Ireland, and Jinyi Shangguan for valuable discussions, and Garrett Somers and Greg Feiden for suggestions on the manuscript. QF and GJH are supported by general grant 11473005 awarded by the National Science Foundation of China and by a course development grant awarded by Peking University. This work has made use of data from many facilities, including the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration; the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation; the European Space Agency (ESA) mission Gaia, processed by the Gaia Data Processing and Analysis Consortium, funded by the institutions participating in the Gaia Multilateral Agreement, and the VizieR catalogue access tool, CDS, Strasbourg, France.