Massive stars in the Small Magellanic Cloud : Evolution, rotation and surface abundances

We study the evolutionary and physical properties of evolved O stars in the Small Magellanic Cloud (SMC), with a special focus on their surface abundances to investigate the efficiency of rotational mixing as a function of age, rotation and global metallicity. We analyse the UV + optical spectra of thirteen SMC O-type giants and supergiants, using the stellar atmosphere code CMFGEN to derive photospheric and wind properties. We compare the inferred properties to theoretical predictions from evolution models. For a more comprehensive analysis, we interpret the results together with those we obtained for O-type dwarfs in a former study. Most dwarfs lie in the early phases of the main-sequence. For a given initial mass, giants are farther along the evolutionary tracks, confirming that they are more evolved than dwarfs. Supergiants have higher initial masses and are located past the terminal age main-sequence. We find no clear trend of a mass discrepancy, regardless of the diagram that was used to estimate the evolutionary mass. CNO abundances are consistent with nucleosynthesis from the CNO cycle. Comparisons to theoretical predictions reveal that the initial mixture is important when the observed trends in the N/C versus N/O diagram are to be reproduced. A trend for stronger chemical evolution for more evolved objects is observed. More massive stars, are on average, more chemically enriched at a given evolutionary phase, qualitatively consistent with evolutionary models. Abundance ratios supports the theoretical prediction that massive stars at low metallicity are more chemically processed than their Galactic counterparts. Finally, models including rotation generally reproduce the surface abundances and rotation rates when different initial rotational velocities are considered. Nevertheless, there are objects for which a stronger braking and/or more efficient mixing is required.


