<figure><img src="image_1.png" /></figure>

# HD 35502: a hierarchical triple system with a magnetic B5IVpe primary

J. Sikora,<fnref n="1" />,<fnref n="2" />* G. A. Wade,<fnref n="2" /> D. A. Bohlender,<fnref n="3" /> M. Shultz,<fnref n="1" />,<fnref n="2" /> S. J. Adelman,<fnref n="4" /> E. Alecian,<fnref n="5" />,<fnref n="6" /> D. Hanes,<fnref n="1" /> D. Monin,<fnref n="3" /> C. Neiner<fnref n="7" /> and the MiMeS and BinaMIcS Collaborations

Accepted 2016 May 3. Received 2016 May 2; in original form 2016 January 5

## ABSTRACT

We present our analysis of HD 35502 based on high- and medium-resolution spectropolarimetric observations. Our results indicate that the magnetic B5IVsnp star is the primary component of a spectroscopic triple system and that it has an effective temperature of 18.4 ± 0.6 kK, a mass of 5.7 ± 0.6 M<sub>☉</sub>, and a polar radius of 3.0<sup>+1.1</sup><sub>-0.5</sub> R<sub>☉</sub>. The two secondary components are found to be essentially identical A-type stars for which we derive effective temperatures (8.9 ± 0.3 kK), masses (2.1 ± 0.2 M<sub>☉</sub>), and radii (2.1 ± 0.4 R<sub>☉</sub>). We infer a hierarchical orbital configuration for the system in which the secondary components form a tight binary with an orbital period of 5.668 66(6) d that orbits the primary component with a period of over 40 yr. Least-Squares Deconvolution profiles reveal Zeeman signatures in Stokes V indicative of a longitudinal magnetic field produced by the B star ranging from approximately -4 to 0 kG with a median uncertainty of 0.4 kG. These measurements, along with the line variability produced by strong emission in Hα, are used to derive a rotational period of 0.853 807(3) d. We find that the measured v sin i = 75 ± 5 km s<sup>-1</sup> of the B star then implies an inclination angle of the star's rotation axis to the line of sight of 24<sup>+6</sup><sub>-10</sub>°. Assuming the Oblique Rotator Model, we derive the magnetic field strength of the B star's dipolar component (14<sup>+9</sup><sub>-3</sub> kG) and its obliquity (63 ± 13 deg). Furthermore, we demonstrate that the calculated Alfvén radius (41<sup>+17</sup><sub>-6</sub> R<sub>*</sub>) and Kepler radius (2.1<sup>+0.4</sup><sub>-0.7</sub> R<sub>*</sub>) place HD 35502's central B star well within the regime of centrifugal magnetosphere-hosting stars.

**Key words:** stars: early-type – stars: individual: HD 35502 – stars: magnetic field.

## 1 INTRODUCTION

Magnetic B-type stars exhibiting strong emission (e.g. σ Ori E, HD 142184, HD 182180; Landstreet & Borra 1978; Grunhut et al. 2012; Rivinius et al. 2013) serve as important testbeds for understanding how stellar winds interact with magnetic fields. Models such as the Rigidly Rotating Magnetosphere (RRM) model (Townsend & Owocki 2005) provide a qualitative description of these systems (e.g. Townsend, Owocki & Groote 2005; Krtička et al. 2009; Oksala et al. 2010); however, detailed comparisons with observations of σ Ori E have uncovered important discrepancies

which require explanations (Oksala et al. 2012, 2015). By relaxing the RRM model's condition that the magnetic field remains undistorted, Townsend & Owocki (2005) proposed the centrifugal breakout scenario in which the field loops episodically break and reconnect in response to an accumulating magnetospheric mass. Although magnetohydrodynamic simulations support this hypothesis (ud-Doula, Townsend & Owocki 2006), no observational evidence of the breakout events (e.g. optical flares) has yet been reported (Townsend et al. 2013).

While these tests of the current theoretical framework provide useful information, their conclusions are based on a relatively small number of case studies. Lately, this number has been increasing as demonstrated by the recent confirmation of HD 23478's centrifugal magnetosphere (CM; Sikora et al. 2015), as well as the discovery

<footnote n="1">¹Department of Physics, Engineering Physics & Astronomy, Queen's University, Kingston, ON K7L 3N6, Canada</footnote>

<footnote n="2">²Department of Physics, Royal Military College of Canada, PO Box 17000 Kingston, ON K7K 7B4, Canada</footnote>

<footnote n="3">³Herzberg Astronomy and Astrophysics, National Research Council of Canada, 5071 West Saanich Road, Victoria, BC V9E 2E7, Canada</footnote>

<footnote n="4">⁴Department of Physics, The Citadel, 171 Moultrie Street, Charleston, SC 29409, USA</footnote>

<footnote n="5">⁵LESIA, Observatoire de Paris, PSL Research University, CNRS, Sorbonne Universités, UPMC Univ. Paris 06, Univ. Paris Diderot, Sorbonne Paris Cité, 5 place Jules Janssen, F-92195 Meudon, France</footnote>

<footnote n="6">⁶UJF-Grenoble 1/CNRS-INSU, Institut de Planétologie et d'Astrophysique de Grenoble (IPAG) UMR 5274, F-38041 Grenoble, France</footnote>

<footnote n="7">⁷LESIA, Observatoire de Paris, CNRS UMR 8109, UPMC, Université Paris Diderot, 5 place Jules Janssen, F-92195 Meudon Cedex, France</footnote>

<footnote>* E-mail: james.sikora@queensu.ca</footnote>
of the candidate CM-host, HD 345439 (Hubrig et al. 2015). The latest addition to this particular subset of magnetic B-type stars, HD 35502, is the focus of this paper.

Over the past 60 yr, the nature of HD 35502 has been redefined in various ways. Located within the Orion OB1 association (likely within the OB1a subgroup; Landstreet et al. 2007), it was initially identified as a B5V star (Sharpless 1952; Crawford 1958). Higher resolution spectra later obtained by Abt & Hunter (1962) revealed both narrow and broad spectral lines, the latter of which being characterized with a $v \sin i$ of 290 km s$^{-1}$. Moreover, Hα lines were reported to be relatively weak; an analysis of early-type stars within Ori OB1 carried out by Nissen (1976) demonstrated that HD 35502's He abundance was approximately half that of the nearby chemically normal field stars. These results motivated its eventual reclassification as a B5IVsnp star (Abt & Levato 1977).

HD 35502's magnetic field was first detected by Borra (1981) and later confirmed by subsequent studies (Bychkov, Bychkova & Madej 2005; Glagolevskij et al. 2010). Following the initial detection, it had been suggested that some of the unusual features apparent in its spectrum may be related to this strong field. In this paper, we use high-resolution spectra to provide a new interpretation of HD 35502 as a spectroscopic triple system whose primary component is a magnetic B-type star hosting a centrifugally supported magnetosphere.

In Section 2, we discuss both the polarized and unpolarized spectroscopic observations used in this study. Section 3 focuses on our derivation of some of the physical parameters of the system including its orbital configuration, along with the effective temperatures, surface gravities, masses, radii, and projected rotational velocities of the three stellar components. The various analytical methods used to derive these parameters, such as the modelling of spectroscopic and photometric data, are also described. In Section 4, we discuss the evidence of rotational modulation from which we derive the rotational period of HD 35502's primary component. In Section 5, the magnetic field measurements of this component are derived along with the field geometry and strength. In Section 6, we discuss and characterize the magnetic B star's magnetosphere. Finally, our conclusions along with our recommendations for further analytical work to be carried out are summarized in Section 8.

# 2 OBSERVATIONS

## 2.1 ESPaDOnS & Narval spectropolarimetry

Spectropolarimetric observations of HD 35502 were obtained over the course of 5 yr (2008 August 23 to 2013 September 24) in the context of the MiMeS (Wade et al. 2015) and BinaMiCS (Aleician et al. 2015) surveys. 19 Stokes V observations were obtained using the high-resolution ($R \simeq 65\,000$) spectropolarimeter Narval installed at the Télescope Bernard Lyot (TBL) over a wavelength range of approximately 3600–10 000 Å. 10 Stokes V spectra were also obtained using the twin instrument ESPaDOnS installed at the Canada-France-Hawaii Telescope (CFHT). Three of these observations exhibited signal-to-noise ratios (SNRs) $\lesssim 100$ and were removed from the analysis. A median SNR of 522 was obtained from the 26 observations. Both the ESPaDOnS and Narval observations were reduced using the LIBRE-ESPRIT pipeline (Donati et al. 1997) yielding final Stokes I and V spectra (for a detailed description of the reduction procedure, see e.g. Silvester et al. 2012). The Heliocentric Julian Dates (HJDs), total exposure times, and SNRs are listed in Table 1.

## 2.2 dimaPol spectropolarimetry

24 medium-resolution spectropolarimetric observations were obtained with dimaPol ($R \simeq 10\,000$) installed at the Dominion Astrophysical Observatory (DAO; Monin et al. 2012) from 2009 February 7 to 2012 February 15. Two of these observations had SNRs $\lesssim 100$ and were removed from the analysis. The remaining 22 Stokes V observations of $H\beta$ were used to derive longitudinal field measurements; the HJDs, exposure times, SNRs, and longitudinal field measurements are listed in Table 2.

## 2.3 FEROS spectroscopy

32 unpolarized spectra were acquired from 2013 December 30 to 2014 January 3 using the spectrograph FEROS mounted on the 2.2 m MPG/ESO telescope located at La Silla Observatory. The instrument has a resolving power of $R \simeq 48\,000$ across a wavelength range of 3600–9200 Å (Kaufer et al. 1999). The spectra were reduced using the FEROS Data Reduction System. The pipeline automatically carries out bias subtraction, flat-fielding, and extraction of the spectral orders; wavelength calibration is carried out using ThAr and ThArNe lamps. Uncertainties in the measured intensities were estimated from the root mean square (rms) of the continuum intensity at multiple points throughout each spectrum (e.g. Wade et al. 2012). The HJDs, total exposure times, and SNRs are listed in Table 3.

## 2.4 Hα spectroscopy

A total of 131 spectroscopic observations of Hα covering various wavelength ranges from approximately 6300–6800 Å are also used in this study. 113 of these observations were obtained at the DAO from 1991 November 26 to 2012 February 4. Seven of the spectra were removed from the analysis on account of their SNRs being $\lesssim 50$. Both the McKellar spectrograph installed at the 1.2 m Plaskett telescope and the spectrograph mounted at the Cassegrain focus of DAO's 1.8 m telescope were used to acquire the spectra.

The remaining 11 observations were obtained at CFHT from 1991 November 21 to 1995 October 3 using the now decommissioned Coudé f/8.2 spectrograph.

## 2.5 uvby photometry

149 uvby photometric measurements were obtained from 1992 January 27 to 1994 March 13 using the 0.75 m Four College Automated Photoelectric Telescope (FCAPT) on Mt. Hopkins, AZ. The dark count was first measured and then in each filter the sky-ch-c-v-c-v-c-v-c-ch-sky counts were obtained, where sky is a reading of the sky, ch that of the check star, c that of the comparison star, and v that of the variable star. No corrections have been made for neutral density filter differences among each group of variable, comparison, and check stars. HD 35575 was the comparison and HD 35008 the check (i.e. second comparison) star. The standard deviations of the ch-c values were 0.006 mag, except for u for which it was 0.008 mag. We adopted uncertainties of 0.005 mag for each measurement based on the highest precision typically achieved with FCAPT. Table 4 contains the complete list of photometry.

# 3 PHYSICAL PARAMETERS

