Search
Close
Search
Advanced Search
Search Menu
AI Discovery Assistant

ABSTRACT

We use a 3D Monte Carlo radiative transfer code to study the polarimetric and photometric variability from stationary corotating interaction regions (CIR) in the wind of massive stars. Our CIRs are approximated by Archimedean spirals of higher (or lower) density formed in a spherical wind originating from the star and we also made allowance for a bright Gaussian spot at the base of the CIR. Comparing results from our code to previous analytical calculations in the optically thin case, we find differences which we attribute mainly to a better estimation of the total unpolarized flux reaching the observer. In the optically thick case, the differences with the analytical calculations are much larger, as multiple scattering introduces additional complexities including occultation effects. The addition of a Gaussian spot does not alter the shape of the polarization curve significantly but does create a small excess in polarization. On the other hand, the effect can be larger on the light curve and can become dominant over the resulting CIR, depending on the spot parameters and density of the wind.

polarization, radiative transfer, methods: numerical, stars: winds, outflow, stars: Wolf–Rayet

1 INTRODUCTION

Radiation-driven winds are a defining feature of massive stars (Puls, Vink & Najarro 2008). High mass-loss rates cause a significant fraction of their envelope to be lost to the interstellar medium in all phases of their evolution, which combined with their hot and intense radiation flux, contributes to the enrichment, ionization, and excitation of the gas and dust in the surrounding medium and ultimately, to the evolution of the stellar populations in galaxies (Kennicutt & Evans 2012; Langer 2012). It is well known that these outflows are inhomogeneous on a small scale (e.g. Moffat et al. 1988), but spectral line variability has revealed that, in some cases, large-scale structures also form in the winds of these stars. These large-scale structures were first discovered in the solar wind and were later generalized to a broader context of stellar winds (e.g. Mullan 1984). The most common evidence for the presence of these large-scale asymmetries are the discrete absorption components (DACs) that are observed in the absorption troughs of ultraviolet (UV) P Cygni profiles of O stars (e.g. Howarth & Prinja 1989; Kaper et al. 1996, 1999), which have been shown by Massa & Prinja (2015) to originate at or very close to the stellar surface. However, they are also revealed as large-amplitude variations in the strong optical emission lines of Wolf–Rayet (WR) stars (e.g. Morel, St-Louis & Marchenko 1997; Morel et al. 1999; Chené & St-Louis 2010). These variations are found to be epoch-dependant, meaning that the periodic changes are generally found at various epochs but their characteristics can evolve.

Hydrodynamic simulations have shown that a perturbation at the base of an optically thin wind, represented by a bright (or dim) spot on the surface of the star, generates large-scale corotating structures to form out of the interaction between high and low-velocity flows as the star rotates, both in the 2D (Cranmer & Owocki 1996) and 3D (Dessart 2004) cases. Brown et al. (2004) presented an analytical model to deduce the kinematics of these structures, from optical depth profiles obtained from spectroscopic observations. These large spiral-like structures have been appropriately named ‘Corotating Interaction regions’ or ’CIRs’. The epoch-dependant nature of the above-mentioned periodic spectroscopic changes could then be attributed to the dissipation and regeneration of the CIRs.

A major effort to study CIRs came from the IUE MEGA campaign (Massa et al. 1995) in which different types of massive stars were monitored in UV spectroscopy; the WN5 star, WR6 (St-Louis et al. 1995), the B0.5 Ib star, HD64760 (Prinja, Massa & Fullerton 1995), and the O4I(f)n star, ζ Puppis (Howarth, Prinja & Massa 1995). More recent efforts include characterizing the wind structure of WR1 using spectropolarization (St-Louis 2013) as well as a detailed study of ζ Puppis from extensive time-dependant photometry and spectroscopy (Ramiaramanantsoa et al. 2018).

In two previous publications, we have presented a simplified analytical model for polarimetric variability from such CIRs by describing them parametrically as a spiral-like density enhancement in an otherwise unperturbed spherical wind. In Ignace, St-Louis & Proulx-Giraldeau (2015), we developed the model, expanding on Ignace, Hubrig & Schöller (2009), in the optically thin electron scattering limit and allowing multiple CIRs to be placed on the star, at arbitrary latitudes and azimuth. Our model polarization curves present clear phase-dependant signatures for one or two CIRs, but for multiple CIRs create more complex behaviours. In St-Louis, Tremblay & Ignace (2018) we extend our model to optically thick winds by accounting for multiple scattering using a ‘core-halo’ approach. This approach defines a pseudo-photosphere beyond the radius of the star, above which the wind can be treated as optically thin. This way, we treated this pseudo-photosphere as the source from which the light emerges, while still having treated the wind and the CIR as initiating from the actual stellar radius. That model was then applied to polarimetric observations obtained from the literature of the WR star WR6, well known to show consistent periodic (P = 3.77 d) but epoch-dependant photometric, polarimetric, and spectroscopic variability without convincing evidence for the presence of a companion. A Levenberg–Marquardt (LM) nonlinear least-squares minimization algorithm was developed to fit 13 different data sets obtained over a time span of about 5 yr (as well as two older data sets from 20 yr before). Two CIRs were used, and a number of parameters related to the stellar wind were adopted. The algorithm was able to fit all observations with consistent stellar parameters and found a stellar inclination of 166° and an orientation of the stellar axis on the plane of the sky of 63°. In all cases, the CIRs were found to be located close to the stellar equator and separated by approximately 90° in longitude. Only their specific locations on the stellar surface were found to differ from one epoch to the next.

In this paper, we expand upon these two papers using Monte Carlo radiative transfer (MCRT). Section 2 describes the MCRT model, and Section 3 presents validation tests to determine the limit between the optically thin and optically thick cases, as well as the number of photons required to obtain significant results. We also present error estimates for different wind densities introduced by varying the random number generator seed. In Section 4 we compare our polarimetric calculations to those obtained both for the optically thin (Ignace et al. 2015) and thick (St-Louis et al. 2018) limits. Finally, in Section 5 we present a parameter study for polarization and light-curve calculations including a CIR and an associated bright Gaussian spot. For these calculations, we have used stellar parameters typical of Wolf–Rayet stars. We conclude in Section 6.

2 THE MONTE CARLO RADIATIVE TRANSFER MODEL

The principle behind MCRT simulations is to follow a large number of monochromatic energy packets, each containing an ensemble of individual photons. So the packet, hereafter referred to as a ‘photon’ may be partially polarized. These packets are randomly emitted from the stellar surface and travel through some medium, in our case a hot stellar wind, until they escape the system. Collectively, these photons represent the luminosity emitted by the star. The scattering during radiation transport is determined by randomly sampling the optical depth while the change in direction (and polarization) of the packet is determined by randomly sampling the Rayleigh phase function. Photons are followed until they leave the envelope where they are placed into the appropriate latitude, azimuth, and frequency bins. To improve our signal-to-noise ratio, we have also implemented a source function sampling procedure. This procedure samples every photon interaction, which measures the scattering source function, and then emits virtual photons in the direction of every individual viewpoint. The result is then weighted by the probability of the virtual photon to scatter and subsequently escape in that direction (often called photon peeling, e.g. Yusef-Zadeh, Morris & White 1984; Whitney 2011).

Our goal with these Monte Carlo simulations was to study the polarimetric and photometric changes caused by the presence of a density perturbation representing the CIRs. We employed a 3D approach for the wind model by adopting a time-independent spherical wind from a hot star, threaded by a CIR as a density perturbation in the shape of a 3D Archimedean spiral. This simplified structure is not quite equivalent to that predicted by hydrodynamic calculations (e.g. Cranmer & Owocki 1996) but our goal with these initial numerical models is to compare our results with our analytical calculations. In future work, we intend to complexify the density and velocity structures of our CIRs. Rotation was simulated by having observers view the star from different azimuths along the same ‘inclination’ or latitude, which act as different rotation phases for the star. As electron scattering is the dominant source of polarization in hot stellar winds, we neglected all other types of scattering. Since electron scattering is wavelength independent, we used a monochromatic approach for the purposes of this paper.

