Exploring the stellar rotation of early-type stars in the LAMOST Medium-Resolution Survey. II. Statistics

Stellar rotation profoundly affects almost every aspect of stellar evolution (Maeder & Meynet, 2000). It plays a crucial role in dynamo-driven magnetic activity (Spruit, 2002), stellar winds (Matt & Pudritz, 2005), surface abundances (Preston, 1974), chemical yields (Prantzos et al., 2018), internal structure (Aerts et al., 2019), and external structure (Rivinius et al., 2013). Consequently, understanding stellar rotation is essential for our understanding of the origin and evolution of stellar angular momentum. However, the connection between models and observations, although much improved in the last few decades, still suffers from several shortcomings.

From a theoretical perspective, stellar rotation on the main sequence (MS) is the result of the long-term evolution of angular momentum loss and redistribution, regulated by intricate processes early in the star formation phases, e.g. the fragmentation of rotating clouds (Boss & Bodenheimer, 1979), magnetic braking (Mouschovias & Morton, 1985), bipolar overflows (Pudritz, 1985), and magnetic locking with the accretion disk (Lin et al., 2011). Once the initial rotation rate at the zero-age MS (ZAMS) has been established, the rotation evolution can be described by a stellar rotation model that encompasses the physical prescriptions of meridional circulation, angular momentum transport, and (radiative, supra-Eddington, and mechanical) mass loss (Ekström et al., 2012). From an observational perspective, the observed stellar rotation rates are usually entangled with projection effects, which further complicates the use of observed rotation rates as benchmark for stellar rotation models.

The distribution of stellar rotation rates can serve as a crucial probe of the formation and evolution processes of stellar populations. A sharp decline in rotation velocities of MS stars is found around $1.2\,\mathrm{M_{\odot}}$, with more massive stars having mean rotation velocities of hundreds of $\mathrm{km\,s^{-1}}$, while lower-mass stars have velocities on the order of a few $\mathrm{km\,s^{-1}}$ (Kraft, 1970). This is mainly so, because more massive stars lack deep convective envelopes or strong magnetic fields that slow down stellar rotation. Many studies (e.g., Abt & Morrell, 1995; Penny, 1996; Huang & Gies, 2006; Huang et al., 2010) have been dedicated to deconvolve the observational distribution and derive a true rotational distribution by means of large population statistics. A bimodal distribution of rotation rates is observed among young solar-mass stars in the Orion nebula (Attridge & Herbst, 1992; Stassun et al., 1999; Barnes, 2003). The origin of such a bimodality is not yet fully understood, but it is thought to be a manifestation of star–disk interactions during pre-MS evolution (Bouvier et al., 1993).

A similar bimodal distribution is also observed among massive stars ($M\gtrapprox 2.5\,\mathrm{M_{\odot}}$), although the underlying mechanism might not be the same. Abt & Morrell (1995) found a bimodal distribution of rotation velocities among A-type MS stars in the field, where the slowly rotating population was mainly ascribed to chemically peculiar (CP) stars. A large-scale study by Zorec & Royer (2012), who reported $v\sin i$ measurements for 2014 B6- to F2-type field stars, found a substantial difference between the rotation distributions of stars with masses below and above $2.5\,\mathrm{M_{\odot}}$. The less massive population exhibits a unimodal profile, while the more massive population is bimodal. They also revealed the complex behavior of stellar rotation, with strong acceleration in the first one-third of the MS evolutionary phase, followed by a mild spin-down. This calls for more efforts to improve our understanding of stellar rotation, in terms of both theory and observations.

Significant progress in stellar rotation statistics is now achievable with multi-object spectrographs acquiring data for large numbers of stars. Large spectroscopic surveys such as the William Herschel Telescope Enhanced Area Velocity Explorer (WEAVE, Dalton et al., 2012), the Galactic Archaeology with High Efficiency and Resolution Multi-Element Spectrograph (GALAH, Buder et al., 2021), the Radial Velocity Experiment (RAVE, Kunder et al., 2017), the Apache Point Observatory Galactic Evolution Experiment (APOGEE, Majewski et al., 2017), and the Large Sky Area Multi-Object Fiber Spectroscopic Telescope (LAMOST, Cui et al., 2012; Liu et al., 2020) will continue to provide millions of stellar spectra and, hence, allow to unveil an unprecedentedly detailed view of stellar features. In Paper I (Sun et al., 2021), we compiled a catalog containing more than 40,000 stars with effective temperatures of $7000\,\mathrm{K}\leq T_{\mathrm{eff}}\leq 14,500\,\mathrm{K}$ from LAMOST data release (DR) 7 medium-resolution spectra. With such a large and homogeneous sample, we aim to review the rotation distribution of early-type stars as function of their masses, ages, and metallicities.