Based on the high-resolution spectra obtained of HD 35502, three distinct sets of spectral lines are apparent: the strong and broad lines
<table><caption>**Table 1.** ESPaDOnS and Narval spectropolarimetric observations. The SNRs per 1.8 km s⁻¹ pixel are reported at 5400 Å. The fifth, sixth, and seventh columns list the radial velocities of the three stellar components (see Section 3.1). The two right-most columns indicate the longitudinal magnetic field derived from Hβ (see Section 5) along with the associated detection status: definite detection (DD), marginal detection (MD), and no detection (ND) as outlined by Donati et al. (1997).</caption><thead><tr><th>HJD</th><th>Total exp. time (s)</th><th>SNR (pix<sup>-1</sup>)</th><th>Instrument</th><th>v<sub>r,B</sub> (km s<sup>-1</sup>)</th><th>v<sub>r,A<sub>1</sub></sub> (km s<sup>-1</sup>)</th><th>v<sub>r,A<sub>2</sub></sub> (km s<sup>-1</sup>)</th><th>⟨B<sub>z</sub>⟩<sub>Hβ</sub> (kG)</th><th>Detection status</th></tr></thead><tbody><tr><td>2454702.138</td><td>1800</td><td>677</td><td>ESPaDOnS</td><td>19.4 ± 1.3</td><td>55.9 ± 1.4</td><td>4.8 ± 2.7</td><td>-0.23 ± 0.19</td><td>DD</td></tr><tr><td>2455849.677</td><td>3600</td><td>662</td><td>Narval</td><td>20.0 ± 1.3</td><td>-15.3 ± 2.6</td><td>71.3 ± 1.7</td><td>-0.03 ± 0.16</td><td>DD</td></tr><tr><td>2455893.623</td><td>3600</td><td>522</td><td>Narval</td><td>21.2 ± 1.4</td><td>-5.1 ± 1.5</td><td>61.3 ± 2.9</td><td>-2.05 ± 0.20</td><td>DD</td></tr><tr><td>2455910.518</td><td>3600</td><td>263</td><td>Narval</td><td>22.7 ± 2.0</td><td>0.7 ± 2.0</td><td>55.9 ± 5.3</td><td>-1.30 ± 0.50</td><td>DD</td></tr><tr><td>2455934.528</td><td>3600</td><td>587</td><td>Narval</td><td>19.6 ± 1.5</td><td>-20.9 ± 1.1</td><td>76.4 ± 2.0</td><td>-2.12 ± 0.19</td><td>DD</td></tr><tr><td>2455936.534</td><td>3600</td><td>472</td><td>Narval</td><td>18.6 ± 1.4</td><td>76.5 ± 1.6</td><td>-19.4 ± 2.8</td><td>-2.32 ± 0.24</td><td>DD</td></tr><tr><td>2455938.525</td><td>3600</td><td>564</td><td>Narval</td><td>20.4 ± 1.5</td><td>18.1 ± 5.5</td><td>37.0 ± 5.9</td><td>0.08 ± 0.19</td><td>ND</td></tr><tr><td>2455944.500</td><td>3600</td><td>494</td><td>Narval</td><td>20.9 ± 1.4</td><td>1.9 ± 2.5</td><td>54.0 ± 2.9</td><td>0.16 ± 0.23</td><td>ND</td></tr><tr><td>2455949.429</td><td>3600</td><td>545</td><td>Narval</td><td>20.5 ± 1.3</td><td>45.4 ± 1.1</td><td>11.2 ± 2.5</td><td>-0.95 ± 0.20</td><td>DD</td></tr><tr><td>2455950.472</td><td>3600</td><td>478</td><td>Narval</td><td>21.4 ± 1.8</td><td>-12.5 ± 2.4</td><td>68.1 ± 2.9</td><td>-0.18 ± 0.23</td><td>ND</td></tr><tr><td>2455951.471</td><td>3600</td><td>431</td><td>Narval</td><td>22.8 ± 1.5</td><td>-22.6 ± 1.0</td><td>77.7 ± 2.0</td><td>-1.00 ± 0.27</td><td>DD</td></tr><tr><td>2455966.376</td><td>3600</td><td>550</td><td>Narval</td><td>19.8 ± 1.5</td><td>48.6 ± 1.5</td><td>7.5 ± 3.2</td><td>-1.87 ± 0.20</td><td>DD</td></tr><tr><td>2455998.332</td><td>3600</td><td>397</td><td>Narval</td><td>20.9 ± 1.4</td><td>53.6 ± 2.6</td><td>3.4 ± 3.4</td><td>0.28 ± 0.28</td><td>ND</td></tr><tr><td>2455999.362</td><td>3600</td><td>402</td><td>Narval</td><td>20.7 ± 1.9</td><td>81.8 ± 2.9</td><td>-25.2 ± 1.0</td><td>-1.64 ± 0.28</td><td>DD</td></tr><tr><td>2456001.309</td><td>3600</td><td>523</td><td>Narval</td><td>20.0 ± 1.5</td><td>-4.6 ± 2.4</td><td>60.0 ± 1.6</td><td>-2.39 ± 0.23</td><td>DD</td></tr><tr><td>2456003.329</td><td>3600</td><td>528</td><td>Narval</td><td>21.0 ± 1.6</td><td>13.5 ± 2.2</td><td>42.8 ± 1.9</td><td>0.51 ± 0.21</td><td>DD</td></tr><tr><td>2456202.665</td><td>3600</td><td>604</td><td>Narval</td><td>21.4 ± 1.3</td><td>65.4 ± 3.4</td><td>-9.4 ± 1.9</td><td>-2.30 ± 0.17</td><td>DD</td></tr><tr><td>2456205.618</td><td>3600</td><td>494</td><td>Narval</td><td>21.1 ± 1.5</td><td>-15.1 ± 2.3</td><td>69.5 ± 2.8</td><td>-0.34 ± 0.22</td><td>DD</td></tr><tr><td>2456224.646</td><td>3600</td><td>450</td><td>Narval</td><td>23.0 ± 1.3</td><td>27.7 ± 1.4</td><td>27.7 ± 1.4</td><td>-0.52 ± 0.24</td><td>ND</td></tr><tr><td>2456246.505</td><td>3600</td><td>505</td><td>Narval</td><td>20.5 ± 1.3</td><td>-15.6 ± 2.6</td><td>69.8 ± 1.5</td><td>-1.72 ± 0.22</td><td>DD</td></tr><tr><td>2456293.881</td><td>1600</td><td>710</td><td>ESPaDOnS</td><td>22.6 ± 1.4</td><td>79.4 ± 2.5</td><td>-24.9 ± 0.9</td><td>-1.39 ± 0.20</td><td>DD</td></tr><tr><td>2456295.787</td><td>1600</td><td>190</td><td>ESPaDOnS</td><td>21.0 ± 2.6</td><td>10.4 ± 1.6</td><td>44.1 ± 1.4</td><td>-3.06 ± 0.71</td><td>ND</td></tr><tr><td>2456295.808</td><td>1600</td><td>231</td><td>ESPaDOnS</td><td>22.8 ± 2.3</td><td>9.3 ± 1.3</td><td>45.2 ± 1.9</td><td>-2.83 ± 0.57</td><td>MD</td></tr><tr><td>2456556.002</td><td>1600</td><td>582</td><td>ESPaDOnS</td><td>24.7 ± 1.9</td><td>43.0 ± 1.9</td><td>10.8 ± 1.8</td><td>-1.37 ± 0.28</td><td>DD</td></tr><tr><td>2456557.140</td><td>1600</td><td>670</td><td>ESPaDOnS</td><td>19.2 ± 1.5</td><td>-18.2 ± 1.9</td><td>70.9 ± 1.6</td><td>-3.27 ± 0.20</td><td>DD</td></tr><tr><td>2456560.077</td><td>1600</td><td>612</td><td>ESPaDOnS</td><td>20.8 ± 1.4</td><td>74.5 ± 1.9</td><td>-20.5 ± 2.1</td><td>-0.05 ± 0.25</td><td>ND</td></tr></tbody></table>

<table><caption>**Table 2.** Spectropolarimetric observations of HD 35502 obtained with dimaPol. Columns 1–3 list the HJDs, exposure times, and SNRs. Column 4 lists the longitudinal field measurements derived from the Hβ Stokes V profiles.</caption><thead><tr><th>HJD</th><th>Total exp. time (s)</th><th>SNR (pix<sup>-1</sup>)</th><th>⟨B<sub>z</sub>⟩<sub>Hβ</sub> (kG)</th></tr></thead><tbody><tr><td>2454869.764</td><td>4800</td><td>410</td><td>-1.95 ± 0.27</td></tr><tr><td>2454872.773</td><td>3600</td><td>240</td><td>-0.64 ± 0.25</td></tr><tr><td>2455109.047</td><td>5400</td><td>240</td><td>-2.88 ± 0.51</td></tr><tr><td>2455110.964</td><td>5400</td><td>260</td><td>-1.60 ± 0.25</td></tr><tr><td>2455167.855</td><td>6000</td><td>290</td><td>-2.27 ± 0.21</td></tr><tr><td>2455168.829</td><td>6000</td><td>320</td><td>-2.21 ± 0.38</td></tr><tr><td>2455169.881</td><td>7200</td><td>280</td><td>-1.72 ± 0.43</td></tr><tr><td>2455170.867</td><td>7200</td><td>190</td><td>0.31 ± 0.31</td></tr><tr><td>2455190.766</td><td>7200</td><td>230</td><td>-2.57 ± 0.33</td></tr><tr><td>2455191.781</td><td>7200</td><td>260</td><td>-2.29 ± 0.30</td></tr><tr><td>2455192.790</td><td>7200</td><td>130</td><td>-3.95 ± 1.09</td></tr><tr><td>2455193.742</td><td>7200</td><td>210</td><td>-1.61 ± 0.33</td></tr><tr><td>2455261.655</td><td>6000</td><td>290</td><td>-2.04 ± 0.41</td></tr><tr><td>2455262.685</td><td>6000</td><td>300</td><td>-2.31 ± 0.26</td></tr><tr><td>2455264.661</td><td>6000</td><td>270</td><td>-0.89 ± 0.37</td></tr><tr><td>2455580.803</td><td>6000</td><td>270</td><td>-0.07 ± 0.38</td></tr><tr><td>2455583.800</td><td>4800</td><td>140</td><td>-2.38 ± 0.94</td></tr><tr><td>2455594.711</td><td>7200</td><td>230</td><td>-2.39 ± 0.60</td></tr><tr><td>2455611.691</td><td>5700</td><td>210</td><td>-2.07 ± 0.47</td></tr><tr><td>2455904.864</td><td>6600</td><td>170</td><td>-3.60 ± 0.56</td></tr><tr><td>2455964.666</td><td>5400</td><td>190</td><td>-2.48 ± 0.34</td></tr><tr><td>2455972.644</td><td>7200</td><td>260</td><td>0.58 ± 0.44</td></tr></tbody></table>

associated with a hot star and two nearly identical components attributable to two cooler stars which are observed to change positions significantly. Based on HD 35502's reported spectral type, the bright, dominant component is presumed to be a hot B5 star (Abt & Levato 1977); the weaker components are inferred to be two cooler A-type stars based on the presence of Fe II lines and the absence of Fe III lines. As will be shown in the next section, the lines of the A stars show velocity variations consistent with a binary system. Hence we conclude that HD 35502 is an SB3 system.

Some of the B star lines (the He I lines, in particular) appear to exhibit intrinsic variability. Such features are commonly found in magnetic He peculiar stars (e.g. Borra, Landstreet & Thompson 1983; Bolton et al. 1998; Shultz et al. 2015).

## 3.1 Orbital solution