2.1 Polarization

We used Stokes vectors to describe the polarization in our simulation, as defined by Chandrasekhar (1960). Stokes vectors are separated into the four parameters I, Q, U, and V by (e.g. Whitney 2011):
(1)
where I describes the intensity of the incoming light, Q and U, the linear polarization components and V, the circular polarization. These components depend on the spherical coordinates θ and ϕ in the star’s frame of reference, and therefore obviously differ depending on the observer’s frame of reference. In our model, V(θ, ϕ) is always 0 as we assume that there are no sources of circular polarization in the wind. Photons emitted by the star have an initial Stokes vector of (1,0,0,0) and travel a distance characterized by a random optical depth of τ = −ln (1 − ϵ), where ϵ is a random number between 0 and 1, which is chosen from a Poisson distribution with unit mean. After travelling this distance, they scatter and are given a new direction (θ′, ϕ′) using the Mueller and phase function matrices which rotate the frame of reference from the previous scattering frame to the next (e.g. Code & Whitney 1995). The resulting Stokes vector is expressed as
(2)
where S′ is the previous Stokes vector. L(i) is the Mueller rotation matrix for an angle i
(3)
and R is the Rayleigh phase function matrix for a scattering angle χ, given by
(4)
Once the photon has escaped from the wind, the Stokes vector reference direction is realigned to the stellar rotation axis, from which we binned the photons into the appropriate observer bins depending on the inclination. Once the simulation is finished, the final Stokes vector is normalized by the total flux for each bin so that
(5)
We then correct for the rotation angle of the stellar axis on the plane of the sky ψ by rotating the final Stokes vectors using L(ψ).

2.2 Wind and grid properties

Our spherical wind was described by a 3D adaptive mesh spherical grid in (r, θ, ϕ), with each cell containing information on the density within it. The density depends on a number of input parameters, as well as on the cell’s distance from the star and whether or not a CIR crosses its path. Following our analytical approach, the density of the wind itself is given by
(6)
where n0 is the number density scaling factor given by
(7)
μe, the mean molecular weight per free electron, mH, the mass of a hydrogen atom, |$\tilde{r}$|⁠, the radial distance from the centre of the star normalized to the stellar hydrostatic radius R*, and |$w(\tilde{r})$| is the normalized velocity of the wind at a location |$\tilde{r}$| given by
(8)
where w0 = v0/v∞, with v0 being the velocity at the base of the wind and v∞, the terminal velocity of the wind. Note that for simplicity, we adopted a standard beta velocity law with an exponent of 1, which is why this exponent does not appear explicitly in the above equation. The number density scaling factor is related to the optical depth parameter τ0 by
(9)
where σT is the Thomson scattering cross-section. We defined grid cells to become radially larger further out in the wind due to the fact that the density varies as 1/r2. However, we have added an extra layer of fine structure for the grid cells containing the CIR in the form of sub-grids for these particular cells. That way, we can resolve the CIR with a much better precision. For our simulations, the coarse grid contains 43 r cells (going from R* to 100R*), 18 θ cells, and 24 ϕ cells. The fine grid for the cells containing a sub-grid have an extra set of (r, θ, ϕ) cells depending on the CIR half-opening angle β0 and the winding rate of the CIR, which is characterized by the ratio of the rotation and wind terminal velocities, vrot/v∞. The star itself contains 1 r cell and the same number of θ and ϕ cells as the wind.
Each CIR that crosses the centre point of a grid cell increments its density. The density in a given cell in our grid is therefore given by
(10)
where η is the density contrast with respect to the spherical wind in a CIR crossing the centre point of the cell given by
(11)
with nCIR the number density in the CIR and nwind that of the spherical wind. η = 0 for a spherical wind. Note that we have adopted the approach, described in Ignace et al. (2015), in which the CIR shape is determined from the radial streamline flow in the rotating frame of reference of the star, corresponding to a spiral pattern as seen by an observer. The centre position ϕs of this spiral is defined by
(12)
where ϕ0 is the azimuth of the CIR at the base of the wind, θCIR is the colatitude of the CIR, and vrot is the equatorial rotation speed of the star.
Note that we characterize the amount of winding of our CIRs with the winding radius defined by
(13)
P being the period of rotation. This parametrization defines the degree of curvature of the spiral. Higher values of v∞/vrot correspond to less curved CIRs.

2.3 Spot model

Our model has the option to add Gaussian spots at the base of the CIRs. The luminosity contribution of a spot at a given point on the surface of our star is given by
(14)
where
(15)
A is the intensity of the spot relative to the star, |$\hat{n}$| is the unit surface vector at a given location on the star, |$\hat{n}_{\rm CIR}$| is that of the centre of the base of the CIR related to the spot and
(16)
with β0 the half width at half-maximum of the spot, which also corresponds to the half-opening angle of the CIR. The luminosity of a given point on the surface of the star is then given by the sum of the contributions from each spot plus the contribution from the base photosphere. This gives a luminosity map of the whole star, which in turn provides the probability per surface area for a photon to be emitted at this position, normalized by the total luminosity.

2.4 Input parameters

Our model requires several input parameters to describe both the spherical wind and the CIRs that we define below. First, there are a number of stellar and wind parameters:

  • Nphot, the number of photons in the Monte Carlo simulation,

  • R*, the stellar radius,

  • vrot, the equatorial rotation speed,

  • τ0, the optical depth scaling factor (see equation 9),

  • v0, the speed of the wind at its base,

  • v∞, the terminal velocity of the wind.

There are also a number of CIR input parameters:

  • NCIR, the number of CIRs,

  • ϕ0, the azimuth of a given CIR,

  • θCIR, the colatitude of a given CIR,

  • β0, the CIR half-opening angle,

  • η, the density contrast of the CIR (see equation 11),

  • Lspot/Lphot, the ratio of the luminosity of the spot at the base of the CIR to the stellar luminosity.

We also require a number of parameters associated with the different viewpoints, such as their total number, the inclination and azimuth relative to the star of each individual viewpoint. To simulate rotation, we take data from a series of viewpoints at the same inclination but different azimuths. These act as different rotation phases of the star.

Unless stated otherwise we have adopted these parameters that are thought to be appropriate for the star WR6:

  • R* = 2.65R⊙ (Hamann, Gräfener & Liermann 2006),

  • v0 = 57 km s−1 (obtained by assuming |$v(r) = v_\infty \left(1-\frac{bR_*}{r}\right)$| with b = 0.97),

  • v∞ = 1900 km s−1 (see St-Louis et al. 1995).

As for default CIR parameters, unless it is explicitly mentioned, CIRs always have:

  • η = 1,

  • θCIR = 90° and ϕ0 = 0°,

  • β0 = 15°.

Also note that Lspot/Lphot = 1, unless explicitly noted that spots were used in the model.

One last input parameter we must set is the seed, a number which initializes a sequence of (pseudo) random numbers that the MCRT code uses. This is useful to reproduce simulations for debugging purposes, but it can also be useful to vary for estimating statistical error, which is what we do in Section 3.3.

3 MODEL VALIDATION

To verify the validity of our model, we have performed tests described in this section. First, we determine the transition between an optically thin and thick wind. Secondly, we explore how polarization values change as a function of the total number of photons. Finally we characterize the effect of the chosen seed for our MCRT simulations.

3.1 Optical thickness