This article is organized as follows. In Section 2, we analyze the uncertainties in $v\sin i$ in more detail. The effects of possible contamination and its influence on the derived rotation rates are discussed in Section 3. Section 4 discusses our statistical approach to derive the true rotational velocity distribution and we further apply this method to the rotation rate distribution in Section 5. We also discuss the dependence of the stellar rotation rates on metallicity, in Section 6. Finally, our conclusions are summarized in Section 7.

Our stellar sample was defined in Paper I, where we used a data-driven method to estimate the stellar parameters and abundances (‘stellar labels’) of early-type candidates selected from the H$\alpha$ and Mg i b line indices provided by LAMOST DR7 medium-resolution spectra. Having adopted a selection based mainly on effective temperatures ($7000\,\mathrm{K}\leq T_{\mathrm{eff}}\leq 14,500\,\mathrm{K}$) and data quality, we derived a final catalog composed of 40,034 stars. The precisions of our estimates are $\sim 75\,\mathrm{K}$, 0.06 dex, 0.05 dex, and $\sim 3.5\,\mathrm{km\,s^{-1}}$ for $T_{\mathrm{eff}}$, $\log g$, [M/H], and $v\sin i$, respectively, for signal-to-noise ratios, $\mathrm{SNR}>60$. However, this only represents the average behavior of our estimates; the actual performance may vary (Sun et al., 2021). As regards the observational uncertainties, they may originate from the line profiles of a single exposure and/or from minor variations among multiple-epoch observations.

To investigate the dependence of the uncertainty on $v\sin i$, we used the difference between the individual $v\sin i$ values from the single-epoch observations and the associated mean value for the co-added spectra. There ought to be only minor differences among the rotation velocities of single-epoch exposures; their difference represents the uncertainty in our estimates. Royer et al. (2002a) found that the error in $\ln(v\sin i)$ can be approximated by a normal distribution, based on multi-epoch (48), high-SNR ($>250$) observations of Sirius. This is because fast rotation blends the absorption lines and affects the accuracy of radial and rotation velocity measurements. Therefore, the variation was calculated in logarithmic space (see Figure 1). For $\ln(v\sin i)<4$ ($v\sin i<50\,\mathrm{km\,s^{-1}}$), the standard deviation (orange squares) exhibits a dependence on the $v\sin i$ value. This is at odds with the result of Royer et al. (2002b, 2007), who found a uniform error value on a logarithmic scale. This is due to the differences between the methods used to estimate rotation velocities.

Because those earlier authors adopted the first zero of the Fourier transform of the line profiles (Carroll, 1933; Gray, 1976) based on high-resolution spectra ($R>20,000$) to determine $v\sin i$, their estimates are less affected by resolution limitations. However, our data from the LAMOST MRS are unavoidably further broadened due to the finite resolution power of the instrument, which is known as instrumental broadening. The observed line profiles are a convolution of the true profiles and the instrumental profile, which makes it hard to estimate very slow rotation velocities and, hence, results in larger uncertainties for slow compared with fast rotators.

Based on our estimates using a spectrograph with $R\sim 4000$ (Sun et al., 2019), the detection limit for the MRS data ($R\sim 7500$) is around $30\,\mathrm{km\,s^{-1}}$. For these slow rotators, the relative uncertainties are comparable to the absolute values, and so such data should be used with caution. We used an inverse proportional function to fit the standard deviations. The median values of $\Delta\ln v\sin i$ are also presented as red triangles in Figure 1. For stars with rotation velocities above the detection limit, the relative bias—$\mathrm{Bias}(v\sin i)/v\sin i$—is smaller than 3%.

