Rapid polarization variations in the O4 supergiant $\zeta$ Puppis

We present time-series linear-polarization observations of the bright O4 supergiant $\zeta$ Puppis. The star is found to show polarization variation on timescales of around an hour and longer. Many of the observations were obtained contemporaneously with Transiting Exoplanet Survey Satellite (TESS) photometry. We find that the polarization varies on similar timescales to those seen in the TESS light-curve. The previously reported 1.78-day photometric periodicity is seen in both the TESS and polarization data. The amplitude ratio of photometry to polarization is ~9 for the periodic component and the polarization variation is oriented along position angle ~70 deg-160 deg. Higher-frequency stochastic variability is also seen in both datasets with an amplitude ratio of ~19 and no preferred direction. We model the polarization expected for a rotating star with bright photospheric spots and find that models that fit the photometric variation produce too little polarization variation to explain the observations. We suggest that the variable polarization is more likely the result of scattering from the wind, with corotating interaction regions producing the periodic variation and a clumpy outflow producing the stochastic component. The H$\alpha$ emission line strength was seen to increase by 10% in 2021 with subsequent observations showing a return to the pre-2018 level.


INTRODUCTION
There is no star in the sky that is both hotter and brighter than  Puppis 1 (HD 66811, Naos); some observational basics are listed in Table 1.On the basis of its O4 I(n)fp classification we will refer to it as a supergiant, although it should be noted that this is a spectroscopic designation; in evolutionary terms,  Pup appears to be a product of binary evolution, and may be core-hydrogen burning (Howarth & van Leeuwen 2019).Additional spectral-type qualifiers indicate moderately broad lines ('(n)'), He ii 4686 and strong N iii 4640 emission ('f'), and unspecified peculiarities ('p').
As well as being bright,  Pup is also luminous (log(/L ⊙ ) ≃ 5.6; Howarth & van Leeuwen 2019), with a strong radiation-driven mass outflow.These characteristics have made it a touchstone for stellarwind studies, starting the discovery of UV P-Cygni profiles (Morton et al. 1969) and continuing to an extensive range of observational, modelling, and theoretical studies (at the time of writing, the Simbad ★ E-mail: j.bailey@unsw.edu.au 1 It is the brightest star in Puppis, with the - assignments being distributed among other components of the former constellation of Argo Navis.Sources: (1) Walborn et al. (2010), Sota et al. (2011).
database identifies 465 papers mentioning  Pup published from 2000 onwards).
The star is known to be variable across the observable electromagnetic spectrum, on a variety of timescales.Ground-based photometry obtained in 1986 and 1989 by Balona (1992) showed irregular microvariability (albeit with possible periods of 5.21 or 1.21 days, the former having been identified in H spectroscopy by Moffat & Michaud 1981), while Marchenko et al. (1998) found a period of 2.56 days in Hipparcos photometry (1989)(1990)(1991)(1992)(1993).Using ∼1000 days of photometry from the Solar Mass Ejection Imager (SMEI) instrument on the Coriolis satellite (2003)(2004)(2005)(2006), Howarth & Stevens (2014) discovered a 1.78-d period, with an amplitude on the order of 10 mmag, varying by a factor of ∼2 on 10-to 100-day timescales.This period persisted in BRIght Target Explorer-Constellation (BRITE-Constellation) data (Ramiaramanantsoa et al. 2018(Ramiaramanantsoa et al. , observations 2014(Ramiaramanantsoa et al. -2015)), and is present in Transiting Exoplanet Survey Satellite (TESS) observations taken in 2019 and 2021 (Burssens et al. 2020, and below).The higher-quality space photometry also confirms the presence of stochastic variability, of comparable amplitude to the periodic (but variable) signal but with shorter timescales (∼hours; e.g., Ramiaramanantsoa et al. 2018).
We began obtaining linear photopolarimetry of  Pup in April 2020.The observations were taken as part of a survey of  ′band polarizations of the 135 stars in the Hipparcos catalogue that are brighter than 3. m 0 and south of declination +30 • ; it extends our earlier study of 50 southern stars within 100 pc of the Sun (Cotton et al. 2016).In the course of this survey several previously unknown polarization variables have been identified, including  Pup, for which results are presented in this paper.

Polarization observations
There have been few previous polarimetric observations of  Pup.Circular spectropolarimetry by Barker et al. (1981), David-Uraz et al. (2014), and Hubrig et al. (2016) sets a limit on any dipole magnetic-field strength of ≲100 G. Linear-polarization structure through the H line was found by Harries & Howarth (1996; see also Harries 2000), who measured an ∼-band continuum polarization of ∼0.04%.Serkowski (1970) reports a -band polarization of (0.09±0.02)%, while Heiles (2000) lists a value of (0.04±0.10)%, in an unspecified passband, in his agglomeration of stellar-polarization catalogues. 2  We obtained 255 linear-polarization observations of  Pup between 2020 Apr 5 and 2021 Mar 8. Most were made with the Miniature High-Precision Polarimetric Instrument (Mini-HIPPI, Bailey et al. 2017) mounted on a 23.5-cm Schmidt-Cassegrain telescope (Celestron 9.25-inch) at Pindari Observatory in Sydney.A smaller number of observations were made with the HIPPI-2 instrument (Bailey et al. 2020a) on the 3.9-m Anglo-Australian Telescope (AAT) at Siding Spring Observatory.Both Mini-HIPPI and HIPPI-2 work in the same way, using a ferro-electric liquid crystal (FLC) modulator operating at 500 Hz and compact photomultiplier tubes as detectors.
All the Mini-HIPPI observations reported here were made using an SDSS g ′ filter.The HIPPI-2 observations were taken through three filters: a 425-nm short-pass filter (425SP), and the SDSS g ′ and r ′ filters.Transmission curves for the filters and for other instrument components are given in Bailey et al. (2020a).The typical polarimetric precision achieved at g ′ on  Pup was ∼30 ppm (parts-per-million) with Mini-HIPPI and ∼5 ppm with HIPPI-2. 2 Heiles cites the Mathewson et al. (1978) compilation (CDS catalogue II/34) for this number; they in turn record their source as "unidentified".We have been unable to locate the primary source of the quoted value; our considered speculation is that it may be an otherwise unpublished observation obtained as part of the programme described by Klare et al. (1972), in which case the formal polarization uncertainty would be perhaps only ∼±0.02%, at  eff ≃ 420 nm.(2000), 4. Donati (2003), 5. Kaufer et al. (1999), 6. Semel, Donati & Rees (1993).† Average of three spectra taken over 10 nights * Average of spectra taken over Dec 5-12.
Full details of the observations and calibration methods, and tables of polarization results are given in appendix A.

Space photometry
The Transiting Exoplanet Survey Satellite (TESS, Ricker et al. 2015) has observed  Pup in 2-minute cadence; here we use data taken in sectors 7 & 8 (2019, Jan 8 to Feb 27) and 34 & 35 (2021, Jan 14 to Mar 6).TESS observes in a broad red band covering 600-1000 nm (a longer wavelength than any of the polarimetric observations).
The PDCSAP (Pre-Search Data Conditioning Simple Aperture Photometry) light-curves were accessed using the Python Lightkurve and Astropy packages (Lightkurve Collaboration, 2018;Astropy Collaboration, 2018) and are shown in Fig. 1 as normalized light-curves with periodograms.

Spectroscopy
Motivated by reports of a possible increase in  Pup's mass-loss rate (Cohen et al. 2020;cf. Section 3.1), we obtained a new spectrum of  Pup on 2021 Feb 2 (during the TESS sector 34 observations) using the HERMES spectrograph at the 3.9-m Anglo-Australian Telescope (AAT).HERMES has four optical bandpasses with a resolving power  ≃ 28 000.For this work, we use only the CCD 3 data, which includes the H feature.Exposure times were 0.1, 1, and 10 seconds, and the signal:noise in the combined spectrum is ∼90 per pixel.Three further spectra were obtained in Oct/Nov 2021 using a Shelyak eShel spectrograph on a PlaneWave CDK 1-m telescope at Mardella Observatory, Western Australia ( ≃ 10 4 , fully resolving stellar spectral features).Further spectra were taken in March 2023 using the CDK 1-m telescope and a 0.35-m Ritchey-Chrétien telescope, also with a Shelyak eShel spectrograph, in Perth, Western Australia.
We also recovered the mean spectrum from AAT/UCLES observations obtained in 2000 as part of an unsuccessful spectropolarimetric search for a magnetic field (Donati & Howarth, unpublished); and 13 further high-resolution spectra taken between 2005 and 2016, from the ESO and CFHT archives.These echelle spectra all have resolving powers in the region of  ≃ 50 000, and signal:noise ratios in excess of 100.A representative subset of the data, spanning the duration of the available observations, is summarized in Table 2, with the H profiles shown in Fig. 2.

Spectroscopy: long-term wind changes?
Cohen et al. (2020) reported X-ray spectroscopy of  Pup obtained using Chandra in 2018-2019.They interpreted the observations as indicating a 30-40% increase in the mass-loss rate, , compared with a similar observation from 2000.Because Balmer emission normally arises through recombination, the H line strength is expected to vary as density squared (in the optically thin limit, which is an adequate approximation in this case).The suggested increase in  would therefore give rise to an easily observable factor ∼2 increase in wind emission.
We examined published H spectra of  Pup spanning a halfcentury to review possible longer-term changes.Sources in addition to those already cited are Ebbets (1980); Bohannan et al. (1986); Harries & Howarth (1996); Hillier et al. (2012), andRamiaramanantsoa et al. (2018).As far as we can ascertain (from often small-scale plots), the peak of the profile is within 2% of 1.12× continuum in all previously published spectra, as is also true for the archival spectra we have examined; the profile was 'normal' as late as September 2016 (Fig. 2).
In that context, the 2021 HERMES and eShel spectra do appear to be exceptional, with peak intensities up to ∼1.18× continuum; compared to the mean of the archival spectra, the emission flux is ∼10% greater, measured over the velocity range ±1000 km s −1 .Although the increase is much less than expected on the basis of the proposed change in  (and could arise at fixed  from changes in the wind's velocity law, (), or in its clumping), this does seem to offer some support for the suggestion of a change in the nature of the outflow around the time of the 2018/19 Chandra observations, with an H signature still present 1 1 /2-2 years later.Note that we have not been able to find any H observations for 2018/2019.However, the morphology appears to have returned to normal in our most recent spectrum (March 2023).

Light-curves
The TESS light-curves and periodograms are given in Fig. 1, with numerical results in Table 3.Previous measurements of the 1.78-d  (Zechmeister & Kürster 2009;Ferraz-Mello 1981), with errors estimated from Monte-Carlo simulations using a residual-permutation algorithm.
The BRITE-Constellation results were obtained in blue and red filters (390-460 nm and 545-695 nm), with an amplitude 7±3% greater in the blue (Ramiaramanantsoa et al. 2018).The single TESS passband has an effective wavelength  eff ≃ 800 nm, and so observed amplitudes may be ∼10% smaller than the BRITE-b values.SMEI also had a very broad response, with  eff ≃ 600 nm for blue stars, roughly similar to the BRITE-r passband.
From the periodograms it is apparent that the previously reported 1.78-day periodicity is present in the TESS data, but is highly variable in amplitude.The signal is strongest in the sector 34 data, at an intermediate level in the sector 7 and 35 data, and is essentially absent in sector 8.The first harmonic is also apparent in sectors 7 and 35 but is variable in strength, indicating that the light-curves are variable in shape and can be non-sinusoidal.The first harmonic was also seen in the SMEI and BRITE-Constellation data reported by Howarth & Stevens (2014) and Ramiaramanantsoa et al. (2018), and shows more strongly in the latter dataset than in the TESS observations.
It is clear from the light-curves that there is also substantial highfrequency variability that is not associated with the 1.78-day periodicity.This is most clearly seen in the sector-8 light-curve where the 1.78-day period is not apparent.However, similar variability is also seen in the other sectors in addition to the periodic component.There are no distinct periodicities apparent in the periodograms associated with this component.This stochastic variability was first discussed in detail by Ramiaramanantsoa et al. (2018), although Balona's (1992) discovery of "irregular microvariability" at a similar level is evidently related.

Polarization variability
The polarization variability of  Pup was apparent after the first few observations, made beginning in April 2020.It was soon established that substantial variation could be seen over a few hours.For example, observations on 2020 May 16 (MJD 58985) show  varying from −687 to +54 ppm in ∼4 hours.
We made intensive polarization observations during the period of TESS observations for sectors 34 and 35 in early 2021.Sequences of up to 12 observations per night were made over this period.Figure 3   shows the g ′ polarization data for the normalized Stokes parameters  and , all given in ppm.
The polarization data are overlaid on the TESS photometery with the mean level subtracted and divided by 15.With this scaling the polarization and photometry amplitudes are similar, and it can be seen that the polarization and photometry vary on similar timescales.The plots are not intended to suggest that the photometry and polarization are correlated; in some cases they vary in opposite directions (see discussion in Section 3.5), but it is clear that the same timescales of variability are seen in both photometry and polarimetry and that the amplitude in photometry is about 15 times that seen in each Stokes parameter.
Fig. 4 shows the periodograms of the polarization data in  and , using the full set of  ′ -band polarization data described in Section A. The highest peaks in both  and  correspond to the 1.78 day period.Most of the other strong peaks correspond to ±1or ±2-day aliases of the 1.78 day period, as expected for the irregular spacing of these ground-based data.Most of the polarization data were obtained only in the  ′ -band, but a limited amount of data were obtained at multiple wavelengths and is shown in Fig. 5.These show that variability in the  ′ and  ′ bands are very similar; there is perhaps a slightly increased amplitude in the blue (425SP) filter, although there are very few points.
In Fig. 6 we show the TESS light-curve for sectors 34 and 35 and the corresponding polarization data plotted against phase.The epoch and period used for the phase determinations were, E = MJD 59230.0,P = 1.78 days.The periodicity is clearly seen in the light-curve and both Stokes parameters.The amplitudes of the binned polarization variation are 440 ppm in  and 426 ppm in .The amplitude in , measured as half the vector distance in the  plane between the furthest two bins, is 598 ppm.The photometric amplitude from the binned phase curve is 5184 ppm (4.8 mmag), intermediate between previously reported results (Table 3).The ratio of photometric to polarimetric amplitude is therefore ∼9 for the periodic component.Fig. 7 shows the phase-folded polarization data taken in 2020.There is no TESS photometry at this time.There is still evidence of periodic polarization variability but the amplitude is smaller than in the 2021 data, indicating that, as for the photometry, the amplitude of 1.78-day polarization signal changes.The phasing of the polarization maximum and the form of the phase curve also look different to the 2021 data, which is again consistent with the changes seen in photometry.
The lower panels of Figs. 6 and 7 show the polarization variation over the 1.78-day period in the  plane.While there is quite a lot of scatter, the data points in Fig. 6 (2021 data) mostly lie along a diagonal line from top left to bottom right.This corresponds to polarization position angles of ∼70 • and ∼160 • .Figure 7 (2020 data) shows a similar pattern but rotated a little to higher angles.

Timescale of stochastic variability
As already noted stochastic variability is present in both the photometry and polarimetry.In Fig. 3 it can be seen that the variability timescales are similar by comparing the typical rise or fall time of observations on timescales of a few hours.To make this analysis more quantitative we have plotted histograms of timescales determined from the local inverse slope of the photometric or polarimetric data.These are given in Fig. 8, with the method further illustrated in Fig. 9.
For the TESS data, each point included in the histogram is derived by measuring the local slope of the photometric data from a pair of points 40 minutes apart.This is then converted to a time (Δt), which is the time corresponding to a change of 1 at that slope, where  is the standard deviation of the flux values in the full TESS datasets (as in Fig. 1).These histograms are similar for all four TESS sectors.They have a peak at about 0.06 days (1.44 hours), a minimum value of about 0.03 days (corresponding to the steepest changes seen in the light curves) and a long tail of larger values which correspond to pairs of points on the turnarounds between the rising and falling sections.
It is apparent that these histograms characterize the stochastic component of the variability and are not very sensitive to the 1.78 day periodic components, since the histograms have very similar shape in Sector 8, which has no periodic variability, and in Sector 34, which has the strongest periodic component.
The bottom panel of Fig. 8 shows the histogram derived in the same way from the polarimetry data.This is measured from pairs of consecutive polarization measurements made on the same night.The typical spacing of these data points is similar to the 40 minute value used with the TESS data.The slope was measured separately for the q and u data, and both sets were combined in the histogram.The polarimetry histogram is noisy due to the availability of far fewer  pairs of observations, and the larger relative errors on the points.However, it can be seen that the shape is generally similar to that seen in the TESS photometry.The peak in the histogram appears in roughly the same place but is broader.The broadening can be understood as due to the errors on the polarimetry which can be comparable to the differences between adjacent points, and thus have a larger effect than in the case of photometry.

Correlation analysis
In Fig. 10 we show plots of the correlation between the TESS photometry and the polarization data (, , and  = √︁  2 +  2 ) for the 104 polarization observations taken during TESS coverage.The TESS photometry has been averaged over 30-minute windows centred on the mid-point time of each polarization observation to provide the corresponding values plotted in the figures.The three top panels in Fig. 10 show the original data that include the periodic as well as the stochastic variability.In the lower panels of Fig. 10 the phase binned curves for the periodic components (as shown by the red dots in Fig. 6) have been subtracted from the data points for both photometry and polarimetry to leave only the stochastic component.
Table 4 presents the standard deviations of the data points plotted in Fig. 10 and the Pearson correlation coefficients, , corresponding to each of the six panels in the figure.Based on these standard deviations the ratio of variability amplitudes in photometry and polarization is 12.2 for the total variability and 19.3 for the stochastic component alone.In Section 3.3 we found an amplitude ratio of 9 for the periodic component alone.The periodic component of the variability thus shows up more strongly in polarization.Nevertheless, it is clear from the lower panels of Fig. 10 that the stochastic component is clearly seen in polarization.The scatter in the data points is many times larger than the statistical errors on the data points.This is shown in particular by the AAT HIPPI-2 data (blue points on the plot) for which the errors are typically ∼5 ppm, whereas the points scatter over ±200 ppm.
The correlation coefficients are also listed in Table 4.For our sample size,  values greater than 0.193 [0.251] are significant at the 5% [1%] level.For the original data, which includes the periodic component, we expect to see some correlation given that the periodic signal is present in all datasets (as shown in Fig. 6).The actual measured correlations between polarization and photometry are not strong (| | ≃ 0.19-0.25).This can be understood as a result of the different phasing and shapes of the phase curves.
For the data sets that have the periodic component subtracted (lower panel of Fig. 10) we find no significant correlation between , , and the TESS photometry (correlation coefficients of −0.05 and 0.11).However, there is a larger correlation ( = 0.25, significant at the 1% level) between  and photometry.The corresponding regression line is shown in the lower-right panel of Fig. 10, and shows that, typically,  increases from ∼100 ppm for the faintest points to ∼220 ppm for the brightest, although there remains a large scatter around this line.Fig. 11, which shows the correlation between  and  for the same datasets as in Fig. 10, helps to explain what we are seeing.The lower panel in this figure shows that there is no preferred direction for the stochastic polarization variations.The correlation coefficient is 0.05.The data seem consistent with the stochastic variability being due to a series of events each of which increases the brightness and polarization but with random orientations, leaving no correlations with  and .
In contrast the upper panel of Fig. 11 for the full dataset shows the clear preferred orientation from top left to bottom right, as already noted in Section 3.3 and shown in Fig. 6.The correlation coefficient between  and  here is −0.42, the largest of any of the correlation plots.

DISCUSSION
The new observations reported here, and in particular the detection of polarization variability associated with both the periodic and stochastic components of the variation, provides some additional constraints on the causes of the variability.Two mechanisms that have been discussed as explanations for the periodic variability are pulsation and rotational modulation arising from surface inhomogeneities ('starspots').

Polarization modelling
We model the polarization produced in stellar photospheres using a version of the synspec spectral synthesis code (Hubeny et al. 1985;Hubeny 2012) modified to do a fully polarized radiative-transfer calculation, using the vlidort code of Spurr (2006).In a hot star the polarization is due to scattering from electrons.For a spherical star the radial symmetry means that the polarization will average to zero.Net polarization will arise only when there is a departure from spherical symmetry.In previous work we have observed and modelled the polarization in hot stars due to departure from spherical symmetry as a result of rotational distortion (Cotton et al. 2017a;Bailey et al. 2020b;Lewis et al. 2022;Howarth et al. 2023), reflection in binary systems (Bailey et al. 2019;Cotton et al. 2020) and nonradial pulsation (Cotton et al. 2022).
The input required to synspec/vlidort is one or more stellaratmosphere model structures.In past work we have generally used atlas9 models, which can be easily calculated for the required combinations of  eff and log .The parameters of  Pup,  eff ≃ 40 kK, log  ≃ 3.5, put it outside the range of atlas9 grids (Castelli & Kurucz 2003;Howarth 2011), as the star is too close to the Eddington limit for stable models to be obtained.
Here, we instead use the OSTAR2002 grid of non-LTE models by Lanz & Hubeny (2003), which includes models close to the Eddington limit, including the range relevant to  Pup.The validity of using a hydrostatic model to represent hot, low-gravity stars is discussed by Lanz & Hubeny (2003); in a real star radiation pressure results in a strong stellar wind.However, to test the hypothesis that surface features (spots, or pulsation) give rise to the photometric and polarization variability these photospheric models should suffice (as argued by Ramiaramanantsoa et al. 2018).
Examples of the polarization predicted by some of these models are given in Figs 12 and 13.Here the polarization and specific intensity are plotted as a function of  = cos  where  is the surface-normal viewing angle.The polarization is positive (perpendicular to the limb of the star) over most of the range, but becomes large and negative (parallel to the limb) for small values of  (close to the limb of the star, Collins 1970).Polarization for some of these model atmospheres has also been calculated using a different method by Harrington (2015).
The dashed line in Fig. 12 is the Harrington model for  eff = 40 kK and log  = 3.5.Harrington included only continuum opacities in his calculations, whereas our models also include lines, which may account for the very slightly higher intensities and polarization in the Harrington (2015) model.
From the models plotted in Fig. 12 and 13 it is possible to see the dependence of polarization on  eff and gravity.Polarization increases with higher temperatures and lower gravities, a trend also seen in cooler-star models (Cotton et al. 2017a;Bailey et al. 2020b).Typically polarization at 460 nm increases by factors of 2.5-3 going from 30kK to 40kK.

Pulsation
Pulsation was suggested by Baade (1986) and Reid & Howarth (1996) as a possible explanation for an ∼8.5 hour period seen in absorption line profiles.Subsequent spectroscopy does not show this period (Baade 1991;Reid & Howarth 1996) and it is not seen in the more-recent space photometry (Howarth & Stevens 2014;Ramiaramanantsoa et al. 2018).Pulsation in low-order (ℓ = 1,2) modes was suggested by Howarth & Stevens (2014) as the likely explanation for the 1.78-day variation.Ramiaramanantsoa et al. (2018) argue against a pulsation mechanism based on the non-sinusoidal phase variation, and the changes in shape of the light-curve.They point out that though some radial pulsators can show non-sinusoidal light-curves, the behaviour seen in  Pup, sometimes showing double-peaked light-curves, is incompatible with pulsation.
Our detection of large-amplitude (∼400 ppm) polarization variations over the 1.78-day period provides additional constraints.Polarization variations in hot stars can be produced by the distortion of the stellar photosphere due to non-radial pulsation.Such effects were suggested and modelled more than 40 years ago in the context of  Cephei pulsators (Odell 1979;Watson 1983), but have only very recently been detected observationally (Cotton et al. 2022), in the bright star  Crucis.The polarization arises from electron scattering in the stellar atmosphere, together with the departure from spherical symmetry due to the pulsations.Only non-radial modes of ℓ = 2 or higher can result in polarization.Radial (ℓ = 0) and dipole (ℓ = 1) modes do not produce polarization variations (Watson 1983) and so can be ruled out as the source of the 1.78-day periodicity in  Pup.
The polarization amplitude seen in  Pup is about 50 times larger than that seen in  Cru by Cotton et al. (2022).However, analytic modelling such as that of Watson (1983) predicts only the relative amplitudes in polarization and photometry.In the case of  Cru the polarization amplitude (in g ′ ) was 35× smaller than the photometric amplitude seen by TESS.The corresponding ratio for  Pup is 9× as described in Section 3.3.For a given mode, this ratio is determined by the ratio of the quantities   and   as defined by Watson (1983).These quantities can be derived from a stellar-atmosphere model and Table 5 shows these quantities and their ratio for  Pup compared with the ratio calculated in the same way for  Cru by Cotton et al. (2022).The stellar-atmosphere model used for  Pup was taken from the OSTAR2002 grid (Lanz & Hubeny 2003) for  eff = 40 kK, log  = 3.5.The calculations were performed for wavelengths appropriate to our observations (  at 460 nm for g ′ polarimetry;   at 800 nm for TESS photometry).
The values for the ratio   /  are greater by factors of 4.3 at ℓ = 2 and 4.1 at ℓ = 3 for  Pup compared to  Cru.For equivalent modes and inclinations, this means the amplitude ratio in polarization relative to photometry will be greater by the same amounts (Watson 1983;Cotton et al. 2022).The observed amplitude ratio of 9 for  Pup is therefore plausible for non-radial pulsation in a similar mode to that observed in  Cru.
However, we nevertheless consider pulsation to be unlikely as the correct explanation for the 1.78-day periodicity.It would require a single mode with strong polarization to be the only period seen.In  Cru, 11 frequencies were detected, with only two being seen in polarization and two others having much higher photometric amplitudes than the modes seen in polarization (Cotton et al. 2022).Additionally, the non-sinusoidal nature and changes in shape of the phase curve remain strong arguments against a pulsation origin.

Rotational modulation by photospheric spots
Rotational modulation has been suggested as the cause of (different) periodicities seen in  Pup in spectroscopy (Moffat & Michaud 1981) and photometry (Marchenko et al. 1998).In their analysis of the BRITE-Constellation photometry Ramiaramanantsoa et al. (2018) argued for the 1.78-day periodicity being due to rotational modulation.They presented models involving two bright photospheric spots that could reproduce the observed double-peaked light-curve (although, as pointed out by Howarth & van Leeuwen (2019), any low-amplitude periodic photometric signal can be reproduced by a spot model).
We here model the polarization produced by a rotating star with photospheric spots.We use a spherical model for the star, as did Ramiaramanantsoa et al. (2018), ignoring the rotational flattening, which we would expect to introduce a fixed polarization offset, but not to change substantially the phase dependence (see Section 4.6.2).We used OSTAR2002 model atmospheres (Lanz & Hubeny 2003), with  eff = 40 kK, log  = 3.75 for the star (corresponding to the green line in Fig. 13), and  eff = 45 kK, log  = 3.75 for the spot.
To determine the integrated polarization we use a similar approach to that used by Bailey et al. (2020b), overlaying a rectangular grid of pixels over the observed view of the star, with spacing 0.01 of the star radius, and calculating the specific intensity and polarization for each pixel using synspec/vlidort as described in section 4.1.This produces a map of the intensity and polarization distribution over the star such as that shown in Fig. 14.Summing the data over all pixels gives the integrated intensity and polarization; repeating the analysis for different rotation angles of the star enables the phase curves to be determined.The integrated intensity was calculated at 800 nm to match the TESS photometry, while the integrated polarization was calculated at 460 nm for the g ′ band polarimetry.
We used an inclination for the star's rotation axis of 33 • .This is the required value if the 1.78-day period is the rotation period, and is constrained by the measured distance, observed flux and  e sin  as described by Howarth & van Leeuwen (2019).The adjustable spot parameters are the number, size, and location (for fixed temperature).We tried models with single spots, and with two identical spots spaced in longitude.We were not able to reproduce the shape of the TESS light-curve with a one-spot model.A single spot near the equator produces too small a width for the bright section of the light-curve.The width can be increased by moving the spot to higher latitude, but then the slopes of the rising and falling branches are not well matched.
A model with two spots near the equator, similar to that used by Ramiaramanantsoa et al. (2018), gives a better match to the lightcurve.Our TESS light-curve does not show the clear double peak seen in the BRITE-Constellation data, so the two spots need to be moved closer together in longitude.A good fit to the light-curve is obtained using a model with two circular spots at latitude 5 • , spaced by 110 • in longitude, with radii of 20.6 • measured from the centre of the star.This model is shown by the blue lines in Fig. 6.
While this model reproduces the TESS light-curve quite well, the resulting polarization variations do not fit the observations at all.The amplitude of the  and  variations falls far short of that observed, and the form of the curves is quite different.Although the light-curve is single peaked, the modelled polarization curves have multiple peaks, and more structure than is seen in the observations.It should be noted that the position angle of the star's rotation axis is unknown, so the polarization model can be rotated arbitrarily in the QU plane.However, it is clear that no such rotation improves the fit to the observations.This is most obvious from the QU plot (Fig. 6, bottom panel), where the circular pattern produced in the model is quite unlike the extended distribution of the data points along the plot diagonal.
While our model is not unique and there are likely to be other spot configurations that can fit the light curve, there is no reason to expect such changes to result in significantly larger polarizations.
The polarization variability therefore does not seem to be consistent with an origin in photospheric spots.If the spot model is the correct interpretation for the periodic photometric variation, then the polarization variation must be produced in some other way.

Polarization due to corotating interaction regions
The 1.78-day period is also seen in He ii 4686 emission (Ramiaramanantsoa et al. 2018) and in X-ray emission (Nichols et al. 2021) both of which indicate that the periodicity is not confined to the photosphere, but extends at least some way into the wind.The polarization could therefore arise from scattering in material just above the photosphere, that is still corotating with the star.
Ignace, St-Louis & Proulx-Giraldeau (2015) and Carlos-Leblanc et al. ( 2019) have presented models of polarization due to corotating interaction regions (CIRs) in the wind from a hot star.CIRs are known to occur in the solar wind (Rouillard et al. 2008) and have been invoked to explain variability in P-Cygni profiles in hot-star winds (Cranmer & Owocki 1996;Morel et al. 1997).Cranmer & Owocki use a simple 'bright spot' model to induce azimuthal wind structure, forming spiral-like density enhancements in the wind.Scattering from the gas in these non-spherically-symmetric structures produces phase-dependent polarization.
As a mechanism for the polarization variability, these models have the advantage that it is easier to explain the amplitude of the polarization variations, as the polarization of light scattered from electrons in optically thin gas can be very high.The mechanism also naturally results in variations repeating with the rotation period of the star.
The models for the polarization due to a single CIR do not provide a good match to the behaviour of  Pup in the QU plane.However, more flexibility can be obtained by including two or more CIRs (Ignace et al. 2015).Using two CIRs in the wind, St-Louis et al. (2018) were able to fit a range of different observed polarization curves for the Wolf-Rayet star WR6 that could not be modelled with a single CIR.If the CIRs are generated by the photospheric spots used to explain the light curve, then two CIRs are to be expected for  Pup.Alternatively, if both the photometric and polarimetric periodicities are due to CIRs, then the sometimes double peaked light curve again suggests the presence of two CIRs.
The B1Iab supergiant  Leo (HD 91316) shows variability with a period of 26.8 days (Aerts et al. 2018), interpreted as "rotational modulation by a dynamic aspherical wind".This may be another example of the same mechanism as in  Pup, albeit with a much longer rotation period.

The stochastic variability
The analysis in Sections 3.3 and 3.5 shows that as well as the periodic component of variability,  Pup shows stochastic variability of polarization on a similar timescale and weakly correlated with the photometric variability.We note that similar variability is present in other hot stars with winds.Polarimetry of Wolf-Rayet (WR) stars (St.-Louis et al. 1987;Drissen et al. 1987;Robert et al. 1989) has shown many of them to vary in a stochastic way with amplitudes similar to or larger than that seen in  Pup.The range of variability (measured as   ) is from 160 ppm (described as essentially instrumental) up to 1550 ppm in WR40 (Moffat & Robert 1991;Robert et al. 1989).The equivalent value for  Pup is 228 ppm (Table 4).Robert et al. (1989) found an anticorrelation between the polarization amplitude of these variations and the terminal wind velocity  ∞ .The amplitude of stochastic photometric variability in WR stars is also correlated with wind terminal velocity (Lenoir-Craig et al. 2022).The typical ratio of photometric to polarimetric amplitude for WR stars is ∼20 (Robert 1992), very similar to what we see for the stochastic component of  Pup.Such observations suggest that the stochastic variations originate in the clumpy winds of these objects.
Models of the polarization variability produced by a wind with optically thin clumps have been given by Richardson et al. (1996); Li et al. (2000); Davies et al. (2007) and Li et al. (2009).The timescale of expected variability is determined in part by the wind flow time given by  ★ / ∞ where  ★ is the stellar radius and  ∞ is the wind terminal velocity.The time that a single clump is sufficiently close to the star to contribute to the variable polarization will be a few times the wind flow time with the factor depending on the velocity law for the wind (Davies et al. 2007).
For  Pup the wind terminal velocity is  ∞ = 2250 km s −1 (Puls et al. 1996) and the stellar radius  ★ ∼ 13R ⊙ (Howarth & van Leeuwen 2019).The resulting wind flow time is 1.1 hours.The variability time scales we see in photometry and polarimetry (∼1.4 hours, see Section 3.4) are consistent with this.
These models also make predictions about the amplitude of polarization variability.The variability is expected to be proportional to the mass loss rate ( ) per wind flow time.This determines the amount of scattering gas close to the star.The amplitude also depends on the clump rate () the number of clumps ejected from the star per wind flow time.A large  reduces the polarization variability due to statistical averaging of the random effects from many clumps.We can use the models in Li et al. (2000) and scale for the mass loss rate of  Pup (  = 3.5 × 10 −6 M ⊙  −1 , Cohen et al. 2010) and the wind parameters given above, to find that the   of 228 ppm would be obtained for  ∼ 20.However, modelling by Davies et al. (2007) produces polarizations about 5 times larger, and requires  > 400 (the largest value modelled) to match the same observed variability.
The models have difficulty matching the ratio of photometric to polarimetric amplitude seen in WR stars (Richardson et al. 1996) which is observed to be ∼20 as also seen in  Pup.Optically thin models for the clumps produce lower ratios than this.Richardson et al. (1996) suggest that the clumps may be optically thick which will reduce the polarization levels due to multiple scattering, and enhance the photometric variability due to contributions from emission (in addition to scattering) from the clumps.However, no detailed modelling of winds with optically thick clumps has been performed.
One of the most variable WR stars is WR40 (HD 96548).Ramiaramanantsoa et al. ( 2019) have shown that the stochastic photometric variations of this star can be explained by a clumpy wind with the clumps scattering light from the star.This star has also been studied with simultaneous photometry and polarimetry (Ignace et al. 2023).The correlations between ,  and photometry for WR 40 show some similarities to those we find for the stochastic component in  Pup as discussed in Section 3.5.The ,  plots in both cases show no preferred angle.WR 40 shows no correlation between polarization and photometry.We find a small positive correlation for  Pup.The ratio of polarization to photometric amplitude is similar in both stars.
Ignace et al. ( 2023) explain the observations of WR 40 in terms of a clumpy wind model in which clumps are ejected from the star in random directions.The lack of correlation between photometry and polarization arises from the different ways in which brightness and polarization vary as the clump moves away from the star.The ,  distribution of the stochastic variation, seen in Fig. 11, and lack of any significant ,  correlation, indicates that the clumps are ejected with random directions, as was also the case in WR 40.
Hot supergiants have not been as well studied for polarization variability as WR stars. Cephei (HD 210839), a supergiant of spectral type O6.5I(n)fp (Sota et al. 2011), is an example of a star that shows stochastic variation similar to  Pup in its TESS light-curves, and has polarization variations (Hayes 1978) of similar amplitude to those in  Pup.Krtička & Feldmeier (2021) use the TESS observations of this star as an example of how stochastic variability can be generated by wind instability.Polarization variations on short timescales have also been observed in a number of OB supergiants (Hayes 1984(Hayes , 1986;;Lupie & Nordsieck 1987).The latter authors describe "random polarimetric fluctuations" in seven out of 10 objects studied which they attribute to "electron scattering off blobs embedded in the stellar wind".
Stochastic photometric variability3 , similar to that observed in  Pup, has been found from recent space photometry, to be common in the light curves of many OB supergiants (e.g.Pedersen et al. 2019;Bowman et al. 2019;Burssens et al. 2020).The cause of this variability is the subject of debate, with possible stellar mechanisms being internal gravity waves (Bowman et al. 2019) or subsurface convection (Cantiello et al. 2021).However, such processes, where the light variation originates at the stellar photosphere, are not expected to result in significant polarization.For  Pup and other cases where the stochastic variability is seen in both photometry and polarimetry the clumpy wind model seems more plausible as the direct cause of the variability.However, this does not rule out other processes in the star, such as those just described, contributing indirectly by driving the formation of clumps in the wind.

Mean polarization
The error weighted mean of all of the g ′ observations is  = −38.6 ± 0.9 ppm,  = −4.6 ± 0.9 ppm (or  = 38.9 ± 0.9 ppm).[The mean is higher for MJD 58790 -59031 but similar for MJD > 59220].The historic polarization observations described in Section 1 support such a small constant component for  Pup.This is surprisingly small for a star at a distance of 332±11 pc.

Interstellar polarization
In general, interstellar polarization increases at a rate of 0.2 to 2 ppm/pc within ∼100 pc of the Sun (Bailey et al. 2010;Cotton et al. 2016), and at ten times that rate beyond that (Behr 1959).Indeed, the simple polarization with distance plot of Gontcharov & Mosenkov (2019) shows a median value of ≈ 3000 ppm, with a minimum of ≈ 400 ppm for this distance; their Figure 9 suggests the position angle, , is likely close to either 90 degrees or 45 degrees (for  Pup this is 93.4 ± 0.7 degrees).
To get a more specific understanding of the space around  Pup, in Fig. 15 we have plotted observations of nearby stars from the agglomerated catalogue of Heiles (2000).The, largely historical, observations have large uncertainties4 , but taken together, the right hand panel indicates a floor in polarization with distance of around 2 ppm/pc, consistent with the value found for Southern hemisphere stars within the Local Hot Bubble by Cotton et al. (2016Cotton et al. ( , 2017b)).However there is a region centred around approximately [116 • , −38 • ] where polarization is seen to increase more steeply -approaching the 20 ppm/pc found by Behr (1959) -beginning from a distance of around 200 pc.The mean polarization of  Pup is well below that which would be expected from these trends.
The discrepancy with the expected interstellar polarization implies either multiple misaligned clouds, whose contributions cancel along the line of sight, or a large constant intrinsic polarization component for  Pup -on the order of 650 ppm -that cancels the interstellar component.Either scenario is plausible.The former is supported by the bifurcated distribution of polarization position angles for the region indicated by Gontcharov & Mosenkov (2019).
Pup seemingly lies on the Sunward side of the centre of the Gum Nebula (centred at [120, −43]5 ) but within its expanding cloud (Woermann et al. 2001).It is a similar distance to the Vela OB2 Association (Choudhury & Bhatt 2009).The Gum Nebula is supposed to be associated with the supernova explosion of  Pup's past companion (Choudhury & Bhatt 2009).The Gum Nebula spans a few hundred parsecs, and reaches 60 pc past  Pup toward the Sun (Woermann et al. 2001), and is thus a good candidate for a contrary interstellar component.

Rotational polarization
A possible intrinsic mechanism for constant polarization is rapid rotation which results in a net stellar polarization due to the departure of the star from spherical symmetry.Recently this phenomenon has been observed and modelled in a number of stars (Cotton et al. 2017a;Bailey et al. 2020a;Lewis et al. 2022;Howarth et al. 2023).To cancel a presumed interstellar polarization of thousands of ppm, would require a larger effect than has hitherto been observed.
The methods previously used for modelling the rotational polarization, described in detail in Bailey et al. (2020a), require a set of stellar atmosphere models covering the variation in local temperature and gravity from the equator to the pole of the rotating star.As explained in Section 4.1 the atlas9 models we have used previously do not extend to the temperatures and gravities needed for  Pup.In order to allow such modelling to extend to higher temperatures we have obtained the required set of models by interpolating between models in the OSTAR2002 grid as described in Section 4.1.However, even this method is not sufficient to reach the gravities required for the equatorial regions of a rapidly rotating  Pup model.Such models require some extrapolation beyond the coverage of the OSTAR2002 grid.It should be noted that the grid is limited to values of Γ rad < 1 where Γ rad is the ratio of radiative to gravitational acceleration (see Section 4 of Lanz & Hubeny 2003).The need to extrapolate the grid therefore means that the Eddington limit is exceeded locally at the equator.
The factors that lead to high polarization are a high rotation rate (specified as /  where  is the equatorial angular velocity of the star, and   is the critical angular velocity) and a high inclination.High temperature and low gravity also favour high polarization as they increase the relative importance of scattering in the atmosphere.The range of parameters possible for  Pup are discussed by Howarth & van Leeuwen (2019).Some of the possible models have slow rotation and are not likely to result in significant rotational polarization.As described in Section 4.3, models that are consistent with    = 1.78 days are constrained to relatively low inclinations ( ∼ 33 • , Howarth & van Leeuwen 2019) and the /  depends on the adopted mass.In Fig. 16 we show the rotational polarization predicted for the Howarth & van Leeuwen (2019)  this configuration as the low equatorial gravity would require substantial extrapolation beyond the range of the OSTAR2002 grid.As explained above this means that the results of our hydrostatic modelling are unlikely to be meaningful, and such stellar parameters may not be realistic.
The second model plotted in Fig. 16 is an example that shows the conditions needed to get a rotational polarization ∼1000 ppm.Both a large /  and a high inclination are needed and this is not compatible with what we expect for  Pup (Howarth & van Leeuwen 2019).

Wind asymmetry
Net polarization could also result if the wind has an asymmetric shape due to rotation.This possibility was investigated by Harries & Howarth (1996) in their analysis of spectropolarimetry of  Pup.These authors had the same problem we have described above, that the interstellar polarization is unknown and so the intrinsic polarization cannot be determined.However they determined a lower limit on the intrinsic polarization of 0.08% (800 ppm) based on the difference between line and continuum polarization and an upper limit of 0.44% (4400 ppm) based on their estimate of the maximum likely interstellar polarization.They then used a model of an asymmetric wind and determined that the lower limit corresponded to an equator-to-pole density ratio of 1.3, and the upper limit to a ratio of 3.
Given our different interpretation of the variable polarization in  Pup their determined lower limit is no longer valid.However, the value is similar to the ∼650 ppm we estimate as a likely interstellar value and so the 1.3 equator-to-pole value is about what might be needed as an asymmetry to cancel such an interstellar polarization.However, it should also be noted that Harries & Howarth (1996) used an inclination of 90 • in their modelling, whereas we now think the inclination is ∼33 • which will result in smaller polarizations.
If there was a substantial wind asymmetry we would also expect to see an asymmetric distribution in the QU plane for the stochastic variability (lower panel of Fig. 11) since this effectively maps the directions at which clumps are ejected.As discussed in Sections 3.5 and 4.5 there is no such effect.

CONCLUSIONS
The discovery of previously unobserved polarization variability in a much studied star like  Pup is an indication that polarimetry of even the brightest stars is a neglected field.Efficient polarimeters on small telescopes such as Mini-HIPPI on the 23.5-cm telescope used here, or that described by Bailey et al. (2023), can make important contributions.
We have made 255 linear polarization observations of  Pup including many made at the same time as TESS observations in early 2021.Spectroscopic observations obtained at the same time show somewhat stronger H emission than seen at other times (2000 -2016, and 2023).This increase may be related to the increased mass-loss rate from  Pup reportedly seen in 2018/2019 Chandra observations (Cohen et al. 2020).
The polarization is found to show rapid variations on similar timescales to those seen in the photometry.The polarization varies over the photometric 1.78-day period, and also shows more-rapid variability corresponding to the high-frequency stochastic component seen in photometry.
The polarization amplitude ratio (photometric amplitude divided by polarimetric amplitude) is ∼12 for the variability as a whole, ∼9 for the periodic component and ∼19 for the stochastic component.The periodic component shows variation along a preferred position angle ∼70 • -160 • .The stochastic variation shows a weak correlation between photometry and polarization and has no preferred direction.
We have tried to fit the 1.78-day variability with a model of a rotating star with bright photospheric spots like that proposed by Ramiaramanantsoa et al. (2018).However, models that fit the lightcurve, produce polarization variations with far too low an amplitude and a quite different form of phase curve to those observed.
The presence of polarization variations rules out pulsation in radial (ℓ = 0) or dipole (ℓ = 1) modes as the origin of the 1.78-day periodicity.Non-radial pulsations in ℓ = 2 or higher modes could, in ideal circumstances, produce polarization amplitudes as high as that observed.However, it seems unlikely that such a polarization-favourable mode should be seen in the absence of any other modes.The nonsinusoidal nature and changes of shape of the phase curve also seems inconsistent with pulsation as noted by Ramiaramanantsoa et al. (2018).
We suggest that a more likely explanation for the polarization variation is scattering from gas in the outflowing wind.This mechanism can more easily produce polarization at the levels observed.Polarization due to a clumpy wind has usually been invoked to explain short timescale polarization variability seen in other hot stars such as Wolf-Rayet stars and OB supergiants (see discussion in Section 4.5).Periodic variations in polarization can be produced by scattering from corotating interaction regions in the wind (Ignace et al. 2015;St-Louis et al. 2018).The stochastic variability could be explained by a model of randomly ejected clumps like that used for the similar variability observed in WR40 (Ramiaramanantsoa et al. 2019;Ignace et al. 2023).
The mean polarization level of  Pup is close to zero, which is surprising, considering that we expect significant interstellar polarization at its 332 pc distance.This is presumably the result of fortuitous cancellation of different polarization components with different position angles.The major contributions are most likely interstellar.the addition of a six-position filter wheel and a new compact electronics unit based on the same design as used for HIPPI-2 (Bailey et al. 2020a).The data were reduced using the methods described by Bailey et al. (2020a).Mini-HIPPI was previously used, on a different telescope, to discover the phase-dependent polarization variability of the binary systems Spica (Bailey et al. 2019) and  1 Sco (Cotton et al. 2020).
Each observation requires measurements of the star at four instrument position angles (0 • , 45 • , 90 • and 135 • ), together with sky measurements at the same four angles.For the Mini-HIPPI observations the total exposure time on star was 800 seconds, but due to overheads in the process of moving to sky and recentering the star, which was inefficient on the small telescope, a typical observation lasted from 30 to 35 minutes.For the HIPPI-2 observations exposure times were 160 seconds in the  ′ filter, 320 seconds in  ′ , and 240 seconds in 425SP.
Table A1 summarizes the polarization observing runs.The telescope polarization (TP) was calibrated using observations of lowpolarization standard stars as described in Bailey et al. (2020a).The low-polarization standards used for the Mini-HIPPI observations were Sirius and  Cen.The measured telescope-polarization values for each run are listed in Table A1 (as  TP ,  TP ).
The polarization position-angle for HIPPI-2 was calibrated using observations of the polarized standard stars listed in Table 4 of Bailey et al. (2020a).For Mini-HIPPI observations, we used two bright standards from that list (HD 84810 and HD 187929) and an additional standard, HD 161471 ( 1 Sco,  = 3.0) for which we adopt the parameters  max = 2.28%,  max = 0.56 m,  = 2.8 • (Serkowski et al. 1975).
The polarization data for the Mini-HIPPI observations are listed in Table A2.Polarization values are given as normalized Stokes parameters  = / and  = /, measured in parts per million (ppm).Observing times are given as the Modified Julian Date (MJD = JD − 2400000.5)for the mid-point of the observation, corrected to the Solar-system barycentre.Table A3 lists, in a similar format, the data obtained with HIPPI-2 on the AAT.The quoted uncertainties include the instrumental positioning error (14.0 ppm for the Mini-HIPPI observations), as discussed in Bailey et al. (2020a).
In an upcoming work (Cotton et al., in prep.)we redetermine the parameters that characterise the efficiency of the HIPPI-2 and Mini-HIPPI modulators based on a decade of on-sky measurements, rather than a mix of laboratory and on-sky determinations as presented in Bailey et al. (2017Bailey et al. ( , 2020a)).The main impact of this is to shift the maximum efficiency  max of the different units by several percent.The difference is significant with respect to the reported errors for objects with large polarizations, as presented here.Consequently, here we adopt interim efficiency parameters, based on the work in progress; these are given in a footnote to Table A1.For the MT modulator used with Mini-HIPPI, the reduced efficiency of the Glan-Taylor prism compared to the Wollaston prism, used with HIPPI-2, is wrapped into these figures.

Figure 1 .
Figure 1.Light-curves (left) and periodograms (right) for  Pup from the four TESS sectors considered here.Light-curves are normalized PDCSAP data plotted in parts-per-million (ppm) with the mean value subtracted.The Lomb-Scargle periodograms are plotted as (semi-)amplitudes in ppm.Dashed and dotted lines are at a period of 1.78 days (0.56 c/d) and its first harmonic respectively.

Figure 3 .
Figure 3. Polarization ( and ) observations for times during and around the TESS sector 34 and 35 observations.The grey curves are the TESS light-curves with the mean subtracted, and divided by 15, and show that the polarization and flux vary on similar timescales.

Figure 4 .
Figure 4. Periodograms of the full set of polarization ( and ) observations (compare with those in Fig.1).The highest peaks in both  and  correspond to the 1.78 day period as marked by the dashed orange line.Other strong peaks correspond to 1 or 2 day aliases (marked by black dotted lines)

Figure 5 .
Figure 5. Polarization wavelength dependence as seen in the HIPPI-2 observations.Lines join data points obtained at closely spaced times (typically less than 30 minutes).

Figure 6 .
Figure 6.TESS light-curve and polarization data folded over the 1.78-day period.The top panel is the TESS light-curve for sectors 34 and 35 folded into 100 phase bins.The middle panels are the corresponding polarization data (for MJD > 59220) as individual points (grey) and averaged into 20 phase bins (red).The lower panel is the same data plotted as a QU diagram with the arrows showing the order from low to high phase values.The blue lines are the predictions of a spot model as described in section 4.3.

Figure 7 .
Figure 7.As Fig. 6 but for polarization data taken in 2020 (MJD 58790 -59031).There are no TESS observations for these dates.Variability over the 1.78-day period is less obvious for these data.

Figure 8 .
Figure 8. Timescale of stochastic variability in the four TESS sectors, and in the polarimetry.The histograms of the inverse slope of the data normalized to the standard deviations () of the data points are plotted.See text for further details, and Fig. 9 for an illustration of how the histogram data are obtained.

Figure 9 .
Figure 9. Illustration of the "inverse slope" method for obtaining the data points used in Fig 8.The slope is measured from two points in the TESS data set 40 minutes apart (the solid orange line).The time (Δt) corresponding to a change of 1 is determined.The histogram of all such pairs of points is then plotted.For the polarization data, pairs of consecutive observations on the same night are used in the same way.

Figure 10 .
Figure 10.Correlations between polarization and photometry for polarization observations with simultaneous TESS data.Upper panels are the original observations.Lower panels have the mean light-curve and polarization curves for the periodic component (the red dots in Fig. 6) subtracted and show only the stochastic variations.

Figure 11 .
Figure 11.Correlation between normalized Stokes parameters ( and ) for the same data sets as in Fig. 10.Upper panel is the original observations and shows similar structure to that seen in Fig. 6 for the phase folded data.The lower panel has the periodic component (the red dots in Fig. 6) subtracted to show the stochastic variations.

Figure 12 .
Figure12.Predicted polarization and specific intensity at 460 nm as a function of  for OSTAR2002 stellar-atmosphere models at log  = 3.5, as described in section 4.1.Positive polarization is perpendicular to the limb of the star, negative polarization is parallel to the limb.The solid lines are calculated with synspec/vlidort.The dashed line is a polarization calculation byHarrington (2015) for one of the same model atmospheres.

Figure 14 .
Figure 14.Example of polarization spot modelling.The distribution of specific intensity and overlaid polarization vectors are shown for the two-spot model described in section 4.3.Without the spots the radial pattern of polarization would integrate to zero over the whole disk.The spots break this symmetry and lead to a net polarization that changes as the star rotates.

Figure 15 .
Figure15.Observations of stars nearby to  Pup fromHeiles (2000) agglomerated polarization catalogue.The left hand panel shows the polarization position angle (North over East) for each star as a headless vector on a sky map.The right hand panel shows polarization as a function of distance, with dashed guidelines drawn at / =0.2, 2 and 20 ppm/pc.Plotted are all observations of stars fall within 600 pc of the Sun and 5 degrees separation of  Pup -which is shown in black.The stars are colour coded according to / and numbered in order of angular separation from  Pup 1: HD 65925, 2: HD 68553, 3: HD 64316, 4: HD 64503, 5: HD 63868, 6: HD 63465, 7: HD 62753, 8: HD 64287, 9: HD 62974, 10: HD 63032, 11: HD 62991, 12: HD 62876, 13: HD 71286, 14: HD 64802, 15: HD 71459, 16: HD 71302, 17: HD 61899, 18: HD 70556.The distance information for each star has been updated from the catalogue value using SIMBAD -resulting in the use of either Gaia Collaboration, (2022) or vanLeeuwen (2007) parallaxes.Not shown on the right hand plot are stars 10 and 7, which have polarizations of 12400 ± 350 ppm and 7200 ± 1000 ppm respectively.These two stars are known, by us, to have large intrinsic variable polarizations; the latter because it is a Be star, the former is the subject of an ongoing program of study.A floor in polarization with distance of around 2 ppm/pc, is indicated.However, polarization increases more steeply -approaching the 20 ppm/pc found byBehr (1959) -beginning from a distance of around 200 pc in a region centred around [116 • , −38 • ].The mean polarization of  Pup is well below that which would be expected from these trends.

Figure 16 .
Figure 16.Predicted polarization level due to the rotational distortion of rapidly rotating stars.The lower curve is for the parameters of  Pup according to the /M ⊙ = 25 model of Howarth & van Leeuwen (2019).The upper plot shows the parameters needed to generate a high polarization, in particular a more rapid rotation and a higher inclination, but these parameters are not consistent with what we expect for  Pup.

Table 1 .
Basic observed properties

Table 2 .
Selected spectroscopic observations of  Pup.

Table 3 .
Results of time-series analyses of space photometry.Semiamplitudes are listed for the fundamental frequencies (and so may differ slightly from actual photometric amplitudes in the presence of significant harmonic content).

Table 4 .
Standard deviations and correlation coefficients of TESS photometry and polarimetry (data points as in Fig.10).

Table 5 .
Pulsation model comparison with  Cru (see Section 4.2)

Table A1 .
Bailey et al. (2020a)runs and telescope-polarization (TP) calibrations.Indicates use of a 2× negative achromatic lens, effectively making the focal ratio f/16.A full description, along with transmission curves for all the components can be found inBailey et al. (2020a).The following parameters were used to calculate modulation efficiency as a function of wavelength; MT:  0 = 504.5, = 1.726,  max = 0.916; MLBLE1:  0 = 455.1,  = 1.969,  max = 0.926. Mean values are given as representative of the observations made of  Pup, and  is the number of observations.

Table A2 .
Polarization observations of  Pup in g ′ obtained with the Mini-HIPPI polarimeter on the Pindari Observatory 23.5 cm telescope.

Table A3 .
Multi-filter polarization observations of  Pup with HIPPI-2 on the AAT.