To demarcate thin and thick limits, we carried out a series of linearity tests to determine how q and u behave as τ0 increases. The linear polarization for various values of a CIR’s ‘winding’ radius, r0, defined in equation (13) were calculated for different inclinations of the stellar rotation axis. Our assumption is that a departure from linearity indicates that multiple scattering effects due to optical depth start to become significant. For these tests, we calculated only the linear polarization for a single viewpoint at (θ, ϕ) = (0, 0) and determined how it varies with τ0 and r0/R*. In Fig. 1 we present the calculated q and u values as filled symbols for the edge-on (top) and pole-on (middle) views when the winding radius is varied from r0/R* = 5 to r0/R* = 100, the latter corresponding to an essentially straight CIR. For technical reasons (the reference direction of the observed Stokes vectors is the projection of the stellar rotation axis on the plane of the sky, which is undefined for i = 0), we cannot strictly use i = 0° so our pole-on calculations are in reality for i = 1°. In the bottom panels, we present the polarization for r0/R* = 100 as viewed from different inclinations. To determine where the linear polarization starts to deviate for a linear behaviour, we fitted a straight line to the first few values at small τ0. The fits appear as solid lines. From this figure, we conclude that departures from a linear behaviour begin roughly around τ 0 = 0.03, although for some of the curves a non-linear behaviour only begins for higher values of τ0. Note that the total electron-scattering optical depth integrated along the line-of-sight to the observer is given by
(17)
with ne(r) = ρwind(r)/μemH. Using equation (6), this can be expressed as
(18)
which results in
(19)
For our adopted value of b = 0.97, this yields τe = 3.6 τ0. In order to encompass as much as possible all different winding radii and inclinations, we have adopted τ0 = 0.03 as our optically thin limit; anything higher will be considered as an optically thick calculation. This limit corresponds to τ e of the order of 0.1, a value much smaller than the value usually considered to be the optically thick limit, i.e. τe ≈ 1. Note that the above results are for a density contrast in the CIR of η = 1.
Linear polarization as a function of τ0 for an edge-on (top row) and pole-on view (middle row) for one CIR located in the stellar equator for various values of the winding radius from 5 to 100R*. In the bottom row, we present the linear polarization for an essentially straight CIR (r0/R* = 100) from different viewing angles. The solid lines represent a linear fit of the first few values showing up to what point the change in polarization remains linear.
Figure 1.

Linear polarization as a function of τ0 for an edge-on (top row) and pole-on view (middle row) for one CIR located in the stellar equator for various values of the winding radius from 5 to 100R*. In the bottom row, we present the linear polarization for an essentially straight CIR (r0/R* = 100) from different viewing angles. The solid lines represent a linear fit of the first few values showing up to what point the change in polarization remains linear.

Open in new tabDownload slide

3.2 Photon numbers

For a spherical wind, the polarization should be zero. In Fig. 2, we present our calculated polarization as a function of phase (left-hand panels) and in the q − u plane (right-hand panels) for a spherical wind only, with an intermediate optical depth of τ0 = 0.1 for a pole-on (top) and edge-on (bottom) view using different numbers of photons for each run, from 1 × 106 photons to 1 × 109 photons. For comparison, we also present the polarization from a CIR at the equator with r0 = 100R* for both viewing angles when using the highest number of photons.

Polarization contribution from a spherical wind for runs with different numbers of photons in our simulation. We have also included our calculations of the polarization from a spherical wind with a CIR with τ0 = 0.1 for comparison.
Figure 2.

Polarization contribution from a spherical wind for runs with different numbers of photons in our simulation. We have also included our calculations of the polarization from a spherical wind with a CIR with τ0 = 0.1 for comparison.

Open in new tabDownload slide

To minimize the random noise while minimizing the numerical error associated with the spherical wind, as well as the computation time, the number of photons used in the simulation needs to be optimized. From Fig. 2, with only a million photons, a residual polarization for our spherical wind of about |$0.03{{\ \rm per\ cent}}$| is obtained. As the number of photons is increased, the residual polarization gradually decreases. Finally, very little difference is apparent between simulations with Nphot = 1 × 108 and Nphot = 1 × 109 for both the pole-on and edge-on views. We therefore decided to use Nphot = 1 × 108 for all simulations presented in this work. We also notice qualitative differences in the shape of the pole-on and edge-on q-u noise curves. Indeed, a wind that is viewed nearly pole-on leads to a circular pattern in this q-u plane, while a wind viewed edge-on does not show such a clear pattern. Instead, the polarization values seem to be aligned along a preferred axis in the q-u plane. This behaviour is readily explained by our source function sampling algorithm. Each emitted photon and photon interaction emits a virtual photon in the direction of each viewpoint weighted by the probability it has to scatter towards them. This results in a correlation between viewpoints due to each of them receiving the same virtual photons (more precisely, source function sample events), which are simply weighted differently. For this particular case (τ = 0.1) the residual polarization of the spherical wind for 1 × 108 photons is of the order of |$0.003{{\ \rm per\ cent}}$| for both pole-on and edge-on views, which is essentially negligible compared to that of a CIR.

3.3 Statistical error from seed values

To determine an approximate numerical error bar for the calculated values from our simulations, we performed a series of simulations with identical input parameters, but for different initial seeds for both the pole-on and edge-on view for one equatorial CIR with a winding radius of r0 = 5R* for different values of τ0. For each τ0 value, we performed 20 simulations and calculated the mean polarization and the standard deviation. Fig. 3 shows our calculated mean polarization values and associated standard deviations as a function of phase for different values of τ0. In the top panel we show a pole-on view for τ0 from 0.01 to 3 while in the bottom panel we show an edge-on view from τ0 = 0.01 to 0.3. Note that computations with large values of τ 0 are very expensive in computing time. Therefore, for those cases, we only calculated the polarization values for a small number of phases for the pole-on view. As expected, both the polarization and its standard deviation generally increases with the number of scatterings, characterized by larger values of τ0.

Pole-on and edge-on polarization values as a function of phase for a CIR with r0 = 5R* located at the equator for different values of τ0 from 0.01 to 3.0. The error bars represent the standard deviation from the mean given by running the same simulation with 20 different seeds.
Figure 3.

Pole-on and edge-on polarization values as a function of phase for a CIR with r0 = 5R* located at the equator for different values of τ0 from 0.01 to 3.0. The error bars represent the standard deviation from the mean given by running the same simulation with 20 different seeds.

Open in new tabDownload slide
To show the behaviour between the standard deviation, σp, and τ0 we use the pole-on case for which the polarization is nearly constant with phase. In Fig. 4, we plot the relative error, defined by the mean of the ratio between the standard deviation and the polarization at each phase as a function of τ0 for the pole-on view. The vertical error bars in this plot were calculated using error propagation. The error on the mean value of the polarization, P, was set to σp/N, where N is the number of simulations (in this case, N = 20), while the error on σp was set by calculating the confidence interval for the variance at the 95 per cent level, assuming a normal distribution. The latter uncertainty dominates over the error on the mean of the polarization leading to asymmetric error bars. The behaviour of the relative error as a function of τ0 can be explained roughly as follows. Assuming that the relative error can be drawn from a Poisson distribution, it should be inversely proportional to the square root of the number of photons that scatter at least once given by
(20)
If we carry out N such calculations, the standard deviation of the mean should be proportional to |$1/\sqrt{N}$| and therefore the measured error should vary as
(21)
At small τ0, δp/p should therefore be proportional to |$\frac{1}{N \cdot N_{\rm phot}\tau _0}$| which seems to be compatible with what we observe. At large τ0, δp/p should be constant and again this seems to be compatible with what we observe. In Fig. 5 we plot the mean of the standard deviations as a function of τ0 for the edge-on view. Note that the scatter on δp is higher at larger τ0, which is simply due to the error bars being dependent on phase, as we can see in the bottom graph of Fig. 3. However the mean standard deviation on the points still seems to rise as we go to denser winds.
Mean relative errors a function of τ0 for the different curves in the top graph of Fig. 3.
Figure 4.