The rotation behavior of our sample stars might not be solely the result of their own evolution but it could also be affected by several external processes. Therefore, possible polluters need to be discarded first. In this section, we present two possible contaminants: binary and CP stars. We check for the incidence of periodic variables and star cluster members in Section 3.3 and Section 3.4, respectively. We then discuss the different rotation properties of these contaminants and the effects of possibly having missed contaminants on the rotation distribution of ‘normal’ (single, non-CP, non-variable, non-cluster member) stars. Following removal of all identifiable contaminants, we are left with 33,655 normal stars (see Section 4.2).

Binaries, particularly close binaries, can drastically affect the rotation properties of both of their components through tides, mass transfer, or even mergers (de Mink et al., 2013). Tidal interactions tend to synchronize stellar and orbital rotation, transferring the stellar spin angular momentum to the orbital angular momentum of the binary system. Mass transfer from the donor star can spin up the accreting star and cause it to approach the break-up limit. Mergers may also change stellar rotation, although this is the least well-understood process.

CP stars represent another possible contaminant type in our sample. They are early-type MS stars exhibiting anomalous chemical abundances. According to Preston (1974)’s classification, they can be divided into four subgroups: metallic line (Am), magnetically peculiar (mAp), stars with enhanced Hg ii and Mn ii (HgMn), and He-weak stars. Most CP stars exhibit slow rotation, but the distribution may differ among subgroups. Abt & Morrell (1995) found that virtually all Am and mAp stars have equatorial rotation velocities below $120\,\mathrm{km\,s^{-1}}$.

To obtain the equatorial rotation velocities, $v$, we need to correct the projected rotation velocities $v\sin i$ for the effects of inclination, $i$. For individual stars, this is only feasible in rare cases, e.g. for eclipsing binaries (where strong tidal effects ensure that the rotation axis is parallel to the orbital axis) and for stars with a rotational modulation signal in their light curves, among others. However, we can correct for projection effects statistically under the assumption of randomly oriented rotation axes (Abt & Morrell, 1995; Royer et al., 2007; Zorec & Royer, 2012).

To estimate the distribution of the equatorial rotation velocities, we followed Royer et al. (2007) and Zorec & Royer (2012), whose method can be summarized in three steps:

1. 1.

   We first smooth the distribution using a Gaussian kernel. The selection of the bandwidth of the kernel-smoothing window is done using Silverman’s rule of thumb (Silverman, 1986) which minimises the mean integrated squared error

   $$h=\left(\frac{4\hat{\sigma}}{3n}\right)^{\frac{1}{5}}$$

   (2)

   where $h$ is the bandwidth of the smoothing parameter, $\hat{\sigma}$ is the standard deviation of the samples, and $n$ is the sample size.

2. 2.

   We then recover the PDF of the true projected rotation velocities, $\theta\equiv\widehat{v\sin i}$, through rectification of the error law. As our measurements are subject to random uncertainties, the observed $v\sin i$ represents the maximum of the probability distribution, whereas $\sigma_{v\sin i}$ is the width of the distribution. Royer et al. (2007) assumed that the error in $\ln(v\sin i)$ was normally distributed, with a mean of 0 and a variance of $\left(\ln(1+\alpha)\right)^{2}$, where the coefficient $\alpha=0.1$ is a constant. However, as shown in Figure 1, our estimates of $v\sin i$ are susceptible to instrumental broadening for small $v\sin i$. Therefore, $\alpha$ depends on the value of $v\sin i$, which must hence be implemented into our calculation. The deconvolution is carried out using the Lucy (1974) iterative technique.

   Although our $v\sin i$ estimates are consistent with literature values up to $\sim 400\,\mathrm{km\,s^{-1}}$, they could be affected by temperature and gravity inhomogeneities introduced by rapid rotation. Gravity darkening suppresses the equatorial region’s contribution to the rotational broadening, which has been discussed extensively in the context of Be stars, i.e. B-type stars exhibiting Balmer-line emission (e.g. Townsend et al., 2004). Nonetheless, a correction for this underestimation is only necessary if the rapid rotators represent a statistically significant fraction (Zorec & Royer, 2012).