Introduction
Although the prime properties of a (single) massive star are set by its initial mass and metallicity, its ultimate destiny and actual path in the Hertzsprung-Russell diagram (HRD) critically depend on mass loss and rotation (Puls et al. 2008;Langer 2012). Furthermore, rotation, the amount of mixing, the metal content, Send offprint requests to: J.-C. Bouret This research is based on observations made with the NASA/ESA Hubble Space Telescope obtained from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. These observations are associated with programmes GO 7437, GO 9434, and GO 11625 Based on data products from observations made with ESO Telescopes at the La Silla Paranal Observatory under programmes ID 67. D-0238, 70.D-0164, 074.D-0109, 079.D-0073, and 079.D-0562 and the mass-loss rate form a feedback loop where more rotation leads to more mixing and changes in the mass-loss rate (e.g. Gagnier et al. 2019), which then affect the rotation rate, and so on. However, the relative roles that these factors play in every phase of single massive-star evolution remain unclear, and evolutionary sequences linking various types of massive stars are poorly defined.
The effects of rotation include an increase in mass-loss rate, a change in the evolutionary tracks in the HRD, a lowering of the effective gravity, an extension of the main-sequence (MS) phase, and the mixing of CNO-cycle processed material up to the stellar surface (Maeder & Meynet 2000;Heger & Langer 2000). A direct consequence of including rotation in stellar evolution calculations is that for a given metallicity, each point in the HRD cannot be uniquely associated with a unique zero-age main-sequence (ZAMS) mass. In stellar evolution codes different settings and prescriptions for rotation-related quantities, such as the transport of chemicals and angular momentum (purely diffusive or advecto-diffusive), instabilities (e.g.Tayler-Spruit dynamo for the internal magnetic fields), shear instabilities, and meridional circulation affect the global structure of a massive star, alter the stellar evolution and the path of a star in the HRD. (Brott et al. 2011;Ekström et al. 2012;Langer 2012).
In the past decade, several studies have investigated the link between surface abundance patterns and rotation, yielding controversial results. The V LT -FLAMES 1 survey (Evans et al. 2005 and V LT -FLAMES Tarantula Survey (VFTS, Evans et al. 2011) found that the majority of the massive OB stars they studied in the Large Magellanic Cloud (LMC) follow the predictions of models well, including rotation. However, the remaining 30-40% of their sample stars present surface nitrogen abundances that cannot be understood in the current framework of rotational mixing for single stars, for instance populations of slowly rotating nitrogen-enriched stars, or rapidly rotating OB-type stars that show little or no nitrogen enhancement (e.g. Hunter et al. 2007Hunter et al. , 2008Hunter et al. , 2009Rivero González et al. 2012;Grin et al. 2017).
At higher metallicity, subsequent studies of surface abundances of OB stars in the Galaxy lead to somewhat different results. Hunter et al. (2009) did not identify clear signs of chemical enrichment at the surface of their B stars, in contrast to Morel et al. (2006Morel et al. ( , 2008. This highlighted the need for samples that include stars that are relatively evolved off the ZAMS to test the effects of rotation. Martins et al. (2015) analysed a large sample of 74 (single) Galactic O-type stars of all spectral types and luminosity classes (20-50 M mass range) that were observed in the context of the MiMeS survey of massive stars (Wade et al. 2016). In this sample, the surface abundances could be explained by models with rotation in a large majority of stars (see also Martins et al. 2017).
To understand the effect of rotational mixing on the chemical evolution of massive stars, studies of stars in the SMC are particularly important because the metallicity in the SMC is only 0.2 Z and half that of the LMC. Rotational mixing is indeed expected to be more efficient at lower metallicities (e.g. Meynet & Maeder 2005;Georgy et al. 2013), because of i) steeper internal gradient of the angular velocity and ii) reduced loss of angular momentum caused by lower mass-loss rates. Furthermore, the interplay between stellar rotation and the stellar wind, which blurs the initial rotational velocity properties of the stars within only a few million years in the case of Galactic massive stars, is expected to be reduced with decreasing metal content. The initial conditions of rotational velocity remain better preserved during the MS life of O-type stars at low metallicity.
Studies of massive star properties at the lower metallicity of the SMC are still relatively rare. Although the fraction of slowly rotating SMC stars with strong nitrogen surface enrichment is similar to the LMC (Hunter et al. 2009), only upper limits on the nitrogen abundance of rapidly rotating B-type stars (with little or no nitrogen enhancement) could be obtained. Dufton et al. (2020) also reported a similar number of nitrogen-rich slowly rotating stars in NGC 346 and in the LMC. Bouret et al. (2013) analysed ultraviolet (UV) and optical spectra of 23 O dwarfs in the SMC and showed that a majority of stars, regardless of their mass, have abundance ratios (N/C versus N/O) that match the predictions of stellar evolution models. 1 Very Large Telescope-Fibre Large Array Multi-Element Spectrograph The fraction of stars close to the ZAMS that show unexplained differences to the standard evolution is 9%, which is far smaller than the 40% outliers of OB-type core-hydrogen burning stars found in the LMC (Hunter et al. 2009;Grin et al. 2017).
Following this initial work, we now focus on more evolved objects. We adopt the same observational strategy and modelling method for optimal consistency between both studies.
The paper is organised as follows. In Sect. 2 we present the observational data sets used in this work. In Sect. 3 we introduce our modelling strategy. We discuss the evolutionary status of the sample stars in Sect. 4, their masses in Sect. 5, and their rotation in Sect. 6. The properties of their surface abundances are discussed in Sect. 7. We then give a summary in Sect. 8.

Observations
Following Bouret et al. (2013), the present work primarily relies on the analysis of UV spectra. Here we use data obtained for our Hubble Space Telescope (HST) COS programme GO 11625 (I. Hubeny, PI), completed with spectra obtained with HST/STIS and available at the Mikulski Archive for Space Telescopes (MAST 2 ). These spectra were supplemented with optical spectra whenever available (see Table 1).
For the purpose of the analysis, we would preferably only work with single stars. The main reason is that predictions of the surface chemical patterns resulting from binary evolution are still very uncertain (see Langer 2012, for a review). On the other hand, we lacked information, a priori, about the multiplicity of the sample stars.
For programme GO 11625, we verified for each selected target that there was no other bright nearby star that might be in the COS aperture, which is 2.5 arcsec in diameter and subtends ∼ 0.8pc at the distance of SMC. For stars from STIS programme GO 7437, Heap et al. (2006) concluded that the targets in the acquisition images did not appear to be multiple. However, they acknowledged that despite the STIS spatial resolution of 0.1", which is much better than that of COS, systems whose separations are smaller than 6300 AU cannot be resolved. However, a detailed analysis of the spectra and of the spectral energy distribution (SED) can reveal a companion star if the companion star is not too faint.

UV and optical data
The sample for our HST/COS programme GO 11625 was extracted from the Massey (2002) UBVR CCD survey of the Magellanic Clouds (see Bouret et al. 2013). We observed each target in a sequence of four science exposures, using gratings G130M and G160M with central wavelengths (λ1291, λ1327) and (λ1577, λ1623) for the adopted settings, respectively. With this sequence, the full UV spectrum from 1132 to 1798 Å is obtained at a resolving power of R ≈ 20,000. Co-addition of overlapping segments improves the final signal-to-noise ratio (S/N), which typically is ≈ 20 -30 per resolution element for a single HST orbit.
For HST/STIS programme GO 7437, stars were observed through the 0.2 -0.2 aperture using the far-UV MAMA detector in the E140M mode. The spectral interval ranges from 1150 Å to 1700 Åand is covered in a single exposure. The effective spectral resolving power of R ≈ 46,000, while the typical S/N per binned data point ranges from 55 to 110 at 1300 Å (the Notes. (a) The primary identification is the AV number (Azzopardi et al. 1975;Azzopardi & Vigneau 1982) . (b) Cross-ID for this star is Sk 80 (Sanduleak 1968), MPG 789 (Massey et al. 1989), or NGC346-001  wavelength with the highest sensitivity), and from 25 to 50 at 1600 Å. We refer to Walborn & Howarth (2000) and Heap et al. (2006) for a more detailed description of the observations and reduction.
The HST/STIS spectrum of AV 232 was retrieved from the MAST archives. This spectrum was obtained for programme GO 9434 (PI: J. Lauroesch), designed for the study of the highly ionised hot component of the interstellar medium (see Jenkins & Tripp 2006). The same detector, gratings, and aperture settings as above were used. The exposure times were set such that an S/N above 30 per resolution element was achieved at 1240, 1393, and 1550 Å.
In the early days of STIS (1999), order 86, which contains N iv λ1718, did not fall on the detector, and stars such as those observed for programme GO 7437 lack this important spectral region. Fortunately, the specific location of different echelle orders on the UV MAMA detector has changed significantly over time, and for spectra after the fourth service mission on HST, order 86 was recorded on the detector. The spectrum of AV 232 therefore does show N iv λ1718. More recently, the standard STIS reduction pipeline stopped extracting order 86 because the regions used to estimate the inter-order background fall off the top of the detector.
Whenever possible, we supplemented the UV spectra with optical spectra that are available from various sources. For all stars from HST/STIS programme 7437, optical coverage exist, either obtained with ESO/CASPEC or with AAT/UCLES. They were extensively presented in , and Heap et al. (2006), to which we refer for more details. For four stars of this programme (AV 15, AV 47, AV 75, and AV 95), ESO/UVES (Dekker et al. 2000) spectra are also available from the ESO Phase 3 data products archive, and we chose to use them as they achieve better S/N for a higher spectral resolving power (typically about 45,000).
AV 232 was observed with ESO/UVES for our programme 079.D.0073 (PI: E. Depagne). The spectrum covers the spectral range ∼ 4720 Å -6830 Å. The reduction process was performed using the UVES context within MIDAS (see Ballester et al. 2000), and it included flat-fielding, bias, and sky-subtraction, and a relative wavelength calibration. Cosmic-ray removal was performed with an optimal extraction method for each spectrum. The S/N is ≈ 130. Another ESO/UVES spectrum of AV 232 obtained in programme 067.D-0238 (PI: P. Crowther) was used to extend the optical coverage down to 3900 Å (the full spectral coverage of this spectrum is ≈ 3730 Å -5000 Å). The S/N of this spectrum is ≈ 125. More details about the instrumental settings and reduction process can be found in Crowther et al. (2002). This star was also part of the ESO VLT/FLAMES survey of massive stars in the SMC , to which we refer for a more detailed description of the properties of the FLAMES data and their reduction). This optical spectrum was modelled and presented in Mokiem et al. (2006Mokiem et al. ( , 2007. The optical spectrum for AV 43 has a lower spectral resolution and was obtained with the AAT/2dF by Evans et al. (2004). These authors provided details about the data and their reduction. Three stars with HST/COS spectra have no optical spectra available (see Table 1), hence they are more prone to uncertainties in the analysis.
The optical spectra obtained with ESO/UVES were normalised interactively using a cubic-spline fit to the pre-defined continuum regions. We examined the regions around important He and CNO lines (see Sect 3.2) and the wings of the Balmer lines with particular care. Details about the normalisation of the spectra obtained with ESO/CASPEC or AAT/UCLES can be obtained in Walborn & Howarth (2000) and in Evans et al. (2004) for the AAT/2dF spectrum.
The non-contemporaneity of the observations in the UV and in the optical may affect the results of the spectroscopic analysis because the significance of a parameter set derived from data obtained several years apart may be questioned (assuming a single set of parameters can indeed be obtained at all). Intrinsic variability for instance could jeopardise our ability to derive such a single set of parameters representing the stellar properties. Variability like this is unknown for our targets. An analogy with Galactic cases indicates that the changes in the stellar and wind properties are not enough to prevent the determination of a single parameter set, where the changes are bracketed within the uncertainties on each physical parameter (e.g. Bouret et al. 2012). Intrinsic variability aside, binarity could also be a cause for confusion for systems that were observed years apart in different spectral bands, as in our case. We have no indication for such binarity in our target stars, except for AV 77, for which the UV spectrum clearly shows the presence of two components (we have no optical spectrum to further confirm this). Confirmation that other stars of our sample are in fact binary systems would require more observations either in the UV or optical or in both, in particular, to search for radial velocity or spectral variations. The ULLYSES programme on HST (the Hubble UV Legacy Library of Young Stars as Essential Standards) or/and its ground-based optical to near-infrared (NIR) range counterpart on ESO/XSHOOTER, which are currently performed, should provide the necessary additional data that we lack at this point.

Model atmosphere
The spectroscopic analysis was performed using model atmospheres and synthetic spectra calculated with the code CMFGEN (Hillier & Miller 1998;Hillier et al. 2003). CMFGEN computes non-local thermodynamic equilibrium (NLTE) line-blanketed model atmospheres, solving the coupled radiative transfer and statistical equilibrium equations in the comoving frame of the fluid in a spherically symmetric outflow. To facilitate the inclusion of extensive line blanketing, a formalism of super-levels is adopted, allowing the incorporation of many energy levels from ions of many different species in the model atmosphere calculations. More specifically, in this work we included ions of H, He, C, N, O, Ne, Mg, Si, P, S, Ar, Ca, Fe, and Ni (see Table 2).
The radiative acceleration was calculated from the solution of the level populations and was used to compute iteratively the hydrodynamical density structure of the inner atmosphere. The inferred velocity law in this quasi-hydrostatic region was then smoothly connected to a β velocity law in the wind. The wind mass-loss rate, density, and velocity are related by the continuity equation. Wind clumping is implemented in CMFGEN assuming the optically thin clumps formalism (see Hillier & Miller 1999, for more details).
After the atmosphere model converged, a formal solution of the radiative transfer equation was computed in the observer's frame (Busche & Hillier 2005), thus providing the synthetic spectrum for a comparison to observations. At this step, we used a radially dependent turbulent velocity (see Hillier & Miller 1998).
Finally, we accounted for the effects of shock-generated Xray emission in the model atmospheres. Because no measured Xray fluxes are available for our targets, we adopted the standard Galactic luminosity ratios logL X /L bol ≈ −7 (Sana et al. 2006), with the bolometric luminosities determined in this study. This scaling ratio might not be valid at low metallicity, depending on the physical nature of the shock generating the X-rays (radiative or adiabatic) and the wind density (Owocki et al. 2013). Depending on the actual relation between the metallicity and the mass-loss rate, logL X /L bol could change accordingly.

Stellar parameters
The photospheric and wind parameters were determined using the same procedure and diagnostics as in Bouret et al. (2013Bouret et al. ( , 2015, to which we refer for more details. Below, we simply outline some details specific to the present work. -Stellar luminosity: We used the flux-calibrated COS and STIS spectra together with optical and NIR photometry (Bonanos et al. 2010, and references therein) to constrain the luminosity, and interstellar extinction as a side-product (the values we derived for E(B − V) are listed in Table 3). In general, it is possible to constrain the intrinsic stellar luminosity by comparing the theoretical SED predicted by a model for a set of fundamental parameters (mostly T eff , log g, R * ,Ṁ, Z) to the observed spectro-photometry when extinction amounts and laws, as well as a distance, are taken into account. We adopted a distance modulus to the SMC of 18.91±0.02 (Harries et al. 2003).
After accounting for a correction for reddening caused by the Galactic foreground (Cardelli et al. 1989), we find that the synthetic SEDs match the observed very well, from the near-IR range (JHK) down to optical (UBVI). This is expected because overall, the SMC extinction laws are quite similar to the Galactic laws in this wavelength range. However, in the UV range, an additional contribution, specific to SMC conditions, was required to account for the flux distribution, which we did using a correction for reddening caused by the SMC interstellar medium (Gordon et al. 2003). While this correction is valid for the SMC bar, some of our targets are located in the outskirts (AV 439 and AV 307) of this galaxy, which might explain why a residual mismatch between the modelled and observed fluxes shortward of ≈ 1250 Å is stronger in these objects. Figure 1 shows a typical fit to the reconstructed SED using UV flux-calibrated spectrum and photometric data, with a theoretical SED computed with CMF-GEN.
For every star of the sample, we adopted a ratio of total-toselective extinction R V = 3.1 (see also Massey et al. 2009), consistent with our work on SMC dwarfs (see Bouret et al. 2013, for a discussion on adopting a different R V ). We did not find obvious indications for deviation from this standard R V with our procedure, even though the extinction is not expected to be accurately described by a single-parameter functional form in the UV (Fitzpatrick & Massa 2007).
-Effective temperature and surface gravities: When optical data were available, we used the ionisation balance of helium lines (equivalent width and line profiles) to derive T eff . Two lines, He i λ5876 and He ii λ4686, were discarded because they often suffer from non-negligible wind contribution in evolved stars. Additionally, we used spectral diagnostics in the UV range, such as ionisation ratios of iron ions (Fe iii to Fe vi), to estimate the effective temperature. For several stars, and in contrast to Bouret et al. (2013), we were unable to use the ratios C iv λ1169 to C iii λ1176 because these lines showed strong signs of wind contamination (the two are sometimes blended).
For stars with both UV and optical spectra, the large number of diagnostics together with the good S/N of the data means that the typical uncertainty on T eff can be as good as 1000 K. How- Notes. The ions listed here are used (or not) according to the effective temperature of the models.  Cutri et al. (2003), 3. Massey (2002), and 4. Evans et al. (2006). The colour excess, E(B − V), is divided into two contributions from the Galactic foreground and the SMC field (see above). ever, for stars with no optical spectra, the uncertainty on T eff can be as large as 1500 K. During the T eff determination process, the relative strength of optical helium lines was used to constrain the helium abundance. For AV 77, AV 307, and AV 439 (UV spectra only), the helium abundance was set to adjust the strength of He ii 1640 Å.
In Table 4, T eff refers to UV-based T eff only for these three stars. For the other stars, the quoted T eff corresponds to that providing the best fit over the whole UV-optical spectral range.
When optical spectra are available, the principal source of uncertainty in the determination of surface gravities is rooted in the use of échelle spectra whose rectification is uncertain, especially in the vicinity of (broad) Balmer lines. This translates into typical uncertainty of 0.1 -0.15 dex (including AV 43 with a 2dF spectrum, i.e. with lower spectral resolution). For stars with UV spectra only, we adopted typical values for the corresponding spectral types (Massey et al. 2009) because photospheric features in the UV range show little sensitivity to changes in surface gravities (but see Heap et al. 2006, for a discussion of using lines from different iron ions in the UV range). The uncertainty was also adopted to be 0.2 dex in this case. Note that UV wind lines do show some sensitivity to log g but we cannot use it for A&A proofs: manuscript no. JCB_AA39890 quantitative purposes until the wind driving is fully understood from first principles. For consistency with our previous study (Bouret et al. 2013), we accounted for the effect of centrifugal forces caused by rotation (see below), and corrected the effective gravity derived from the spectroscopic analysis following the approach outlined in Repolust et al. (2004). These corrected log g values (Table 4) were used to derive the spectroscopic masses.
-Surface abundance: We used the numerous photospheric lines of iron ions (Fe iii to Fe vi) in the UV to constrain the iron content of the stars. The sensitivity of iron line strengths and line ratios to temperature, gravity, iron abundance, and/or microturbulent velocity (see also Heap et al. 2006) has been discussed previously by Bouret et al. (2015). For the photospheric parameters listed in Table 4, the iron lines in the UV spectra of our sample stars are well reproduced with models having Fe/Fe = 0.2, further confirming earlier metallicity determinations in the SMC (Venn 1999). For consistency with Bouret et al. (2013), we used abundances for Mg and Si from Hunter et al. (2007), while the standard solar abundances from Asplund et al. (2009), scaled down by a factor of 5 (or equivalently, by −0.7 dex), were adopted for all other elements included in the models (Ne, P, S, Ar,Ca, Ni), but excluding CNO.
CNO surface abundances and associated errors were determined following the same procedure as in Martins et al. (2015). When the fundamental parameters were constrained, we ran several models in which we changed only the CNO abundances. For a set of selected lines of given element, a χ 2 analysis was performed in which these selected lines were combined. The computed χ 2 values were then renormalised to a minimum value of 1.0, which defines the best-fit abundance. 1σ uncertainties are defined by χ 2 = 2.0.
-Carbon: C iv λ1169 and C iii λ1176 serve as prime indicators for carbon abundances. For stars of the latest type in the sample, where C iii is the dominant ionisation stage of carbon, C iii λ1247 was given more weight than C iii λ1176 because the latter shows some sensitivity to the wind density (see also Bouret et al. 2015). Powerful diagnostic lines in the optical for stars later than O 7.5 are C iii 4070, and C iii 4153, C iii 4157, and C iii 4163 (Martins et al. 2015(Martins et al. , 2017. However, they are at best very weak or even absent in our optical spectra, and could seldom be used. Overall, in the optical the several diagnostic lines that are available for carbon (Martins et al. 2015) are weaker than the UV lines. While C iii λλλ4647, 4650, and 4651 are often detected as the strongest carbon lines in the optical, they are notoriously sensitive to the details of the atomic physics and atmosphere model (Martins & Hillier 2012). Although the effect of Fe metal lines found by Martins & Hillier (2012) should be lower than for galactic O stars because the Fe metal abundance is down a factor of five, we used them only as consistency checks.
-Nitrogen: In the UV, the photospheric lines N iii λλ1183, 1185 and N iii λλ1748, 1752 (for stars observed with COS) were used as primary diagnostics. For stars with COS spectra, N iv λ1718 is also detected but was mostly used as consistency checks because it is potentially affected by the stellar wind (and hence the mass-loss rate and clumping filling-factor). In the optical, N iv 4058 is rarely detected although it is quite frequently seen in Galactic evolved O-type stars with similar spectral types . Other lines that are often used as main diagnostics of nitrogen in the optical include N iii 4511, N iii 4515, N iii 4518, and N iii 4524. These are often weak in several stars, such that we opted not to rely on them except for a few cases such as AV 95, AV 232, and AV 327. We opted to use N iii λλ4634, 4642 as the main diagnostic for the nitrogen abundance in the optical (but see Rivero González et al. 2011González et al. , 2012, for a thorough discussion of their formation processes, in particular, their sensitivity to the background metallicity, nitrogen abundance, and wind strengths).
-Oxygen: In the UV, the O iv λλ1338, 1343 lines generally serve as a diagnostic for the oxygen abundance, with O vλ1371 as a secondary indicator (e.g. Bouret et al. 2013). However, they are affected by wind contribution for stars with significant massloss, such as several of the giants and supergiants of our sample. O iii λλ1150,1154 could be used instead, but these lines are on the shortest side of the UV spectra where the flux level is more uncertain. Only one line is strong enough to be detected (depending on the T eff ) in the optical spectra (O iii λ5592). In all cases, the oxygen content was tuned to improve the fit quality to the lines listed above. However, because of the paucity of sensitive and reliable diagnostics, the values we derive suffer from larger uncertainties than those for carbon and nitrogen.
-Rotation and macro-turbulence: For comparison to observations, the synthetic spectra were convolved with appropriate instrumental profiles, as well as rotational and macroturbulent profiles that contribute to line-broadening. A pure isotropic Gaussian profile was used to mimic the effect of macroturbulence (but see Aerts et al. 2009, for a warning about adopting a Gaussian formulation of macro-turbulence).
To derive the projected rotational velocity (v sin i) of the stars, we first used optical spectra, when available. We used the IACOB-BROAD (Simón-Díaz & Herrero 2014) analysis tool to derive the projected rotational velocity (v sin i) and the contribution of macro-turbulence to the total line-broadening (v mac ). Depending on the spectral type and spectrum quality of the target, we relied on diagnostic lines such as He i λ4471, He i λ4713, He i λ4922, Si ii λ4552, and He ii λ4541. In several instances, however, the optical spectra suffered from significant nebular contamination and/or an insufficient S/N such that unambiguous and accurate determination of v sin i and v mac with this method was not guaranteed. We refer to Ramírez-Agudelo et al. (2013) for a thorough discussion of the different biases and problems with this method in the case of the VFTS.
We verified that the derived v sin i and v mac yielded consistent UV spectral profiles. Many UV spectra show partial or fully resolved components of the C iii λ1176 multiplet, indicating moderate (v sin i < 100 km s −1 ) rotational velocities, with a typical precision of 20 km s −1 . Other useful constraints in the UV are provided by the O ivλλ1338-1343 and numerous iron lines.
For three stars with no optical spectra, we simply convolved the synthetic spectra with rotational profiles (plus instrumental) and varied v sin i until we achieved a good fit to the observed photospheric profiles, adopting a standard limb-darkening law for the convolution. We refer to Hillier et al. (2012) for a discussion about the effect of limb darkening on some spectral lines of fast rotators. For these stars the line-broadening parameter therefore is v sin i alone, which means that the values we quote are upper limits on the true v sin i.
-Wind parameters: For stars with saturated (UV) P Cygni profiles, the wind terminal velocities, v ∞ , were measured from the maximum velocity associated with the saturated black portion of the wind profiles. In cases of unsaturated P Cygni profiles, v ∞ was estimated from the blue-most edge of the absorption component of UV P Cygni profiles. The typical uncertainty for this determination of v ∞ is 100 km s −1 (depending on the maximum microturbulent velocity we adopt). For stars with underdeveloped P Cygni profiles, we confirmed our values by scaling down the wind terminal velocities found in Galactic stars of the same spectral types. More precisely, we used extracts from Table 1 in Kudritzki & Puls (2000), and applied a scaling v ∞ ∝ Z n with n = 0.13 (Leitherer et al. 1992).
Mass-loss rates were derived from the analysis of UV P Cygni profiles of C iv and sometimes S iv resonance doublets, plus the Hα line for stars with optical spectra. For these stars, a unique value ofṀ allows a good fit to both the UV lines and Hα. In all cases, less weight was given to the N v λ1238-1242 doublet as it is blended with the adjacent Lyα line, and is also notoriously sensitive to the interclump medium and to the X-ray flux, which is not constrained for our sample stars (e.g. Zsargó et al. 2008).
The β exponent of the wind velocity law was derived from the fit of the shape of the P Cygni profile. Clumping parameters, f and v cl , were derived in the UV domain, following Bouret et al. (2005). We found no need to tune v cl to improve the fit to observed lines, therefore we adopted a canonical v cl = 30 km s −1 (Bouret et al. 2005). Some photospheric lines in the optical also show some sensitivity to the adopted filling factor (and scaleḋ M). For photospheric H i and He lines, for instance, this is essentially caused by a weaker wind contribution (emission) in clumped models, thus producing deeper absorption than smoothwind models. The wind properties for our giant and supergiant SMC stars will be extensively discussed in a forth-coming paper. Table 4 summarises the fundamental parameters we derived for the giants and supergiants of programmes HST GO 7437 and 11625 introduced in Sect. 2.1. In the following sections, we merge these 13 stars with the sample of 23 dwarfs presented in Bouret et al. (2013) and discuss the properties of the whole (unevolved + evolved stars) sample. Table 6 presents the surface abundances and ratios measured in the present analysis together with those measured for the dwarf sample. Details about the other parameters concerning the unevolved sample can be found in Bouret et al. (2013).
The SEDs and spectra of the evolved sample stars are satisfactorily fit by single-star models, with the exception of AV 77. For this star alone, the spectra distinctively indicate binarity. Its UV spectrum clearly shows the presence of two components. Nevertheless, we cannot rule out that more members of the evolved sample are in fact binary systems, until more observations are obtained, in particular to search for radial velocity variations. The sample of (23) dwarf stars likely contains up to six binaries (Bouret et al. 2013). Appendix B provides the bestmodel fits for each object in the present study.

Evolutionary properties
In Fig. 2 we show the location of the individual stars in the classical HRD (top) and in the Kiel diagram (KD, bottom). The latter depends only on T eff and log g, which quantities are independent of distance and reddening, and can be directly determined from spectroscopy. To build the KD, we corrected the measured surface gravities for the centrifugal acceleration (see Repolust et al. 2004, for a derivation of a first-order analytical form for this correction).
To build the HRD, we used the flux-calibrated UV spectra together with UVBIJHK photometry to derive log L/L for each star. This quantity is well constrained because of the very moderate extinction towards the stars in our SMC sample and because of the accurate distance modulus (see Sect. 3.2). In this sense, the physical properties and evolutionary status of our sample stars, derived from location in the HRD, are more secure than for Galactic (field) stars, where distances are both uncertain and varied. Therefore we expect that the HRD on one hand and a distance-independent plot such as the KD (or rather equivalently, a spectroscopic HR diagram as introduced by Langer & Kudritzki 2014) on the other hand lead to the same conclusions concerning the evolutionary properties of the sample stars.
We compared the locations of the stars in the HRD and KD to evolutionary tracks and isochrones for models with rotation predicted by the Geneva code for metallicity typical of SMC (Georgy et al. 2013). These models start on the ZAMS with an equatorial surface velocity of v init = 0.4v crit , (v crit being the critical surface velocity), which for the mass range we considered is close to ≈ 300 km s −1 .

Star
Sp. Type For stars with UV + optical spectra, uncertainties on T eff are ± 1000 K, ± 0.1 dex for log g, ± 0.1 dex on log L/L . As for the wind quantities, an uncertainty of ± 0.2 dex was estimated for the mass-loss rates (given in M yr −1 ), while v ∞ was measured within ± 100 km s −1 . For stars without optical data, log g is accurate within ± 0.15 -0.2, while uncertainties on T eff are ± 1500 K. For both v sin i and v mac , uncertainties are ± 10 -20 km s −1 .
This initial rotation rate corresponds to velocities of between 110 and 220 km s −1 at the middle of the MS (Ekström et al. 2012). It is therefore a good match to the average v sin i ≈ 110 km s −1 of the full sample, although it must be kept in mind that the initial rotational velocity for evolutionary models is without projection effects. This may be of importance for stars such as the fastest two rotators of our sample (v sin i ≥ 300 km s −1 ), which are also close to the ZAMS, and it can therefore be expected that the typical velocity assumed in this grid of models is probably lower than the actual initial rotational velocities for these dwarfs.
The error bars on the position of stars in the HRD and KD reflect the uncertainties on the effective temperature, surface gravities, and luminosities from the analysis above. The bulk of the sample stars has masses ranging from about 17 to 55 M in both diagrams. The isochrones indicate that the star ages range from 1 Myr up to 10 Myr. We refer to Bouret et al. (2013) for a full discussion of the evolutionary status and ages of the dwarfs, in particular, concerning the stars in NGC 346.
As a whole, the sample stars cover the regions of the HRD and KD where O-type stars are expected. Although the two observing programmes were designed to cover the full range in spectral type and luminosity, there is a significant gap in the distribution, in particular, between the hottest two stars of the sample. The different luminosity classes are relatively well segregated in both diagrams. The observed overlap between stars of different luminosity classes may be the consequence of a misclassification. We recall that we collected the spectral types and luminosity classes from Massey (2002) and Walborn & Howarth (2000), which did not systematically use the same criteria to assign a spectral type and luminosity class, nor did they use the same spectroscopic and photometric data. These uncertainties are possibly related to the sensitivity to abundance and metallicity effects of the luminosity classification criteria for O-type stars (see Massey et al. 2009;Walborn et al. 2014, for a careful presentation of these effects).
Another possible cause for the observed overlap between stars of different luminosity classes is binarity. If binarity is important, luminosities may be overestimated for some stars, lead-ing to a greater luminosity spread. However, only one star in the evolved sample (AV 77) shows unambiguous signs of binarity.
In the HRD, all dwarfs are located in the first part of the Hburning phase, well before the reversal of T eff on the MS (the terminal-age main-sequence, TAMS). The ages indicated by the isochrones are shifted towards higher values than those in Bouret et al. (2013), where we relied on the tracks and isochrones of Brott et al. (2011) with an average v rot = 180 km s −1 . This is due to different input parameters such as the overshooting. Higher initial rotation (about 300 km s −1 in our models) indeed produces a longer MS phase, with a shift of the TAMS to lower temperatures (e.g. Maeder & Meynet 2001). For the three dwarfs that passed the 7 Myr isochrone, a timescale argument excludes that they are at a later evolutionary stage (looping back to the blue for a lower initial mass track). This is further supported by the moderate surface nitrogen enrichments we measure (cf. Sect 7). Models with rotation (Brott et al. 2011;Ekström et al. 2012;Georgy et al. 2013) predict that such enrichment should indeed occur during the MS phase, and is not expected to be delayed to post-MS phases.
The supergiants are clearly separated into two mass bins, one above 45 M and one around 30 M (three and two stars, respectively). All but one (the hottest) are beyond the turn-off, in the part where the tracks evolve to the red at constant luminosities. A closer look shows, however, that this apparently least evolved supergiant, in the second half of the MS next to the M * = 50 M , could also be located on the M * = 44 M and have passed the TAMS. This reading of the HRD is supported by the location of the star in a log g -log(N/C) diagram, showing that the chemical abundances at the surface of the star underwent significant processing compatible with a post-MS phase (see Sect. 7.3). The two ∼30 M stars are located farther along the evolutionary track than the more massive supergiants, indicating that they are more evolved, as further indicated by the nearby 7 Myr isochrone.
The KD also shows a dichotomy for the supergiants, although with different mass distributions. The KD indicates that the five supergiants have masses in the range 40 -60 M . The two supergiants with the highest log g are less evolved than their fellows with lower log g, which are all well beyond the TAMS. Supergiants are the only luminosity class for which the KD indicates systematically higher masses than the HRD, by more than 8 M on average (see next section). Finding some of the most massive stars of our sample among the supergiants could be expected because the most evolved O stars are usually biased toward higher masses (the advanced evolutionary stages of lower mass stars occur at cooler temperatures). Overall, our sample contains six to eight stars that are more massive than 40 M (dwarfs and supergiants), regardless of whether the HRD or the KD is used to determine the stellar mass. This is noteworthy because this mass limit seems to be connected to the properties of stellar evolution for a metallicity typical of the SMC (see e.g. Ramachandran et al. 2019;Dufton et al. 2020).
A direct reading of the HRD indicates that the giants of our sample gather along two evolutionary tracks, with initial masses about 20 M (for the stars O9 IIII and later) and 32 M (the O7-O8 III stars), respectively. Those with M ≈ 32 M lie half-way between dwarfs and supergiants, that is, they lie in the second part of the MS and after the TAMS. The less massive giants (M ≈ 20 M ), however, are well past the TAMS, and as such should be as evolved as supergiants. Analysis of their surface abundance ratios confirms this advanced evolutionary status (Sect. 7.3).
In the KD diagram, the dichotomy between the O7-8 and the O9-B0.5 is still present, although the dispersion is higher. The first group of stars lies in the mass range bracketed by the 25 M and 40 M tracks, in the second part of the MS, reminiscent of the findings by Martins et al. (2017) for Galactic O7-8 giant stars. The later type giants have masses between 20 M and 25 M , and appear to be clearly more evolved (past the TAMS for the 20 M track).
Most dwarfs (slightly less than 90%) have higher log g values than giants, while the lowest log g are found for supergiants. In the KD, supergiants are therefore farther along the evolutionary tracks than the giants and the dwarfs for a given initial mass, which confirms that they are more evolved (see e.g. the 40 M track in the bottom plot of Fig. 2). However, there is a grey zone in which some dwarfs have log g similar to giants (log g = 3.65 and log g = 3.75). One of these two dwarfs is probably a binary (Bouret et al. 2013;Dufton et al. 2019), while the other is ascribed a luminosity class IV by Walborn & Howarth (2000). Furthermore, half of the stars in the giants sample have surface gravities 3.25 ≤ log g ≤ 3.50, which is expected to be typical of supergiants in most calibrations (although these are for Galactic stars, e.g. Martins et al. 2005). The giants located past the TAMS are thus evolved objects with lower masses and luminosities than supergiants. That they have weaker winds explains their luminosity class classification Martins & Palacios (2021) The three stars with the highest surface gravities lie between the 2 and 4 Myr isochrones in the KD and are therefore among the youngest in our sample. This status is qualitatively supported by their location in the HRD, although an age determination based on their location in the HRD is inaccurate because the lower part of the HRD shows notoriously poor sensitivity to age as a result of the small separation between isochrones, and the position of these stars overlaps with different isochrones within their uncertainties. Two of these objects are fast rotators, with v sin i ≥ 300 km s −1 , and the third is classified as a O7 Vz star (see Bouret et al. 2013). It is misleading to link their higher surface gravity to a younger evolutionary stage than the other Otype dwarfs in our sample, however, because three other, hotter dwarfs share the same age range while displaying lower log g (see Sabín-Sanjulián et al. 2014, for a thorough discussion of the O Vz phenomenon, and related classification criteria).

Stellar masses
A comparison of the evolutionary masses derived from the classical HRD and the KD is presented in Fig. 3 by the ratio M t HRD to M t KD as a function of M t HRD . A majority of stars present a ratio M t HRD / M t KD within 1σ of the mass ratio centered at unity, three stars show a ratio M t HRD / M t KD = 1. Twelve objects show discrepancies higher than 1σ (i.e. 25% of the sample), which include the five supergiants, three dwarfs, and two giants. The outlier beyond the 2σ limit in Fig. 3 is MPG 355. This star is found to be significantly less massive in the HRD (85 M ) than in the KD (120 M ), although the masses may marginally agree within the error bars. The optical spectrum of MPG 355 has also been modelled by Massey et al. (2009). Interestingly, their HRD indicates a mass of 75 M , while a KD would still yield a mass of 120 M , although Fig. 3. Ratio of evolutionary masses as derived from the HRD, and KD, respectively. The symbols and colour-coding are the same as in Fig. 2. The full line indicates the one-to-one relation, the dotted line indicates the mean value of the mass ratio. The dashed lines and the shaded area correspond to a difference of ± 1σ, and 2σ from the mean value, respectively. Black circle indicate stars with v sin i ≥ 200 km s −1 , while pink crosses indicate stars whose analysis relied on UV spectra only (no optical data available). A typical error bar for the mass ratio is indicated. the star is marginally farther along its evolutionary sequence in their case. Figure 3 also shows that the majority of stars with v sin i ≤ 200 km s −1 , as well as those for which no optical analysis could be performed to derive the surface gravity, only show very moderate differences of M t HRD over M t KD , regardless of the evolutionary status, as indicated by the luminosity class. Further confirmation requires a larger sample, however, because we are hampered here by small number statistics.
The only systematic property we observe is that the five supergiants have lower M t HRD than M t KD , and that the difference is higher than 1σ. Uncertainties in either the surface gravities or the distance plus reddening correction are not expected to be larger for these five stars, for which we have a full spectral coverage with very good S/N in most cases. Overall, a small majority of stars tend to exhibit M t HRD lower than M t KD (20 stars out of 36 for the full sample). Although very marginal, this is reminiscent of the results by Sabín-Sanjulián et al. (2017) for the LMC, but this was obtained with different evolutionary tracks (i.e. Brott et al. 2011).
In contrast to the results of Markova et al. (2018), we do not observe a systematic trend for higher ratios M t HRD / M t KD as the HRD mass increases. We recall that although Markova et al. (2018) used the HRD and the spectroscopic HRD (sHRD, Langer & Kudritzki 2014) to derive evolutionary masses for a sample of Galactic stars, neither KD nor sHRD require knowledge of stellar distances, and masses are almost independent of stellar radii. We now focus on a comparison of spectroscopic masses to evolutionary masses derived from the HRD and KDs (Fig. 4). As usual, the spectroscopic masses were obtained by combining the measured log g, corrected for centrifugal acceleration (Repolust et al. 2004), with the stellar radius (from T eff from the spectroscopic analysis, and the luminosity from the modelling the SED).
No obvious trend emerges from one plot or the other; the mass ratios are evenly distributed about the 1:1 line. It is noteworthy that fast rotators do not present any specific peculiarities in the mass discrepancy, which concurs with the conclusion of Markova et al. (2018) for Galactic stars.
There might be a weak trend for a stronger mass discrepancy with lower spectroscopic masses (typically below 30 M ) for the evolutionary masses from the HRD and KD. This is reminiscent, although not as clear, of the results of Markova et al. (2018), who found a trend for a stronger mass discrepancy with lower spectroscopic masses for M evol (HRD) and for masses derived from the sHRD (Langer & Kudritzki 2014) in their sample of Galactic stars.
The spectroscopic versus evolutionary mass discrepancy is shown in the left plot in Fig. 4. All stars fall within 2σ of the 1:1 M evol (HRD)/M spec mass ratio, although for 15 stars (i.e. ≈ 41% of the sample) M evol (HRD) differs from M spec by more than 1σ. For M evol (KD)/M spec , the right plot in Fig. 4 shows that for 10 stars (≈ 28% of the sample), this ratio differs by more than 1σ, and in one supergiant (and MPG 355), this ratio differs by more than 2σ. In both cases, our values are therefore fairly close to the normal distributions, according to which 68% of the objects should lie within ±1σ.
The stellar distance is known to a high accuracy for SMC stars, therefore the error in log g is the dominant source of error in the spectroscopic mass, hence in the ratio of the evolutionary mass to the spectroscopic mass. These uncertainties are partially rooted in the difficulties associated with rectifying echelle spectra, especially in the vicinity of the broad H lines. On the other hand, it is striking that the stars for which no optical data are available, hence whose log g should be most affected by uncertainties, do not show a conspicuous trend of a larger mass discrepancy. In addition, in the specific cases of the M evol (KD) versus M spec mass discrepancy, the intrinsically weak dependence of the KD on the distance suggests that the observed mass discrepancy is not a consequence of poorly constrained distances (or alternatively but somewhat counter-intuitively, that the error on the distance and reddening correction that goes into the calculation of the luminosity used in the classical HRD is smaller than the weak distance dependence of KD on the distance).
A possible explanation for the larger M t HRD (or M evol (KD)) compared to M spec is that the spectroscopic masses are underestimated because our CMFGEN models do not account for turbulent pressure when the photospheric stratification is calculated, the neglect of which leads to lower log g than the actual value (see e.g. Massey et al. 2013;Markova et al. 2018). This effect is expected to be stronger for higher mass stars and/or for evolved stars with larger atmospheric scale height. The evidence for this trend is lacking in our sample.

Projected rotational velocities
The top panel of Fig. 5 shows the distribution of v sin i in our sample. A wide range of values is covered, from 20 to more than 350 km s −1 , with a clear bias toward low to moderate values (e.g. 25 stars have v sin i ≤ 120 km s −1 ). The general shape of this histogram is likely a consequence of the way in which the sample is constructed. Although rotation velocity was not a selection criterion of the COS sample, it was meant to ensure that we would have a broad distribution of rotation velocities, especially compared to the the STIS sample alone, which was biased toward stars with narrow-line spectra, and did not include stars rotating faster than 120 km s −1 (Heap et al. 2006). The actual distribution of true rotational velocities is possibly different owing to Fig. 4. Ratio of evolutionary masses (as derived from the HRD and KD, left and right, respectively) to spectroscopic mass. The symbols and colour-coding are the same as in Fig. 3. The dashed lines and the shaded area correspond to a difference of ± 1σ and 2σ from the one-to-one value, respectively. Typical error bars for the mass ratios are indicated. the projection factor. Based on considerations of their mass discrepancies (spectroscopic versus evolutionary from the HRD), Heap et al. (2006) argued that the stars in the STIS sample were preferentially viewed pole-on (see also Hillier et al. 2003). This is important when measured surface abundances as a function of rotation are compared to theoretical predictions (see Sect. 7.4).
Using projected values rather than true rotational velocities may for instance introduce a bias when polar and equatorial regions of a rapidly rotating stars are or are not homogeneously enriched.
To further examine the rotation rates, we split the full sample into two groups, 15 stars with log g ≥ 3.95 (i.e. stars on the MS), and 21 stars for which 3.1 < log g < 3.95 (i.e. post-MS stars (see e.g. Hunter et al. 2009;Dufton et al. 2018)). The bottom panel in Fig. 5 shows the cumulative distribution functions (CDFs) of the unevolved and evolved sub-samples. This plot shows that there is a larger fraction of unevolved stars for projected rotation velocities below 100 km s −1 . The group of unevolved stars also has a fraction of ∼ 26% of stars with v sin i ≥ 200 km s −1 , while the fraction is ∼ 15 % for the evolved stars. This is reminiscent of the findings by Mokiem et al. (2006), who found that their sample of SMC stars contains relatively more rapidly rotating unevolved stars than evolved stars. Furthermore, the difference in our case is significantly higher than what was found by Penny & Gies (2009) in their analysis of SMC stars observed with FUSE, and taken at face value, this seems to support the predictions of stellar evolution that stars slow down as they evolve on the MS (Brott et al. 2011;Ekström et al. 2012). However, a Kolmogorov-Smirnof (KS) test in the combined sub-samples does not allow us to firmly conclude that the distribution of v sin i for the unevolved stars is different from that of the evolved stars (significance level of 0.05 for a p-value = 0.016). This is caused for a large part by the large spread of v sin i in SMC dwarfs.
We are interested in the comparison of the surface abundances as a function of metallicity, therefore we further compared the distribution of projected rotation velocities of our sample stars to that of the sample of Galactic stars analysed by Martins et al. (2015), who used the same tools and method to measure v sin i and CNO abundances. Using the same criteria as above to define unevolved and evolved stars, we built two subsamples of 19 unevolved and 45 evolved stars.  ing faster than 110 km s −1 . We considered merging the sample of Martins et al. (2015) with that of Cazorla et al. (2017), for instance, as these authors analysed a sample of Galactic O-type with fast rotation using CMFGEN and the same method used here and in Martins et al. (2015). However, only one star in their sample would be considered unevolved (log g > 3.95), which would not change the shape of the CDF. In the end, we refrained from doing so to avoid biasing the resulting sample towards high rotation. The left panel of Fig. 6 also shows that 80% of the sample of unevolved Galactic stars has v sin i < 100 km s −1 , while this fraction is smaller than 60% for the SMC sample.
The right panel focuses on evolved stars and shows the opposite: the distribution of v sin i extends to higher values for SMC stars than for their Galactic counterparts. Quantitatively, the CDFs in this panel indicate that the fraction of evolved Galactic stars with v sin i above 200 km s −1 is 8% compared to 14% for the SMC stars.
We then performed a KS test on the sample sets of Galactic versus SMC. We first tested for the hypothesis that the v sin i distribution for SMC dwarf stars is shifted towards faster rotators than in the Milky Way (MW) dwarf stars; the null hypothesis for both distributions having the same mean can be rejected at a significance level of 0.025 (p-value = 0.015). In a similar test for the evolved samples in the SMC compared to the Galaxy, the null hypothesis can be rejected at a significance level of 0.01 (pvalue =0.006); as for the dwarfs, the v sin i distribution for SMC evolved stars are shifted towards faster rotators than the MW stars. We therefore conclude that there is significant evidence that the v sin i distributions for SMC and MW stars are different.
These results should not be viewed as evidence that SMC stars rotate faster than Galactic stars, although Martayan et al. (2007) have shown that this was the case for B stars (see also Penny 1996;Mokiem et al. 2006;Penny & Gies 2009, for O stars). Rather, this indicates that the velocity distributions are different in the Galactic and SMC sample by the construction of these samples. We remind that the sample used by Martins et al. (2015) was biased towards low rotational velocity to maximise magnetic field detection (Wade et al. 2016). This should be kept in mind when the effects of rotation on stellar evolution are discussed. In view of these different conclusions, our result must be considered with caution given the small number of stars in each sub-sample as well.

Surface abundances
In this section we first discuss the nature of the surface abundance patterns and their relation to the CNO cycle. We then focus on the metallicity dependence of surface chemical enrichment. Finally, we discuss the effects of rotation and age on the the degree of chemical processing.

CNO ratios
In Fig. 7 we show the N/C as a function of N/O to probe the chemical evolution of massive stars because these ratios vary substantially throughout the various phases of CNO-cycle (see e.g. Maeder et al. 2014). The measured ratios are correlated, showing a trend of higher N/C with higher N/O. This is consistent with the expectation that surface abundances are affected by products of the CNO cycle operating in the interior of massive stars (e.g. Maeder et al. 2014). More precisely, a trend for stronger chemical evolution for more evolved objects is observed from dwarfs and giants on one hand to supergiants on the other. However, this trend is not as clearly defined as the trend that was observed for Galactic O-type stars by Martins et al. (2015). Most dwarfs and giants populate the same region in the plot in Fig. 7, and a few dwarfs exhibits N/C as high as the supergiants.
Of the seven stars with the highest (log positive) N/Os, five are supergiants and two are dwarfs. Furthermore, AV 15 and AV 75, the two supergiants with the lowest value of log(N/O), are also the least evolved. Of the two dwarfs, one (MPG 355) has the highest initial mass of the sample, and the other (AV 177) is more massive than any of the giants of our sample and is also a fast rotator, with v sin i = 220 km s −1 .
For further information, we overplot in Fig. 7 the expectations of nucleosynthesis through the CNO cycle, more precisely, the limiting cases of the complete CNO, and partial CN cycles derived by Maeder et al. (2014). These curves express the degree of dilution of the CNO abundances at equilibrium mixed with the initial CNO abundances. These limiting cases above are expected to bracket the (more realistic) solutions obtained from numerical models. In numerical models, predicted surface N/C and N/O will not only depend on nucleosynthesis, but are also affected by the adopted rotation rates and the assumptions used to describe mixing processes at different depths in the stars. We adopted a solar initial CNO mixture from Asplund et al. (2009), as in the models we use below (Georgy et al. 2013). Figure 7 shows a clear trend for a majority of dwarfs and giants to lie on or even above the CNO-equilibrium line. Only the supergiants all fall between the two limiting case. The observed dispersion of N/C for a given N/O is a consequence of the dispersion of the evolutionary stage of each individual object as much as the result of various internal conditions favouring either the CN or the CNO cycle. Despite this scatter, Fig. 7 provides evidence that internal mixing produces N enrichment at the expense of carbon and oxygen.
An explanation for some of the N/Cs lying outside the theoretical predictions is that the initial mixture adopted to compute the analytical solutions of Maeder et al. (2014) is not relevant for SMC stars. To illustrate this, we collected different values for the SMC CNO baseline available in the literature. Figure 8 present a comparison of the effect of these initial CNO compositions on the expected trends for the CNO cycle. The CNO initial abundances have been measured from various sites, from B-type stars (Hunter et al. 2009;Rolleston et al. 2003, upper left and right, respectively) to H II regions (Kurt & Dufour 1998, lower left) and using supernova remnant radiative shocks (Dopita et al. 2019). We list them in Table 3. All of them lead to N/O and N/Cs that are significantly different from the solar ratios that were adopted for Fig. 7. The four panels clearly display the effect of the initial N/O and N/Cs on the theoretical expectations for the CN and the CNO equilibrium limiting cases. They also emphasise that these initial ratios may lead to different interpretations of the chemical evolution of our sample stars.
The upper left panel, using initial CNO abundances from Hunter et al. (2009), shows that several stars of our sample have N/O or/and N/Cs that are lower than the baseline ratios. The latter are clearly affected by the high initial value of the mean nitrogen abundance for the SMC derived by Hunter et al. (2009). A recent study by Dufton et al. (2020) based on B-type stars in NGC 346, however, found that N , the SMC baseline abundance for nitrogen is around 6.5 (where N = log [N/H]+12), that is, about 0.8dex lower than Hunter et al. (2009). This offset between initial CNO from Hunter et al. (2009) has been described in Martins & Palacios (2013) for Galactic stars. On the other hand, a majority of our sample stars fall nicely between these two extreme lines, homogeneously covering the parameter space bracketed by the theoretical CN and CNO lines. This suggests that they are still evolving towards the full CNO equilibrium.
The upper right (Rolleston et al. 2003) and lower left panels (Kurt & Dufour 1998) present similar trends: a good fraction of our sample stars fall below the CN curve. Although the measured ratios can be reconciled within the error bars with the expectations for the CN and the CNO equilibrium for the adopted CNO baselines, it is possible that these SMC initial abundances do not represent the initial mixture of the O-type stars of our sample. More quantitatively, this is a direct consequence of too high an N/C but too low an N/O in Rolleston et al. (2003) and Kurt & Dufour (1998). Kurt & Dufour (1998) measured the CNO baseline abundances from H II regions, which are expected to best represent the composition of recently formed stars. However, except for N66, which is located in NGC 346 where several of our sample stars are located, the other H II regions used by Kurt & Dufour (1998) lie well outside the main bar of the SMC, and thus might yield abundances that are not representative of O-type field stars in the SMC.
The CN and CNO equilibrium trends relying on the SMC baseline CNO composition derived by Dopita et al. (2019) (lower right panel) best bracket the measured CNO ratios. These SMC initial abundances therefore likely best represent the CNO mixture out of which our sample stars formed. In terms of chemical evolution, the location of a large majority of the sample stars between the two extreme lines indicates that they are indeed evolving away from the CN-equilibirum towards the full CNO-equilibrium. The individual values of CNO elements differ mildly from those used (scaled solar) in the models of Georgy et al. (2013) that we used in this study. More importantly, the N/O and N/Cs in Dopita et al. (2019) are lower by a factor of 2 at most than in the Sun. For these abundances the number of outliers is the smallest (only three, and only one when the error bars are taken into account). Remarkably, the scaled solar abundances come very close in terms of the smallest number of outliers (cf. Table 5).
There is now ample evidence that nitrogen in the SMC is depleted relative to the other elements (cf. e.g. the references for the discussion above). Regardless of the true value is, this calls for caution when the observed surface abundances of our sample stars are compared to predictions of stellar evolution models built on (scaled) solar abundances. Variations in initial abundances induce change in the zero-points of the nuclear paths. As a consequence, the relative variations of abundances throughout the evolution must be estimated taking the differences in the initial ratios into account. Georgy et al. (2013) noted for instance that when the initial composition measured for the SMC by Hunter et al. (2009) is used, the nitrogen enrichment for a 15 M star in the middle and at the end of the MS differs by slightly more than a factor of two with respect to the enrichment obtained at the same stage when scaled solar initial abundances for the heavy elements are used.
Finally, regardless of the initial mixture we considered, one giant lies well below the baseline N/O and N/Cs. This star is AV 69, which is classified as an OC7.5 giant (Walborn & Howarth 2000). Hillier et al. (2003) analysed AV 69 in detail and found its nitrogen abundance to be very low, more precisely, lower than the SMC baseline. The analysis of Hillier et al. (2003) also indicated that AV 69 is also non-standard in terms of C and O (see Table 6) with respect to the SMC baselines (see Table 5). Hillier Article number, page 13 of 30 A&A proofs: manuscript no. JCB_AA39890  Notes. The scaled Asplund et al. (2009) abundances for CNO are the baseline that we consider in the rest of the analysis as it is used in the Geneva models (but see text). et al. (2003) further argued that this star is a slow rotator that did not experience significant mixing. Slow rotation is indeed an obvious possibility to explain why the star does not show signatures of strong enough mixing that would take the N/C to higher levels. A direct reading of the CNO abundances suggests that they would have to be be primordial in AV 69. Martins et al. (2016) similarly concluded that Galactic OC supergiants likely experienced very little processing through the CN or CNO cycles. They further argued that because the OC phenomenon in Galactic stars is related to the chemically unprocessed surface abundances, it is best seen in more evolved objects for which chemical processing is already strong in morphologically normal stars. As a consequence, for Galactic stars, the OC phenomenon does not exist among dwarfs that are only slightly chemically evolved. When we extrapolate this to environments with lower metallicities such as the MCs, where stronger and faster chemical processing is expected, the OC phenomenon may appear in earlier evolutionary stages and spectral type such as giants and even dwarfs. This is supported by the identification of two LMC early-type stars as OC stars by Evans et al. (2015), or by MPG 113, which is an SMC OC dwarf (Bouret et al. 2003(Bouret et al. , 2013.
In the following sections, we investigate whether the level of surface chemical processing we measured are consistent with expectations from stellar evolution theory. We focus more specifically on the metallicity, rotation, and the evolutionary status of our sample stars.

Dependence on metallicity
In the past decade, several studies have been published on the surface abundances of massive stars in the Galaxy (e.g. Martins et al. 2015;Mahy et al. 2015;Martins et al. 2016Martins et al. , 2017Cazorla et al. 2017;Markova et al. 2018;Carneiro et al. 2019) or in the Magellanic Clouds (e.g. Rivero González et al. 2012;Bouret et al. 2013;Grin et al. 2017). Evolutionary computations predict that for stars of similar masses and evolutionary stages, the chemical enrichment is stronger at lower metallicities for identical initial v init /v crit (e.g. Maeder et al. 2009;Brott et al. 2011;Georgy et al. 2013). Hence the relations between the observed N/C and N/Os may differ as a function of mass and evolutionary stage. Observationally, this prediction still lacks confirmation.
Here we compare our results to those of Martins et al. (2015), which were obtained with similar tools and methods. Figure 9 presents the two samples in a log(N/O) -log(N/C) diagram. No clear trend with metallicity is visible. Although a few SMC dwarfs appear as evolved as some Galactic O supergiants, the opposite is also true: the O supergiants of our sample appear to be somewhat less evolved than some Galactic O dwarfs.
To better isolate any potential metallicity effect between the Galactic and SMC stars, we restricted the samples to stars that cover the same region of the KD (see Fig. 10, top panels). Their surface abundances are shown in the bottom panels of Fig. 10. No conspicuous trend is observed, but a higher N/C (for a given N/O) in SMC stars compared to Galactic stars may be spotted.
A potential pitfall here is that this visual indication might be biased by initial v init /v crit being higher in low metallicity stars than in Galactic stars. In this case, the (potentially) stronger surface enrichment in metal-poor stars might indeed be the consequence of their higher initial rotational velocities rather than be due to steeper internal gradients of angular velocities. This word of caution is supported by the shape of the cumulative distribution functions of v sin i for both populations (cf. Sect. 6, see also the discussion in Penny & Gies 2009, of this topic).
For a more robust view, we applied a principal component analysis (PCA) to the log(N/O) -log(N/C) data to verify whether the SMC stars are on average more enriched in log(N/C) than their MW counterparts for a given log(N/O) value. For a twodimensional data sample, the PCA analysis is equal to a single rotation of the axes, yielding two new components, PC1 and PC2, the first of which points to the direction of maximum data variance, whereas the second is orthogonal to that direction. In our case, PC1 should in particular coincide with the ordinary least-squares regression line, and PC2 would correspond to orthogonal deviations from this regression.
To define the principal components, we joined the two samples (MW + SMC) because we wish to find the directions of maximum variance in the abundance space (see e.g.   Figure 11 shows the abundance space log (N/O) -log(N/C) in this new coordinate system. The dashed black line is the ordinary least-squares regression line when the two data samples are joined. A quick glance at this plot shows that the SMC data have higher PC2 scores than MW data. This is reinforced by the disposition of the SMC points with respect to the grey area, which corresponds to the 2σ confidence band for a locally estimated scatterplot smoothing (LOESS) regression (e.g. Cappellari et al. 2013) of the Milky Way log(N/O) -log(N/C) relation. The LOESS curve is found by a non-parametric regression that fits simple models to localised subsets of the data and smooths the relations. It can be compared with the blue line, which corresponds to the ordinary least-squares regression line for the Milky Way data alone.
We used a two-sample Kolmogorov-Smirnoff test to confirm the null hypothesis that the distribution of PC1 and PC2 scores for the MW and SMC comes from the same distribution. While the data samples do not differ statistically in the first principal component, the null hypothesis can be discarded for the PC2 scores with a probability (p-value) ∼ 4.9 × 10 −5 .
By joining the two samples before applying the PCA, our analysis took a conservative approach. Had we defined the first principal component by using only the MW log(N/O) -log(N/C) data, the PC1 axis would have coincided with the blue line in Fig.  11, and the difference between the PC2 scores for the MW and SMC stars would have been larger. From this we conclude that the SMC stars are on average more enriched in log(N/C) than the MW for a given log(N/O) value.
The PCA of our sample and Martins et al. (2015) supports stellar evolution predictions that low metallicity indeed favours stronger chemical processing. We acknowledge that a more definitive conclusion still requires the analysis of a larger sample of SMC stars (at least commensurate to the available sample of Galactic stars) to define more statistically significant bins of the same initial masses. These should be analysed as homogeneously as possible using a larger variety of tools (atmosphere codes) and methods (e.g. line fitting, equivalent widths, and curve of growth).  Figure 12 shows log (N/C) plotted against the (log of) surface gravity. The N/C indicates the evolutionary level of a star in terms of nucleosynthesis, while the surface gravity is a proxy for the evolutionary point reached by a star on its mass track. While there is clear separation in log g between dwarfs on one hand and the giants and supergiants sample on the other, there is a significant overlap of these two categories of evolved stars that correspond to stars in similar evolutionary states at or just past the TAMS (cf. our comments on Sect. 4, see also Fig. 2). Dwarfs display a wide range of N/C values, while supergiants show a re-markably homogeneous set of N/C values. Giants have a somewhat intermediate behaviour. The two dwarfs with the lowest log g, typical of what we find for some giants, display log(N/C) values similar to those of their more evolved (giants) counterparts. This property further suggests that these two stars are indeed more evolved than the remaining stars in the dwarf sample. We also note that while high N/Cs are not systematically associated with fast rotation, all the fast rotators in our sample display significant signs of chemical processing (cf. Fig. 12). This is further discussed in the next section.

Surface abundances and evolutionary status
To better isolate the effect of evolution, we again plot the N/Cs as a function of log g, but we sort the stars by bins of (ini- Fig. 11. PCA analysis of the full merged SMC + Galactic sample. Squares show SMC stars (this study), and circles correspond to Galactic stars (Martins et al. 2015). See text for details. tial) masses in the two panels of Fig. 13. We used the average mass calculated for the whole sample, that is, ≈ 30 M , to separate the sample. Also illustrated are the evolutionary tracks for 15 M and 40 M with v init = 0.1v crit , v init = 0.2v crit , and v init = 0.4v crit .
The models of Georgy et al. (2013) that we used for SMC stars predict a non-monotonic behaviour of the maximum surface enrichment with mass, unlike the models of Ekström et al. (2012) for Galactic stars (see e.g. Martins et al. 2015): from 15 to 40 M , this maximum value decreases before it increases for masses between 40 and 85 M . The reason for this behaviour in the low-metallicity models of Georgy et al. (2013) is manyfold and discussed in Appendix A.
In the mass range 15 M to 30 M (left panel), all but three stars have log(N/C) falling in the part of the diagram covered by the 15M and 32 M tracks for various initial rotation rates, demonstrating the ability of models with rotation to account for the surface abundances of a large majority of stars in that mass range.
Thirteen out of the sixteen dwarfs are reasonably well reproduced by the models with various degrees of rotation. Most of these objects are little evolved in terms of surface abundances. This is expected given their position in the HRD (or KD), where they lie in the first half of the MS. The N/Cs of the giants in the 15 -30 M mass bin are correctly reproduced by the models. They are also higher than those of the dwarfs in the same mass bin, supporting the prediction from theoretical models that surface chemical processing increases with time.
The three outliers not accounted for by any model are all dwarfs. Two have N/Cs similar to those we find for supergiants, as shown in Fig. 12. The location of these stars in the diagram argues for a much earlier and stronger processing and mixing than displayed in the tracks available in the public grid of models 3 . The third outlier has log(N/C) = -1, which is 0.4 dex below the (solar-like) baseline adopted for the models of Georgy et al. (2013). The measured log(N/C) for AV 468 could be reconciled with a baseline for the SMC if values such as Kurt & Dufour (1998) were adopted instead, so that at this point we cannot consider this object incompatible with model predictions.
The evolved status of the dwarf with the lowest surface gravity (MPG 12, log g = 3.65) is suggested by the luminosity class IV ascribed to this star. It is further supported by a high enhancement in surface nitrogen (see also Bouret et al. 2013, and references therein for more comments). The location of this star in the diagram is closer to the location of the giants of this mass sub-sample, both in terms of log g and N/Cs. The initial mass of MPG 12 is also very similar (if not identical within the error bars) with those of these four giants.
Stars with masses above 30 M , that is six dwarfs, four giants and the five supergiants of the original sample, are shown in the right panel of Fig. 13. All but one star in this mass bin have M < 60 M , the exception being MPG 355 with M = 85 M . The five supergiants show the highest N/C of the whole sample, in agreement with qualitative expectations for their evolutionary status. The average log(N/C) for these five stars is 0.85, with a small dispersion around this value (0.1). Similarly, dwarfs with M * ≥ 30 M have an average log(N/C) ≈ 0.55, again with a relatively small dispersion, which is far smaller that what we observed for the lower mass bin (left panel).
For the dwarfs and supergiants the N/Cs are typically higher than predicted by models for the three values of v init we considered. Interpreting the behaviour of our measurements with respect to model predictions is made difficult by the nonmonotonic evolution of N/C with mass, that is, the reversal above 40 M . Nevertheless, the systematic increase of the N/C with the initial rotation rate in the theoretical models suggests that the enhanced abundances more in line with the observed abundance might be obtained if a larger v init /v crit were used in the evolutionary models.
The four giants in the right panel of Fig. 13 display very different N/Cs, from a significantly processed (AV 95 with log(N/C) = 0.3) to an essentially unprocessed ratio (AV 47 with log(N/C) = -0.61). These three stars have very similar physical parameters and masses (≈ 32 ± 2 M ), and the HRD (Fig. 2) indicates that they are located in the second half of the MS, before the TAMS, hence they have similar evolutionary status. The conspicuous outlier among giants is AV 69 (see Sect. 7.1). In the high-mass bin, above 30 M , the models therefore quantitatively account less well for the observed surface abundances. Qualitatively, observations confirm that more evolved objects show on average more surface chemical processing, in agreement with expectations from mixing due to rotation. For each luminosity class, there is of course a dispersion in each mass bin due to different rotational velocities (cf. next section), but the average trend is relatively clear: in our sample, more massive O-type stars are on average more chemically enriched.
We emphasise that the same conclusion holds qualitatively when initial masses from the KD are used instead of those from the H-R diagram. Fewer than a handful of stars shift from one mass bin to another, and this does not affect the general conclusion we just reached. For a more quantitative test of the ability of theoretical models to reproduce the surface chemistry of SMC O stars, it is obvious that 1) a larger sample and 2) dedicated population synthesis models are necessary.

Surface abundances and rotation
In the specific case of the SMC, and to the best of our knowledge, systematic studies on the effect of rotation on surface abundances have been limited to B-type stars (Hunter et al. 2007(Hunter et al. , 2008(Hunter et al. , 2009Dufton et al. 2018Dufton et al. , 2020. These studies showed that evolutionary models based on the code of Yoon et al. (2006) (see also Brott et al. 2011) failed to reproduce the position of up to 40% of core-hydrogen burning stars in a log (N/H) -v sin i diagram (e.g. a population of slowly rotating nitrogen-rich stars). Our analysis extends the study of Bouret et al. (2013) towards more evolved evolutionary phases, and allows a full comparison to the studies above.

O versus B stars
In Fig. 14 we show the N/Cs of our sample as a function of the (projected) rotation velocity, compared to the location of the B stars from Hunter et al. (2009). With this comparison we determine whether there is a trend of more chemical processing for more massive stars, as predicted by models including rota-tion (Brott et al. 2011;Ekström et al. 2012) and observed in the Galaxy (Martins et al. 2017).
We separated the stars according to their surface gravity. In a first group we gathered stars located on the MS, while in the second group we collected stars near the TAMS and the early post-MS phases. In doing so, we tried to reduce the effect of age on surface abundances. In each sub-sample, N/C should thus depend mainly on rotation and initial mass. The bottom panels of Fig. 14 display the log(N/C) -v sin i diagrams for both subsamples. In the MS sample, O stars with small v sin i have N/C surface values that cover a wide range. At high v sin i the maximum N/C values seem to be the same as in the low v sin i sample, but fast-rotating O stars with low N/C are not observed. This is an indication that high v sin i implies on average high N/C. The opposite is not true because there are slow rotators with N/C as high as fast rotators.
The comparison B stars are all slow rotators and thus no conclusion regarding the effect of v sin i on surface abundance can be drawn for this sub-sample. At low v sin i B stars show a range of N/C values that is almost comparable to that of O stars. In the sub-sample covering the TAMS and early post-MS phases, only a few stars have v sin i > 200 km s −1 . It is thus impossible to discuss any trend of chemical enrichment with rotation. At low v sin i (i.e. up to ∼150 km s −1 ), the few O stars have a range of N/C values similar to that of the MS O stars. The B stars seem to have N/C in the same range, although only two objects reach the N/C of O supergiants. The top panel of Fig. 14 shows in that subsample B stars, O giants and O supergiants correspond to groups of stars having increasing masses. The supergiants are significantly more enriched than both O giants and B stars, which share the same range of N/Cs. Surface chemical processing thus appears to increase with mass above ∼30 M , while it is essentially constant below that limit. This qualitative behaviour is reminiscent of what models predict (see Appendix A).
We stressed in the previous section that for the most massive stars the maximum N/C observed value seemed to correlate with mass but does not do so for the least massive stars, although this last point is not robustly demonstrated. This may be a first observational evidence of the theoretical effects just described. We acknowledge, however, that this needs to be further investigated with larger samples because the number of stars in the current study remains small (about 10 to 20 objects at most). If true, and given our previous conclusion that among O stars the enrichment is stronger for more massive stars, this would indicate that surface abundances depend on mass only for the most massive stars.

Comparison to theoretical predictions
In Fig. 15 we compare the N/Cs as a function of v sin i for the sample as a whole. We further compare the positions of the stars to the different parts of the tracks that correspond to their range of surface gravities as deduced from the spectroscopic analysis. In the left panel, the models are those we adopted throughout this work, that is, they have an initial rotation v init = 0.4v crit (Georgy et al. 2013). The non-monotonic behaviour of the maximum surface enrichment with mass predicted by the models of Georgy et al. (2013), discussed in Appendix A, is also visible in Fig. 15.
For most of the dwarfs located in the area covered by the evolutionary tracks, the measured N/C (around 0.5 ± 0.2 for values of v sin i ( 160 km s −1 ) is predicted to be reached at a later evolutionary stage, regardless of the initial value of the mass in the models. Similarly, all the supergiants fall close to the tracks at the correct evolutionary stage, although for masses different from those measured from the evolutionary diagrams. For both cases, the location of the stars on the tracks qualitatively supports the assumptions of rotational mixing in the models, although quantitatively, the measured N/C appear to require stronger mixing than currently implemented in the models.
Another conspicuous feature of Fig. 15 is the population of dwarfs and giants with v sin i ≤ 160 km s −1 , and log (N/C) < 0, which fall in a region of the log (N/C) -v sin i plot that no tracks with v init = 0.4v crit reach. This may indicate that the initial rotation rates selected for the models are not relevant for a fraction of our sample. To further explore this assumption, the right panel of Fig. 15 shows that models with a range of initial rotation (v init = 0.1v crit , and v init = 0.2v crit ) can indeed reproduce the location of the stars of our sample, although they too seem to struggle to do so for the correct evolutionary stage or mass for several objects. The same comments likely hold for stars with v sin i larger than 250 km s −1 , but for v init that would be larger than 0.4v crit .
Again, we acknowledge that this discussion is limited by the size of the sample. Nevertheless, the analysis of the properties of our sample stars provides an indication that most objects can be accounted for by models including rotation, although a small Fig. 15. log (N/C) as a function of projected rotational velocity (v sin i) for the sample stars. Left: The evolutionary tracks are labelled by initial mass (initial rotation v init = 0.4v crit ). We overplot in blue, green, and red the parts of the tracks that correspond to the minimum and maximum value of log g as measured for the dwarfs, giants, and supergiants, respectively. Right: Models with a range of initial rotation v init = 0.1v crit , 0.2v crit , and 0.4v crit , and initial masses of 15 M and 40 M are shown. In both panels, a horizontal scaling corresponding to a factor π/4 to account for an average projection effect has been applied to the theoretical curves. The horizontal dotted line indicates the baseline N/C for the SMC used in the models. fraction clearly escape these predictions. For these objects, faster braking or more efficient mixing may be needed.

Summary
We have analysed the UV and optical spectra of 13 giants (luminosity class III) and supergiants (class I) O stars in the SMC. Using atmosphere models computed with the code CMFGEN, we derived their physical and evolutionary properties and their surface abundances. These properties and abundances were discussed together with those of the 23 SMC dwarf stars (luminosity class V) presented by Bouret et al. (2013) and analysed with the same method. The properties of the full sample as a whole were compared to predictions of stellar evolution models in Georgy et al. (2013). A discussion of the wind properties of this (full) sample will be presented in a coming paper.
Evolutionary masses determined from the HRD and the KD usually agree relatively well given the respective uncertainties. Spectroscopic and evolutionary masses do not show a clear systematic difference as a function of mass, as reported for instance by Markova et al. (2018) for Galactic O stars. The highresolution spectra used in this work to derive surface gravities are echelle spectra. The rectification of these spectra is notoriously difficult and often uncertain in the vicinity of the broad H lines that are the primary diagnostics of log g. This is even more severe when the signal-to-noise ratios are limited. Because improving log g determinations is crucial for determining stellar masses by spectroscopic methods, a major undertaking in the near future would be to design reliable and accurate spectrum rectification methods.
The surface C, N, and O abundances are consistent with nucleosynthesis from the hydrogen-burning CNO cycle. Comparisons to theoretical predictions reveal the importance of the initial mixture for reproducing the observed trends in the N/C versus N/O diagram.
A PCA analysis of the abundance ratios for the whole (dwarfs and evolved stars) sample supports the theoretical prediction of models including rotation that massive stars at low metallicity are more chemically processed than their Galactic counterparts.
A trend for stronger chemical evolution for more evolved objects is observed, from dwarfs to giants and supergiants. Similarly, for a given evolutionary stage, there is evidence that above ∼30 M , more massive O stars are more chemically processed at their surface. However, comparison with literature data indicates that this may not be the case below ∼30 M when B stars are considered. This pattern is expected from stellar evolution at low metallicity, although the limit above which the mass dependence kicks in is quantitatively different from what we observe. The N/C as a function of the projected rotation is overall well accounted for by models when different initial rotation rates are adopted for the models. However, for the most massive stars of our sample, a stronger mixing efficiency is required to reproduce the N/Cs. In addition, a handful of slowly rotating O stars with high N/Cs is not reproduced by evolutionary models. A similar result has been found for slowly rotating B stars (Hunter et al. 2008;Grin et al. 2017;Dufton et al. 2020). For these objects, faster braking may be required.
Extending our work to much larger samples of O stars is needed if we are to fully understand the factors controlling surface abundances and the interplay between surface abundances and rotation. This is the goal of ambitious programmes such as ULLYSES and XShootU, which will target more than 200 OB stars in the MCs, and provide full spectral coverage from the UV to the optical and NIR. The enabled science supported by the ULLYSES -XShootU data spans from constraining the properties of massive stars at low metallicity and their evolution to building spectral templates for stellar population synthesis, and using them to study stellar populations at low metallicity. Note : For CNO abundances, X = 12 + log (X /H ). For comparison, we give the abundances for He, C, N, and O that we consider as the baseline in our the analysis (scaled solar, Asplund et al. 2009).