Mean relative errors a function of τ0 for the different curves in the top graph of Fig. 3.

Open in new tabDownload slide
Mean of the standard deviations for the different curves in the bottom graph of Fig. 3 as a function of τ0.
Figure 5.

Mean of the standard deviations for the different curves in the bottom graph of Fig. 3 as a function of τ0.

Open in new tabDownload slide

4 COMPARISON WITH THE ANALYTICAL RESULTS

4.1 Optically thin

In Ignace et al. (2015), we presented analytical calculations for the periodic polarization variability from CIRs embedded in an otherwise spherical wind for the optically thin case. In our second paper (St-Louis et al. 2018) we applied a similar approach to optically thick winds. For our MCRT calculations, in order to differentiate between the effects of the CIRs and that of the wind, we proceeded in two steps. First, we only included CIRs by imposing that the density in the grid cells not associated with them was 0. As a second step, we used a non-zero density for the spherical wind and added the density to that of the CIRs.

In Fig. 6, we present a series of polarization intensity images (⁠|$P = \sqrt{Q^2+U^2}$|⁠) for the case of an equatorial CIR only (no wind) and for a value of τ0  = 0.03, which corresponds to the limit of an optically thin case. The four columns correspond to different rotational phases. Phase 0 is for the base of the CIR facing the observer and phase 0.5 is for when it is behind the star. The first four rows are for a CIR with a modest winding radius of r0/R* = 5 and for an inclination of the stellar axis of respectively i = 1° (nearly pole-on), i = 30°, i = 60°, and i = 90° (edge-on). The last row present the edge-on view of an almost straight CIR (r0/R* = 100). The polarization intensity images are proportional to the density of the gas and therefore in addition to showing the distribution of the polarization, they also provide a map of the density structure of the CIR. We have added a circle to the maps to indicate the size of the stellar disc.

Polarization intensity images for one CIR placed at the stellar equator with τ0 = 0.03. Each column corresponds to a different rotation phase, with phase 0 being the CIR footprint facing the observer. The first four rows are for r0/R* = 5 with the stellar inclination varying from i = 1° to i = 90° and the last r0/R* = 100. The circle shows the position of the stellar disc.
Figure 6.

Polarization intensity images for one CIR placed at the stellar equator with τ0 = 0.03. Each column corresponds to a different rotation phase, with phase 0 being the CIR footprint facing the observer. The first four rows are for r0/R* = 5 with the stellar inclination varying from i = 1° to i = 90° and the last r0/R* = 100. The circle shows the position of the stellar disc.

Open in new tabDownload slide

4.1.1 CIRs only

In this section we compare our Monte Carlo polarization calculations with the results from our analytical model, first in the optically thin case and then in the optically thick limit. Our goal is to confront both approaches to verify if they agree and to bring to light any differences there may be.

We first discuss our results when assuming that only the grid cells containing the CIR have a density incremented by equation (10) and that all other cells have a nil density. Note that we retain the same contrast for the CIRs with the wind. This way we will be able to compare the cases with and without a spherical wind. In Fig. 7, we show a comparison between analytical (curves) and MCRT (symbols) calculations for one CIR with r0 = 5R* placed at different latitudes (θCIR from 20° to 80°). In the top panels, we show p normalized by τ0 as a function of phase and q versus u for an inclination of 30° and in the bottom panels, the same plots for an inclination of 60°.

Comparison between our MCRT linear polarization values (shapes) and those from the analytical model (lines) when only a CIR with r0 = 5R* is included, meaning that the wind density is set to 0. We vary the CIR latitude θCIR as we did in Ignace et al. (2015) and present results for an inclination of i = 30° for the top row and i = 60° for the bottom row. Calculations are carried out for τ0 = 0.03. Filled circles (black) are for θCIR = 20°, triangles (red) for θCIR = 40°, squares (green) for θCIR = 60° and diamonds (blue) for θCIR = 80°.
Figure 7.

Comparison between our MCRT linear polarization values (shapes) and those from the analytical model (lines) when only a CIR with r0 = 5R* is included, meaning that the wind density is set to 0. We vary the CIR latitude θCIR as we did in Ignace et al. (2015) and present results for an inclination of i = 30° for the top row and i = 60° for the bottom row. Calculations are carried out for τ0 = 0.03. Filled circles (black) are for θCIR = 20°, triangles (red) for θCIR = 40°, squares (green) for θCIR = 60° and diamonds (blue) for θCIR = 80°.

Open in new tabDownload slide

We note the close similarities between the analytical and the Monte Carlo model. The curve is either single or double peaked, depending if the CIR is viewed in a more stationary manner by a given viewpoint. However, there still are some small differences that are larger than the numerical error, discussed in Section 3.3. In general, the MCRT values are below the analytical ones with the largest deviations at phases near 0.5 when the CIR is located behind the star. Strangely, there does not seem to be a coherent pattern in the deviations as a function of θCIR with the best agreement for θCIR = 20° and θCIR = 80° and the worst for θCIR = 40°. This behaviour is most likely specific to this particular viewing configuration.

4.1.2 CIR with wind

Even though the net contribution from the spherical wind should in principle be zero, we have carried out the Monte Carlo simulations with a non-zero density for the wind. We present our results in Fig. 8, superimposed on the same analytical curves as in Fig. 7. One can see immediately that the differences are much more pronounced than for the case without a wind. Although the general form of the curves are the same, the polarization is attenuated for all phases and maybe even slightly shifted (see for example i = 60°, θCIR = 40°). Once again the largest differences are for phases near 0.5.

Comparison between our MCRT linear polarization values (shapes) and those from the analytical model (lines) when both the CIR with r0 = 5R* and the wind are included. We vary the CIR latitude θCIR as we did in Ignace et al. (2015) and we present results for an inclination of i = 30° for the top row and i = 60° for the bottom row. Calculations are carried out for τ0 = 0.03. The symbols and colours are as in Fig. 7.
Figure 8.

Comparison between our MCRT linear polarization values (shapes) and those from the analytical model (lines) when both the CIR with r0 = 5R* and the wind are included. We vary the CIR latitude θCIR as we did in Ignace et al. (2015) and we present results for an inclination of i = 30° for the top row and i = 60° for the bottom row. Calculations are carried out for τ0 = 0.03. The symbols and colours are as in Fig. 7.

Open in new tabDownload slide

4.1.3 What causes the differences?

Part of the difference might be explained in the way our polarization values are normalized in equation (5). For the analytical calculations, since the scattering in the envelope is expected to be small compared to the direct star light, the total intensity was assumed to be that coming directly from the star, i.e. I = I*. For the Monte Carlo calculations, all scattering contribution as well as the pre- and post-scattering attenuation are included by default since each photon run in the simulation contributes to I. In this interpretation, the differences are so much larger when we include the spherical wind because there are simply many more scatterings.

To attempt to verify this hypothesis we ran simulations where the polarization bins were normalized by the flux of the escaping stellar photons I* only, instead of the total number of photons. In Fig. 9, we compare the analytical values for a pole-on view of a equatorial CIR with r0 = 5R* (blue curve) embedded in an optically thin wind with τ0 = 0.03 as well as the Monte Carlo results when our polarization values are normalized using all photons (black curve) or only the photons escaping directly from the star (green curve). Although the difference is not compensated for completely, the resulting polarization is definitely much closer to the analytical results, with a deviation of only about 7 per cent between the two curves, as opposed to a difference of about 20 per cent between our original result and the analytical calculation. Note that this simple test does not remove the pre- and post-scattering attenuation terms which still contributes to the Monte Carlo value but obviously not to the analytical one.

Comparison between a polarization normalization by stellar escaping photons only and a normalization by all photons, for a wind with τ0 = 0.03 and a CIR with r0/R* = 5 in the pole-on view. Note that we still find a small difference with the analytical results, even with the different normalization. However the resulting polarization curve is much closer to the analytical results.
Figure 9.