3. 3.

   Following correction for the error law, the final step is to correct the distribution for inclination effects. Under the assumption of a random distribution of rotation axes, the PDF of the corrected projected rotation velocities, $\Psi(\theta)$, becomes

   $$\Psi(\theta)=\int\gamma(v)P(\theta|v)\mathrm{d}v=\int\gamma(v)\frac{\theta}{v}\frac{H(v-\theta)}{\sqrt{v^{2}-\theta^{2}}}\mathrm{d}v,$$

   (3)

   where $\gamma(v)$ is the true equatorial rotation velocity distribution and $P(\theta|v)$ is the conditional probability. $H(x)$ represents the Heaviside step function, which is 1 if $x\geq 1$ and otherwise 0.

We provide two examples for different mass ranges (from 1.8 to $2.0\,\mathrm{M_{\odot}}$ and from 2.6 to $2.8\,\mathrm{M_{\odot}}$) to illustrate the process of deconvolution: see Figure 4. The original observed distribution of $v\sin i$ is shown as the cyan histograms with a bin width of $20\,\mathrm{km\,s^{-1}}$. Its smoothed $v\sin i$ distribution, error-corrected true projected rotation velocity ($\theta$) distribution, and projection-rectified true equatorial rotation velocity ($v$) distribution are represented by the blue dotted line, the orange dashed line, and the red solid line, respectively.

These two mass ranges were selected to represent the typical situation in the less massive ($M<2.5\,\mathrm{M_{\odot}}$) and more massive ($M>2.5\,\mathrm{M_{\odot}}$) regimes. Zorec & Royer (2012) found that less massive stars have a unimodal equatorial velocity distribution, while more massive stars trace a bimodal distribution. Our results are largely consistent with previous results, except for several significant differences:

• •

  The velocity distribution of the less massive stars (Figure 4, left) shows a peak velocity at $v\approx 110\,\mathrm{km\,s^{-1}}$, which is $\sim 40\,\mathrm{km\,s^{-1}}$ smaller than the peak velocity reported by Zorec & Royer (2012) for the same mass range.

• •

  For the more massive stars, we found an excess of stars with $v\sin i<100\,\mathrm{km\,s^{-1}}$, which results in a broad and rather flat rectified distribution. This is different from the bimodal distribution found by Abt & Morrell (1995) and Zorec & Royer (2012). The distribution in the right panel of Figure 4 could be interpreted as showing two populations, with one group representing slow rotators that share similar rotation velocities, and a second group comprising fast rotators whose velocities are widely spread. If this were correct, we could still claim a bimodal distribution for our more massive sample. Note that the peak velocity of the slow rotators is similar to that found by Zorec & Royer (2012). However, the number fraction of slow rotators is rather distinctive in both papers. Zorec & Royer (2012) found that the maximum PDF for slow rotators is smaller than that for fast rotators, and the number fraction for the former population is less than 20%. However, slow rotators provide a more distinct contribution to our sample.

• •

  Zorec & Royer (2012) claimed that the sum of two Maxwellian functions can fit the rectified distributions. They used the derived Maxwellian functions to quantify the typical velocities and number fractions of slow and fast rotators. This was based on Deutsch (1967), who showed that the distribution of $v$ can be described by Maxwellian functions. However, we note that this practice may not always yield a reliable representation of the $v$ distribution, especially for the rapidly rotating population.

We suggest that the difference in peak velocity for the less massive stars and the excess of slow rotators among more massive stars can be attributed to missed contamination by CP stars. As discussed in Section 3.5, CP stars are mainly slow rotators and we expected to have $\sim 15$% of such stars misclassified as normal stars. Given these missing CP stars, the peak velocity of the low-mass stars could be biased toward the lower end and a possibly bimodal distribution could be hidden. This may also explain an excess of slow rotators compared with previous studies.

However, recent papers have suggested the existence of a population of slowly rotating normal early-type stars (Sikora et al., 2020; Qin et al., 2021). Most of these studies are based on photometric variability caused by the inhomogeneous structure of the stellar surface, which could be used to infer the rotation period. Qin et al. (2021) searched for low-velocity rotating normal (non-CP) A-type stars from the LAMOST–Kepler project and found a significant fraction of slow rotators among normal A-type stars. Such a contradiction is also reported for the 30 Doradus starburst region in the Large Magellanic Cloud (LMC). The O-type stellar sample (Ramírez-Agudelo et al., 2013) shows a similar rotation distribution as the right panel of Figure 4), exhibiting an excess of slow rotators, while the B-type stellar sample (Dufton et al., 2013) present a bimodal population of very slow and fast rotators, with few stars rotating at the low-velocity peak seen for the O-type stars. These controversies may suggest that previous studies of the rotation properties of early-type stars may have missed a significant number of slow rotators.