The radial velocity ($v_r$) of the central B star (B) in each observation was determined using spectral lines for which no significant contribution from the two A stars ($A_1$ and $A_2$) was apparent. C II λ4267 was found to be both relatively strong (with a depth of 10 per cent of the continuum) and only weakly variable. The Hα spectra encompassed a limited range of wavelengths with few lines from which $v_r$ could be accurately determined. We used C II λ6578 and He I λ6678 in order to estimate $v_r$ from all of the spectra (spanning a 22 yr period). However, measurements made from He I λ6678 were subject to systematic errors associated with strong variability (see Section 6). Moreover, the shallower depth of C II λ6578 (<4 per cent
<table><caption>**Table 3.** Unpolarized spectra obtained using FEROS. The SNRs per 2.8 km s⁻¹ pixel listed in column 3 are estimated from the rms of the continuum near λ = 5400 Å. The fourth, fifth, and sixth columns list the radial velocities of the three stellar components (see Section 3.1).</caption><thead><tr><th>HJD</th><th>Total exp.<br>time (s)</th><th>rms<br>SNR</th><th>v<sub>r,B</sub><br>(km s<sup>-1</sup>)</th><th>v<sub>r,A1</sub><br>(km s<sup>-1</sup>)</th><th>v<sub>r,A2</sub><br>(km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>2456656.568</td><td>300</td><td>271</td><td>19.5 ± 1.6</td><td>75.4 ± 2.4</td><td>-25.4 ± 1.5</td></tr><tr><td>2456656.609</td><td>600</td><td>315</td><td>20.9 ± 1.4</td><td>76.0 ± 2.1</td><td>-26.2 ± 1.1</td></tr><tr><td>2456658.535</td><td>600</td><td>262</td><td>19.2 ± 1.9</td><td>10.3 ± 3.1</td><td>39.2 ± 2.9</td></tr><tr><td>2456658.542</td><td>124</td><td>145</td><td>19.2 ± 2.7</td><td>10.1 ± 1.8</td><td>39.2 ± 1.5</td></tr><tr><td>2456658.679</td><td>300</td><td>215</td><td>20.5 ± 1.8</td><td>2.9 ± 1.8</td><td>46.5 ± 2.5</td></tr><tr><td>2456658.683</td><td>300</td><td>200</td><td>21.8 ± 1.8</td><td>2.8 ± 1.7</td><td>46.5 ± 3.0</td></tr><tr><td>2456658.686</td><td>300</td><td>195</td><td>22.3 ± 1.5</td><td>2.8 ± 1.7</td><td>46.8 ± 2.0</td></tr><tr><td>2456658.690</td><td>300</td><td>178</td><td>20.2 ± 2.0</td><td>2.2 ± 2.1</td><td>47.0 ± 2.4</td></tr><tr><td>2456658.694</td><td>300</td><td>215</td><td>20.9 ± 1.8</td><td>2.3 ± 2.3</td><td>47.2 ± 2.7</td></tr><tr><td>2456658.698</td><td>300</td><td>248</td><td>21.4 ± 1.5</td><td>2.1 ± 2.0</td><td>47.4 ± 2.9</td></tr><tr><td>2456658.701</td><td>300</td><td>219</td><td>20.4 ± 1.7</td><td>1.9 ± 1.8</td><td>47.6 ± 2.1</td></tr><tr><td>2456658.702</td><td>300</td><td>181</td><td>21.3 ± 1.6</td><td>1.6 ± 2.3</td><td>47.6 ± 1.6</td></tr><tr><td>2456658.709</td><td>300</td><td>286</td><td>21.8 ± 1.9</td><td>1.3 ± 1.5</td><td>48.1 ± 2.2</td></tr><tr><td>2456658.713</td><td>300</td><td>251</td><td>21.5 ± 3.3</td><td>1.3 ± 2.3</td><td>48.1 ± 2.4</td></tr><tr><td>2456658.763</td><td>600</td><td>252</td><td>19.7 ± 1.6</td><td>-1.5 ± 4.7</td><td>51.1 ± 2.4</td></tr><tr><td>2456659.628</td><td>600</td><td>267</td><td>18.4 ± 1.6</td><td>-29.3 ± 1.3</td><td>77.8 ± 2.3</td></tr><tr><td>2456659.671</td><td>300</td><td>222</td><td>19.1 ± 1.5</td><td>-29.9 ± 1.2</td><td>78.1 ± 1.7</td></tr><tr><td>2456659.677</td><td>400</td><td>284</td><td>19.6 ± 1.6</td><td>-30.1 ± 1.5</td><td>78.3 ± 1.8</td></tr><tr><td>2456659.682</td><td>300</td><td>233</td><td>19.9 ± 1.9</td><td>-30.0 ± 1.6</td><td>78.5 ± 2.3</td></tr><tr><td>2456659.685</td><td>300</td><td>209</td><td>19.7 ± 1.7</td><td>-29.8 ± 1.3</td><td>78.0 ± 1.6</td></tr><tr><td>2456659.689</td><td>300</td><td>299</td><td>18.8 ± 1.8</td><td>-29.9 ± 1.7</td><td>78.0 ± 1.6</td></tr><tr><td>2456659.693</td><td>300</td><td>245</td><td>17.5 ± 2.9</td><td>-30.3 ± 1.9</td><td>78.1 ± 2.0</td></tr><tr><td>2456659.697</td><td>300</td><td>310</td><td>19.9 ± 1.7</td><td>-30.0 ± 1.1</td><td>78.1 ± 1.3</td></tr><tr><td>2456659.700</td><td>300</td><td>240</td><td>20.5 ± 1.7</td><td>-30.1 ± 2.0</td><td>78.2 ± 1.8</td></tr><tr><td>2456659.704</td><td>300</td><td>222</td><td>20.4 ± 1.6</td><td>-30.1 ± 1.6</td><td>78.1 ± 1.9</td></tr><tr><td>2456659.708</td><td>300</td><td>231</td><td>20.7 ± 2.0</td><td>-30.2 ± 1.3</td><td>78.2 ± 1.3</td></tr><tr><td>2456659.746</td><td>600</td><td>300</td><td>20.5 ± 1.5</td><td>-29.8 ± 1.3</td><td>78.0 ± 1.4</td></tr><tr><td>2456660.614</td><td>600</td><td>249</td><td>20.6 ± 1.6</td><td>-4.4 ± 3.6</td><td>53.2 ± 4.3</td></tr><tr><td>2456660.653</td><td>600</td><td>276</td><td>21.6 ± 1.6</td><td>-2.4 ± 3.1</td><td>51.7 ± 4.1</td></tr><tr><td>2456660.724</td><td>600</td><td>265</td><td>21.5 ± 1.6</td><td>1.1 ± 2.7</td><td>48.1 ± 2.9</td></tr><tr><td>2456660.763</td><td>600</td><td>299</td><td>19.7 ± 1.7</td><td>2.8 ± 2.2</td><td>46.5 ± 3.1</td></tr><tr><td>2456660.801</td><td>600</td><td>250</td><td>21.5 ± 2.0</td><td>5.1 ± 1.9</td><td>43.9 ± 2.6</td></tr></tbody></table>

<table><caption>**Table 4.** uvby photometry obtained from 1992 January to 1994 March. Each measurement is given as a magnitude difference between the variable source ('v', HD 35502) and a non-variable comparison star ('c', HD 35575). This table will appear only in the electronic version of the paper.</caption><thead><tr><th>HJD</th><th>u(v-c)<br>(mag)</th><th>v(v-c)<br>(mag)</th><th>b(v-c)<br>(mag)</th><th>y(v-c)<br>(mag)</th><th>HJD</th><th>u(v-c)<br>(mag)</th><th>v(v-c)<br>(mag)</th><th>b(v-c)<br>(mag)</th><th>y(v-c)<br>(mag)</th></tr></thead><tbody><tr><td>2448648.764</td><td>1.25</td><td>1.04</td><td>0.96</td><td>0.89</td><td>2449291.008</td><td>1.27</td><td>1.06</td><td>0.98</td><td>0.91</td></tr><tr><td>2448657.697</td><td>1.27</td><td>1.06</td><td>0.97</td><td>0.89</td><td>2449292.006</td><td>1.28</td><td>1.08</td><td>0.99</td><td>0.91</td></tr><tr><td>2448676.640</td><td>1.27</td><td>1.07</td><td>0.98</td><td>0.91</td><td>2449294.930</td><td>1.27</td><td>1.07</td><td>0.99</td><td>0.92</td></tr><tr><td>2448678.635</td><td>1.26</td><td>1.05</td><td>0.97</td><td>0.90</td><td>2449297.965</td><td>1.28</td><td>1.08</td><td>0.99</td><td>0.91</td></tr><tr><td>2448681.635</td><td>1.27</td><td>1.07</td><td>0.97</td><td>0.90</td><td>2449311.949</td><td>1.27</td><td>1.06</td><td>0.98</td><td>0.91</td></tr><tr><td>2448687.621</td><td>1.28</td><td>1.07</td><td>0.98</td><td>0.90</td><td>2449312.945</td><td>1.26</td><td>1.06</td><td>0.97</td><td>0.90</td></tr><tr><td>2448692.632</td><td>1.27</td><td>1.08</td><td>0.99</td><td>0.91</td><td>2449318.933</td><td>1.26</td><td>1.05</td><td>0.98</td><td>0.90</td></tr><tr><td>2448695.633</td><td>1.27</td><td>1.06</td><td>0.98</td><td>0.91</td><td>2449324.920</td><td>1.26</td><td>1.06</td><td>0.97</td><td>0.90</td></tr><tr><td>2448696.605</td><td>1.26</td><td>1.06</td><td>0.97</td><td>0.90</td><td>2449325.917</td><td>1.26</td><td>1.05</td><td>0.97</td><td>0.90</td></tr><tr><td>2448697.604</td><td>1.26</td><td>1.06</td><td>0.97</td><td>0.90</td><td>2449328.909</td><td>1.27</td><td>1.06</td><td>0.97</td><td>0.91</td></tr></tbody></table>

of the continuum) and its blending with Hα resulted in both a decrease in precision and a larger scatter in v<sub>r</sub> compared with those values derived from C II λ4267.

The radial velocities were calculated by fitting a rotationally broadened Voigt function to the C II λ4267, C II λ6578, and He I λ6678 lines. The uncertainties were estimated through a bootstrapping analysis involving the set of normalized flux measurements (I/I<sub>c</sub>) spanning each line. A random sample of 61 per cent of

the data points was selected to be removed at each iteration. These points were then replaced by another set that was randomly sampled from I/I<sub>c</sub>. The fitting routine was then repeated on this new data set. 1000 iterations of the bootstrapping routine were carried out and a probability distribution was obtained for each fitting parameter. The uncertainties in each of the fitting parameters were then taken as the 3σ standard deviations associated with each probability distribution.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** ESPaDOnS, Narval, and FEROS Stokes I LSD profiles (right) and dynamic Stokes I LSD profiles (left) generated such that the three spectral components are emphasized. The observations are phased by the A star orbital period of 5.668 66(6) d. The vertical dashed black lines indicate the surface of the B star located at $v = \pm v \sin i = \pm 75$ km s$^{-1}$. The dashed black sinusoids correspond to the fits obtained for $v_r$ of the two A stars.</figcaption></figure>

The value of $v_r$ inferred from C II λ4267 was found to exhibit a median uncertainty of 1.5 km s$^{-1}$ and a standard deviation of 1.4 km s$^{-1}$. Larger uncertainties were derived using He I λ6678 and C II λ6578 ranging from 1 to 53 km s$^{-1}$. Similarly, $v_r$ inferred from C II λ6578 and He I λ6678 yielded larger standard deviations of 9 and 18 km s$^{-1}$, respectively. In the case of He I λ6678, the decrease in precision and increase in scatter relative to the more stable C II λ4267 measurements is likely the result of the intrinsic line variability. No significant $v_r$ variability was detected using C II λ4267 ($\langle v_r \rangle = 21 \pm 2$ km s$^{-1}$), C II λ6578 ($\langle v_r \rangle = 30 \pm 15$ km s$^{-1}$), or He I λ6678 ($\langle v_r \rangle = 30 \pm 8$ km s$^{-1}$).

The radial velocities of the two A stars were calculated from Stokes I profiles produced using the Least Squares Deconvolution (LSD) method (Donati et al. 1997; Kochukhov, Makaganiuk & Piskunov 2010). The LSD line mask used to carry out the procedure was compiled using data taken from the Vienna Atomic Line Database (VALD; Kupka et al. 2000). In order to isolate the A stars from the dominant B star component in the Stokes I LSD profiles, we used a line list associated with an 8000 K star having a surface gravity of log g = 4.0 (cgs) and a microturbulence of $v_{mic} = 0$. Fig. 1 shows the LSD profiles generated using a different line mask in which both the A and B star components are apparent. The majority of the radial velocities were then determined by simultaneously fitting two Gaussians to the sharp components of the Stokes I profiles. In the case of the Narval observation obtained at HJD = 2456224.646, the sharp line profiles were completely blended and the radial velocities were estimated by fitting a single Gaussian and adopting the resultant velocity for both components. The $v_r$ errors were estimated using a 1000 iteration bootstrapping analysis. We note that the contribution of the B star to the Stokes I LSD profiles was generally weak; however, in certain observations, small contributions were present which resulted in small deformations in the continuum between the two A star profiles. In these

cases, the A star line appearing closest to these deformations was more affected than the other A star thereby yielding slightly higher uncertainties in the fitting procedure. The values of the two A stars' radial velocities were found to range from -30.4 to 78.6 km s$^{-1}$ with an average uncertainty of 2.4 km s$^{-1}$.

The spectral characteristics of the two A-type components are nearly identical; therefore, it is not possible to unambiguously attribute a particular line profile in each spectrum to a particular star. Nevertheless, the importance of this ambiguity can be reduced by making simplifying assumptions.

First, we assumed that the two A stars are gravitationally bound and therefore orbit a common centre of mass (having a radial velocity $v_{cm}$) with a period $P_{orb}$. Furthermore, we assumed that the orbits are circular implying that the A star $v_r$ variations are purely sinusoidal and described by

$$v_{r,i}(t) = v_{cm} + K_i \sin(2\pi t / P_{orb} + \phi_i), \quad (1)$$

where $K_i$ and $\phi_i$ are the semi-amplitude and phase shift of the ith A-type component, respectively. The fact that the radial velocities are observed to oscillate symmetrically about a constant average radial velocity of $\langle v_{r,A} \rangle = (v_{r,1} + v_{r,2})/2 = 25 \pm 3$ km s$^{-1}$ suggests that (1) $K_1 = K_2$ and (2) $|\phi_1 - \phi_2| = \pi$. With these assumptions, we applied the following procedure:

(i) define a grid of possible orbital periods;

(ii) define an amplitude and phase shift of the radial velocity variations based on the maximum observed $v_r$ separation;

(iii) for every period, determine which sinusoidal model the blue and red shifted spectral lines must be associated with in order to minimize the residuals.

The two components in each observation were then identified using whichever period returned the minimal residual fit. A traditional period fitting routine (e.g. Lomb–Scargle) could then be applied to the $v_r$ time series of each star separately thereby yielding more precise periods, amplitudes, and phase shifts for each model.

We chose a grid of periods ranging from 0.1 to 10 d in increments of $10^{-5}$ d (~1 s). The amplitudes ($K_1 = K_2$) and phase shifts ($\phi_1 = |\phi_2 - \pi|$) were defined by the maximum $v_r$ separation of 109 km s$^{-1}$ (i.e. phase 0.994 where this phase corresponds to the phase of the B star's maximum longitudinal magnetic field derived in Section 5). The analysis then involves assigning radial velocities to each of the A star components based on a best-fitting period of 5.6687 d. An alternative means of identifying the orbital period uses the fact that the quantity $|v_{r,1} - v_{r,2}|$ varies with a period of $2P_{orb}$, as outlined by Hareter et al. (2008). Applying this method yields a similar value of $P_{orb} = 5.6680(6)$ d.

With the radial velocities of the two A stars correctly assigned to each individual component, a more precise analysis of the binary orbital parameters was carried out using ORBITX, a FORTRAN code later adapted to IDL which determines the best-fitting $P_{orb}$, time of periastron passage (T), eccentricity (e), longitude of the periastron (ω), semi-amplitudes of each component's radial velocities ($K_1$ and $K_2$), and the radial velocity of the centre of mass (γ) (Tokovinin 1992). This calculation yielded $P_{orb} = 5.668$ 66(6) d, $T = 2456$ 658.172 ± 1.652, $e = 0.003^{+0.006}_{-0.003}$, $\omega = 82 \pm 105^{\circ}$, $K_1 = 55.5 \pm 0.4$ km s$^{-1}$, $K_2 = 52.7 \pm 0.4$ km s$^{-1}$, and $\gamma = 26.5 \pm 0.2$ km s$^{-1}$. These results imply a mass ratio of $M_1/M_2 = 1.05 \pm 0.02$, a projected total mass of $(M_1 + M_2)\sin^3 i = 0.186 \pm 0.008$, and a projected semimajor axis of $a \sin i = 0.0564 \pm 0.0008$. These values are listed in Table 5. Fig. 2 shows the radial velocities of the A stars phased by the 5.668 66 d
<table><caption>**Table 5.** Orbital parameters of the A+A binary.</caption><tr><td>P<sub>orb</sub> (d)</td><td>5.668 66(6)</td></tr><tr><td>T</td><td>2456 658.172 ± 1.652</td></tr><tr><td>e</td><td>0.003<sup>+0.006</sup><sub>-0.003</sub></td></tr><tr><td>ω (°)</td><td>82 ± 105</td></tr><tr><td>K<sub>1</sub> (km s<sup>-1</sup>)</td><td>55.5 ± 0.4</td></tr><tr><td>K<sub>2</sub> (km s<sup>-1</sup>)</td><td>52.7 ± 0.4</td></tr><tr><td>γ (km s<sup>-1</sup>)</td><td>26.5 ± 0.2</td></tr><tr><td>M<sub>1</sub>/M<sub>2</sub></td><td>1.05 ± 0.02</td></tr><tr><td>(M<sub>1</sub> + M<sub>2</sub>)sin<sup>3</sup>i (M<sub>☉</sub>)</td><td>0.186 ± 0.008</td></tr><tr><td>a sin i (au)</td><td>0.0564 ± 0.0008</td></tr></table>

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Top: the nearly sinusoidal fits to the radial velocities of the two A stars phased by a period of 5.668 66 d. Bottom: the reduced χ² distribution yielded by the period fitting routine applied to one of the A stars. The 5.668 66 d period is indicated by the red arrow.</figcaption></figure>

orbital period and compared with the radial velocities computed from the orbital solution.

Comparing $\langle v_{r,A} \rangle = 25 \pm 3 \text{ km s}^{-1}$ with the average B star $v_r$ of $20.5 \pm 1.6 \text{ km s}^{-1}$ and noting that no significant variability in $\langle v_{r,A} \rangle$ was detected over the 22 yr observing period implies a very long orbital period of the A binary about the B star. A lower limit for this period is derived in Section 3.2.

## 3.2 SED fitting

Photometric fluxes of HD 35502 have been measured throughout the UV, visible, and near-infrared spectral regions thereby allowing the temperatures and radii of the three stellar components to be constrained. Ultraviolet measurements were previously obtained at four wavelengths – 1565, 1965, 2365, and 2749 Å – by the S2/68 instrument on board the TD1 satellite (Thompson et al. 1978). Photometry spanning the visible spectrum were taken from the Geneva Observatory’s catalogue of U, B, V, B₁, B₂, V₁, and G filters (Rufener 1981). Additionally, infrared observations obtained by 2MASS (J, H, and Ks filters; Skrutskie et al. 2006) and WISE (W1 and W2 filters; Wright et al. 2010) were used. The reported Geneva, 2MASS, and WISE magnitudes were converted to the flux units of erg s⁻¹ cm⁻² Å⁻¹ using the zero-points reported by Rufener & Nicolet (1988), Cohen, Wheaton & Megeath (2003), and Wright et al. (2010).

The reported photometric measurements of HD 35502 include the contributions from each of the three stellar components. This renders a spectral energy distribution (SED) fitting analysis partic-

ularly susceptible to degenerate solutions; however, speckle interferometry measurements obtained by Balega et al. (2012) provide additional photometric constraints on the system. They detected magnitude differences of 1.45 ± 0.02 and 1.21 ± 0.02 mag using filters centred on 5500 Å and 8000 Å, respectively. The sources were reported to have angular separations of 69 ± 1 mas and 68 ± 1 mas in the two filters. The speckle companion is also identified in observations obtained by Horch et al. (2001) with a consistent angular separation of $\rho < 59$ mas.

In conjunction with the distance to HD 35502, the angular separations may be used to determine the associated linear separation. A distance of $d = 430 \pm 120$ pc was inferred from the $2.35 \pm 0.68$ mas Hipparcos parallax (van Leeuwen 2007). However, assuming HD 35502 to be a member of the Orion OB1a subassociation, we inferred a moderately more precise value of 400 ± 80 pc based on the subassociation’s reported average distance modulus of $\langle \text{dm} \rangle = 8.00 \pm 0.46$ mag (Brown, de Geus & de Zeeuw 1994). The projected linear separation between the two speckle sources was then found to be $27 \pm 5$ au. The minimum orbital period of the A star binary system around the B star can then be approximated by assuming upper limit masses of 8 and $3 M_{\odot}$ for the B and A stars, respectively (the actual masses are derived in Section 8). This implies an orbital period of $P_{\text{orb}} \gtrsim 40$ yr, which is consistent with the fact that no significant variations were detected in either the B star radial velocities or the A star binary’s systemic radial velocity.

The observed photometry was fit using ATLAS9 synthetic SEDs generated from the atmospheric models of Castelli & Kurucz (2003). The grid consists of models with effective temperatures ranging from 3.5 to 50.0 K and surface gravities spanning log g = 0.5–5.0 (cgs), as described in detail by Howarth (2011). This grid was linearly interpolated in order to produce models with a uniform temperature and surface gravity resolution of 125 K and 0.01 dex for $T_{\text{eff}} = 5–25$ K and log g = 3.0–4.75 (cgs). All of the SEDs were then multiplied by the transmission functions associated with each of the narrow band filters: TD1 UV (Carnochan 1982), Geneva (Rufener & Nicolet 1988), 2MASS (Cohen et al. 2003), and WISE (Wright et al. 2010).

Modelling the photometry of unresolved multistar systems using synthetic SEDs requires a large number of fitting parameters and therefore the solution is expected to be highly degenerate. The contribution to the total flux from each of the three stellar components depends on, among other factors, their effective temperatures, surface gravities, and radii. In order to reduce the number of solutions, we adopted a solar metallicity and a microturbulence velocity of $v_{\text{mic}} = 0 \text{ km s}^{-1}$. As with many Bp stars, HD 35502’s primary exhibits chemical spots on its surface (see Section 6); however, on average, a solar metallicity may be adopted.

The high-resolution spectra of HD 35502 obtained by Narval, ESPaDOnS, and FEROS suggest that the two cooler A star components are approximately identical in terms of their $T_{\text{eff}}$, log g, and line-broadening parameters (see Section 3.3). If we assume that the two A stars contribute identically to the SED, the number of independent models required in the fitting routine is reduced from three to two thereby resulting in a total of six free parameters: $T_{\text{eff}}$, log g, and the stellar radius, R, for both the B star and the (identical) A stars.

The effective temperature of a star inferred from fitting model SEDs to photometry is highly dependent on the assumed colour excess, $E(B - V)$. Given HD 35502’s probable location within the Orion OB1a association (Landstreet et al. 2007), the extinction caused by gas and dust is expected to be significant. Indeed,
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Marginalized posterior probability distributions returned by the MCMC algorithm that was applied to the SED fitting. Each frame demonstrates the correlations that are apparent between various parameters, while the vertical and horizontal blue lines indicate the value of each parameter associated with the maximum likelihood solution (i.e. $T_B$, $T_A$, $R_B$, and $R_A$). The contours approximately correspond to $1\sigma-3\sigma$ confidence regions in increments of $0.5\sigma$; $\sigma_{T_B}$, $\sigma_{T_A}$, $\sigma_{R_B}$, and $\sigma_{R_A}$ indicate the value of each parameter's $1\sigma$ region.</figcaption></figure>

Sharpless (1952) and Lee (1968) report values (without uncertainties) of $E(B-V)$ of 0.13 and 0.14, respectively. We used the method of Cardelli, Clayton & Mathis (1989) with an adopted to selective total extinction ratio of $R(V) = 3.1$ in order to deredden the observed photometry. Small differences in the resulting best-fitting parameters of $<3$ per cent were found by using an $E(B-V)$ of 0.13 or 0.14. Although we investigated how our analysis was affected by varying the colour excess from 0.0 to 0.2, the final effective temperatures are reported after assuming $E(B-V) = 0.14$.

We found that a Markov Chain Monte Carlo (MCMC) fitting routine provided a suitable means of determining the most probable solution while simultaneously revealing any significant degeneracies. This was carried out by evaluating the likelihood function yielded by a set of randomly selected fitting parameters drawn from a prior probability (see e.g. Wall & Jenkins 2003). For each iteration, the derived likelihood is compared with that produced by the previous iteration. If a new solution is found to yield a higher quality of fit (higher likelihood), these parameters are adopted, otherwise, the previous solution is maintained. In order to broadly sample the solution space, the MCMC algorithm is designed to adopt poorer fitting solutions at random intervals thereby preventing a local (but not global) maximum likelihood from being returned.

Uniform prior probability distributions (flat priors) were defined for $T_{\text{eff}}$, log g, and R, where the latter was constrained within 1.0–10.0 R⊙. The two speckle observations (Balega et al. 2012) were then included in the total prior probability as monochromatic flux ratios (i.e. magnitude differences) at 5500 and 8000 Å. We assumed that the reported 0.02 mag uncertainties correspond to $1\sigma$ significance. The marginalized posterior probability distributions produced after $10^6$ iterations in the Markov Chain are shown in Fig. 3.

The most probable effective temperatures for the B and A star models were found to be $18.4 \pm 1.2$ and $8.9 \pm 0.6$ kK, respectively, where the uncertainties correspond to the 93rd percentile (approximately $2\sigma$). The fitting parameters used to derive the stellar radii, $R_B$

<table><caption>**Table 6.** Stellar parameters of HD 35502.</caption><tbody><tr><td>Sp. type<sup>a</sup></td><td>B5IVsnp+A+A</td></tr><tr><td>π (mas)<sup>b</sup></td><td>2.35 ± 0.68</td></tr><tr><td>(dm) (mag)<sup>c</sup></td><td>8.00 ± 0.46</td></tr><tr><td>d (pc)</td><td>400 ± 80</td></tr><tr><td colspan="2"><strong>Photometry</strong></td></tr><tr><td>V (mag)<sup>d</sup></td><td>7.331 ± 0.004</td></tr><tr><td>E(B-V) (mag)<sup>e</sup></td><td>0.14</td></tr><tr><td colspan="2"><strong>B star parameters</strong></td></tr><tr><td>T<sub>eff</sub> (kK)</td><td>18.4 ± 0.6</td></tr><tr><td>log (g) (cgs)</td><td>4.3 ± 0.2</td></tr><tr><td>v sin i (km s<sup>-1</sup>)</td><td>75 ± 5</td></tr><tr><td>log L/L<sub>⊙</sub></td><td>3.0 ± 0.3</td></tr><tr><td>M/M<sub>⊙</sub></td><td>5.7 ± 0.6</td></tr><tr><td>R<sub>p</sub>/R<sub>⊙</sub></td><td>3.0<sup>+1.1</sup><sub>-0.5</sub></td></tr><tr><td>R<sub>eq</sub>/R<sub>⊙</sub></td><td>3.1<sup>+1.8</sup><sub>-0.4</sub></td></tr><tr><td>τ<sub>age</sub> (Myr)</td><td>20 ± 20</td></tr><tr><td>P<sub>rot</sub> (d)</td><td>0.853 807(3)</td></tr><tr><td colspan="2"><strong>A star parameters</strong></td></tr><tr><td>T<sub>eff</sub> (kK)</td><td>8.9 ± 0.3</td></tr><tr><td>log (g) (cgs)</td><td>4.3 ± 0.3</td></tr><tr><td>v sin i (km s<sup>-1</sup>)</td><td>12 ± 2</td></tr><tr><td>log L/L<sub>⊙</sub></td><td>1.4 ± 0.3</td></tr><tr><td>M/M<sub>⊙</sub></td><td>2.1 ± 0.2</td></tr><tr><td>R/R<sub>⊙</sub></td><td>2.1 ± 0.4</td></tr><tr><td>τ<sub>age</sub> (Myr)</td><td>&lt;630</td></tr></tbody></table>

Table references: <sup>a</sup>Abt & Levato (1977), <sup>b</sup>van Leeuwen (2007), <sup>c</sup>Brown et al. (1994), <sup>d</sup>Rufener (1981), <sup>e</sup>Lee (1968).

and $R_A$, depend on the distance to HD 35502 (i.e. as a scaling factor given by $R_*^2/d^2$). Although the posterior probability distributions for $R_B$ and $R_A$ both yield $2\sigma$ uncertainties of $0.2 R_⊙$, the consideration of the relatively large distance uncertainty ($d = 400 \pm 80$ pc) implies larger uncertainties of $\delta R_B = 0.6 R_⊙$ and $\delta R_A = 0.4 R_⊙$. The most probable radii and their uncertainties found from the MCMC analysis are then given by $R_B = 3.0 \pm 0.6 R_⊙$ and $R_A = 2.1 \pm 0.4 R_⊙$. The derived temperatures and stellar radii are listed in Table 6. The analysis was insensitive to changes in log g as indicated by an essentially flat posterior probability distribution; therefore, no definitive surface gravity can be reported. Comparisons between the observed photometry and the best-fitting model are shown in Fig. 4, where we have adopted log g = 4.3 for both the A and B models as derived in Section 3.3.

The model B star flux ($F_B$) and the model binary A star flux ($2F_A$) can be used to verify our initial assumption that the two components detected in the speckle observations do indeed correspond to the central B star and the A star binary system. The models can be compared with the speckle observations by calculating the flux ratios, $2F_A(\lambda)/F_B(\lambda)$, at the speckle observation wavelengths. Both $F_A(\lambda)$ and $F_B(\lambda)$ are integrated over wavelength intervals of 200 and 1000 Å (i.e. the FWHM of the filters used by Balega et al. 2012) centred at 5500 and 8000 Å, respectively. We then obtain magnitude differences of $\Delta m_{\text{syn}}(5500\AA) = 1.46$ and $\Delta m_{\text{syn}}(8000\AA) = 1.23$. These values yield a negligible discrepancy with the speckle observations of 1 per cent at $\lambda = 5500\AA$ and 2 per cent at $\lambda = 8000\AA$.

## 3.3 Spectral line fitting

Several properties of HD 35502's three stellar components may be estimated through comparisons with synthetic spectra (e.g. the
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Comparisons between the observed photometry (red points) and the model SED (solid black curve). The dashed blue and dot-dashed black curves correspond to the model SEDs of the composite A star components (i.e. $2F_A(\lambda)$) and the B star, respectively. The black crosses indicate the flux obtained by multiplying the model SED by the transmission function of the associated filter.</figcaption></figure>

surface gravity, line broadening characteristics, etc.). We carried this out using local thermodynamic equilibrium (LTE) models generated with SYNTH3 (Kochukhov 2007). The code computes disc-integrated spectra using spectral line data provided by VALD (Kupka et al. 2000) obtained using an extract stellar request for a specified effective temperature, surface gravity, and microturbulence velocity in conjunction with ATLAS9 atmospheric models (Kurucz 1993). The synthetic spectra can then be convolved with the appropriate functions in order to account for instrumental and rotational broadening effects.

The ESPaDOnS, Narval, and FEROS observations were normalized using a series of polynomial fits to the continuum. The relatively shallow (≈5 per cent of the continuum) and narrow lines produced by the two A stars made the spectral line modelling inherently uncertain. For instance, the typical rms of the continuum near the A stars' Mg I λ4703 lines was found to be approximately 14 per cent of the line depth. Thus, the SNRs of the majority of the A star lines were relatively low. This was mitigated to some extent by binning the observed spectra with a bin width of ≈0.03 Å (i.e. 2 pixels). In order to account for the instrumental profile of the ESPaDOnS and Narval observations, the synthetic spectra were convolved with a Gaussian function assuming a resolving power of R = 65 000; similarly, the FEROS spectra were fit after convolving the synthetic spectra assuming R = 48 000.

The quality of fit yielded by the total normalized synthetic spectrum ($F_{\text{tot}}$) depends not only on $T_{\text{eff}}$ of the three models but also on their (relative) luminosities: $F_{\text{tot}} = (\sum_i L_i F_i)/\sum_i L_i$, where $L_i$ and $F_i$ are the luminosities and synthetic spectra of the ith component. We adopted the 18.4 and 8.9 kK values associated with the B and two A stars obtained from the SED fitting (Section 3.2). Moreover, we assumed that the luminosities of the two A stars are equal, thereby reducing the number of degrees of freedom in the spectral line fitting analysis.

With $T_{\text{eff}}$ specified, the stellar luminosities can be estimated through various methods. We found that the best results were obtained by letting the luminosity ratio of the B and A star models, $L_B/L_A$, be a free parameter and subsequently finding the minimum $\chi^2$ fit for a given surface gravity (log g) and rotational broadening ($v \sin i$). This was carried out using the observed spectra for

which the two A stars were most widely separated in wavelength (phase 0.994 in Fig. 2) in the wavelength range of 4200–4300 Å. This region was chosen because of the presence of the strong and essentially non-variable C II line produced by the B star along with many A star lines of various elements (e.g. Fe, Ti, Cr, Mn). Most importantly, this wavelength range is free of B star lines exhibiting obvious chemical abundance anomalies and variability such as those observed from He. A subsample of this region containing C II λ4267 is shown in Fig. 6. This technique yielded $L_B/L_A = 5.2$, which is consistent with the median value implied by the SED fitting analysis – within 4000 Å ≤ λ ≤ 6000 Å – of $8.4_{-2.4}^{+4.1}$. Ultimately, using an $L_B/L_A$ of 5.2 instead of 8.4 produced a moderate increase of $1.1\sigma$ in the best-fitting solution's overall quality of fit.

log g and $v \sin i$ of HD 35502's three components were then fit using various B and A star lines while recalculating the best-fitting $L_B/L_A$ for every change in the parameters. As a result of the presumed chemical peculiarities and line variability of the B star (see Section 6) and the limited number of lines, log g$_B$ could not be reliably constrained using He or metal lines (e.g. Mg I and Mg II). Instead, we relied upon the wings of the strong and broad Balmer lines. In particular, the observations of Hβ, Hγ, and H δ obtained at a phase of 0.49 in the B star's rotational period (the phase of minimum emission) were used in order to minimize the effects of emission. Several metal lines were used to constrain $v_{\text{B}} \sin i$ such as C II λ4267, S II λ5640, and Fe II λ5780. log g$_B$ and $v_{\text{B}} \sin i$ were found to be $4.3 \pm 0.2$ (cgs) and $75 \pm 5$ km s$^{-1}$, respectively. Examples of the best-fitting model spectra are shown in Fig. 5; the adopted range in $v_{\text{B}} \sin i$ is shown in Fig. 6.

The surface gravity and rotational broadening of the two A stars were fit simultaneously using several Fe II and Mg I lines. Their best-fitting log g and $v \sin i$ values were found to be $4.3 \pm 0.3$ (cgs) and $12 \pm 2$ km s$^{-1}$, respectively. Two examples of the modelled A star lines are shown in Fig. 7.

## 3.4 Hertzsprung–Russell diagram

The masses, ages, and polar radii of HD 35502's three stellar components may be estimated by comparing their positions on the Hertzsprung–Russell diagram (HRD) with theoretical isochrones.

In order to determine the B star's luminosity, we used the $18.4 \pm 0.6$ kK effective temperature and $3.0 \pm 0.6$ R$_{\odot}$ radius derived from the three star SED fit discussed in Section 3.2. The Stefan–Boltzmann law then yields a luminosity of $\log L/L_{\odot} = 3.0_{-0.5}^{+0.4}$. Similarly, the position of the A stars on the HRD can be identified using $T_{\text{eff}} = 8.9 \pm 0.3$ kK and $R_A = 2.1 \pm 0.4$ R$_{\odot}$. We calculate an A star luminosity of $\log L/L_{\odot} = 1.4 \pm 0.3$. The HRD positions of the B and A stars are shown in Fig. 8.

The masses (M) and polar radii (R$_p$) associated with a given $T_{\text{eff}}$ and log L/L$_{\odot}$ were determined using a grid of Geneva model isochrones generated by Ekström et al. (2012). The grid is calculated for the evolutionary time-scale beginning with the zero-age main sequence up until the core carbon-burning phase for masses of 0.8–120 M$_{\odot}$. The microturbulence velocity was fixed at $v_{\text{mic}} = 0.0$ km s$^{-1}$ and a solar metallicity of Z = 0.014 was assumed. In the case of HD 35502's central B star, the ratio of the angular velocity to the critical angular velocity, $\Omega/\Omega_c$, is known to be significant based on the 0.853 807(3) d rotational period (see Section 4). Its position on the HRD was therefore compared against several additional grids calculated using $\Omega/\Omega_c = 0.4–0.9$ in increments of 0.1 (Georgy et al. 2013). While no significant difference
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Top: comparisons between the best-fitting synthetic (red) and observed (black) H lines, Hγ (left) and Hβ (right). The filled blue region indicates the total uncertainty associated with log g<sub>B</sub> and log g<sub>A</sub>. The observations occur at phase 0.49 in the B star's rotational period (see Section 4) when the weakest emission is visible in the wings of Hβ. Bottom: the residuals associated with the model spectrum (red).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 6.** A subsample of the wavelength region used to determine the luminosity ratio, $L_B/L_A$. The red curve corresponds to the best-fitting model spectrum; the filled blue region indicates the total uncertainty associated with a $v_B \sin i$ of ±5 km s⁻¹ as determined by fitting various metal lines. The black curve shows the observed spectrum. The Fe I λ4273 lines of the two A star components are labelled as 'A₁' and 'A₂'.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 7.** Comparisons between best-fitting model (red) and observed (black points) Fe II (left) and Mg I (right) lines used to constrain the surface gravity and rotational broadening of the two A stars; the filled blue region indicates the total uncertainties associated with $v_A \sin i$ and log g<sub>A</sub>. The phase corresponds to the maximum observed separation between the two binary components.</figcaption></figure>

in the inferred *M* was apparent (i.e. <4 per cent), *R*<sub>p</sub> was found to decrease by as much as 15 per cent.

In order to select the most accurate grid of isochrones and thus, the most accurate $R_p$, $\Omega/\Omega_c$ must first be estimated. Since $\Omega_c$ depends on both the mass and polar radius, it was calculated using the parameters derived from each grid of isochrones. Using $P_{\text{rot}}$ inferred in Section 4 to determine $\Omega$, a range of $\Omega/\Omega_c$ values were found. The appropriate grid was then chosen based on whichever $\Omega/\Omega_c$ most closely agreed with the $\Omega/\Omega_c$ associated with the isochrone grid. A calculated $\Omega/\Omega_c$ of 0.53 yielded the best agreement; we found $R_{p,B} = 3.0_{-0.5}^{+1.1} R_{\odot}$ and $M_B = 5.7 \pm 0.6 M_{\odot}$ using the $\Omega/\Omega_c = 0.5$ isochrones. These results imply an equatorial radius of $R_{\text{eq,B}} = 3.2_{-0.6}^{+1.6} R_{\odot}$. Using von Zeipel's law (von Zeipel 1924), we estimate that the ratio of $T_{\text{eff}}$ between the pole and the equator is approximately 1.02.

The two A stars' *M* and $R_p$ were inferred using the $\Omega/\Omega_c = 0.0$ isochrone grid. The $T_{\text{eff}}$ and $\log L/L_{\odot}$ derived from the three star SED fit then yielded $R_{p,A} = 2.0_{-0.5}^{+0.8} R_{\odot}$ and $M_A = 2.1 \pm 0.2 M_{\odot}$.

We note that both $R_{p,B}$ and $R_{p,A}$ are consistent with $R_B$ and $R_A$ derived in Section 3.2.

## 4 ROTATIONAL PERIOD

Several observed properties of HD 35502 exhibit periodic variability with varying significance. In order to correctly interpret the origin of these variations, it is crucial to identify the periods, the phases at which the maxima and minima occur, and amplitudes with which they occur. This was carried out using the same procedure discussed in Section 3.1 in which the orbital solution of the A star binary was derived.

We first assumed a sinusoidal fit to the data, $f(t)$, given by $f(t) = C_0 + C_1 \sin(2\pi[t - t_0]/P + C_2)$, where P is the period of variability, t is the observation's HJD, $t_0$ is the HJD corresponding to phase 0, and $C_0, C_1$, and $C_2$ are fitting parameters. A $\chi^2$ distribution was then generated using periods ranging from 0.1 to 10.0 d in increments ~1 s. The best-fitting period was inferred from the minimal $\chi^2$ solution and the $3\sigma$ $\chi^2$ interval was taken as the
<figure><img src="image_9.png" /><figcaption>Figure 8. The positions of HD 35502's A and B star components are indicated by the blue diamond and black square. The evolutionary tracks (black dot-dashed lines) and isochrones (red dotted lines) assume a non-rotating star of solar metallicity (Ekström et al. 2012).</figcaption></figure>

associated uncertainty. The uncertainties in the three fitting parameters were estimated using a 1000 iteration bootstrapping analysis. The statistical significance of each derived period was evaluated by comparing the quality of the sinusoidal fit to that yielded by a constant fitting function given by $f(t) = C_0$, where $C_0$ is a time-independent fitting parameter. The difference between the minimal $\chi^2$ values associated with the constant fit ($\chi^2_{\text{const}}$) and sinusoidal fit ($\chi^2_{\text{sin}}$) were then calculated. Any sinusoidal fit having $\chi^2_{\text{const}} - \chi^2_{\text{sin}} \ge 3\sigma$ was considered to be statistically significant.

Various periods were found when this procedure was applied to the longitudinal field measurements ($\langle B_z \rangle$, see Section 5) along with the multiple photometry and equivalent width (EW) measurements (see Section 6). The analyses of nearly all data sets yielded statistically significant variability, with over half corresponding to a unique period near 0.85 d. They were found to be equal to one another within ≈10 s with typical uncertainties ≲10 s and were therefore averaged to obtain a period of 0.853 82(7) d. However, when the Hα EWs were phased with this period, the oldest measurements showed an ≈0.1 phase offset relative to the more recent measurements. This discrepancy was resolved by adopting the best-fitting ephemeris derived from the Hα EWs of

$$JD = 2456295.812\,850 \pm 0.853807(3) \times E, \quad (2)$$

where the reference JD (2456295.812850(3)) corresponds to the epoch of $\langle B_z \rangle$ maximum magnitude. Therefore, while the general accuracy of the rotational period is established by the diverse photometric, spectroscopic, and magnetic data sets, the adopted value and its precision correspond to those implied by the Hα EWs.

The periodic variability of $\langle B_z \rangle$ can be explained, in part, as a consequence of a stable oblique magnetic field configuration that is modulated by the star's rotation. Similarly, rotationally modulated variations exhibited by the EWs of various spectral lines can be produced by at least two mechanisms: (1) non-uniform distributions of chemicals on the stellar surface and (2) hot plasma accumulating in the star's magnetosphere resulting in emission and absorption. All

of these phenomena are commonly exhibited by magnetic B-type stars (e.g. Landstreet & Borra 1978; Leone et al. 2010; Bohlender & Monin 2011). Therefore, we conclude that the ephemeris given by equation (2) is the B star's rotational period.

## 5 MAGNETIC FIELD

Zeeman signatures produced by a magnetic star in the HD 35502 system were detected in circularly polarized (Stokes V) ESPaDOnS and Narval observations. They were found to be coincident with the B star's spectral lines regardless of the inferred velocities of the two A stars. We therefore assumed the detected field to be entirely produced by the B star.

The Stokes V Zeeman signature associated with Hβ yielded 18 definite detections (DDs), 1 marginal detection (MD), and 7 non-detections (NDs) based on the detection criterion outlined by Donati et al. (1997). The SNRs of the observed signatures were optimized using the LSD procedure (Donati et al. 1997; Kochukhov et al. 2010) introduced in Section 3.1. A master line mask containing He and metal lines was generated using data obtained from VALD (Kupka et al. 2000) with a specified $T_{\text{eff}}$, log g, and microturbulence velocity ($v_{\text{mic}}$). All Balmer lines were also removed along with any regions affected by atmospheric absorption (i.e. telluric lines). Several single element line masks were subsequently generated from the He+metal mask by retaining only specific chemical elements including He, C, Si, Fe, and Mg. Clearly, the magnitudes of $\langle B_z \rangle$ derived using different elements will be affected by any non-uniform distribution of chemicals across the star's surface. Therefore, our analysis also includes measurements obtained using H lines (both from LSD profiles and Hβ) which do not typically exhibit non-solar abundances or non-homogeneous surface distributions (i.e. chemical spots). A H line mask was generated in which the H lines exhibiting moderate emission (e.g. Hβ and Hα) were removed. The resultant mask contained three H I lines: H I λ3970, H I λ4102, and H I λ4340.

Two approaches were used to isolate the B star lines from the A star lines using LSD. The first method used a mask generated with $T_{\text{eff}} = 15$ kK, log g = 4.0 (cgs), and $v_{\text{mic}} = 0$ km s⁻¹ which yielded LSD profiles with clear Stokes I contributions from all three stellar components. The narrow line components associated with the two A stars were then fit by Gaussian functions which were subsequently subtracted from the Stokes I profiles. We found that this method could not be consistently applied to all observations. Moreover, the quality of the Gaussian fits was dramatically reduced when applied to the C, Si, and Fe line masks.

The second method used the same log g and $v_{\text{mic}}$ with a significantly higher temperature of $T_{\text{eff}} = 26$ kK. This yielded LSD profiles with minimal contributions from the two A stars. The Stokes I and V LSD profiles generated using the $T_{\text{eff}} = 26$ kK line mask for H, He+metal, He, C, Si, and Fe are shown in Fig. 9, phased according to equation (2). Aside from the Si and Fe LSD profiles, no strong contributions from the A star Stokes I profiles can be discerned.

Each ESPaDOnS and Narval spectropolarimetric observation includes a diagnostic null which may be used to evaluate the significance of any polarized signal. No spurious signals were detected in any of the diagnostic null profiles.

$\langle B_z \rangle$ was inferred from each of the Stokes I and V LSD profiles, as well as from Hβ, using equation 1 of Wade et al. (2000). We used a wavelength of 500 nm with a Landé factor of 1.2 for the He and metal mask measurements and a Landé factor of unity for the H mask and Hβ measurements. The Doppler shift produced by the B star's radial velocity of ≈20 km s⁻¹ was subtracted from each LSD
<figure><img src="image_10.png" /><figcaption>**Figure 9.** H, He+metal, He, C, Si, and Fe Stokes I and V LSD profiles generated from a $T_{\text{eff}} = 26$ kK, $\log g = 4.0$ (cgs) line mask. The profiles are all phased by the B star's rotational period of 0.853 807 d.</figcaption></figure>

profile. The Stokes I and V profiles were then normalized to the continuum intensity at a velocity of $v = -125$ km s$^{-1}$, where the average Stokes V intensity is approximately zero. An integration range of $v \in [-110, 110]$ km s$^{-1}$ was then used in the calculation of $\langle B_z \rangle$ for each of the LSD profiles (i.e. H, He+metal, He, C, Si, and Fe); for the H LSD profiles and Hβ line, this integration range corresponds to the width of the Doppler core. The values of $\langle B_z \rangle_{H\beta}$ inferred from the Narval, ESPaDOnS, and dimaPol observations, along with the status of their detections, are listed in Tables 1 and 2. $\langle B_z \rangle$ derived from the H, He, and metal LSD profiles are listed in Table 7.

High-resolution spectropolarimetry is essentially insensitive to the polarization in the wings of the Balmer lines. As an example of how the Doppler cores of Balmer lines may be used to infer $\langle B_z \rangle$, see fig. 2 of Landstreet et al. (2015), who explain the method in some detail. Our measurements obtained in the con-

text of the current paper, as well as those obtained by Sikora et al. (2015), demonstrate that this method results in longitudinal field intensity and variability in good agreement with other approaches.

All of the $\langle B_z \rangle$ measurements were found to exhibit statistically significant variations with best-fitting periods ranging from 0.853 80 to 0.853 89 d. Only the $\langle B_z \rangle$ values obtained using the Fe LSD profiles, ($\langle B_z \rangle_{\text{Fe}}$) yielded more than one period. Figs 10 and 11 show the $\langle B_z \rangle$ measurements inferred from Hβ and the LSD profiles, respectively, phased by the B star's rotational period (equation 2). It is clear that the scatter of the $\langle B_z \rangle_{\text{He}}$, $\langle B_z \rangle_{\text{Si}}$, and $\langle B_z \rangle_{\text{Fe}}$ measurements is significantly larger than that yielded by $\langle B_z \rangle_{\text{H}}$, $\langle B_z \rangle_{\text{He+metal}}$, and $\langle B_z \rangle_{\text{C}}$. This is likely caused by the presence of He, Si, and Fe chemical spots which are commonly observed on the surfaces of Bp stars. As discussed in Section 6, we find strong evidence for He and Si spots.
<table><caption>**Table 7.** Longitudinal field measurements derived using the LSD profiles that were generated from the ESPaDOnS and Narval spectra. Columns 3–6 list $\langle B_z \rangle$ derived from H, He+metal, He, C, Si, and Fe LSD profiles. This table will appear only in the electronic version of the paper.</caption><thead><tr><th>HJD</th><th>Instrument</th><th>⟨B<sub>z</sub>⟩<sub>H</sub><br>(kG)</th><th>⟨B<sub>z</sub>⟩<sub>He + metal</sub><br>(kG)</th><th>⟨B<sub>z</sub>⟩<sub>He</sub><br>(kG)</th><th>⟨B<sub>z</sub>⟩<sub>c</sub><br>(kG)</th><th>⟨B<sub>z</sub>⟩<sub>Si</sub><br>(kG)</th><th>⟨B<sub>z</sub>⟩<sub>Fe</sub><br>(kG)</th></tr></thead><tbody><tr><td>2454702.138</td><td>ESPaDOnS</td><td>0.39 ± 0.30</td><td>-0.27 ± 0.19</td><td>-0.53 ± 0.35</td><td>-0.52 ± 0.60</td><td>0.32 ± 0.38</td><td>0.09 ± 0.57</td></tr><tr><td>2455849.677</td><td>Narval</td><td>-0.14 ± 0.26</td><td>0.12 ± 0.22</td><td>-0.96 ± 0.41</td><td>0.88 ± 0.60</td><td>1.91 ± 0.93</td><td>1.44 ± 0.65</td></tr><tr><td>2455893.623</td><td>Narval</td><td>-2.75 ± 0.32</td><td>-2.96 ± 0.27</td><td>-3.47 ± 0.48</td><td>-2.29 ± 0.80</td><td>-2.53 ± 0.61</td><td>-0.56 ± 0.88</td></tr><tr><td>2455910.518</td><td>Narval</td><td>-1.26 ± 0.74</td><td>-1.50 ± 0.55</td><td>-1.63 ± 1.04</td><td>-2.27 ± 1.60</td><td>0.17 ± 0.88</td><td>-1.15 ± 2.26</td></tr><tr><td>2455934.528</td><td>Narval</td><td>-2.27 ± 0.28</td><td>-2.86 ± 0.29</td><td>-3.51 ± 0.56</td><td>-2.64 ± 0.73</td><td>-6.63 ± 1.57</td><td>-0.57 ± 0.94</td></tr><tr><td>2455936.534</td><td>Narval</td><td>-2.04 ± 0.35</td><td>-2.21 ± 0.29</td><td>-2.93 ± 0.51</td><td>-1.27 ± 0.90</td><td>-3.25 ± 1.13</td><td>-2.82 ± 0.77</td></tr><tr><td>2455938.525</td><td>Narval</td><td>0.48 ± 0.32</td><td>-0.08 ± 0.25</td><td>-0.42 ± 0.45</td><td>0.11 ± 0.73</td><td>0.90 ± 0.54</td><td>1.83 ± 0.85</td></tr><tr><td>2455944.500</td><td>Narval</td><td>0.25 ± 0.37</td><td>0.30 ± 0.28</td><td>0.25 ± 0.53</td><td>1.05 ± 0.77</td><td>1.30 ± 0.54</td><td>1.94 ± 1.01</td></tr><tr><td>2455949.429</td><td>Narval</td><td>-0.89 ± 0.28</td><td>-1.09 ± 0.24</td><td>-1.69 ± 0.44</td><td>-1.36 ± 0.81</td><td>-1.08 ± 0.51</td><td>-1.19 ± 0.81</td></tr><tr><td>2455950.472</td><td>Narval</td><td>0.41 ± 0.40</td><td>-0.10 ± 0.31</td><td>-0.86 ± 0.60</td><td>0.03 ± 0.96</td><td>1.57 ± 0.69</td><td>0.51 ± 0.80</td></tr></tbody></table>

<figure><img src="image_11.png" /><figcaption>**Figure 10.** Hβ longitudinal magnetic field measurements phased by the rotational period (equation 2). The measurements obtained from Narval (open circles), ESPaDOnS (filled circles), and dimaPol (filled triangles) are shown along with the best-fitting sinusoid (solid black).</figcaption></figure>

The mean and amplitude of the phased $\langle B_z \rangle$ measurements are defined by the fitting parameters $B_0$ and $B_1$ associated with the sinusoidal fitting function $\langle B_z \rangle = B_0 + B_1 \sin(2\pi\theta + \phi)$, where $\theta$ is the phase calculated using equation (2) and $\phi$ is the phase shift. The most precise $B_0$ and $B_1$ values – as indicated by the uncertainties estimated using a 1000 iteration bootstrapping analysis – were derived using Hβ, and the H and He+metal LSD profiles. They were found to be consistent within $2\sigma$. The lowest uncertainties were obtained from the $\langle B_z \rangle_H$ measurements, which yielded a mean and amplitude of $B_0 = -1.41 \pm 0.11$ kG and $B_1 = 1.64 \pm 0.16$ kG where the uncertainties correspond to $3\sigma$.

If we assume that the field is characterized by an important dipole component, the sinusoidal variations in $\langle B_z \rangle$ imply that the dipole's axis of symmetry is inclined (i.e. has an obliquity angle $\beta$) with respect to the star's rotational axis. This interpretation, first described by Stibbs (1950), is known as the Oblique Rotator Model (ORM). Under the assumptions of the ORM, $\beta$ can be calculated from equation 3 of Preston (1967) which depends on $r \equiv |\langle B_z \rangle|_{\text{min}}/|\langle B_z \rangle|_{\text{max}}$ and the inclination angle, $i$, of the star's axis of rotation. The

<figure><img src="image_12.png" /><figcaption>**Figure 11.** Longitudinal magnetic field measurements phased by the rotational period (equation 2). The measurements obtained from the H, He+metal, He, C, Si, and Fe line masks are shown along with their best-fitting sinusoids (solid black). The fit to $\langle B_z \rangle_H$ is shown as dashed black curves in the other panels.</figcaption></figure>
value of *i* can be determined using $P_{\text{rot}}$ given by equation (2), $v \sin i = 75 \pm 5$ km s⁻¹ derived in Section 3.3, and $R_{\text{eq}} = 3.2_{-0.6}^{+1.6} R_{\odot}$ listed in Table 6. We obtained a value of $i = 24^{+8}_{-9}^\circ$. The value of *r* was determined from $B_0$ and $B_1$. Using the values inferred from the $\langle B_z \rangle_{\text{H}}$ measurements, we obtained $r = 0.08_{-0.07}^{+0.09}$. Finally, the obliquity angle was found to be $\beta = 63^\circ \pm 13^\circ$ using equation 3 of Preston (1967).

In addition to the obliquity, the strength of the magnetic field's dipole component, $B_p$, can be calculated by inverting equation 1 of Preston (1967) and letting $t=0$ correspond to $\langle B_z \rangle_{\text{max}}$. We used a linear limb darkening constant that was averaged over the values derived by van Hamme (1993) for the U, B, V, R, and I bandpasses. These specific filters were selected because of the approximate correspondence with the ESPaDOnS and Narval wavelength range. A value of $u = 0.265$ was obtained after interpolating the published table for an effective temperature and surface gravity of 18.4 kG and 4.3 (cgs). *i*, *β*, and *u* then yield $B_p = 14_{-3}^{+9}$ kG. Similar obliquity angles and dipolar field strengths are derived using $H\beta$ along with the He+metal, He, C, and Si LSD profiles. $\langle B_z \rangle_{\text{Fe}}$ exhibits significantly weaker and more uncertain values of $\beta = 33_{-27}^{+36}^\circ$ and $B_p = 8_{-3}^{+6}$.

# 6 EMISSION AND VARIABILITY

Hot magnetic B-type stars are commonly found to exhibit spectral line variability either as a result of chemical spots (e.g. Kochukhov et al. 2015; Yakunin et al. 2015) or from the presence of a hot plasma beyond the stellar surface (e.g. Landstreet & Borra 1978). Furthermore, photometric variability correlated with both of these phenomena, as well as with strong, coherent magnetic fields has been previously reported (e.g. Shore et al. 1990; Oksala et al. 2010).

Along with the *uvby* photometric measurements listed in Table 4, we also analysed *Hipparcos* Epoch Photometry for variability. The catalogue (Perryman et al. 1997) contains 98 observations of HD 35502 which were obtained over a period of 3.1 yr. Three of these measurements have multiple quality flags reported and were therefore removed from our analysis. The remaining measurements have an average of 7.331 mag, a standard deviation of 0.011 mag, and an average uncertainty of 0.009 mag.

The period searching routine described in Section 4 was applied to both the *uvby* and *Hipparcos* data sets. All of the *uvby* measurements were found to exhibit statistically significant variability; however, only *u*(v-c) and *v*(v-c) yielded unique periods of 0.8537(3) d. The analysis of the *Hipparcos* magnitudes ($H_p$) resulted in a best-fitting period of 0.8630(2) d along with five other statistically significant periods ranging from 0.46 to 1.7 d. Fig. 12 shows the sinusoidal fits to *u*(v-c), *v*(v-c), *b*(v-c), *y*(v-c), and $H_p$ obtained when phased by the B star's 0.853 807 d rotational period given by equation (2).

The variability of the spectral lines associated with HD 35502's central B star is most easily detected by calculating EWs. We carried this out for a number of lines for which no significant absorption produced by the two A stars was evident. This included He I, C II, and Si III lines. The EWs of the He and metal lines were calculated using integration ranges of [−100, 100] km s⁻¹ and were normalized to the continuum just outside these limits. The Balmer line EWs ($H\alpha$, $H\beta$, and $H\gamma$) were measured by normalizing to the flux at $|v| \gtrsim 700$ km s⁻¹ and integrating over a velocity range of [−600, 600] km s⁻¹. The uncertainties in the EW measurements were then estimated using a bootstrapping analysis with 1000 iterations. All of the calculated EWs and uncertainties are listed in and Tables 8 and 9.

<figure><img src="image_13.png" /><figcaption>**Figure 12.** Photometric observations obtained using various filters and phased by the B star's rotational period of 0.853 807 d. The top four panels show the *u*(v-c), *v*(v-c), *b*(v-c), and *y*(v-c) magnitude differences between HD 35502 ('v') and the non-variable comparison star, HD 35575 ('c'). The bottom panel shows the *Hipparcos* Epoch Photometry measurements. The horizontal dotted lines indicate the best constant fit to the data.</figcaption></figure>
<table><caption>**Table 8.** EWs of Hα calculated using all of the DAO, CFHT f/8.2 spectrograph, Narval, ESPaDOnS, and FEROS spectra. The uncertainties correspond to 3σ. This table will appear only in the electronic version of the paper.</caption><thead><tr><th>HJD</th><th>Instrument</th><th>EW<sub>Hα</sub><br>(Å)</th><th>HJD</th><th>Instrument</th><th>EW<sub>Hα</sub><br>(Å)</th></tr></thead><tbody><tr><td>2448529.029</td><td>DAO</td><td>2.77 ± 0.10</td><td>2455910.518</td><td>Narval</td><td>2.66 ± 0.08</td></tr><tr><td>2448582.005</td><td>CFHT f/8.2</td><td>3.11 ± 0.09</td><td>2455934.528</td><td>Narval</td><td>2.28 ± 0.10</td></tr><tr><td>2448582.103</td><td>CFHT f/8.2</td><td>3.09 ± 0.09</td><td>2455936.534</td><td>Narval</td><td>2.49 ± 0.08</td></tr><tr><td>2448582.901</td><td>CFHT f/8.2</td><td>3.09 ± 0.09</td><td>2455938.525</td><td>Narval</td><td>3.14 ± 0.09</td></tr><tr><td>2448582.999</td><td>CFHT f/8.2</td><td>3.05 ± 0.09</td><td>2455944.500</td><td>Narval</td><td>3.16 ± 0.09</td></tr><tr><td>2448583.112</td><td>CFHT f/8.2</td><td>2.64 ± 0.08</td><td>2455949.429</td><td>Narval</td><td>2.84 ± 0.10</td></tr><tr><td>2448597.819</td><td>DAO</td><td>1.81 ± 0.10</td><td>2455950.472</td><td>Narval</td><td>3.11 ± 0.10</td></tr><tr><td>2449592.079</td><td>CFHT f/8.2</td><td>3.10 ± 0.10</td><td>2455951.471</td><td>Narval</td><td>2.77 ± 0.08</td></tr><tr><td>2449592.086</td><td>CFHT f/8.2</td><td>3.04 ± 0.10</td><td>2455953.653</td><td>DAO</td><td>2.66 ± 0.09</td></tr><tr><td>2449593.086</td><td>CFHT f/8.2</td><td>3.14 ± 0.09</td><td>2455953.663</td><td>DAO</td><td>2.74 ± 0.09</td></tr></tbody></table>

<table><caption>**Table 9.** EWs of various spectral lines calculated using the Narval, ESPaDOnS, and FEROS spectra. The uncertainties correspond to 3σ. This table will appear only in the electronic version of the paper.</caption><thead><tr><th>HJD</th><th>Instrument</th><th>EW<sub>CIIλ4267</sub><br>(× 10<sup>−6</sup>Å)</th><th>EW<sub>HeIλ4713</sub><br>(× 10<sup>−6</sup>Å)</th><th>EW<sub>SiIIIλ4553</sub><br>(× 10<sup>−6</sup>Å)</th><th>EW<sub>NIIλ4643</sub><br>(× 10<sup>−6</sup>Å)</th><th>EW<sub>Hβ</sub><br>(× 10<sup>−6</sup>Å)</th><th>EW<sub>Hγ</sub><br>(× 10<sup>−6</sup>Å)</th></tr></thead><tbody><tr><td>2454702.138</td><td>ESPaDOnS</td><td>1.53 ± 0.07</td><td>1.23 ± 0.05</td><td>0.89 ± 0.02</td><td>0.25 ± 0.01</td><td>59.38 ± 0.70</td><td>56.40 ± 0.67</td></tr><tr><td>2455849.677</td><td>Narval</td><td>1.49 ± 0.07</td><td>1.30 ± 0.05</td><td>0.84 ± 0.03</td><td>0.26 ± 0.01</td><td>58.36 ± 0.71</td><td>56.13 ± 0.69</td></tr><tr><td>2455893.623</td><td>Narval</td><td>1.47 ± 0.07</td><td>1.98 ± 0.06</td><td>0.78 ± 0.04</td><td>0.28 ± 0.02</td><td>56.84 ± 0.72</td><td>56.01 ± 0.69</td></tr><tr><td>2455910.518</td><td>Narval</td><td>1.46 ± 0.07</td><td>1.49 ± 0.05</td><td>1.00 ± 0.04</td><td>0.31 ± 0.02</td><td>57.59 ± 0.69</td><td>55.84 ± 0.67</td></tr><tr><td>2455934.528</td><td>Narval</td><td>1.49 ± 0.07</td><td>1.78 ± 0.06</td><td>0.87 ± 0.04</td><td>0.24 ± 0.02</td><td>57.09 ± 0.69</td><td>55.69 ± 0.65</td></tr><tr><td>2455936.534</td><td>Narval</td><td>1.50 ± 0.07</td><td>1.63 ± 0.06</td><td>0.78 ± 0.03</td><td>0.34 ± 0.02</td><td>57.33 ± 0.66</td><td>55.97 ± 0.65</td></tr><tr><td>2455938.525</td><td>Narval</td><td>1.50 ± 0.07</td><td>1.24 ± 0.04</td><td>0.81 ± 0.03</td><td>0.30 ± 0.01</td><td>58.16 ± 0.70</td><td>55.76 ± 0.67</td></tr><tr><td>2455944.500</td><td>Narval</td><td>1.50 ± 0.07</td><td>1.18 ± 0.05</td><td>0.84 ± 0.03</td><td>0.31 ± 0.01</td><td>58.23 ± 0.70</td><td>55.93 ± 0.67</td></tr><tr><td>2455949.429</td><td>Narval</td><td>1.55 ± 0.07</td><td>1.34 ± 0.05</td><td>0.85 ± 0.03</td><td>0.32 ± 0.02</td><td>58.42 ± 0.70</td><td>55.74 ± 0.67</td></tr><tr><td>2455950.472</td><td>Narval</td><td>1.51 ± 0.07</td><td>1.16 ± 0.05</td><td>0.93 ± 0.03</td><td>0.31 ± 0.02</td><td>58.64 ± 0.70</td><td>55.89 ± 0.67</td></tr></tbody></table>

The contributions of the two A stars to the total measured Balmer line EWs were approximated by comparing the synthetic EWs associated with the SYNTH3 models discussed in Section 3.3. We found that the total synthetic spectrum (including the B star and the two A stars) yielded EWs of 5.6 Å averaged over Hα, Hβ, Hγ, and Hδ. A similar calculation applied to the single B star model yielded average Balmer line EWs of 4.7 Å suggesting that the presence of the two A stars increase the EW measurements by a factor of ≈1.2.

Statistically significant variations were detected from EW measurements of Hα, Hβ, Hγ, He Iλ4713, and Si IIIλ4553. We note that telluric absorption lines were not removed or minimized in the calculation of these EWs. A range of best-fitting periods were derived; however, only the Hα and He Iλ4713 EW measurements yielded unique periods of 0.853 807(3) d and 0.853 77(3) d, respectively. The strongest variability was measured from Hα for which an amplitude of 0.51 ± 0.03 Å was derived. Similar variability – both in terms of the phase of maximum emission and the best-fitting period – was also detected in Hβ although, at a much lower amplitude of 0.08 ± 0.01 Å. The phased Hα and Hβ EWs exhibit a maximum emission at a phase of 0.99 ± 0.04 and 0.0 ± 0.3, respectively, and are therefore in phase with the ⟨Bz⟩ measurements. The He Iλ4713 EWs are approximately in antiphase with respect to the ⟨Bz⟩ variation with minimum absorption occurring at a phase of 0.5 ± 0.1.

Along with EWs, dynamic spectra were also computed by comparing various spectral lines with their respective average normalized intensity (⟨I/I⟩). Both the dynamic spectra and the EWs of C IIλ4267, He Iλ4713, Si IIIλ4553, Hβ, and Hγ are shown in

Fig. 13. It is evident that He Iλ4713, Si IIIλ4553, and to a lesser extent, C IIλ4267, show absorption features crossing from negative to positive velocities. These features suggest the presence of chemical spots on the B star’s surface. The most obvious spot is associated with He I, which exhibits a maximum absorption at a phase of 0.0 ± 0.1 and is therefore coincident with the epoch of maximum ⟨Bz⟩ magnitude. Assuming that the star’s magnetic field consists of a strong dipole component as discussed in Section 5, this result suggests that He is more concentrated near the field’s negative pole. Enhanced He abundances on the surfaces of magnetic Bp stars have been commonly reported to coincide with either the magnetic equator or magnetic poles (e.g. Neiner et al. 2003; Bohlender & Monin 2011; Grunhut et al. 2012; Rivinius et al. 2013).

A similar plot of the dynamic spectrum and EWs of Hα is shown in Fig. 14, where the additional DAO and CFHT spectra are also included. In order to reduce normalization errors, all of the Hα spectra were consistently normalized using a linear fit to the measured flux at velocities of ±600 km s⁻¹. The observed spectra are compared with the synthetic spectrum (Isyn) discussed in Section 3.3 rather than the average observed spectrum. Isyn includes the contributions from the A stars which move (in velocity space) relative to the B star throughout the B star’s rotational period. Therefore, this method results in a greater contrast between the emission and absorption features associated only with the B star. Strong, nearly symmetrical emission peaks are observed at a distance of ≈4R* at a phase of 0.0. The intensity of this emission is observed to decrease by a factor of ≈2 at a phase of 0.5. Similarly, the ratio between the maximum core emission (at phase 0.25) and minimum core emission (at phase 0.5) is also found to be ≈2.
<figure><img src="image_14.png" /><figcaption>**Figure 13.** Left: dynamic spectra of various H I, He I, C II, and Si III lines. A low-pass filter has been applied to the spectrum in order to minimize the A star contributions and tellurics. Each set of observations are compared with the average spectrum (dashed red). Right: measured EWs associated with the lines shown in the dynamic spectra plots. All of the measurements are phased by the B star's rotational period of 0.853 807 d.</figcaption></figure>

The standard interpretation of the broad Hα emission peaks that are associated with a small number of magnetic B-type stars is that they are produced by two dense clouds of hot plasma, trapped in the magnetic field above the stellar surface, which corotate with the star (e.g. Walborn & Hesser 1976; Landstreet & Borra 1978). Under the assumption that the cloud is optically thin, one would expect the same blue shifted emission feature to be observed half a rotational cycle later shifted towards redder wavelengths (e.g. in Fig. 14, the blue emission peak occurring at phase 0.0 should reappear red shifted at phase 0.5). The fact that the strength of both the blue and red shifted emission peaks decrease between phase 0.0 and phase 0.5 suggests that the plasma clouds are, to an extent, optically thick. The relatively large decrease in emission is currently unprecedented amongst the known CM hosting stars; however, a more moderate decrease in HR 5907's Hα emission is shown in fig. 15 of Grunhut et al. (2012).

Adopting the standard interpretation, the trajectories of the Hα-emitting clouds may be approximately inferred by fitting the velocities at which the peak emission is found on either side of the Hα core as a function of rotational phase. The resulting fits suggest that the plasma clouds follow nearly circular trajectories as indicated by

<figure><img src="image_15.png" /><figcaption>**Figure 14.** Top: observed Hα spectrum ($I/I_c$, solid black curve) compared with the average synthetic spectrum (dotted red). Middle: dynamic spectra of $I/I_c - I_{syn}$ where $I_{syn}$ is the synthetic spectrum. A low-pass filter has been applied to the spectrum in order to minimize the A star contributions and tellurics. Bottom: measured EWs of Hα using ESPaDOnS (open red triangles), Narval (filled green triangles), FEROS (yellow squares), CFHT f/8.2 spectrograph (filled blue circles), and DAO (open black circles) observations. The black curve indicates the best-fitting sinusoid. Both the EWs and dynamic spectra are phased by the B star's rotational period of 0.853 807 d.</figcaption></figure>

the two dashed curves shown in Fig. 14. The mechanism by which this plasma is confined is discussed in the following section.

# 7 MAGNETOSPHERE

As described by ud-Doula, Owocki & Townsend (2008), various characteristics of a star's magnetosphere may be inferred by comparing two parameters: the Kepler radius, $R_K$, and the Alfvén radius, $R_{Alf}$. $R_K$ is the radius at which the gravitational force is balanced by the centrifugal force in a reference frame that is corotating with the star. $R_{Alf}$ characterizes the point within which the magnetic field dominates over the wind and approximately corresponds to the extent of the closed field loops (ud-Doula & Owocki 2002; ud-Doula et al. 2008). Their ratio, $R_{Alf}/R_K$, can therefore be used to define a magnetosphere as either dynamical ($R_{Alf}/R_K < 1$) or centrifugal ($R_{Alf}/R_K > 1$; Petit et al. 2013). It also serves as an indicator of the volume of the magnetosphere: those stars having comparatively larger $R_{Alf}/R_K$ will be capable of confining the emitted wind at larger radii. Furthermore, since a stronger field would be capable of confining more mass, a correlation between the Alfvén radius and the magnetosphere's density may be expected.
<table><caption>**Table 10.** Magnetospheric parameters derived from $\langle B_z \rangle_H$ and $\langle B_z \rangle_{He+metal}$ measurements.</caption><thead><tr><th></th><th>H<sub>LSD</sub></th><th>He+metal<sub>LSD</sub></th></tr></thead><tbody><tr><td>i (°)</td><td>24<sup>+8</sup><sub>-9</sub></td><td>24<sup>+8</sup><sub>-9</sub></td></tr><tr><td>β (°)</td><td>63 ± 13</td><td>66<sup>+9</sup><sub>-13</sub></td></tr><tr><td>B<sub>p</sub> (kG)</td><td>14<sup>+9</sup><sub>-3</sub></td><td>15<sup>+10</sup><sub>-4</sub></td></tr><tr><td>M&#x0302;<sub>B=0</sub> (M<sub>☉</sub> yr<sup>-1</sup>)</td><td>(1.3<sup>+6.0</sup><sub>-1.0</sub>) × 10<sup>-10</sup></td><td>(1.3<sup>+6.0</sup><sub>-1.0</sub>) × 10<sup>-10</sup></td></tr><tr><td>V<sub>∞</sub> (km s<sup>-1</sup>)</td><td>1100<sup>+40</sup><sub>-110</sub></td><td>1100<sup>+40</sup><sub>-110</sub></td></tr><tr><td>η<sub>*</sub></td><td>(2.6<sup>+7.9</sup><sub>-1.3</sub>) × 10<sup>6</sup></td><td>(3.0<sup>+6.5</sup><sub>-1.5</sub>) × 10<sup>6</sup></td></tr><tr><td>R<sub>Alf</sub> (R<sub>*</sub>)</td><td>41<sup>+17</sup><sub>-6</sub></td><td>42<sup>+14</sup><sub>-7</sub></td></tr><tr><td>R<sub>K</sub> (R<sub>*</sub>)</td><td>2.1<sup>+0.4</sup><sub>-0.7</sub></td><td>2.1<sup>+0.4</sup><sub>-0.7</sub></td></tr><tr><td>R<sub>Alf</sub>/R<sub>K</sub></td><td>19<sup>+20</sup><sub>-5</sub></td><td>20<sup>+18</sup><sub>-6</sub></td></tr></tbody></table>

Using the mass and rotational period of HD 35502's B star, we find a Kepler radius of $R_K = 2.1_{-0.7}^{+0.4} R_*$, where $R_*$ is the stellar radius at the magnetic equator. We approximate $R_*$ using $R_{eq}$ since this corresponds to the stellar radius at the latitude where the plasma is expected to accumulate. The Alfvén radius is estimated using equation 9 of ud-Doula et al. (2008) for a dipole magnetic field. This expression requires the calculation of the wind confinement parameter, $\eta_*$, which in turn depends on the dipole magnetic field strength, the equatorial radius, the terminal wind speed ($V_\infty$), and the wind mass-loss rate in the absence of a magnetic field ($\dot{M}_{B=0}$). Following the recipe outlined by Vink, de Koter & Lamers (2000), $\dot{M}_{B=0}$ and $V_\infty$ are derived for a B star having $12.5 < T_{eff} \le 22.5$ kK using $V_\infty/V_{esc} = 1.3$, where $V_{esc}$ is the escape velocity. We obtain $\dot{M}_{B=0} = (1.3_{-1.0}^{+6.0}) \times 10^{-10} M_{\odot} \text{yr}^{-1}$ and $V_\infty = 1100_{-110}^{+40} \text{ km s}^{-1}$. Finally, $\eta_*$ is found to be $(2.6_{-1.3}^{+7.9}) \times 10^6$ using the value of $B_p$ derived from $\langle B \rangle_H$, which then yields $R_{Alf} = 41_{-6}^{+17} R_*$. The magnetospheric parameters associated with both the H and He+ metal longitudinal field measurements derived in Section 5 are listed in Table 10.

Given that the hot plasma surrounding HD 35502's B star is corotating with the star at a distance of $\approx 4R_*$, i.e. between $R_K$ and $R_{Alf}$, it is likely that the plasma is being confined by the strong magnetic field. Similar examples of magnetic B-type stars producing H emission well beyond the stellar radius (at distances of $\approx 2-4R_*$) have been previously reported (e.g. Bohlender & Monin 2011; Grunhut et al. 2012; Oksala et al. 2012). In each of these cases, the star's Alfvén radius exceeds its Kepler radius by approximately an order of magnitude (Petit et al. 2013; Shultz et al. 2014). Using the $R_{Alf}$ value obtained from the $\langle B \rangle_H$ measurements, we derived an $R_{Alf}$ to $R_K$ ratio of $19_{-5}^{+20}$. Therefore, the fact that we observe strong Hα emission is in agreement with this $R_{Alf}/R_K \gtrsim 10$ empirical limit.

The magnetic confinement-rotation diagram compiled by Petit et al. (2013) allows $R_K$ and $R_{Alf}$ to be understood within the broader context of all known O and B stars that host magnetospheres. Our characterization of the magnetosphere hosted by HD 35502's B star suggests that it is well within the CM regime. Only two other stars have been discovered exhibiting similar $R_{Alf}$ and $R_K$ values within the derived $R_{Alf} > 40R_*$ and $1.5 \le R_K \le 2.5R_*$. Although approximately four other stars have lower limits of $R_{Alf}$ and $R_K$ that are consistent with HD 35502's, only HD 182180 (Rivinius et al. 2013) and HD 142184 (Grunhut et al. 2012) have reported upper and lower uncertainties. These two examples have similar effective temperatures, surface gravities, radii, and masses to HD 35502's magnetic B star. However, HD 182180 and HD 142184 are slightly faster rotators ($P_{rot} \approx 0.5$ d) and host magnetic fields with weaker dipolar components ($B_p \approx 10$ kG).

## 8 CONCLUSION

The analysis presented here demonstrates a number of new discoveries regarding the nature of HD 35502. The high-resolution spectroscopic and spectropolarimetric observations obtained using ESPaDOnS, Narval, and FEROS indicate that it is an SB3 system containing a central magnetic B-type star and two cooler A-type stars, all of which lie on the main sequence. We confirm that HD 35502's speckle companion reported by Balega et al. (2012) is indeed the A star binary system. Our analysis indicates that both the A stars are physically nearly identical with a mass ratio of $1.05 \pm 0.02$, masses of $2.1 M_{\odot}$, and effective temperatures of 8.9 kK.

Based on radial velocity measurements, we find that the two A stars form a binary system with an orbital period of 5.668 66(6) d. No radial velocity variations of the B star were detected over the 22 yr observing period, which can be explained by the inferred orbital period of $P_{orb} \gtrsim 40$ yr. However, two other explanations can account for the lack of detected radial velocity variations: (1) the inclination angle associated with the A star binary's orbit about the B star may be $\sim 0^{\circ}$ or (2) the A star binary may lie along the line of sight but not be gravitationally bound to the B star. A number of factors favour the triple system description such as the consistent flux ratios between the three components derived here. Specifically, if the two A stars are significantly closer or further than HD 35502's $400 \pm 80$ pc distance, they would no longer lie on the main sequence. Furthermore, the radial velocity of the Orion OB1a subassociation, of which HD 35502 is most likely a member, has a reported velocity of $\approx 24 \text{ km s}^{-1}$ (Morrell & Levato 1991) and is therefore consistent with both the B star's average radial velocity of $21 \pm 2 \text{ km s}^{-1}$ and the radial velocity of the A star binary's centre of mass ($27 \pm 3 \text{ km s}^{-1}$).

Our analysis of HD 35502's central B star revealed the following:

(i) it has an effective temperature of 18.4 kK, a mass of $5.7 M_{\odot}$, and a polar radius of $3.0 R_{\odot}$;

(ii) it rotates relatively rapidly with a rotational period of 0.853 807(3) d;

(iii) we detect a strong magnetic field and derive the magnitude of its dipolar component ($B_p \approx 14$ kG) and its obliquity ($\beta = 63^{\circ}$);

(iv) it exhibits significant line variability in the form of emission and chemical spots. The emission is predominantly observed in Hα at a distance of approximately four times the stellar radius. Strong He abundance variations indicate a higher concentration near the negative pole of the magnetic field's dipole component;

(v) we derive an Alfvén and Kepler radii of $R_{Alf} \approx 41 R_*$ and $R_K \approx 2.1 R_*$, unambiguously indicating that HD 35502's B star hosts a large centrifugally supported magnetosphere.

Our analysis indicates that the 'sn' classification appearing in HD 35502's historical B5IVsnp spectral type (Abt & Levato 1977) is most likely related to the presence of the sharp-lined binary companion along with the strong H emission. We therefore propose that this system be reclassified as a B5IVpe+A+A system.

Stars hosting strongly emitting CM can provide useful insights towards our understanding of both stellar winds and stellar magnetism. Therefore, it is important that these rare systems be studied in detail. While this particular class of magnetic stars is certainly growing, the number of confirmed examples are still insufficient to solve various outstanding issues. For instance, the inferred magnetospheric material densities of CM-hosting stars are largely inconsistent with currently predicted values (e.g. Rivinius et al. 2013; Townsend et al. 2013; Shultz et al. 2014). Moreover, testing the
validity of theoretical models describing the physical nature of these magnetospheres (e.g. the RRM model derived by Townsend & Owocki 2005) requires detailed comparisons with a diversity of observations such as those recently carried out by Oksala et al. (2012, 2015).

We conclude that the strong variability – both in terms of $(B_z)$ and the observed line variability – makes HD 35502 a favourable subject of a magnetic Doppler imaging analysis (Piskunov & Kochukhov 2002). However, any attempt would require the three spectral components to be disentangled. Our orbital solution provides the necessary first step towards accomplishing this task.

## ACKNOWLEDGEMENTS

GAW acknowledges support in the form of a Discovery Grant from the Natural Science and Engineering Research Council (NSERC) of Canada. SIA thanks the United States National Science Foundation for support on his differential photometric studies. We thank Dr Jason Grunhut for the helpful discussion regarding the interpretation of the magnetospheric emission. This work is based on spectropolarimetric observations obtained at the CFHT which is operated by the National Research Council of Canada, the Institut National des Sciences de l'Univers (INSU) of the Centre National de la Recherche Scientifique of France, and the University of Hawaii, observations obtained using the Narval spectropolarimeter at the Observatoire du Pic du Midi (France), which is operated by the INSU, and observations obtained at the DAO, NRC Herzberg, Programmes in Astronomy and Astrophysics, National Research Council of Canada.

## REFERENCES

* Abt H. A., Hunter J. H. J., 1962, ApJ, 136, 381

* Abt H. A., Levato H., 1977, PASP, 89, 797

* Alecian E. et al., 2015, in Meynet G., Georgy C., Groh J., Stee P., eds, Proc. IAU Symp. 307, New Windows on Massive Stars: Asteroseismology, Interferometry, and Spectropolarimetry. Cambridge Univ. Press, Cambridge, p. 330

* Balega Y. Y., Dyachenko V. V., Maksimov A. F., Malogolovets E. V., Rastegaev D. A., Romanyuk I. I., 2012, Astrophys. Bull., 67, 44

* Bohlender D. A., Monin D., 2011, AJ, 141, 169

* Bolton C. T., Harmanec P., Lyons R. W., Odell A. P., Pyper D. M., 1998, A&A, 337, 183

* Borra E. F., 1981, ApJ, 249, L39

* Borra E. F., Landstreet J. D., Thompson I., 1983, ApJS, 53, 151

* Brown A. G. A., de Geus E. J., de Zeeuw P. T., 1994, A&A, 289, 101

* Bychkov V. D., Bychkova L. V., Madej J., 2005, A&A, 430, 1143

* Cardelli J. A., Clayton G. C., Mathis J. S., 1989, ApJ, 345, 245

* Carnochan D. J., 1982, MNRAS, 201, 1139

* Castelli F., Kurucz R. L., 2003, in Piskunov N., Weiss W. W., Gray D. F., eds, Proc. IAU Symp. 210. Astron. Soc. Pac., Uppsala

* Cohen M., Wheaton W. A., Megeath S. T., 2003, AJ, 126, 1090

* Crawford D. L., 1958, ApJ, 128, 185

* Donati J.-F., Semel M., Carter B. D., Rees D. E., Cameron A. C., 1997, MNRAS, 291, 658

* Ekström S. et al., 2012, A&A, 537, 1

* Georgy C., Ekström S., Granada A., Meynet G., Mowlavi N., Eggenberger P., Maeder A., 2013, A&A, 553, 1

* Glagolevskij Y. V., Chountonov G. A., Shavrina A. V., Pavlenko Y. V., 2010, Astrophysics, 53, 133

* Grunhut J. H. et al., 2012, MNRAS, 419, 1610

* Hareter M. et al., 2008, A&A, 492, 185

* Horch E., van Altena W. F. M. G. T., Franz O. G., López C. E., Timothy J. G., 2001, AJ, 121, 1597

* Howarth I. D., 2011, MNRAS, 413, 1515

* Hubrig S. et al., 2015, A&A, 3, 1

* Kaufer A., Stahl O., Tubbesing S., Nørregaard P., Avila G., Francois P., Pasquini L., Pizzella A., 1999, The Messenger, 95, 8

* Kochukhov O., 2007, in Romanyuk I. I., Kudryavtsev D. O., eds, Physics of Magnetic Stars. Special Astrophysical Observatory of the Russian Academy of Science, Russia, p. 109

* Kochukhov O., Makaganiuk V., Piskunov N., 2010, A&A, 524, 1

* Kochukhov O., Rusomarov N., Valenti J. A., Stempels H. C., Snik F., Rodenhuis M., Piskunov N., 2015, A&A, 574, 1

* Krtička J., Mikulášek Z., Henry G. W., Zverko J., Žižovský J., Skalický J., Zvěřina P., 2009, A&A, 499, 567

* Kupka F. G., Ryabchikova T. A., Piskunov N. E., Stempels H. C., Weiss W. W., 2000, Balt. Astron., 9, 590

* Kurucz R., 1993, ATLAS9 Stellar Atmosphere Programs and 2 km/s grid, Kurucz CD-ROM No. 13. Smithsonian Astrophysical Observatory, Cambridge

* Landstreet J. D., Borra E. F., 1978, ApJ, 224, L5

* Landstreet J. D., Bagnulo S., Andretta V., Fossati L., Mason E., Silaj J., Wade G. A., 2007, A&A, 698, 685

* Landstreet J. D., Bagnulo S., Valyavin G. G., Gadelshin D., Martin A. J., Galazutdinov G., Semenko E., 2015, A&A, 580, A120

* Lee T. A., 1968, ApJ, 152, 913

* Leone F., Bohlender D. A., Bolton C. T., Buemi C., Catanzaro G., Hill G. M., Stift M. J., 2010, MNRAS, 401, 2739

* Monin D., Bohlender D., Hardy T., Saddlemyer L., Fletcher M., 2012, PASP, 124, 329

* Morrell N., Levato H., 1991, ApJS, 75, 965

* Neiner C. et al., 2003, A&A, 411, 565

* Nissen P. E., 1976, A&A, 50, 343

* Oksala M. E., Wade G. A., Marcolino W. L. F., Grunhut J., Bohlender D., Manset N., Townsend R. H. D., MiMe S. C., 2010, MNRAS, 405, L51

* Oksala M. E., Wade G. A., Townsend R. H. D., Owocki S. P., Kochukhov O., Neiner C., Alecian E., Grunhut J., 2012, MNRAS, 419, 959

* Oksala M. E. et al., 2015, MNRAS, 451, 2015

* Perryman M. A. C. et al., 1997, A&A, 323, L49

* Petit V. et al., 2013, MNRAS, 429, 398

* Piskunov N., Kochukhov O., 2002, A&A, 381, 736

* Preston G. W., 1967, ApJ, 150, 547

* Rivinius T., Townsend R. H. D., Kochukhov O., Stefl S., Baade D., Barrera L., Szeifert T., 2013, MNRAS, 429, 177

* Rufener F., 1981, A&AS, 45, 207

* Rufener F., Nicolet B., 1988, A&A, 206, 357

* Sharpless S., 1952, ApJ, 116, 251

* Shore S. N., Brown D. N., Sonneborn G., Landstreet J. D., Bohlender D. A., 1990, ApJ, 348, 242

* Shultz M., Wade G., Rivinius T., Townsend R., 2014, in Sigut A., Jones C., eds, ASP Conf. Ser., Bright Emissaries: Be Stars as Messengers of Star Disk Physics. Astron. Soc. Pac., ON, Canada

* Shultz M., Rivinius T., Folsom C. P., Wade G. A., Townsend R. H. D., Sikora J., Grunhut J., Stahl O., 2015, MNRAS, 449, 3945

* Sikora J. et al., 2015, MNRAS, 451, 1928

* Silvester J., Wade G. A., Kochukhov O., Bagnulo S., Folsom C. P., Hanes D., 2012, MNRAS, 426, 1003

* Skrutskie M. F. et al., 2006, AJ, 131, 1163

* Stibbs D. W. N., 1950, MNRAS, 110, 395

* Thompson G. I., Nandy K., Jamar C., Monfils A., Houziaux L., Carnochan D. J., Wilson R., 1978, Catalogue of Stellar Ultraviolet Fluxes (TD-1). Science Research Council, London

* Tokovinin A., 1992, in McAlister H. A., Hartkopf W. I., eds, ASP Conf. Ser. Vol. 32, IAU Colloq. 135: Complementary Approaches to Double and Multiple Star Research. Astron. Soc. Pac., San Francisco, p. 573

* Townsend R. H. D., Owocki S. P., 2005, MNRAS, 357, 251

* Townsend R. H. D., Owocki S. P., Groote D., 2005, ApJ, 630, 81

* Townsend R. H. D. et al., 2013, ApJ, 769, 33

* ud-Doula A., Owocki S. P., 2002, ApJ, 576, 413

* ud-Doula A., Townsend R. H. D., Owocki S. P., 2006, ApJ, 640, L191

* ud-Doula A., Owocki S. P., Townsend R. H. D., 2008, MNRAS, 385, 97
* van Hamme W., 1993, AJ, 106, 2096
van Leeuwen F., 2007, A&A, 474, 653

Vink J. S., de Koter A., Lamers H. J. G. L. M., 2000, A&A, 362, 295

von Zeipel H., 1924, MNRAS, 84, 684

Wade G. A., Donati J.-F., Landstreet J. D., Shorlin S. L. S., 2000, MNRAS, 313, 851

Wade G. A. et al., 2012, MNRAS, 419, 2459

Wade G. A. et al., 2015, MNRAS, 456, 2

Walborn N. R., Hesser J. E., 1976, ApJ, 205, L87

Wall J., Jenkins C., 2003, Practical Statistics for Astronomers. Cambridge Univ. Press, Cambridge

Wright E. L. et al., 2010, AJ, 140, 1868

Yakunin I. et al., 2015, MNRAS, 447, 1418

## SUPPORTING INFORMATION

Additional Supporting Information may be found in the online version of this article:

**Table 4.** uvby photometry obtained from 1992 January to 1994 March. Each measurement is given as a magnitude difference between the variable source ('v', HD 35502) and a non-variable comparison star ('c', HD 35575).

**Table 7.** Longitudinal field measurements derived using the LSD profiles that were generated from the ESPaDOnS and Narval spectra. Columns 3–6 list $\langle B_z \rangle$ derived from H, He+metal, He, C, Si, and Fe LSD profiles.

**Table 8.** EWs of $H\alpha$ calculated using all of the DAO, CFHT f/8.2 spectrograph, Narval, ESPaDOnS, and FEROS spectra. The uncertainties correspond to $3\sigma$.

**Table 9.** EWs of various spectral lines calculated using the Narval, ESPaDOnS, and FEROS spectra. The uncertainties correspond to $3\sigma$ (http://www.mnras.oxfordjournals.org/lookup/suppl/doi:10.1093/mnras/stw1077/-/DC1).

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the article.

This paper has been typeset from a TeX/LATeX file prepared by the author.