Comparison between a polarization normalization by stellar escaping photons only and a normalization by all photons, for a wind with τ0 = 0.03 and a CIR with r0/R* = 5 in the pole-on view. Note that we still find a small difference with the analytical results, even with the different normalization. However the resulting polarization curve is much closer to the analytical results.

Open in new tabDownload slide

4.2 Optically thick

In this section we compare our numerical calculations with the analytical ones presented in St-Louis et al. (2018), which treats optically thick winds in an approximate way using a ‘core-halo’ approach. Since we already noted significant differences between the MCRT and the optically thin analytical calculations, we elected to adopt the simplest possible configuration. Therefore the calculations presented in this section are for an essentially straight CIR (r0 = 100R*) placed at the stellar equator and we consider only the pole-on and edge-on views. In Fig. 10 we compare analytical (solid curves) and MCRT (symbols) results for three different values of τ 0 : for an optically thin wind (τ0 = 0.03), for a moderate optical depth (τ0= 0.5), and a strongly optically thick wind (τ0 = 2.0). The top row shows the total linear polarization p (left), and Stokes parameter q (right) as a function of phase for the pole-on view while the bottom row shows these same parameters for an edge-on view. For the optically thin calculations (black curve and circles), the difference between the MCRT and analytical values are the same order of magnitude as those shown in Figs 7 and 8. For the pole-on view, for example, the difference in p is around |$0.008{{\ \rm per\ cent}}$|⁠. Note that unlike in the optically thin case, here we do not normalize our polarization values by τ0 because in the optically thick cases, the polarization does not scale linearly with τ0, as can readily be seen. The differences between the analytical and MCRT models are considerable for both optically thick cases.

Comparison between our MCRT linear polarization values (shapes) and those from the analytical solution (lines) for a wind containing a straight CIR with r0 = 100R* at different τ0. Top row has a pole-on view, i = 1°, and the bottom row has an edge-on view, i = 90°.
Figure 10.

Comparison between our MCRT linear polarization values (shapes) and those from the analytical solution (lines) for a wind containing a straight CIR with r0 = 100R* at different τ0. Top row has a pole-on view, i = 1°, and the bottom row has an edge-on view, i = 90°.

Open in new tabDownload slide

For the pole-on view, the amplitude of the q MCRT polarization curves and the phase-independant values of p are wildly different from the analytical values, with deviations of about 0.55 per cent for the τ0 = 0.5 curve (blue) and about 0.70 per cent for the τ0 = 2.0 curve ( red). The relative attenuation of the red curve compared to the blue curve in the MCRT calculations is also much larger than in the analytical model. Finally, the shift between the blue and red curves in the statistical approach is much more pronounced than it is in the analytical one, reaching an amount of almost 0.1 in phase, compared to a 0.01 change in phase in the analytical case. The general trend is that, as the wind becomes increasingly optically thick, the differences between the MCRT and analytical models increase. These differences reach extremely large values.

The q versus phase curve for the pole-on view presents an additional intriguing characteristic. Indeed, as τ0 becomes larger, the maximum of the curve gradually shifts from the value of ϕ = 0.25 expected analytically to phases that are increasingly larger. For the edge-on case, the q curve also presents this shift in its maximum but in addition, two dips to negative q values appear on either side of phase 0 (CIR between the star and the observer). The u curve (not shown here) shows values all close to 0, as expected, as the polarization vector is horizontal on the plane of the sky. These two dips are not predicted by our analytical model and renders p curves (⁠|$\sqrt{q^2+u^2}$|⁠) complex looking.

4.3 Interpretation

The behaviour exhibited in these optically thick MCRT simulations are complex and the difference with the analytical calculations are large. In this section, we will present our interpretation of these results for both the pole-on and the edge-on views.

4.3.1 Pole-on view

Fig. 11 shows a linear polarization image of a spherical wind and CIR at phase 0.25 viewed from the pole. Superposed on the images are lines of various colours corresponding to isocontours on a linear scale. Moving radially from the centre of the star, the polarization rises, reaches a maximum (region in white) and decreases again. This is a well-known behaviour for extended atmospheres of early-type stars (e.g. Brown & McLean 1977; Cassinelli, Nordsieck & Murison 1987). As with the analytical model, we see that the CIR causes an excess in polarization throughout most of the wind, except in the interior near the maximum wind polarization. For the analytical calculations and our optically thin model this produces a double-wave q curve with maxima at positive values at phases 0.25 and 0.75 and minima at negative values at phases 0 and 0.5 (see Fig. 10). However, in the optically thick cases, there is a deficit in the region where the polarization peaks at the location of the CIR. This can readily be seen as a break in the dark blue isocontour. This deficit introduces an important new contribution to the polarization of the wind as it breaks the previously axisymmetric polarization of the wind. The resulting curve has maxima at positive q values at phase 0 and 0.5 and minima at negative q values at phases 0 and 0.25. Both contributions (CIR and deficit) therefore vary in antiphase, which greatly reduces the amplitude of the resulting polarization.

24R* by 24R* linear polarization map (p) observed from the pole from a star with a spherical wind containing an essentially straight CIR (r0 = 100R*) for τ0 = 2.0, at phase 0.25. Isophotes have been drawn, representing surfaces of constant polarization.
Figure 11.

24R* by 24R* linear polarization map (p) observed from the pole from a star with a spherical wind containing an essentially straight CIR (r0 = 100R*) for τ0 = 2.0, at phase 0.25. Isophotes have been drawn, representing surfaces of constant polarization.

Open in new tabDownload slide

There is one final ingredient that explains the apparent gradual shift in the q curve with increasing τ0 that is seen in Fig. 10. As can be seen in Fig. 11, even though our chosen CIR is essentially straight (r0 = 100R*), a slight curvature is still present. This can be readily seen by measuring the position of the centre of the CIR on the yellow and red isophotes. While the centre of the deficit in the wind at r = 1.7R* is on the horizontal, the centre of the yellow isophote at the position of the CIR is clearly below and the centre of the CIR on the red isophote is even lower. This will produce a q polarization curve that is slightly shifted from the one from an optically thin wind. As the wind becomes increasingly thick, the CIR will emerge at higher and higher radii shifting the curve accordingly to higher phases. When the two polarization contributions are added (deficit plus slightly curved CIR) we obtain a curve with a greatly reduced amplitude with maxima that gradually shift towards higher phases as the optical thickness of the wind increases, as seen in Fig. 10.

4.3.2 Edge-on view

For the edge-on view, there is an extra level of complexity because there are now occultation effects. In the bottom right panel of Fig. 10, two dips to negative values can be seen in the q polarization curve, on either side of phase 0 for optically thick calculations (τ0 = 0.5 and τ0 = 2.0). As the optical thickness of the wind increases, the dips become deeper and they gradually move away from phase 0. The q values remain positive in the other parts of the curve. Our interpretation of these dips is that we are seeing the occultation of sections of the wind polarization by the dense CIR on either side of phase 0.

In Fig. 12 we present a 2D sketch of an edge-on view of the wind and CIR at a phase around ϕ = 0.1. It can be seen that a large fraction of the wind polarization in the horizontal direction is occulted by the optically thick CIR. As a consequence, the balance between the horizontal and vertical components, previously leading to a nil polarization, is now broken and produces a net vertical polarization, i.e. negative q values. This leads to the two dips in the q curve. These two negative dips are superposed on a curve that is identical to the one that can be seen for the pole-on view and plotted in the top right panel of Fig. 10. Indeed, the contributions from the spherical wind are the same whether they are viewed pole-on or edge-on.