For a single star that is not affected by binary interactions, we can generalize the formalism of the observed rotation velocity, $v$, by writing out its dependence on mass and stellar age:

where the stellar age is expressed as the relative age and the initial rotation rate is given by $\omega_{\mathrm{init}}=\Omega/\Omega_{\mathrm{crit}}$ (Georgy et al., 2013). The parameter $\Omega$ is the angular velocity and $\Omega_{\mathrm{crit}}$ is the critical velocity where the surface gravity can no longer maintain equilibrium with the centrifugal motion. Since $v$ changes as a function of time, and also depends on the stellar mass, the parameter $v$ itself is not suitable for comparison among samples with different properties.

We suggest this dependence on mass and age can be alleviated by introducing $v/v_{\mathrm{crit}}$ (e.g. Netopil et al., 2017). For slow rotators, this ratio is almost constant, particularly for rigid rotators on the MS (Zorec & Royer, 2012). We illustrate this by exploring the evolutionary models of the Geneva group (Ekström et al., 2012; Georgy et al., 2013). In Figure 5, we show the $v/v_{\mathrm{crit}}$ evolution for different initial rotation rates. The ratio is almost constant for $0.1\leq t/t_{\mathrm{MS}}\leq 0.9$ if $\omega<0.5$. For fast rotators, the ratio exhibits a mild increase with time. This is to account for the enhanced ‘mechanical’ mass loss in the equatorial region as a star reaches the critical velocity. Even though the $v/v_{\mathrm{crit}}$ ratio unavoidably depends on evolution, it is still a good proxy for the initial rotation rates and only weakly dependent on mass, age, or metallicity. The earliest phase, $t/t_{\mathrm{MS}}\leq 0.1$, corresponds to the time taken to reach a quasi-equilibrium state, while the late phase at $t/t_{\mathrm{MS}}\geq 0.9$ represents the star’s deceleration as it evolves off the MS.

Therefore, we apply a similar analysis as in Section 4 to the $v/v_{\mathrm{crit}}$ ratio. Considering the significant $v/v_{\mathrm{crit}}$ variations that occur during the early and late phases of MS evolution, we limit our analysis to stars characterized by $0.1\leq t/t_{\mathrm{MS}}\leq 0.9$, which covers 86% of our sample. The critical velocity was adopted from Georgy et al. (2013) by interpolating the grid of masses and relative ages for a given metallicity. Note that rapid rotation can increase the MS lifetime (Georgy et al., 2013) by transporting fresh hydrogen into the core. Depending on the initial rotation rate, the stellar lifetime could be extended by up to 30%. This was not considered in estimating the ages for our sample based on the non-rotating stellar evolution models. Nevertheless, estimates of the relative age, which represents the fraction of a star’s lifetime spent on the MS, is not affected by this effect.

An overview of the distribution of $v/v_{\mathrm{crit}}$ as a function of mass for normal stars is displayed in Figure 6. Using the same color scale as before, the figure is color-coded by the relevant number density pertaining to the one-dimensional normalized distribution.

For less massive stars ($M<2.5\,\mathrm{M_{\odot}}$), the distribution of $v/v_{\mathrm{crit}}$ is unimodal, with a significant concentration around the peak velocity, particularly in the mass range $1.8\,\mathrm{M_{\odot}}\leq M\leq 2.0\,\mathrm{M_{\odot}}$. The velocity ratio exhibits a clear dependence on stellar mass, in the sense that it increases from $\sim 0.3$ to $\sim 0.55$ as mass increases, which is consistent with Zorec & Royer (2012). However, a lack of slow rotators ($v\leq 100\,\mathrm{km\,s^{-1}}$) previously reported for this mass range is not seen in our data. Since Zorec & Royer (2012) found that slow rotators are underrepresented even when they are not removed from the sample, they claimed that this lack is unlikely due to contamination by CP stars. Instead, they attributed it to sampling effects in the observed volume. However, in our sample we note that there are fewer slow rotators following our selection of a subsample with SNR $>50$ in the g band and overlapping with LAMOST DR5 to reproduce the parent sample of Qin et al. (2019), Hümmerich et al. (2020), and Paunzen et al. (2021). This way, we can reduce the impact of contamination by CP stars by generating a cleaner sample than that used for Figure 6. However, this relative absence of slow rotators only affects stars with $v/v_{\mathrm{crit}}<0.2$ ($v\leq 70\,\mathrm{km\,s^{-1}}$); a significant fraction of stars with $70\,\mathrm{km\,s^{-1}}\leq v\leq 100\,\mathrm{km\,s^{-1}}$ remains. Therefore, we conclude that, even if the reported relative lack of stars is real, it must be affected by contamination by CP stars but does not extend to $v\sim 100\,\mathrm{km\,s^{-1}}$.

This trend between mass and $v/v_{\mathrm{crit}}$ continues to develop as stars become more massive than $2.5\,\mathrm{M_{\odot}}$ and may continue up to $\sim 3.0\,\mathrm{M_{\odot}}$, where a mild turnover/upper limit is seen. The branch of slow rotators becomes prominent for $M>3.0\,\mathrm{M_{\odot}}$, but the bimodality is already present for $M>2.5\,\mathrm{M_{\odot}}$. The reason why the slowly rotating branch is not prominent in this figure is that it covers an extended distribution across $v/v_{\mathrm{crit}}$. In fact, Zorec & Royer (2012) found that a large fraction of stars have intermediate rotation velocities compared with the slow or fast rotators, and it is not until $M>3.0\,\mathrm{M_{\odot}}$ that a gap between the slowly and rapidly rotating branches becomes apparent. The connecting structure is better visible in Figure 7, where the slowly rotating branch has bifurcated from the rapidly rotating branch at $M\sim 2.5\,\mathrm{M_{\odot}}$. For masses greater than $\sim 3.25\,\mathrm{M_{\odot}}$, the gap between the slowly and rapidly rotating branches is well-established.

As stellar masses increase to beyond $\sim 3.6\,\mathrm{M_{\odot}}$, our analysis suffers from poor sampling in the $\log T$ vs $\log L$ parameter space (Figure 5 Paper I). Massive MS stars in early evolutionary phases are missing in our sample. Nevertheless, we can make an educated guess of their properties based on the available information. If we assume that younger stars behave similarly as their older counterparts, it appears that more stars reside on the slowly rotating branch compared with stars with masses in the range $2.5\,\mathrm{M_{\odot}}\leq M\leq 3.25\,\mathrm{M_{\odot}}$.

A similar bimodal distribution in rotation period has been observed for late-type field dwarfs (e.g. McQuillan et al., 2013), mainly by the Kepler and K2 missions (Borucki et al., 2010). Gordon et al. (2021) detected a gap in the rotation periods for $0.57\,\mathrm{M_{\odot}}<M<0.76\,\mathrm{M_{\odot}}$ based on a sample of $\sim$70,000 K2 targets and attributed it to a departure from the Skumanich (1972) spin-down law. As a young star with its envelope initially decoupled from its core loses angular momentum due to magnetic braking, it could slow or even halt its spin-down if the core and envelope begin to recouple at a later time. Despite the differences between our and their samples in terms of the mass regime covered, such a mechanism might also occur in early-type stars (e.g. Zorec & Royer, 2012).

Other than on mass, stellar rotation also depends on metallicity (e.g., Ekström et al., 2008). The meridional circulation cells that carry angular momentum from the inner regions to the surface become less efficient at lower metallicity. In addition, stellar winds also cause angular momentum losses. Wind strengths are greater for higher metallicities. The models of Georgy et al. (2013) and Szécsi et al. (2015) suggest that normally evolving stars with low mass loss (e.g., low metallicities) may increase their rotation during MS evolution and that the rotation distribution of early-type MS stars in lower-metallicity environments might be different from that in higher-metallicity environments. However, thus far, such differences have not been observed conclusively. Studies of rotational velocities of O- and B-type stars in the LMC and the Small Magellanic Cloud (SMC) (e.g., Penny et al., 2004; Penny & Gies, 2009; Dufton et al., 2019) found little compelling evidence for differences with metallicity. Hunter et al. (2008) found that massive stars (with masses greater than $8\,\mathrm{M_{\odot}}$) in the SMC rotate faster than those in the solar neighborhood. (There is no significant difference between the distribution of rotation rates in the Galaxy and the LMC.) However, the analysis was hindered by field star contamination (similar to the situation in Section 3.4) because of the complicated composition of their sample, which includes target stars projected toward young clusters and a fraction of field stars.