2D sketch of an observer (eye) viewing the star (centre circle) with a wind and CIR (2d cone) from an edge-on point of view at phase ϕ = 0.1. The max polarization band in the wind is represented by the region between two disconnected and partially filled circles. The red or lighter band represents the region occulted by the CIR while the black or darker band represents the region occulted by the star. Note that the CIR is not centred on the opening due to its curvature, which is approximated in this sketch by a slightly displaced CIR.
Figure 12.

2D sketch of an observer (eye) viewing the star (centre circle) with a wind and CIR (2d cone) from an edge-on point of view at phase ϕ = 0.1. The max polarization band in the wind is represented by the region between two disconnected and partially filled circles. The red or lighter band represents the region occulted by the CIR while the black or darker band represents the region occulted by the star. Note that the CIR is not centred on the opening due to its curvature, which is approximated in this sketch by a slightly displaced CIR.

Open in new tabDownload slide

The effects of the curvature are also visible in this edge-on view. At τ0 = 0.5, the CIR emerges closer to the star than at τ0 = 2.0. Therefore at ϕ = 0, the polarization reaches almost 0 for the τ0 = 0.5 case as the dominant part of the CIR is then symmetrical in our line of sight, which is not quite the case for the τ0 = 2.0 case. The curvature also manifest itself through the slight asymmetry of the two dips around phase 0, and the two peaks around phase 0.5. The slight curvature inward when the CIR is at phase 0.25 will scatter more photons into the line of sight than the outward curvature at 0.75 would. This also implies that at phase 0.1 more photons will be scattered out of the line of sight than at phase 0.9.

5 GAUSSIAN SPOT MODELS

In this section, we present MCRT calculations of the polarization, now including the total light intensity of a spherical wind and a CIR along with a stellar spot on the surface of the star at the footpoint of the CIR. Different spot models have been studied for a variety of star types, for example, in Al-Malki (1992), whose model generated small variations in polarization due to asymmetries in the photosphere. Here a Gaussian spot model will be used.

5.1 MCRT models including spots on the stellar surface

Here we present results of MCRT models for a spherical wind spanning three values of τ0; and optically thin wind (0.01), a moderately thick wind (0.1), and a thick wind (1.0). We also include a CIR with r0 = 100R* and a density contrast of η = 1 at the stellar equator. Finally, we include a spot on the surface of the star, at the base of the CIR with the same angular extent as the CIR. We will vary the opening angle of the spot and CIR and the intensity of the spot with respect to the rest of the star.

5.1.1 The effect on the polarization curves

In Fig. 13, we present q polarization curves for an edge-on view for three configurations. First, in the top for a spherical wind with only a spot at the surface of the star (no CIR). In the middle panels, we show a spherical wind with a CIR only (no spot). Finally in the bottom panels, we show results for the combination of a spot and a CIR. We also vary the spot and CIR parameters to get a better idea of their effect. In each plot, the black curve is a spot with a luminosity contrast of 1.2 and a half-opening angle of 15°. The blue curve is for the same opening angle but a luminosity contrast of 1.5. Finally, the red curve is for a spot with a 1.2 luminosity contrast but for a wider spot with β0 = 30°.

Contributions to the q polarization from the Gaussian spot (top row), the CIR (middle row), and both combined (bottom row) for three different τ0; τ0 = 0.01 (left column), τ0 = 0.1 (centre column), and τ0 = 1.0 (right column). Three different combinations of spot luminosity Lspot/Lphot and spot opening angular radius β0 were used. Note that the scale on the first row is different to the two other rows in order to see more clearly the amplitude of variation of intensity caused by the spot.
Figure 13.

Contributions to the q polarization from the Gaussian spot (top row), the CIR (middle row), and both combined (bottom row) for three different τ0; τ0 = 0.01 (left column), τ0 = 0.1 (centre column), and τ0 = 1.0 (right column). Three different combinations of spot luminosity Lspot/Lphot and spot opening angular radius β0 were used. Note that the scale on the first row is different to the two other rows in order to see more clearly the amplitude of variation of intensity caused by the spot.

Open in new tabDownload slide

The effects on the q polarization curves of the CIR only are as discussed in the previous section. Here in addition, we can see that increasing the opening angle increases the amplitude of the curve and the depth of the eclipses of the wind and can be explained within the framework of our interpretation.

The effects on the polarization of a spot are illustrated in the top panels. First note that the amplitude of polarization is a factor of ∼10 smaller than in the case of a CIR only. Secondly, as expected, for an optically thin or modestly thick wind, the effect of increasing the brightness ratio of the opening angle is to increase the amplitude of the curve. This curve has two maxima per cycle, one at 0.25 and the other at 0.75 when the scattering angle is 90° and two minima at q = 0 when the spot is in the line of sight of the observer at phase 0 (forward scattering) or behind the star (occulted). The behaviour for the optically thick wind seems more complex. When the spot is in front of the star at ϕ = 0, the q polarization can either be positive if its contrast is higher (1.5) or negative if it is lower (1.2). When the spot is behind the star at ϕ = 0.5, the q polarization is either 0 for a brighter spot (1.5) or negative for a lower luminosity contrast (1.2). This can be accounted for mainly by numerical noise, as the errors on the polarization values at τ0 = 1 in Fig. 15 are quite large (∼0.03) at this scale. As for the CIR, the amplitudes for the optically thick case at phase 0.25 and 0.75 are not quite equal. This can most likely be explained by the fact that the CIR is slightly curved, even for r0/R* = 100 and the leading and tailing edges then cause an asymmetry in the polarization (see Section 4.3.2).

In the bottom panels, we present the combined effects of the spot and CIR. the most important conclusion is that the effect of the spot on the polarization is similar in nature as that of the CIR (excluding the eclipse effects) but that they are of much smaller amplitude. Therefore, they do not affect significantly the shape of the polarization curves.

5.1.2 The effect on the light-curves

In Fig. 14, we present light-curves associated with the polarization curves presented in Fig. 13. For the spot only, the curves are very much as expected with an increasing amplitude, when the spot is brighter and when it has a bigger surface. As the wind becomes thicker, the amplitude becomes smaller and smaller and the eclipse becomes less sharp. This is because as τ0 increases, light from the spot is diffused outwards, making the spot larger and blurrier. For the CIR only, the behaviour is also as expected. When the CIR is in front at phase 0, it eclipses part of the star creating a dip. Of course, if the CIR is wider, the eclipse is also wider but also deeper. At phases 0.25 and 0.75, it scatters the light into the line of sight, generating excess light. For a wider CIR, these excesses are stronger. Finally, when the CIR is behind the star, it is totally invisible and the relative flux is unaffected (=1.0) for thin and moderately thick winds. For thick winds, some flux seems to reach the observer (>1.0) indicating that when it emerges, the CIR is slightly larger than the stellar photosphere. As for the relative amplitudes between the effects of the spot and that of the CIR, when the wind is thin (0.01) the spot dominates, but when the wind is thick (1.0) the CIR dominates. This is true even for moderately thick winds (0.1).

Contributions to the intensity from the Gaussian spot (top row), the CIR (middle row), and both combined (bottom row) for three different values of τ0; τ0 = 0.01 (left column), τ0 = 0.1 (centre column), and τ0 = 1.0 (right column). Three different combinations of spot luminosity Lspot/Lphot and spot opening angle β0 were used. Note that the scale on the first row is different to the two other rows in order to see more clearly the amplitude of variation of intensity caused by the spot.
Figure 14.

Contributions to the intensity from the Gaussian spot (top row), the CIR (middle row), and both combined (bottom row) for three different values of τ0; τ0 = 0.01 (left column), τ0 = 0.1 (centre column), and τ0 = 1.0 (right column). Three different combinations of spot luminosity Lspot/Lphot and spot opening angle β0 were used. Note that the scale on the first row is different to the two other rows in order to see more clearly the amplitude of variation of intensity caused by the spot.

Open in new tabDownload slide