To address the possible metallicity dependence, we selected metal-poor ($\mathrm{[M/H]}<-0.2\,\mathrm{dex}$) and metal-rich ($\mathrm{[M/H]}>0.2\,\mathrm{dex}$) subsamples. Each comprises more than 5000 normal stars to ensure sufficient sample sizes for statistical analysis. Using the same procedures as discussed in Section 5, we analyzed their rotation distributions as a function of stellar mass.

The top and bottom panels of Figure 8 show the rotational profiles for metal-poor and metal-rich stars, respectively. There is a significant difference between both distributions, in the sense that metal-poor stars only exhibit a single branch of slow rotators, with $v/v_{\mathrm{crit}}$ centered at $0.3$–$0.4$, while metal-rich stars clearly show two branches of the stellar velocity ratio. The rapidly rotating branch is similar to the feature we saw in Figures 6 and 7, but only for $M\geq 2\,\mathrm{M_{\odot}}$. The velocity ratio of the less massive, metal-rich stars peaks at $v/v_{\mathrm{crit}}\sim 0.2$. This is much smaller than previously observed for our entire sample. Moreover, the slowly rotating branch is also located at lower velocity ratios than that of the metal-poor subsample.

The shape and location of the velocity distribution pertaining to the metal-rich subsample could be affected by slowly rotating CP stars. Since most CP stars are metal-rich (e.g. Qin et al., 2019), there is a larger fraction of missing CP stars in this metal-rich subsample compared with both the entire sample and the metal-poor subsample, where contamination is expected to be minimal. Netopil et al. (2017) found a linear correlation between the masses and rotation rates of mAp stars, with $v/v_{\mathrm{crit}}$ increasing from 0.08 to 0.17 as mass increases from $\sim 1.8\,\mathrm{M_{\odot}}$ to $\sim 4.5\,\mathrm{M_{\odot}}$ ($\sim 30\%$ scatter). This could be the reason for the low rotation rate of the slowly rotating branch as seen for the metal-rich subsample. Note that the incidence of CP stars is not large enough to completely explain the rotation characteristics for $M\leq 2\,\mathrm{M_{\odot}}$. Thus, among less massive stars (which are most numerous), the rotation rates of normal metal-rich stars are lower than those for the metal-poor stars. The contribution of the rapidly rotating branch becomes obvious for $M\geq 2.5\,\mathrm{M_{\odot}}$.

The lack of a rapidly rotating branch for the metal-poor subsample is intriguing. One possible explanation is that this feature is determined within the first few Myr of a star’s lifetime through ‘disk locking’ (e.g. Bouvier et al., 1993). Magnetic interactions between a star and its accretion disk can remove angular momentum and form intrinsically slow rotators (Lin et al., 2011), and the longer disk lifetime results in slow rotation velocities on the MS (Bastian et al., 2020). However, Yasui et al. (2009, 2010) found that the near-infrared disk fraction of low-metallicity stars declines rapidly on timescales $<1\,\mathrm{Myr}$, much faster than the $\sim 6\,\mathrm{Myr}$ observed for solar abundance. If this were to also apply to our metal-poor subsample, and the initial rotation rates are determined by this mechanism, the metal-poor subsample should rotate faster than their metal-rich counterparts. This may explain the difference in velocity ratio in the slowly rotating branch between the metal-poor and metal-rich subsamples, but not the absence of a rapidly rotating branch.

Another possible mechanism to slow down stellar rotation at birth is magnetic braking (Meynet et al., 2011). The coupling of wind material and magnetic field lines that extend well beyond the stellar surface exerts a torque on the stellar surface layers and thus slows down rotation. Unlike magnetic fields in late-type stars, which are thought to be generated through dynamo action, magnetic fields in early-type stars are likely of fossil origin (Donati & Landstreet, 2009). Sanyal et al. (2017) proposed that high-mass ($M\geq 25\,\mathrm{M_{\odot}}$) stars develop envelope inflation during their lifetime on the MS, which is favorably correlated with metallicity. Therefore, stronger fossil fields are expected for massive stars with lower metallicity than for higher-metallicity stars under the fundamental physics law of flux conservation, where the surface magnetic field strength is expected to decrease as the stellar radius expands. If confirmed, magnetic braking could be more prominent in the metal-poor subsample. However, the incidence rate of magnetic fields in early-type MS stars is around 10% within the solar neighborhood (Sikora et al., 2019), and thus the lack of a rapidly rotating branch cannot be explained exclusively by magnetic braking.

Instead, the lack of a rapidly rotating branch could also be related to the fact that the metal-poor subsample is generally much older than the full sample (see Figure 9). If the overall spin down of metal-poor stars is stronger than predicted by stellar evolutionary models and more (specific) angular momentum is lost during the MS phase, the rapidly rotating branch of the older metal-poor subsample may have been moderated and may even have merged with the slowly rotating branch. This would make it appear as a single population with slightly higher rotation rates than characteristic of the slowly rotating branch.

Using the catalog of Paper I, in this paper we have studied the statistical characteristics of stellar rotation in the MS evolutionary phase. The distribution of the equatorial rotational velocities has been estimated by rectifying the error distribution and correcting for projection effects. We also introduced the ratio $v/v_{\mathrm{crit}}$ to alleviate the dependence of the rotational velocity on mass and age.

Our main results are the following:

• •

  The velocities of the less massive stars ($M<2.5\,\mathrm{M_{\odot}}$) exhibit a unimodal distribution, with their peak velocity ratio increasing from 0.3 to 0.5 as stellar mass increases.

• •

  There is no clear sign of a lack of slow rotators for $M\sim 1.9\,\mathrm{M_{\odot}}$. Even if contamination by CP stars is better corrected for, the gap is smaller than previously reported.

• •

  A bimodal rotation distribution emerges as stars increase in mass ($M>2.5\,\mathrm{M_{\odot}}$). The gap between slow and fast rotation becomes prominent for $M>3.0\,\mathrm{M_{\odot}}$. Compared with the less massive stars, the rapidly rotating branch exhibits a milder increase with stellar mass and a turnover at $M\sim 3.0\,\mathrm{M_{\odot}}$.

• •

  For stars more massive than $3.6\,\mathrm{M_{\odot}}$, more stars reside on the slowly rotating branch compared with stars with masses in the range $2.5\,\mathrm{M_{\odot}}\leq M\leq 3.25\,\mathrm{M_{\odot}}$.

• •

  Metal-poor ([M/H] $<-0.2$ dex) stars only exhibit a single branch of slow rotators, with $v/v_{\mathrm{crit}}$ centered at 0.3–0.4, while metal-rich ([M/H] $>0.2$ dex) stars clearly show two branches. The difference could be attributed to unexpectedly high spin-down rates in the metal-poor subsample. Strong magnetic fields associated with stars with low metallicity could partly explain the result.

• •

  The bulk of the metal-poor stars are characterized by higher rotation rates than the metal-rich stars.

• •

  Roughly 10% of CP (Am and mAp) stars exhibit abnormally high rotational velocities.

Our understanding of early-type stellar rotation could be improved in two ways: (1) currently, sample decontamination largely depends on external catalogs, whose homogeneity and completeness are difficult to assess. This caveat is partially caused by the insufficiently high SNR, limited number of observation epochs, and narrow wavelength coverage of our current spectroscopic survey. In the future, armed with more observations collected by LAMOST-MRS, combined with LAMOST-LRS, we can improve the catalog and apply a more consistent decontamination approach. (2) The description of stellar rotation ($v/v_{\mathrm{crit}}$) still depends on the stellar evolutionary models adopted. Our sample’s homogeneity and size allow us to quantitatively study the rotation distribution as a function of stellar mass, age, and metallicity. For instance, a variate-controlled sample could be used to test the characteristics of stellar evolution on the MS, which in turn could help refine and determine the free parameters remaining in stellar rotation models. Further detailed analysis of the present stellar sample and models of rotating stars will be pursued in forthcoming papers.