5.2 MCRT CIR polarization curves for a range of densities

In Fig. 15, we present polarization curves for a wind with an essentially straight equatorial CIR (r0/R* = 100) in a pole-on (top row) and edge-on (bottom row) view for p (left-hand column) and q (right-hand column) as a function of phase for different values of τ0. For these models, we have also added a Gaussian spot with Lspot/Lphot = 1.2. Here 1.5 phase cycles are shown to make the shape of the curve more clear. For the pole-on view, we can see that the total linear polarization, p, increases with τ0 until it reaches a maximum value of |$\sim 0.25{{\ \rm per\ cent}}$| at τ0 = 0.3−0.5. Above this value, increasing τ0 gradually decreases the value of p until it reaches a value of ∼0.15 at τ0 = 2.0. Our calculation at τ0 = 3.0 gives a very similar polarization value. These effects can also be seen in the amplitude of the q curves shown in the top right panel (the u curves are in antiphase with the q curves). In addition to these variations in the amplitude of the q curve with τ0, we can also see the gradual shift in the maxima of the curves, already described in Section 4.3.1. This shift begins to become significant after τ0 = 0.5, approximately when the maximum in p is reached. This is consistent with our interpretation that at a certain value of τ0 (0.3–0.5) the optical depth in the CIR becomes important enough to break the symmetry of the wind polarization, hereby generating a new linear polarization source that varies in antiphase with the polarization curve generated by the CIR itself. As τ0 increases, the CIR emerges at larger and larger distances from the star and because even with r0/R* = 100 it still presents a slight curvature, the polarization curve from the CIR becomes gradually shifted to higher phases as τ0 increases.

Polarization values p and q as a function of phase of a spherical wind with a single, essentially straight CIR (r0 = 100R*) located at the equator for different values of τ0 from 0.01 to 3.0, with a spot included for a pole-on view i = 1° (top row) and an edge-on view i = 90° (bottom row). Note that our τ0 = 3.0 simulation has half the amount of phase points, making the curve a bit less resolved.
Figure 15.

Polarization values p and q as a function of phase of a spherical wind with a single, essentially straight CIR (r0 = 100R*) located at the equator for different values of τ0 from 0.01 to 3.0, with a spot included for a pole-on view i = 1° (top row) and an edge-on view i = 90° (bottom row). Note that our τ0 = 3.0 simulation has half the amount of phase points, making the curve a bit less resolved.

Open in new tabDownload slide

For the edge-on view, we can see the gradual appearance of the double dips caused by the eclipse of the wind by the CIR on either side of phase 0, also starting around τ0 = 0.3. These dips become deeper and wider as τ0 increases as the wind polarization becomes larger and the CIR occults a larger and larger fraction of the wind polarization.

One interesting thing to note here is that, in the edge-on view, we can see at what τ0 the two peaks in the polarization curve around phase 0 start appearing, in this case around τ0 = 0.3. With increasing density the peaks become higher.

As for the pole-on view, we can see that at in between τ0 = 0.3 and τ0 = 0.5, the polarization peaks reach a maximum and start decreasing. There also seems to be slight phase shifting in q becoming most noticeable around τ0 = 1.0. Note that the decrease in polarization seems to have stopped in between τ0 = 2.0 and τ0 = 3.0, since the amplitude has stayed the same; however the shifting in q still continues.

In view of the fact that the polarization curves we obtain using MCRT for CIRs in a spherical wind are quite different from those obtained using analytical models, the fits of the observations of the WR star WR6 presented in St-Louis et al. (2018) need to be re-done. The fact that our MCRT curves have a much smaller amplitude and, depending on the optical depth of the wind, have peaks that are shifted compared to those obtained with the analytical models will certainly result in different output parameters, such as the density contrast or the opening angle of the CIR and perhaps even in a different orientation of the stellar axis with respect to our line-of-sight.

However, performing such new fits is beyond the scope of this paper. First, we will need to calculate a grid of MCRT models by varying the many parameters of our model. Calculations, particularly at higher values of τ0 are very expensive in computing time. Then, we would need to adjust the model curves to the observations using a robust fitting method such as, for example, a Monte Carlo Markov Chain technique such as the one used in the emcee package. We intend to pursue such fits in an upcoming paper.

5.3 MCRT light curves for a range of densities

Our MCRT calculations also include monochromatic light curves for scattering of star light from both the wind and CIR. Total intensities include light from the star and scattered light, I = 1 corresponding to the intensity from stellar light only. In this section, we present results showing how these vary for various wind densities.

In Fig. 16, we illustrate how the light curves evolve as a function of τ0 for a spherical wind containing a single straight CIR (r0 = 100R*) at the equator with a spot on the surface of the star for pole-on (left) and edge-on (right) views. As expected, the pole-on case yields relatively constant values of I with an increase in amplitude for higher values of τ0. The slight systematic variations at this scale are due to the slight inclination (∼1°). The edge-on case shows curves with the same characteristics as those presented in Section 5.1.2, where the CIR contribution becomes gradually more important as τ0 increases. For low τ0 values, we can see a broad contribution from scattered light centred on phase 0 for this essentially straight CIR. As τ0 increases, this excess becomes more and more reduced by the narrow dip generated by the eclipse of the wind of the star by the CIR. Around phases 0.25 and 0.75 however, as the CIR exits the line of sight of the star, we notice two bumps in the light curve, indicating an excess of photons scattered into the line of sight.

Intensity values as a function of phase of a spherical wind with a single straight CIR (r0 = 100R*) located at the equator for different values of τ0 from 0.01 to 3.0, with a spot included for a pole-on view i = 1° (left) and an edge-on view i = 90° (right). Note the scale for the pole-on intensity is different from the edge-on intensity to highlight the difference between the different curves.
Figure 16.

Intensity values as a function of phase of a spherical wind with a single straight CIR (r0 = 100R*) located at the equator for different values of τ0 from 0.01 to 3.0, with a spot included for a pole-on view i = 1° (left) and an edge-on view i = 90° (right). Note the scale for the pole-on intensity is different from the edge-on intensity to highlight the difference between the different curves.

Open in new tabDownload slide

Fig. 17 presents the light curves from a slightly curved CIR (r0 = 5R*) instead of a straight one, for different values of τ0 for a pole-on view (top) and edge-on view (bottom). For these curves, we have repeated the calculation 20 times and present the mean values in the figure. The error bars correspond to the standard deviation of these means. It can readily be seen that these intensity error bars are much smaller than their polarization counterparts in Fig. 3, as the errors are in general of the order of 0.0001.

Pole-on (top graph) and edge-on (bottom graph) intensity values as a function of phase for a spherical wind containing a single CIR with r0 = 5R* located at the equator for different values of τ0 from 0.01 to 3.0 for the pole-on view and from 0.01 to 0.3 for the edge-on view, with a spot included. The error bars represent the standard deviation given by running the same simulation with 20 different seeds.
Figure 17.

Pole-on (top graph) and edge-on (bottom graph) intensity values as a function of phase for a spherical wind containing a single CIR with r0 = 5R* located at the equator for different values of τ0 from 0.01 to 3.0 for the pole-on view and from 0.01 to 0.3 for the edge-on view, with a spot included. The error bars represent the standard deviation given by running the same simulation with 20 different seeds.

Open in new tabDownload slide

6 CONCLUSIONS

In this paper, we have shown that the results from our Monte Carlo statistical approach for treating CIRs differ significantly from the analytical models. While the Monte Carlo model fits relatively well with the optically thin results of Ignace et al. (2015), only with some minor differences when we do not include the spherical wind with the CIR, the polarization becomes much more attenuated compared to the analytical model when we do include the wind in the Monte Carlo simulations. We interpret this as indicating that the scattered light and/or the pre- and post-scattering attenuation have a much more important impact than previously envisaged on the polarization with a decrease of about |$\sim 20{{\ \rm per\ cent}}$|⁠. When we compare our results with those of St-Louis et al. (2018) for the optically thick limit, the differences become even more important, as multiple scattering adds complexity to the polarization curves. First, the scattering of the photons by the CIR towards the line of site at phases 0.25 and 0.75 for an edge-on view is increasingly reduced as the optical depth becomes higher because of the eclipse by the CIR of the polarized spherical wind on either side of phase 0 (CIR in front). Secondly, because of multiple scattering, the optically thick CIR introduces a deficit in the region of maximum wind polarization, yielding a polarization contribution almost completely in antiphase with the polarization generated by the CIR further out in the wind where the density is smaller. These two contributions in almost complete antiphase greatly reduce the amplitude of the resulting polarization. Of course, the fits to observations presented in our previous paper using analytical curves need to be revisited and will certainly yield different stellar and CIR parameters.

Adding spots on the surface of the star at the base of the CIRs has a small effect on the polarization curve, where a slight excess can be observed. However in the total light curves, three cases can be distinguished depending on the spot parameters and the optical thickness of the wind: The first is when the spot dominates in the optically thin limit, the second is when the CIR dominates in the strongly optically thick limit, and the last is when both contributions are significant, in the moderately optically thick limit.

Although this statistical model presented in this paper is relatively simple, we believe it provides a base on which we will be able to build upon. In the future we plan to treat more complex wind and CIR geometries and kinematic structures such as those that result from hydrodynamical simulations.

ACKNOWLEDGEMENTS

NSL acknowledges financial support from the Natural Sciences and Engineering Research Council (NSERC) of Canada. RI acknowledges support by the National Science Foundation under Grant No. AST-1747658. Computations were made on the supercomputer Briarée from the Université de Montréal, managed by Calcul Québec and Compute Canada. The operation of this supercomputer is funded by the Canada Foundation for Innovation (CFI), the Ministère de l’économie, de la science et de l’innovation du Québec (MESI), and the Fonds de recherche du Québec - Nature et technologies (FRQ-NT).

REFERENCES

Al-Malki
M. B.
,
1992
,
PhD thesis
, University of Glasgow

Brown
J. C.
,
McLean
I. S.
,
1977
,
A&A
,
57
,
141

 

Brown
J. C.
,
Barrett
R. K.
,
Oskinova
L. M.
,
Owocki
S. P.
,
Hamann
W.-R.
,
de Jong
J. A.
,
Kaper
L.
,
Henrichs
H. F.
,
2004
,
A&A
,
413
,
959
10.1051/0004-6361:20031557

Crossref
Search ADS
 

Cassinelli
J. P.
,
Nordsieck
K. H.
,
Murison
M. A.
,
1987
,
ApJ
,
317
,
290
10.1086/165277

Crossref
Search ADS
 

Chandrasekhar
S.
,
1960
,
Radiative Transfer
. D
over Publications
,
New York
.

Google Scholar

OpenURL Placeholder Text

 

Chené
A.-N.
,
St-Louis
N.
,
2010
,
ApJ
,
716
,
929
10.1088/0004-637X/716/2/929

Crossref
Search ADS
 

Code
A. D.
,
Whitney
B. A.
,
1995
,
ApJ
,
441
,
400
10.1086/175363

Crossref
Search ADS
 

Cranmer
S. R.
,
Owocki
S. P.
,
1996
,
ApJ
,
462
,
469
10.1086/177166

Crossref
Search ADS
 

Dessart
L.
,
2004
,
A&A
,
423
,
693
10.1051/0004-6361:20040543

Crossref
Search ADS
 

Hamann
W.-R.
,
Gräfener
G.
,
Liermann
A.
,
2006
,
A&A
,
457
,
1015
10.1051/0004-6361:20065052

Crossref
Search ADS
 

Howarth
I. D.
,
Prinja
R. K.
,
1989
,
ApJS
,
69
,
527
10.1086/191321

Crossref
Search ADS
 

Howarth
I. D.
,
Prinja
R. K.
,
Massa
D.
,
1995
,
ApJ
,
452
,
L65
10.1086/309710

Crossref
Search ADS
 

Ignace
R.
,
Hubrig
S.
,
Schöller
M.
,
2009
,
AJ
,
137
,
3339
10.1088/0004-6256/137/2/3339

Crossref
Search ADS
 

Ignace
R.
,
St-Louis
N.
,
Proulx-Giraldeau
F.
,
2015
,
A&A
,
575
,
A129
10.1051/0004-6361/201424806

Crossref
Search ADS
 

Kaper
L.
,
Henrichs
H. F.
,
Nichols
J. S.
,
Snoek
L. C.
,
Volten
H.
,
Zwarthoed
G. A. A.
,
1996
,
A&AS
,
116
,
257

Crossref
Search ADS
 

Kaper
L.
,
Henrichs
H. F.
,
Nichols
J. S.
,
Telting
J. H.
,
1999
,
A&A
,
344
,
231

 

Kennicutt
R. C.
,
Evans
N. J.
,
2012
,
ARA&A
,
50
,
531
10.1146/annurev-astro-081811-125610

Crossref
Search ADS
 

Langer
N.
,
2012
,
ARA&A
,
50
,
107
10.1146/annurev-astro-081811-125534

Crossref
Search ADS
 

Massa
D.
et al. .,
1995
,
ApJ
,
452
,
L53
10.1086/309707

Crossref

 

Massa
D.
,
Prinja
R. K.
,
2015
,
ApJ
,
809
,
12
10.1088/0004-637X/809/1/12

Crossref
Search ADS
 

Moffat
A. F. J.
,
Drissen
L.
,
Lamontagne
R.
,
Robert
C.
,
1988
,
ApJ
,
334
,
1038
10.1086/166895

Crossref
Search ADS
 

Morel
T.
,
St-Louis
N.
,
Marchenko
S. V.
,
1997
,
ApJ
,
482
,
470
10.1086/304122

Crossref
Search ADS
 

Morel
T.
et al. .,
1999
,
ApJ
,
518
,
428
10.1086/307250

Crossref
Search ADS
 

Mullan
D. J.
,
1984
,
ApJ
,
283
,
303
10.1086/162307

Crossref
Search ADS
 

Prinja
R. K.
,
Massa
D.
,
Fullerton
A. W.
,
1995
,
ApJ
,
452
,
L61
10.1086/309708

Crossref
Search ADS
 

Puls
J.
,
Vink
J. S.
,
Najarro
F.
,
2008
,
A&A Rev.
,
16
,
209
10.1007/s00159-008-0015-8

Crossref
Search ADS
 

Ramiaramanantsoa
T.
et al. .,
2018
,
MNRAS
,
473
,
5532
10.1093/mnras/stx2671

Crossref
Search ADS
 

St-Louis
N.
,
2013
,
ApJ
,
777
,
9
10.1088/0004-637X/777/1/9

Crossref
Search ADS
 

St-Louis
N.
,
Dalton
M. J.
,
Marchenko
S. V.
,
Moffat
A. F. J.
,
Willis
A. J.
,
1995
,
ApJ
,
452
,
L57
10.1086/309706

Crossref
Search ADS
 

St-Louis
N.
,
Tremblay
P.
,
Ignace
R.
,
2018
,
MNRAS
,
474
,
1886
10.1093/mnras/stx2813

Crossref
Search ADS
 

Whitney
B. A.
,
2011
,
Bull. Astron. Soc. India
,
39
,
101

 

Yusef-Zadeh
F.
,
Morris
M.
,
White
R. L.
,
1984
,
ApJ
,
278
,
186
10.1086/161780

Crossref
Search ADS
 
© 2019 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society
This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://dbpia.nl.go.kr/journals/pages/open_access/funder_policies/chorus/standard_publication_model)
Issue Section:
Article
Download all slides