Multiepoch VLT-FORS spectro-polarimetric observations of supernova 2012aw reveal an asymmetric explosion

We present VLT-FORS spectropolarimetric observations of the type II supernova (SN) 2012aw taken at seven epochs during the photospheric phase, from 16 to 120d after explosion. We correct for the interstellar polarization by postulating that the SN polarization is naught near the rest wavelength of the strongest lines - this is later confirmed by our modeling. SN2012aw exhibits intrinsic polarization, with strong variations across lines, and with a magnitude that grows in the 7000A line-free region from 0.1% at 16d up to 1.2% at 120d. This behavior is qualitatively similar to observations gathered for other type II SNe. A suitable rotation of Stokes vectors places the bulk of the polarization in q, suggesting the ejecta of SN2012aw is predominantly axisymmetric. Using an upgraded version of our 2D polarized radiative transfer code, we model the wavelength- and time-dependent polarization of SN2012aw. The key observables may be explained by the presence of a confined region of enhanced 56Ni at ~4000km/s, which boosts the electron density in a cone having an opening angle of ~50deg and an observer's inclination of ~70deg to the axis of symmetry. With this fixed asymmetry in time, the observed evolution of the SN2012aw polarization arises from the evolution of the ejecta optical depth, ionization, and the relative importance of multiple versus single scattering. However, the polarization signatures exhibit numerous degeneracies. Cancellation effects at early times imply that a low polarization may occur for ejecta with a large asymmetry. An axisymmetric ejecta with a latitude-dependent explosion energy can also yield similar polarization signatures as an asymmetric 56Ni distribution. Despite these uncertainties, SN2012aw provides additional evidence for the generic asymmetry of type II SN ejecta, of which VLT-FORS spectropolarimetric observations are a decisive and exquisite probe.


Introduction
Photometry and spectroscopy provide critical information on type II supernova (SN). They characterize, for example, the timescale over which the energy stored in the ejecta is released, the color of the escaping radiation, the level of ionization of the spectrum formation region, and the composition of the ejecta. Such data do not, however, provide direct constraints on the geometry of the explosion. Lacking spatially resolved observations for type II-Plateau (II-P) SNe, spectropolarimetry is the best technique to provide constraints on the morphology of SN ejecta, especially at early times (Shapiro & Sutherland 1982). The polarization is linear, caused by scattering with free electrons, and is typically less than 1% in the continuum (Wang & Wheeler 2008). The current spectropolarimetric dataset of type II SNe is limited to 10 − 20 objects, mostly because of the scarcity of nearby events. Indeed, high signal-to-noise ratio observations are required to reveal such a low-level of polarization.
SN 1987A is a hydrogen-rich, core collapse SN classified as type II-peculiar. Unlike type II-P SNe, which arise from redsupergiant (RSG) star explosions, the progenitor of SN 1987A was a blue supergiant star. It is the first SN for which multiepoch multiband polarimetric observations were obtained (Mendez et al. 1988;Vidmachenko et al. 1988;Barrett 1988;Cropper et al. 1988;Jeffery 1991b;Wang et al. 2002). The 0.5% contin-uum polarization level observed in the optical during the first two months is understood to arise from scattering with free electrons in a prolate or oblate ejecta (Hoflich 1991;Jeffery 1991a). The polarization angle shows variations in time and across lines, suggesting some departures from axial symmetry. Overall, an asymmetric distribution of 56 Ni may also account for the observed polarization (Chugai 1992).
Subsequently, spectropolarimetric observations have been gathered for a number of nearby type II SNe, but now of the plateau type, starting with SN 1999em (Leonard et al. 2001). For SN 2004dj, Leonard et al. (2006) collected spectropolarimetric observations throughout the photospheric and the nebular phase, up to about a year after explosion. The intrinsic polarization of the SN was revealed unambiguously by the strong variation in polarization across lines. Furthermore, the continuum polarization evolved from being small during the photospheric phase to being maximum at the onset of the nebular phase, subsequently decreasing as the inverse of the time squared. Chugai (2006) proposed that the polarization was associated with the excess excitation and ionization originating from the presence of high velocity 56 Ni fingers or blobs. Although intimately associated with ejecta asymmetry, Dessart & Hillier (2011) suggest that the polarization peak at the onset of the nebular phase is caused by a radiative transfer effect as the ejecta turns optically thin, while Article number, page 1 of 20 arXiv:2101.00639v2 [astro-ph.HE] 24 May 2021 A&A proofs: manuscript no. ms its subsequent evolution reflects the drop in ejecta optical depth. Chornock et al. (2010) presented spectropolarimetric observations for a few more type II SNe, confirming the general trend observed for SN 2004dj. It thus appears that type II SNe typically show a temporal increase in continuum polarization through the photospheric phase, reaching a maximum at the onset of the nebular phase, and subsequently dropping.
To supplement this dataset, we initiated a VLT − FORS spectropolarimetric program in 2008. Over a number of years, we gathered multiepoch high quality observations for SNe 2008bk, 2012aw, and 2013ej, as well as a more sparse dataset for a few additional objects. Although this dataset has been presented in a terse format previously (Leonard et al. 2012(Leonard et al. , 2015, we are now in a position to present the observations in more detail, and with modeling results. We start in this paper with SN 2012aw. In the next section, we summarize the results from the polarization modeling of Dessart & Hillier (2011). Section 3 presents the numerical approach used for the present study, which is also presented in detail in Hillier & Dessart (in prep.). As in Dessart & Hillier (2011), we adopt a 2D axisymmetric ejectum but we now model the polarized spectrum for the whole optical range rather than focus on isolated lines and the overlapping continuum. In Section 4 we present the spectropolarimetric observations of SN 2012aw that we collected with VLT − FORS while section 5 presents the results of our modeling of SN 2012aw covering both the photospheric and the nebular phase. In Section 6, we investigate the degeneracy of the polarization results by investigating the impact of our assumptions (magnitude of the asymmetry, adoption of the mirror symmetry with respect to the equatorial plane, and nature of the asymmetry). We present our conclusions in Section 7. Dessart & Hillier (2011) studied the linear polarization in 2D axisymmetric Type II SN ejecta using a long-characteristic code as well as a Monte Carlo code. They explored how the intrinsic continuum and line polarization change with wavelength, bound-bound transition, albedo, SN age, or some properties of the asymmetry. Below, we summarize their results since they are relevant for the present study and the interpretation of polarization signatures in general.

Summary of results from Dessart & Hillier (2011)
Even though electron scattering produces a gray opacity, the absorption opacity varies between photoionization edges and causes a variation of the albedo with wavelength. Consequently, the continuum polarization in type II SNe may vary across the optical range. In the modeling, we take this effect into account by separating the different sources of opacity at each wavelength, namely electron scattering on the one hand and on the other hand the bound-bound and bound-free processes (see also Hillier 1994Hillier , 1996. Scattering in lines (i.e., bound-bound transitions), and especially resonance lines, is assumed to be unpolarized. In reality, line scattering has a polarizing effect but it is much weaker than the polarizing effect caused by free electrons (Jeffery 1991a).
The residual polarized flux arises from a non-cancellation of the local polarization integrated on the plane of the sky. A net polarization can result from an asymmetric distribution of scatterers, but it can also stem from an asymmetric distribution of the flux (for example due to optical-depth effects ;Hillier 1994;Dessart & Hillier 2011;Vlasis et al. 2016;Bulla et al. 2015).
Even for large asphericity, a large polarization may not result because of strong cancellation effects. This situation holds particularly at early times because of the confinement of the spectral formation region, which favors a low residual polarization. The integrated polarized flux is also inhibited by multiple scatterings since these tend to randomize the scattering directions, making the radiation field more isotropic. These effects suggest that SN ejecta should produce a small (or at least a smaller) continuum polarization at earlier times, whatever the level of asymmetry. At nebular times, when the ejecta optical depth is below unity along all sight lines, photons typically scatter once with free electrons. Because of the weaker cancellation effects, this situation produces a greater residual polarization for a given asymmetric configuration. Under optically thin conditions, the continuum polarization scales linearly with continuum optical depth (Brown & McLean 1977). Despite these advantages, the faintness of SNe at nebular times presents a major challenge for spectropolarimetric observations (Leonard et al. 2006).
The effects of optical depth can lead to complex polarization signatures, distinct from what is obtained under optically thin conditions. One such feature is a sign reversal of the polarization (equivalent to a 90 deg change in the polarization angle) in the continuum across the optical range (in part driven by the change in albedo between the Balmer edge and the Paschen edge), as well as sign reversals across line profiles and the adjacent continuum. In observations, sign reversals (especially if they are weak) may be hard to diagnose if the interstellar polarization is not known accurately. Sign reversals have been seen in spectropolarimetric observations of some type II-P SNe (Chornock et al. 2010).
In optically thick lines such as Hα, the polarization is generally zero somewhere between the location of maximum absorption in the P-Cygni trough and the location of maximum flux near line center (see for example Figs. 8 and 16 in Dessart & Hillier 2011, in cases where the zero crossing is associated with a sign reversal of the polarization). This location may thus be useful for constraining the interstellar polarization, especially if other lines can also be used for that purpose (for example the Ca ii near-infrared triplet or Hβ). In contrast, weak lines may not cause much reduction to the continuum polarization (for example Hδ) and thus only induce an inflection in the continuum polarization level. In the opposite regime of strong line overlap, for example in regions of line blanketing, both the total flux and the polarized flux are strongly reduced.
For a fixed ejecta asymmetry, the continuum polarization is expected to reach a maximum when the core is revealed. This is understood as stemming from the progressive growth of the extent of the spectrum formation region during the photospheric phase and the drop in envelope optical depth (favoring single instead of multiple scatterings). It may occur at fixed asymmetry, if the asymmetry increases toward the inner ejecta regions, and may even occur if the asymmetry is confined to the outer, H-rich ejecta layers. Such a strong increase in continuum polarization has been seen in a few type II SNe (Leonard et al. 2006;Chornock et al. 2010), and this feature is also present in the SN that we discuss in this paper. At nebular times, when the total electron scattering optical depth of the ejecta is below unity, the continuum polarization is expected to follow a 1/t 2 evolution for a fixed ionization (since the asymmetry and viewing angle do not change with time; see also Brown & McLean 1977), following the drop in ejecta optical depth. In SN ejecta, however, deviations from this could arise if the asymmetry varies with depth and if the ionization evolves in a complicated fashion -for example following the growth of the γ-ray photon mean free path with time.
All of these conclusions apply to the new models presented here using the upgraded long-characteristic code (Section 3.3;

1D-X1 SN2012aw
Fig. 1. Multiepoch spectral comparison between SN 2012aw and our reference 1D CMFGEN model 1D-X1 (see Section 3.1). The data correspond to the observed total flux for SN 2012aw obtained with VLT − FORS in spectropolarimetric mode. The model is redshifted, reddened, and then normalized to the observed flux at 7000 Å. The V-band brightness offset, calculated from the observed photometry (Bose et al. 2013; corrected for distance and extinction) and the model photometry, is given for each epoch and is typically on the order of 0.1 mag. This model is named x1p5 in Hillier & Dessart (2019) and is our best match model for SN 2012aw. see also Hillier & Dessart, in prep.). However, they are now more directly visible since the upgraded code delivers the full polarized spectrum (across the optical or near-infrared range) rather than just the polarization for a single isolated line and its adjacent continuum, as previously done in Dessart & Hillier (2011).

Modeling approach
The spectropolarimetric modeling we perform in this study requires a few preparatory steps. All simulations are based on 1D nonlocal thermodynamic equilibrium (nonLTE) time-dependent radiative transfer simulations of RSG star explosion models performed with CMFGEN, computed recently for a separate study on the diversity of Type II SNe (Hillier & Dessart 2019). These CMFGEN simulations are used as initial conditions for the 2D polarized radiative transfer modeling. We use an updated version of the long-characteristic code described in Dessart & Hillier (2011). The main improvements are the computation of the full optical polarized flux as well as extra flexibility for setting up the 2D axisymmetric ejecta (see Hillier & Dessart, in prep.).
Here, a 2D axially-symmetric ejecta is built using various approaches. The first approach is to introduce a latitudinal scaling to the density (i.e., mass density, ion density, atom density, and electron density), opacity, and emissivity computed from a 1D CMFGEN model. The second approach is to "stretch" these quantities in radial space, by an amount that depends on latitude. One limitation of these two approaches is that the adopted scaling or radial offset may not reflect well the properties of axisymmetric and asymmetric ejecta. In the third approach, we assign distinct 1D CMFGEN models to a range of latitudes so all latitudes are characterized by a physical, albeit 1D, model. One may use two models (defining the properties of the 2D ejecta along the equator and the pole) or a larger set of models to introduce smaller scale asymmetries covering a smaller opening angle.
In this work, we used a combination of approaches. For the first epoch of spectropolarimetric observations of SN 2012aw, we explored the influence of a latitudinal scaling or a latitudinal radial displacement of a 1D model (options 1 and 2 above). This somewhat artificial approach is not followed further. Instead, the bulk of the modeling is based on a mapping in latitude of two distinct 1D CMFGEN models to produce a large scale asymmetry of the ejecta. This approach is also applied from early to late times and thus allows one to see the evolution of the SN polarization for a given asymmetric ejecta. Below, we review the properties of the 1D CMFGEN models and the various 2D asymmetric (but axisymmetric) ejecta that we build from these, before describing in detail how the polarization calculations are performed. The nomenclature as well as the distinguishing properties of the 1D and 2D simulations performed in this study are summarized in Table 1.
3.1. The reference 1D CMFGEN model SN 2012aw was studied by Hillier & Dessart (2019) as part of a large investigation on the origin of the diversity of type II SNe, in particular in terms of V-band decline rate and spectral diversity during the photospheric phase. The reference model we use for SN 2012aw is named here 1D-X1 (originally called x1p5 in Hillier & Dessart 2019) and corresponds to a 15 M zero age main sequence star, exploded to yield an ejecta kinetic energy of 1.2 × 10 51 erg, an ejecta mass of 12.12 M , and a 56 Ni mass of 0.056 M . Figure 1 illustrates the agreement in the optical range (the data are from our spectropolarimetric VLT − FORS Table 1. Summary of 1D and 2D model properties. The first column gives the dimensionality of the simulation, 1D referring to the sphericallysymmetric nonLTE time-dependent simulations with CMFGEN and 2D to the polarized axially-symmetric steady-state radiative transfer simulations with LONG_POL. The second column gives the model name. The third column provides some characteristics of each simulation. The rightmost column indicates in which section the corresponding simulation is discussed. All 2D simulations adopt mirror symmetry with respect to the equatorial plane unless otherwise stated. ). This agreement suggests that a spherical ejecta model yields a satisfactory match to the photometric and spectroscopic observations of SN 2012aw. The deviations from spherical symmetry that we introduce must therefore remain small enough not to deteriorate the agreement between observations and the results for model 1D-X1.

Application of a latitudinal scaling to the reference model
In this approach, the density ρ at an ejecta location (r, β), where r is the local radius and β is the polar angle (i.e., the angular separation from the axis of symmetry), is given by with B set (in all cases) to 1/(1 + A/3) (chosen so that 1 0 B(1 + A cos 2 β) sin β dβ = 1). With this choice, prolate (oblate) density configurations correspond to a positive (negative but greater than −1) A value. Other choices are possible. For example, by varying the exponent on the cosine one can modulate the latitudinal confinement of the asymmetry. This option is not explored here. With this adopted density scaling, and at a given radius r, opacities and emissivities scale with B 2 (1 + A cos 2 β) 2 , while the electron-scattering opacity scales with B(1 + A cos 2 β).

Application of a radial scaling to the reference model
The second option is to apply a radial scaling of the 1D CMFGEN model and use the same density, opacity, and emissivities as the 1D CMFGEN model on that scaled grid. For example, the density (as well as opacities and emissivities) at r along β corresponds to its counterpart in the original 1D CMFGEN model at r/B(1 + A cos 2 β). If the density declines outward, positive (negative) A values correspond to prolate (oblate) density configurations. For an ejecta characterized by a power-law density with exponent n ρ , the radial scaling causes a maximum density contrast between pole and equator of (1+A) n ρ . For n ρ = 12 (which is representative of the outer layers of the ejecta in model 1D-X1), this maximum density contrast corresponds to 14.6 for A = 0.25 and to 130.0 for A = 0.5 (these values are used in two models presented in Section 5.2). In ejecta regions of constant density, a radial scaling yields no latitudinal variation in density (in realistic ejecta, the density gradient varies with depth so a stretching may correspond to complicated variations in density with both inclination and depth). In contrast, the latitudinal scaling (see Eq. 1) always introduces a variation in density with latitude, and this variation between pole and equator has the same magnitude at any given depth.  and 1D-X2b. Bottom: Mass density, initial 56 Ni mass fraction, electron density, gas temperature, and hydrogen ionization fraction vs. velocity for models 1D-X1 (solid) and 1D-X2b (dashed) at 44 d after explosion. These two models vary in the outer ejecta structure, the inner core structure, and the 56 Ni mixing adopted (strong or weak, presence or absence of a 56 Ni-rich shell at a Lagrangian mass of 12 M , corresponding to a velocity of about 4000 km s −1 in the resulting ejecta). The two models have the same mass-density profile outside of the inner ejecta (with the exception of the regions beyond 10000 km s −1 ), but different electrondensity profiles because of the different 56 Ni profiles.
3.2.3. Combination of the reference 1D CMFGEN model with another model In this approach, we build the axisymmetric ejecta assigning distinct 1D CMFGEN models to the equatorial and to the polar directions, and interpolate between the models at intermediate latitudes. The mapping of a given model over a range of latitudes is flexible. Multiple models could also be used but here we use only two to facilitate the interpretation. We define β 1/2 as the maximum polar angle over which the "polar" model applies.
Hence, besides the reference model 1D-X1, we used an alternate model 1D-X2b. This model deviates from model 1D-X1 both for the progenitor and for the ejecta. The core material within the inner 5 M is made completely homogeneous (this implies strong mixing) and is also given a fixed density (conserving total mass and yields), thereby erasing the jump in density present in model 1D-X1 at the He-core edge. The progenitor is enshrouded by a dense extended atmosphere prior to explosion, which produces a lower maximum ejecta velocity as well as the presence of a dense shell in the outer ejecta. The most important feature is however the presence of a 56 Ni rich shell (in 1D), which we add "by hand" at a Lagrangian mass of 12 M , so that it is located around 4000 km s −1 in the resulting ejecta once homologous expansion is established. In practice, this is done in the radiation-hydrodynamics calculation with V1D (see Livne 1993 andDessart 2019 for discussion), at 1000 s after the explosion trigger (and therefore in mass space since the velocity profile has not settled at that time). Because there is very little mass between 12 M (or the ejecta location where the velocity is about 4000 km s −1 ) and the outermost ejecta layers, an excess in 56 Ni (relative to model 1D-X1) is present throughout the outer ejecta. However, the 56 Ni mass fraction remains very small, typically below 0.01, and the 56 Ni is mixed beyond about 2000 km s −1 with a material rich in H and He (the material that used to be in the H-rich envelope of the progenitor star). Model 1D-X2b corresponds to a (spherically-symmetric) ejecta with a kinetic energy of 1.29 × 10 51 erg, a total mass of 12.43 M , and a 56 Ni mass of 0.047 M . Figure 2 compares the spectral properties of model 1D-X2b with the observations of SN 2012aw in the optical range. Because of the adjustments made, the resulting SN radiation differs sizably from that obtained for model 1D-X1, and thus deviates from SN 2012aw. This was the goal, and is in line with the earlier findings of Chugai (2006). In model 1D-X2b, the stronger mixing and the presence of a 56 Ni shell at 4000 km s −1 (added to give a source of asymmetry in a 2D ejecta; see below) both lead to a stronger and broader Hα line in the second half of the plateau phase. The Hα mismatch suggests that there is no strong 56 Ni mixing along all ejecta-centered directions, but that strong 56 Ni mixing, if present, must be limited to a modest solid angle. The lower 56 Ni mass also causes an underestimate of the optical brightness compared to SN 2012aw. As long as model 1D-X2b covers a small fraction of the whole 2D ejecta volume, the global ejecta and radiation properties should primarily reflect those of model 1D-X1. We compare the ejecta properties and the bolometric light curves of the 1D models 1D-X1 and 1D-X2b in Fig. 3.
In most of the simulations we adopt mirror symmetry with respect to the equatorial plane. Model 1D-X2b is assigned the polar directions and model 1D-X1 is assigned the equatorial direction. The properties at intermediate latitudes is flexible. In the reference setup described in Section 5, model 1D-X2b represents the ejecta for all polar angles between zero and 22.5 deg. Beyond 33.75 deg, we use model 1D-X1 instead. In between these Article number, page 5 of 20 A&A proofs: manuscript no. ms The ejecta asymmetry stems from the asymmetric distribution of 56 Ni, as well as differences in the outer ejecta (through the presence of a dense RSG progenitor atmosphere or not) and differences in the inner ejecta (complete mixing and smoothing of the He core material) -see Fig. 3 and Section 3.2.3 for details. two angles, we linearly interpolate in cos β. This interpolation applies to all relevant quantities for the 2D radiative transfer and in particular the opacities and emissivities. We thus interpolate between two physically consistent 1D CMFGEN models (also derived from physically consistent radiation-hydrodynamics simulations of the explosion). This is superior to using one 1D CMFGEN model and applying some density scaling or radial scaling on its properties (see previous two sections). Figure 4 illustrates some of the properties of the resulting 2D ejecta. The left panel shows the normalized electron density (in the log base 10) at 44 d after explosion. With this setup, we can explore the configuration for an enhanced 56 Ni abundance confined to a narrow range of latitudes. Such a configuration is known to yield a residual polarization (Chugai 2006). This setup is not an exact representation of what may occur in nature. Detailed 3D simulations of neutrino-driven explosions suggest that the 56 Ni material may be asymmetrically distributed in a complicated structure of blobs, shells, fingers, elongated toward large velocities, and made of essentially pure 56 Ni (Wongwathanarat et al. 2015). This material is macroscopically mixed in space during the dynamical phase of the explosion, so that it may coexist at a given velocity with H-rich material but it is not microscopically mixed with it. In our present approach, we instead apply a macroscopic and a microscopic mixing so that both 56 Ni and H-rich material coexist at a given velocity. This approach is perhaps not as problematic as it seems since, with respect to the polarization of light, it is the influence of 56 Ni on the electron density that matters. As time passes, the γ-rays associated with the radioactive decay of 56 Ni and 56 Co are less efficiently trapped locally, and the more so for the 56 Ni advected out to large velocities during the explosion. Thus, the volume influenced by this decay heating extends beyond that occupied by the 56 Ni-rich material. An illustration of this effect is presented in Dessart et al. (2012), in particular their Figs. 1 and 2. Hence, the configuration corresponding to an extended cocoon of enhanced ionization (and thus electron density) around any pure 56 Ni blob is well captured by our setup.

Our approach for polarization modeling
The simulations presented in this study were performed with the long characteristic code LONG_POL first implemented to treat continuum (Hillier 1994) and line (Hillier 1996) polarization in the context of multiscattering axisymmetric envelopes of hot star winds. This code was subsequently extended to treat the case of Type II SN ejecta (Dessart & Hillier 2011) but limited to solving for the polarized flux for an individual line and its overlapping continuum.
The code has since been improved in two important ways (Hillier & Dessart 2020, in prep.). First, the code now computes the full optical polarized spectrum and is thus no longer limited to a single line. This implies that line overlap is automatically treated, irrespective of the number of overlapping lines. Second, the code is no longer limited to 2D ejecta produced by the latitudinal density scaling or radial stretching of a precomputed 1D CMFGEN model. We can now also combine a set of precomputed 1D CMFGEN models, assigning a given model to a specific range of latitudes. The benefit from this latter approach is the improved physical consistency since the opacities and emissivities imported into LONG_POL are those from the precomputed nonLTE 1D CMFGEN models. This is less artificial than prescribing an asymmetry through parametrized scalings. A remaining weakness is that we still interpolate between 1D models at a given latitude. Further, the opacities and emissivities (with the exception of the electron-scattering emissivity) assigned to this 2D ejecta are not computed from a fully consistent 2D nonLTE model. Practically, the 2D radiative transfer is performed as a formal solution (i.e., with the opacities and emissivities fixed).  Note: All observations were obtained at the VLT using the 300 line/mm grism ("GRIS_300V") along with the "GG435" ordersorting filter to prevent second-order contamination. This resulted in a useable spectral range of 4350 Å − 9200 Å. All observations were made through a 1 . 0 slit, which delivered a resolution of ∼ 12 Å, as derived from the FWHM of night-sky lines. The flux standard used for all observations was BD+26 • 2606 (Oke & Gunn 1983), observed on the same nights as the science observations.
All observations were made with the position angle of the spectrograph slit set to 0 deg, which often differed significantly from the parallactic angle (Filippenko 1982); when coupled with the relatively narrow slit used, the potential for differential light loss -and, hence, inaccurate relative flux calibration -is significant. The other benefit is that we can combine any variety of models, differing in explosion energy, composition such as 56 Ni mass or 56 Ni mixing, or clumping.
In LONG_POL, the polarization is produced by electron scattering and described by the Stokes parameters I, Q, U, and V (Chandrasekhar 1960). The polarization is thus linear and V is identically zero. We use the same nomenclature as in Dessart & Hillier (2011) and define the corresponding flux-like polarization quantities F I , F Q , and F U output by the code (the calculations are done in I l , I r , U and V space where I l = 0.5(I + Q), I r = 0.5(I − U); see Hillier 1994 for details) -the definitions of these various quantities are given in Section 2 of Dessart & Hillier (2011). The geometry of the axisymmetric ejecta is described with a right-handed set of unit vectors (ξ X , ξ Y , ξ W ). The origin is the center of the ejecta, the axis of symmetry lies along ξ W , ξ Y is in the plane of the sky, and the observer lies in the XW plane (Hillier 1994(Hillier , 1996. We adopt the same sign convention as in Dessart & Hillier (2011). The flux F Q is positive (negative) when the polarization is parallel (perpendicular) to the symmetry axis. Because of the imposed axial symmetry and the adopted orientation, F U is zero. Hence, the polarization fluxes that we present from our simulations with LONG_POL will be limited to F Q and to P = 100 × |F Q |/F I (and thus corresponding to a percentage), where F I is the total flux. To keep with the convention in spectropolarimetric observations, the normalized Stokes parameters q and u are defined as q = F Q /F I and u = F U /F I . In this study, we only discuss the flux-like polarized quantities obtained with LONG_POL. We refer the reader to Dessart & Hillier (2011) for a discussion of the maps of the polarization intensities on the plane of the sky for various asymmetric configurations and wavelengths. In general, because of the very low level of polarization from asymmetric SN ejecta, these polarization intensity maps are difficult to analyze.
For a LONG_POL calculation, we first construct the 2D ejecta (see Section 3.2). The initial information in this 2D model consists of radius, density, velocity, electron scattering opacity, absorptive opacity as well as emissivity for the selected range of wavelengths. The 2D ejecta is mapped with about one hundred radial points 1 and nine polar angles between pole and equator (equally-spaced; mirror symmetry with respect to the equatorial plane is by default assumed). The radiation field I at each ejecta location (r, β) is discretized in azimuth φ and latitude θ. We use thirteen φ angles while the latitudinal rays are set according to the impact parameter. Fourteen points (or rays) are used to cover uniformly from the ejecta center to the innermost ejecta radius, and the radial grid is used beyond that to define the impact parameter.
In LONG_POL, most of the computing time is taken by a nested β and φ angle loop (in which the 2D radiative transfer equation is solved for). This section is parallelized with MPI and LONG_POL is run with N β N φ processors (hence typically of order a few hundred processors). Convergence is reached after 20 − 30 iterations, although the total and polarized fluxes change little after about ten iterations. This convergence is faster and better for optically thin ejecta since the amount of scattering is in that case much reduced. A fully converged model takes approximately 24h with 200 processors (or 400 processors if mirror symmetry is not adopted).

Observations of SN 2012aw
We obtained seven epochs of spectropolarimetry of the Type II-P supernova SN 2012aw using the European Southern Observatory's Very Large Telescope (VLT) at Cerro Paranal in Chile (Appenzeller et al. 1998   showing the total flux F I at the top, the observed polarization angle θ P (we show the quantity 180 + 1 2 arctan(F U , F Q ), in degrees), and the observed percentage polarization P (we show here the traditional polarization defined as 100 q 2 + u 2 ; see, e.g., Leonard et al. 2001). The quantities θ P and P are displayed with a binning of 15 Å/bin. Each epoch is color coded (see labels and Table 2). To help gauge the reality of spectropolarimetric features, in the bottom plot we show a representative uncertainty per 15Å bin for each epoch (color coded error bars in the upper left corner), taken as the mean 1σ uncertainty in the total polarization calculated from 6900 Å to 7400 Å, a region largely free of line features at all epochs. (the particular uncertainty in any given bin rises and falls with the S/N, becoming larger in the P-Cygni troughs and smaller in the peaks.) We overplot the adopted interstellar polarization based on a Serkowski law (Serkowski et al. 1975;Cikota et al. 2017) derived by assuming that the intrinsic polarization is zero at the regions of the strong, consistent, depolarization seen near the flux emission peaks of Hβ, Na i D, Hα, and the Ca ii near-infrared triplet at epoch 6 (thick black line; see Section 4.2 for details). mounted at the Cassegrain focus of the Antu (UT1) telescope. The data sample days 16 -120 following the date of explosion (16.09 March, 2012 UT) estimated by Bose et al. (2013). Details of the observational epochs, instrumental setup, exposure times, and observational conditions are given in Table 2. The data are presented in Fig. 5.
Each complete observational sequence consisted of four separate exposures with the half-wave retarder plate positioned at each of 0, 45, 22.5, and 67.5 deg; multiple complete sequences were executed at all epochs except for epoch 2, for which only one complete sequence was obtained (see Table 2). All data were bias subtracted and flat-field corrected using dome flats obtained through the polarimetry optics at a waveplate position angle of 0 deg on the same nights as the observations 3 ; reductions carried out with unflatfielded data were found to be virtually identical to those made with the flatfielded data, as expected due to the re-dundant number of half-wave positions obtained (see, e.g., Patat & Romaniello 2006).
One-dimensional sky-subtracted spectra were extracted using both the optimal (Horne 1986) and standard extraction techniques, with an extraction width set to encompass at least three times the effective "seeing" of the particular observations (see Table 2). The wide extraction was used to minimize the effect of spurious polarization features being introduced by interpolating the counts in fractional pixels at the edges of narrow extraction apertures (see Leonard et al. 2001). We subtracted a linear interpolation of the median values of background windows on either side of the object from the object's spectrum; for consistency, the "sky" background region was set for all extractions to be ± 7 . 5 : 2 . 5 from the centroid of the object's spatial profile; an exception was made for the final, seventh epoch, where poor seeing necessitated that a more remote sky region (± 8 . 75 : 3 . 75) be used to assure minimal object light in the sky window. Each spectrum was then wavelength and flux-calibrated, corrected for continuum atmospheric extinction, and rebinned to a 5 Å bin −1 linear scale.
We formed the normalized Stokes q and u parameters and their statistical uncertainties for each complete sequence of observations according to the methods outlined by Miller et al. (1988) and Cohen et al. (1997). For each observational sequence, we derived results using both the optimal and standard extractions separately, and then compared them to each other. In general, the two extraction algorithms yielded virtually identical results, and so the optimal extractions were preferred due to their (slightly) higher signal-to-noise ratio and ability to discount minor cosmic ray strikes. However, in certain spectral regions that contained particularly sharp spectral features that also affected the "sky" region (e.g., regions of strong telluric absorption or, more rarely, severe cosmic-ray strikes across the spatial profile of the object on the CCD), the optimal extraction algorithm produced spurious features in the polarimetry that were not found in the standard extraction results (and, were not consistent from one observational set to the next on the same night). Presumably, such features occur due to a breakdown of the assumption of a "smoothly varying" spatial profile for the object, which underlies the successful application of the optimal weighting algorithm (Horne 1986). This occurred quite frequently at the telluric "A-band" feature near 7600 Å; in one instance, it also occurred in the vicinity of a particularly bad cosmic-ray strike. In these limited regions (typically ∼ < 40 Å), the standard extraction's polarimetric results were inserted; otherwise, the results derived from the optimal algorithm were used.
The fact that multiple, identically obtained, complete observational sequences were taken at nearly all epochs (see Table 2) permitted us to perform valuable tests on the reduced data and calculated uncertainties. To build confidence in our calculated statistical uncertainties, we utilized the data for epoch 3, for which 13 complete, independent, and identically obtained sets were acquired, and subjected them to the following test. First, we computed for each 5Å wavelength bin the 1σ spread of values from all 13 sequences for a Stokes parameter (we took the normalized Stokes q); in principle, this should give a rough estimate of the 1σ uncertainty for each single set of data. (The calculated mean uncertainty derived in this way across the entire spectrum was ∼ 0.25% per bin.) We then compared this with the statistical uncertainty computed through our algorithm for a single observational sequence, and found the two estimates to agree with high precision: The mean difference between the two calculations was less than 0.01% across the spectrum. With confidence in our calculated statistical uncertainties, then, we performed a final check on the veracity of all statistically significant features seen in the polarimetry at all epochs (except epoch 2, for which only one complete set was obtained; see Table 2) by intercomparing the individual sequences' results, carefully examining each one for spurious features not seen in the others.
We then corrected all data for the small ( ∼ < 0.1%) amount of instrumental polarization known to exist at the VLT (Cikota et al. 2017). To correct for the wavelength dependence of the measured polarization angle induced by chromatism of the retarder plates, we applied the chromatic zero angle correction provided by the VLT website 4 to all data. In addition, we checked the absolute zero-point of the derived polarization angle in the sky coordinate system on each night by observing polarized standard stars from the lists of Clemens et al. (1990;HD 127769) and Cikota et al. (2017;Hiltner652;BD-144922;HDE316232;and Vela1). With the exception of the final (seventh) epoch, all polarized standard stars' derived polarization angles (in the Vband; nominal wavelength range of 5050 − 5950 Å) agreed with the cataloged values to within 0.6 deg. On the seventh epoch, the polarization angle derived for HD 127769 differed from the value reported by Clemens et al. (1990) by 2.4 deg; the cause for this discrepancy is not known, although we do note the unusually poor seeing (Table 2) and short exposure time for the star (one second per waveplate position; more typical of other nights was 3 − 5 seconds) on this night. In all cases, including this final epoch, we kept the nominal polarization angle derived directly from the data, and did not apply any overall zero-point offset based on our observations of polarized standard stars. To check for instrumental polarization beyond that already accounted for through application of the linear relations provided by Cikota et al. (2017), we also observed the null polarization standard HD 109055 (Clemens et al. 1990;Berdyugin et al. 1995) at every epoch, and always found it to be null to within 0.05%.
For each epoch, the normalized Stokes q and u parameters derived from all complete observational sequences were combined according to their statistical weights, with new uncertainties calculated. We then removed the interstellar polarization derived in Section 4.2, removed a redshift of 778 km s −1 (de Vaucouleurs et al. 1991, via NED), rebinned to 5 Å bin −1 , and created a final, rest-frame dataset with a wavelength range of 4350 Å -9150 Å. We also computed results binned to 15, 25, 50 and 100 Å bin −1 for presentation purposes. To generate the total flux spectrum for each epoch, a final correction for telluric absorption bands (Wade & Horne 1988) was made.
To illustrate the continuum polarization as a function of time we show in the left panel of Fig. 6 the normalized Stokes q and u measurements for all seven epochs in three spectral regions devoid of strong lines around 6000, 7000, and 8000 Å. All measurements lie in a single quadrant of the (q, u) plane, and to first order lie along a straight line. The biggest departure from a straight line occurs for the data at 120 d. The offset is significant, although the data at this epoch also has the worst signalto-noise ratio. If SN 2012aw were axisymmetric all polarization data would lie along a straight line in the (q, u) plane, and, in the absence of interstellar polarization, would pass through the origin. The presence of a dominant axis in the data sets indicates that SN 2012aw must possess a strong degree of axial symmetry. However, the scatter (and in some cases one could argue for an epoch dependent offset) of the polarization observations about the straight line indicate the structure cannot be simply axisymmetric.
A study of the polarization behavior for Hα (right panel of Figure 6) yields a similar conclusion to that reached using the continuum polarization. Again, all measurements lie in a single quadrant of the (q, u) plane, and to first order lie along a straight line. As for the continuum, it is the observations at 120 d that show the largest offset. Because of the presence of a dominant axis, we rotate the observations so that the bulk of the polarization lies in q (Section 4.3). This has the consequence that sign reversals in q can easily be seen in the same plot, whereas such reversals are not seen in P and imply a 90-deg change in the angle of polarization.

Interstellar polarization
Correcting for interstellar polarization is a delicate problem with no perfect solution. In this study, we attempted two different techniques that rely on quite different theoretical assumptions.
Our first approach was to assume that the intrinsic polarization of the SN is zero at the first epoch; thus, all observed polarization should be due to the ISP. Theoretically, this may be motivated by the fact that instabilities, born in the progenitor core during the explosion or at the H/He interface after shock passage, influence little the outer ejecta. Unfortunately, there is clearly some variation in polarization across lines even at the earliest epoch, incompatible with the expected smooth featureless polarization from the interstellar medium. Nonetheless, we fit a spline to the "continuum" regions (i.e., those away from strong line features) as a first estimate of the ISP.
For our second approach we considered the theoretical prediction that during the optically thick, "photospheric" phase of an SN II's development (up through about day ∼ 115 for SN 2012aw, when its light curve begins to drop off of the plateau; Bose et al. 2013), broad, unblended emission features (e.g., Hα, the Ca ii near-infrared triplet) may possess zero intrinsic SN polarization (e.g., Trammell et al. 1993;Hoeflich et al. 1996;Dessart & Hillier 2011). If complete depolarization of all SN light (line and continuum) has occurred in such spectral regions -a further assumption commonly employed to estimate the total ISP in SN studies (e.g., Howell et al. 2001;Leonard et al. 2006;Nagao et al. 2019) as well as polarization studies of emission-line stars (e.g., Meyer et al. 2002) -then any measured polarization there is due solely to ISP. For multiepoch data, this should present as fixed, unchanging levels of polarization (and, polarization angle) with time in these spectral regions, since ISP does not change appreciably over the timescales involved whereas intrinsic SN polarization may.
Examination of Fig. 5 shows that the expectation is quite well realized empirically for SN 2012aw, where sharp, consistent, depolarizations are seen for several strong lines, including Hα, Hβ, Na i D, and the Ca ii near-infrared triplet. That is, while nearby regions of the spectrum vary by up to ∆P ∼ 1.5 % through the observational sequence, these spectral regions remain consistent to within ∆P ∼ 0.1 %. Motivated by this consistency, we fit a Serkowski law (Serkowski et al. 1975) through the strong depolarizations seen in Hα, Hβ, Na i D, and the Ca ii near-infrared triplet to the epoch 6 data (similar results were found using epochs 3 to 7). Using these four spectral locations as a reference for the interstellar polarization, we searched for the set of values for λ max,ISP , q max,ISP , and u max,ISP that produced the best fitting curve given by a Serkowski law. Each q max,ISP , and u max,ISP correspond to θ ISP and p max,ISP through θ ISP = 0.5 arctan(u max,ISP /q max,ISP ) (2) and p max,ISP = q max,ISP cos(2θ ISP ) + u max,ISP sin(2θ ISP ).
While both estimates of the ISP have weaknesses, we encouragingly found them to agree to better than ∼ 0.2% over the full spectral range covered by our observations. We ultimately adopted the ISP described by our second approach and proceeded to correct all spectropolarimetric observations of SN Article number, page 10 of 20 Luc Dessart et al.: Modeling of the SN 2012aw spectropolarimetric obsvervations 2012aw by it. However, given the difference between our estimates, we expressly recognize that a ∼ 0.2 % systematic uncertainty in the absolute level of all calculated "intrinsic" polarizations exists in all that follows. 5 Fortunately, since the ISP is a smooth function of wavelength, this uncertainty in its absolute value has no effect on interpretation of specific line features, which is the main focus of this work. Further, as the degree of continuum polarization ultimately seen in SN 2012aw approaches ∼ 1.2 % in the latest epochs, the effects of a 0.2 % uncertainty do not qualitatively alter our fundamental conclusions regarding the overall polarization level when comparing to our models.
Adopting this interstellar polarization and focusing on the spectral range between 6900 and 7200 Å (Fig. 5), SN 2012aw yields a nonzero polarization even at the earliest times, although interpretation of this is tempered by the fact that its value, ∼ 0.2 %, is of order our systematic uncertainty in the ISP. SN 2012aw then shows a strong temporal increase in polarization during the second half of the plateau phase, reaching a maximum of about 1.2% as the SN starts plunging in brightness and transitioning to the nebular phase. This has been seen in SN 2004dj (Leonard et al. 2006 (Leonard et al. 2015). SN 2012aw could not be observed at nebular times as it went behind the sun. It was too faint when it later emerged.

Rotation of the Stokes fluxes
For an axially symmetric ejecta the polarization (at all times and at all wavelengths) lies on a diagonal line in the (q, u) plane which passes through the origin. Since the slope of this line simply reflects the orientation of the ejecta symmetry axis on the plane of the sky, one can rotate it so that all the polarization lies in q. This then allows a direct comparison with our model calculations in which all the polarization lies in q. In observations, noise and departures from strict axial symmetry leave a residual signal in u even after such a rotation has been applied.
Hence, when we compare to our polarization results for axisymmetric models, we apply a rotation of the observed Stokes vectors by an angle θ rot = 0.5 arctan(δu/δq), where δu and δq are measured in the (q, u) plane. This is done to put most of the observed polarization in q. The rotated fluxes are then given by q rot = q cos(2θ rot ) + u sin(2θ rot ) and u rot = −q sin(2θ rot ) + u cos(2θ rot ).
The data of SN 2012aw, corrected for interstellar polarization, yield an observed polarization angle of 118 deg in the plane of the sky. However, we choose to rotate the Stokes vectors by the angle orthogonal to this -that is, of 28 deg -in order to produce a sign of the rotated Stokes vectors that agrees with the models, for ease of comparison (rotating through an angle 90 deg different simply changes the sign of the polarization). With this rotation, the bulk of the polarization is now in q, and u is generally close to zero. This works with only modest success at early times when the polarization level is low, but it seems suitable for epochs four and later, as shown in the next section.

Modeling results and comparisons to observations of SN 2012aw
As discussed in Section 3.2, we used three different ways to introduce the asymmetry in our 2D ejecta. In the first two sections below, we describe the results for the first epoch adopting a latitudinal density scaling (Section 3.2.1) or a latitudinal stretching (Section 3.2.2) of a 1D CMFGEN model. In the subsequent section, we use the combination of two models (Section 3.2.3) and compare to all observed epochs of SN 2012aw (Section 4). For clarity, we call these models the scaled model, the stretched model, and the hybrid model. Throughout this section, we present multipanel figures that have the same structure. The top panel shows the total flux F I for both the model and the observations (normalized to unity at 7100 Å). The next two panels show the normalized, rotated Stokes parameters q rot and u rot (we plot the percentage value), also corrected for interstellar polarization (see previous section).
Finally, the bottom panel shows the percentage polarization P for both the observation (black) and the model (colors other than black). The model results are sometimes shown for all inclinations. The inclination of zero degree yields zero polarization. The other inclinations are spaced uniformly every 11.25 deg, colored from blue to red (colorblind rainbow colormap) as the inclination is increased.
Compared with the simulations presented in Dessart & Hillier (2011), the present simulations provide the total and polarized flux spectrum over the entire optical range. We explicitly account for line overlap or blanketing and can directly compare to the optical VLT − FORS spectropolarimetric data that we collected for SN 2012aw.  . 1), corresponding to a density contrast at a given radius of 1.5 and 5.0 between pole and equator (prolate density configurations). The simulations assume mirror symmetry with respect to the equatorial plane. For both models, the asymmetry introduces only a minor change with inclination for the total flux F I (top panel). The polarized flux F Q (or equivalently q) is positive throughout most of the optical range, which means that the electric vector is parallel to the symmetry axis. This is an optical depth effect because in the optically-thin regime, F Q would be negative in this case. Switching to an oblate configuration would result in a sign change of the polarization (rotation of the polarization angle by 90 deg). In absolute terms, the polarized flux is here greater where the flux is greater. The percentage polarization is also greater at shorter wavelengths, albeit modestly. Theoretically, one may expect a greater polarization at shorter wavelength because of the greater albedo toward the Balmer edge (Dessart & Hillier 2011), especially at this young SN age when the material in the spectrum formation region is at least partially ionized. The polarization changes between lines and continuum, as well as within the lines. We see a sign reversal of the polarization in Hα, and for many lines in the case of a stronger asymmetry (right panel of Fig. 7). For a large scaling (A = 4.0), the model P rises to a maximum of 0.6% in the largely line-free region between Fe ii 5169 Å and Hα.

Results with the latitudinal density scaling at early times
The observed polarization of SN 2012aw at the corresponding epoch differs in a number of ways from the model predic-A&A proofs: manuscript no. ms . The ordinate range is the same for q rot and u rot . A binning of 25 Å is used for q rot , u rot , and P, while a binning of 5 Å is used for F I . As described in the caption to Fig. 5, a representative 1σ uncertainty in the polarization was calculated, and is shown in the bottom panels. The polarization displayed in the bottom panels for the SN 2012aw data is the "traditional polarization", as defined in the caption to Fig. 5. For all models u rot is null across the spectrum, and so is not plotted.
tions. It is greater at redder wavelengths, although at this epoch the polarization is quite low, and only about 0.2% at 7000 Å, thus significantly lower that the maximum value obtained with the models discussed above (about 0.6%). Observations also reveal a sizable change in polarization across Hβ and Hα, but not other lines (perhaps because the data are too noisy). While Hβ and Hα in the models also show polarization, the structure does not match the observations. In q the model shows an enhanced polarization in the P Cygni absorption component (not seen in the observations), and the polarization absorption trough is offset from that observed. In P, the models show a very complicated structure which contrasts with the relatively smooth structure seen in the observations. The difficulty in interpreting this epoch is exacerbated by the uncertainties in the interstellar polarization, the low level of polarization, and the significant signal that remains in the rotated u parameter.

Results with the latitudinal stretching at early times
For an asymmetric ejecta resulting from a radial stretching of the 1D CMFGEN model 1D-X1 (see section 3.2.2 for explanations), the results are quite different from models with a latitudinal density scaling (Fig. 8). The changes to the total flux are greater than for the models with the density scaling. The flux F I at 7100 Å increases by 40% (90%) as we vary the inclination from 0 to 90 deg in model 2D-X1-STR0p25 (2D-X1-STR0p50; see Sec-tion 3.2.2 for details on the nomenclature). Once we renormalize the spectra at 7100 Å, the variation with inclination appears mostly in lines in the form of a wavelength shift. This shift is more noticeable in model 2D-X1-STR0p50. The polarization is also overall lower. Varying the magnitude of the radial stretching changes little the resulting polarization, despite the large density contrast between pole and equator (corresponding to values of 14.6 and 130.0 for models 2D-X1-STR0p25 and 2D-X1-STR0p50 since the ejecta follow a power-law density with exponent twelve in the outer regions). In both cases, we see strong sign reversals of the continuum polarization when the inclination is changed. There are also strong reversals of the polarization across line profiles for inclinations corresponding to low latitudes (near equator-on views). Such signatures are likely due to optical depth effects with a complicated cancellation of the local polarization on the plane of the sky. These can occur at early times because asymmetry impacts the distribution on the plane of the sky of both the scatterers and the flux. Hence, for the case of a latitudinal stretching, the model fails to reproduce the greater polarization at longer wavelengths, is unable to reproduce the level of polarization observed, and fails to match the structure across Hα.  Fig. 9. Evolution of the shape factor γ(r) in the hybrid model versus velocity for the eight epochs modeled (the last epoch at 160 d has no observational counterpart). The symbols correspond to the photospheric radii in the equatorial direction (square) and the polar direction (dot).

Results with a bipolar explosion using a combination of 1D CMFGEN models
In this section, we use two distinct 1D CMFGEN models to build an asymmetric ejecta. The ejecta are symmetric with respect to the equatorial plane so it corresponds to a bipolar explosion in which a cone of opening angle 50 deg along the symmetry axis is characterized by distinct properties relative to lower latitudes. In our setup, the difference is the steeper density fall-off beyond about 10,000 km s −1 and the larger 56 Ni content within the H-rich lay-ers of the ejecta (see Fig. 3). In the context of polarization, the main impact of the 56 Ni enhancement is not the change in opacity or abundance (both have a minor influence) but instead the boost it induces to the density of free electrons (the mass density below about 10,000 km s −1 hardly changes with angle, so the impact is on the ionization). Hence, there is at least one source of asymmetry at any given velocity, which implies that this setup may produce a net polarization at any epoch. The initial conditions used in this case are presented in detail in Section 3.2.3, to which we refer the reader for additional information (see also Fig. 3 −4).
One way to characterize the asymmetry of our ejecta is by comparing the photospheric radius along the pole and along the equator. Although our hybrid model is made of the same combination of 1D ejecta at all epochs computed, the gas properties evolve and in a distinct way along the polar and the equatorial directions. The 56 Ni enhancement affects mostly the innermost ejecta layers, so the asymmetry should be weak early on and grow as the inner aspherical ejecta layers are revealed by the receding photosphere. Consequently, in our hybrid model, the ratio of polar and equatorial photospheric radii evolves from 1.02, to 1.08, 1.22, 1.36, 1.59, 1.52 and 1.26 as the SN ages from 16.1, to 46.0, 63.0, 71.9, 91.9, 107.9, and 120.0 d.
Another method for illustrating the structure is with the shape factor defined by Brown & McLean (1977) which (together with the optical depth and the inclination) fully describes the continuum polarization properties of an axially-symmetric ejecta illuminated by a point source. In this paper, we introduce   Fig. 7, but now for the hybrid model 2D-X1-X2b (for an inclination of 67.5 deg) and the observations of SN 2012aw for epochs 1, 2, 3, and 4. a depth-dependent shape factor γ(r) and define it as where N e (r, µ) is the free-electron density at (r, µ), and µ is the cosine of the polar angle β). When r = R min we obtain the expression of the shape factor as in Brown & McLean (1977). Because of the complicated dependence of the emissivities and opacities on depth, no single parameter can fully describe the asymmetric  nature of the ejecta. The shape factor (defined as in Eq. 8) for our hybrid model is shown for all epochs in Fig. 9. Thus revealed, the level of asymmetry of the free-electron density distribution is modest at early times, but continuously grows as time progresses. This arises because the ionization level is naturally higher early on, irrespective of decay heating. 6 However, when the ejecta recombines, the excess decay heating around 4000 km s −1 and within the opening angle β 1/2 generates a local enhancement in free-electron density. The contrast grows in this region at all epochs and should boost the polarization level as time passes. Another striking property is that the asymmetry, as revealed by the shape factor, is really confined within the regions beyond 2000 km s −1 , which are H rich. In contrast, the metal-rich core is moderately asymmetric.
In Figs. 10 − 11, we compare the multiepoch spectropolarimetric observations of SN 2012aw with the results of the 2D polarized radiative transfer code based on this new axisymmetric ejecta configuration (the figure layout is the same as that used in Fig. 7). Figure 10 shows the results for the first four epochs (at 16.1, 46.0, 63.0, and 71.9 d after explosion), and Fig. 11 shows the results for the last three epochs (91.9, 107.9, and 120.0 d after explosion) together with the predictions of the model at 160 d (for which there is no data). The later epochs are better suited for analyzing the SN polarization since the polarization level is higher and the data has a higher signal to noise ratio. The variations across lines are better resolved. For this comparison, we selected the inclination of 67.5 deg since it yielded a good match to the observations at all epochs except the first one.
At the first two epochs, the observed polarization is low, about 0.1 − 0.3 %. This polarization level is roughly matched by the hybrid model but the rise of the polarization at longer wavelength (more obvious at the first epoch) is not predicted (as for the previous calculations described in Sections 5.1 and 5.2). This rise in the red toward the Paschen jump is hard to comprehend as there are no obvious discontinuities in this region. If anything, the albedo decreases smoothly so the polarization is expected to drop as scattering weakens toward the red. This peculiar feature may arise from a cancellation effect due to the competing influence of flux and optical depth (i.e., between the asymmetry of the emission and of the scatterers as seen on the plane of the sky).
At later epochs (63 to 107.9 d), the hybrid model captures the qualitative and the quantitative evolution observed in SN 2012aw, both for the total flux (and thus photometry) and for the polarized flux. This is obtained without any adjustment to the ejecta structure so that the change in polarization arises from the intrinsic evolution of the 2D ejecta as they expand and cool etc.
Across the profile of a strong line, the polarization is expected to show a jump in the blue part of the P-Cygni profile trough. This is observed in Na i D and the Ca ii near-infrared triplet, but rarely in other lines. In contrast, the model shows this as a narrow feature in most lines, including Hα. The difficulty of observing this feature is aggravated by the rebinning of the data, which was used to enhance the signal-to-noise ratio per bin. The feature originates from the deficit of unscattered (direct and thus unpolarized) flux from the SN while the residual flux arises from photons scattered back into the line of sight. In an asymmetric ejecta, such scattered photons yield a net polarization. In the Ca ii near-infrared triplet, the model shows the polarization jump further to the blue than in the observations. Close inspection of the total flux F I shows that the model overestimates the width of the Ca ii near-infrared triplet line at late epochs (more so at epochs 5, 6 and 7), so the mismatch in q or P is probably rooted in the slight overestimate by the model of the ejecta expansion rate.
In strong lines, the polarization is zero or close to zero somewhere between the P-Cygni trough and the location of maximum flux. This feature of line depolarization was used earlier to constrain the interstellar polarization toward the line of sight to SN 2012aw. It is clearly observed at all epochs (more so for epochs three to seven) in Hβ, Na i D, Hα, and the Ca ii nearinfrared triplet. It is also reproduced by the model. This results from the depolarization of line photons. Recombination lines emit photons isotropically and thus cannot produce a net polarization. The large line optical depth also leads to absorption of possibly polarized photons, thereby destroying the polarization they carried. The line photons that scatter with free electrons within the ejecta can however produce a residual polarization. However, this is associated with a redshift (because of the large bulk velocity of the scatterers) so this polarized flux appears in the red wing of each line. 7 As we move to the red wing of lines, the polarized flux is large, and larger for stronger lines. This arises from line photons scattered by free electrons. The effect yields an excess total flux in the red wing of all lines. Because of the asymmetry of the distribution of free electrons, it also yields a clear polarized flux F Q , whose strength scales with F I . However, P is at the same level in the red wing of lines and in the adjacent continuum. This suggests that line photons and continuum photons pick up the same level of polarization when they last scatter in the ejecta (and that statistically they bear the same level of polarization before that last scattering). The observed polarization behavior across lines is quite distinct in type Ia SNe (the metal-rich composition implies a lower number of free electrons per unit mass and a reduced importance of electron scattering for the opacity), which show a maximum polarization within the absorption trough of the strongest lines and a low level of polarization elsewhere (see, e.g., Wang et al. 2003).
In general there is good agreement between the observed and theoretical line profiles, especially after the first two epochs. One discrepancy is that the model profiles exhibit a larger variation in polarization at the blue edge of the P Cygni absorption than do the observations. Particularly striking is the agreement for the polarization across the O i 7774 Å multiplet. In the flux spectrum, O i 7774 Å is an inconspicuous spectral feature. However, in both the observed and model polarization (after the first two epochs) it is conspicuous, and is only slightly less prominent than the feature due to Hα. In the continuum, the level of polarization increases with wavelength in the optical range, although this effect, which is matched by the model, is quite weak. At late times (epochs four to seven), in the region of line blanketing shortward of 5000 Å, the model predicts no clean (line free) continuum region so strong line depolarization occurs. In the red part of the spectrum, and in particular between Hα and O i 7774 Å, the continuum polarization shows a smooth and near-constant magnitude both in the observations and in the model (the model shows fluctuations associated with weak lines, such as the Ca ii 7300 Å doublet). This is the spectral region where it is maximum at any given epoch.
Other discrepancies convey some important information on the source of polarization and departures from what is adopted in the hybrid model. Clearly visible at epochs 5 and 6 is the overestimate of the polarized flux associated with the scattered and red-shifted Hα photons, appearing as a bump in q or P. In contrast, the observations suggest the same level of polarization for these photons and the adjacent continuum-only photons further to the red. This mismatch is, however, only so obvious in Hα and not, for example, in Na i D. The origin of the discrepancy might be an inadequate 56 Ni distribution which causes a too extended  boost to the electron density (for example, a more localized enhancement in 56 Ni at large velocity may be more suitable).
Another important discrepancy, which could be related to the point discussed in the previous paragraph, is the lack of model polarization in the blue part (say < 5500 Å) of the spectrum at the recombination epoch (epochs 3 and later). The observations reveal a lower polarization than around 7000 or around 8000 Å, but only by about half. The model shows low or zero polarization because the polarized photons from the inner asymmetric ejecta are destroyed by interactions with lines (these regions are strongly affected by line blanketing due to Fe ii and Ti ii). One simple way of resolving this issue is to invoke a more distant, that is a more external scattering source, for example in the form a higher velocity 56 Ni enhancement. In this context, this external scatterer would yield a polarization that would scale with the number of incoming photons. The polarization in the blue could then be slightly smaller than at long wavelength because of the residual influence of line opacity at large velocities in the ejecta.
In the comparisons shown in Figs. 10 and 11, the observed polarization (second and third panels from top) was rotated so that the bulk of the polarization lies in the flux F Q . Doing this, the flux F U is close to zero and shows only modest variations with wavelength, in contrast to what is observed in F Q (this is more easily seen at late epochs when the polarization in linefree regions reaches 0.5 to 1.0%). At all epochs except the first one, the flux F Q is negative at all wavelengths in both the model and the observations. Together with the fact that F U is close to zero and much smaller (in magnitude) than F Q , this property suggests that the ejecta is primarily axisymmetric. This property is also evident from the data distribution in the (q, u) plane shown in Fig. 6. Figure 12 summarizes the evolution of the continuum polarization (taken from the averaged value over the range 6900 to 7200 Å) as well as the V-band brightness for both SN 2012aw and for the hybrid model discussed in this section. The model matches well the light curve for the plateau duration and brightness, and for the tail brightness. This occurs because the asymmetry has little effect on the total observed flux -the hybrid model displays the same light curve and spectral evolution as the spherically symmetric model 1D-X1 (with the exception of a slight overestimate of the Hα line strength and width). The hybrid model also matches the low level of polarization at early times and its steep rise in the second half of the plateau. The hybrid model then predicts that the continuum polarization will drop at late times, essentially following the drop in optical depth. Under optically-thin conditions, the polarization scales linearly with τ (Brown & McLean 1977) and in SN τ drops as 1/t 2 at late times (assuming constant ionization). Such an evolution has been observed over a sizeable timeline (Leonard et al. 2006) and explained by polarization modeling (Dessart & Hillier 2011).

Uniqueness of the solution
In the preceding section, we have explored various configurations for the aspherical (2D axisymmetric) ejecta. Here, we study the sensitivity of our results to changes in the 2D ejecta structure in the hybrid model. The hybrid model used model 1D-X2b for the polar latitudes (corresponding to a cone with an opening an-Article number, page 17 of 20 A&A proofs: manuscript no. ms gle of 50 − 60 deg), and model 1D-X1 for other latitudes. In this section, we describe the results when the half opening angle of the cone is increased (Section 6.1), when we drop the assumption of mirror symmetry (equivalent to comparing a bipolar and a unipolar explosion; Section 6.2), and when the explosion energy of the ejecta in the polar regions is increased (Section 6.3).
6.1. Variation with cone opening angle in the hybrid model Figure 13 shows the influence of the half opening angle β 1/2 on the radiative signatures for model 2D-X1-X2b at 107 d. In this particular setup, model X2b extends from the pole up to β 1/2 of 11.25, 22.5, or 33.75 deg. For an inclination of 67.5 deg, we find that the normalized polarized flux (middle panel of Fig. 13) is qualitatively independent of β 1/2 , while quantitatively, the percentage polarization increases with β 1/2 (bottom panel of Fig. 13). For lower inclinations, optical depth effects cause a sign flip (rotation by 90 deg) of the polarization F Q for different values of β 1/2 . Such optical depth effects complicate the analysis of the observed polarization during the photospheric phase. Furthermore, the polarization increase with increasing β 1/2 occurs up to a half opening angle of about 60 deg, beyond which the configuration is somewhat equivalent to swapping model X1 and X2b for a half opening angle of π/2 − β 1/2 . Figure 14 shows the results for the hybrid model 2D-X1-X2b with and without the assumption of mirror symmetry. The angle β 1/2 is 22.5 deg and the inclination is 67.5 deg with respect to the symmetry axis. With mirror symmetry, the ejecta has a bipolar morphology, with the properties of model X2b along the poles. Without mirror symmetry, the ejecta is unipolar so that model X2b only covers polar angles smaller than 22.5 deg. In Fig. 14, we adopt an inclination of 67.5 deg. In this case, the  model with mirror symmetry (i.e., bipolar explosion) shows a percentage polarization in the 7000 Å region that is more than twice greater than when mirror symmetry is ignored (i.e., unipolar explosion). The offset differs from a factor of two because of optical depth effects. For example, for an inclination of 10 deg, the percentage polarization is the same in both cases because the asymmetry from the receding part of the 2D ejecta is obscured. For an inclination of 45 deg, the level of polarization is very different (sign flip, different absolute polarization level) between the two cases, again because of optical depth effects. For an inclination of 90 deg (edge on), the percentage polarization in the 7000 Å region is exactly twice larger if mirror symmetry is assumed rather than ignored. In the mirror-symmetry case, each hemisphere contributes the same residual polarization and the same flux at any wavelength. Without mirror symmetry, only one hemisphere contributes, yielding exactly half the continuum polarization level obtained with mirror symmetry. Finally, whatever the inclination, the total flux F I is essentially identical in both cases. Figure 15 compares the ejecta and radiative properties of the hybrid models 2D-X1-X2b and 2D-X1-X2B, which differ in that the polar regions have an enhanced 56 Ni mass fraction (former model) or an enhanced explosion energy (latter model). Although these two asymmetric ejecta are produced by very different means, they yield very similar total and polarized fluxes. The reason is that in the polar direction, the higher explosion energy leads to a greater density at large velocities. Although the ionization is the same along all latitudes (because of the similar influence of the deeply embedded 56 Ni), the electron density varies with latitude, reflecting the variation in mass density (see left panel of Fig. 15). This polar enhancement in electron density (not caused by an excess in 56 Ni), produces the continuum polarization in model 2D-X1-X2B. In contrast, the mass density in model 2D-X1-X2b is essentially constant with latitude but the greater 56 Ni mass fraction along the pole yields a greater ionization and thus a greater density of free electrons along the poles than at lower latitudes. The results in Fig. 15 indicate that such distinct physical processes (polar enhanced 56 Ni mixing or polar enhanced explosion energy) yield a very similar polarization signal. This emphasizes the degeneracy of polarization signatures.

Conclusion
We have presented VLT − FORS spectropolarimetric observations of the type II SN 2012aw spanning seven epochs of the photospheric phase until the onset of the drop-off from the plateau. Unfortunately, no spectropolarimetric observation of SN 2012aw was made during the nebular phase. SN 2012aw presents polarization characteristics that are qualitatively and quantitatively similar to what has been observed in type II SNe so far, such as the emblematic SN 2004dj (Leonard et al. 2006). A generic property of Type II SNe, corroborated here with 2012aw, is the relatively low polarization early in the plateau phase (although see the counter-example of SN 2013ej; Leonard et al. 2015;Mauerhan et al. 2017Nagao et al. 2019) and the progressive rise to a maximum polarization of about 1 % as the ejecta turns optically thin at the end of the high-brightness phase. By rotation of the Stokes vectors, it is possible to place most of the SN 2012aw polarization along a single, fixed polarization axis, which implies that the ejecta of SN 2012aw has a dominant axis and that the polarization arises from ejecta with an oblate or prolate morphology.
We have modeled the SN 2012aw linear polarization using a long characteristic code that assumes axial symmetry (Hillier 1994(Hillier , 1996Dessart & Hillier 2011;Hillier & Dessart 2021). The code provides the full optical total and polarized flux F Q (for symmetry and geometry reasons, F U is zero) at any epoch during the photospheric or nebular phase. At present, the favored operation mode for the code is to build a 2D axially-symmetric ejecta using two 1D CMFGEN models and assigning them a specific range of latitudes. Such a hybrid 2D model can thus mimic an explosion with a higher energy or a higher 56 Ni mass fraction along specified latitudes.
Our modeling results support the notion that the intrinsic polarization is negligible within strong lines somewhere between the location of maximum absorption and the location of maximum line flux. We use this property to correct for the interstellar polarization. With this choice, the intrinsic polarization of SN 2012aw is small but non zero at the earliest epochs.
The parameter space for producing 2D hybrid models is extended since we may use any 1D CMFGEN from the large set of models that we have calculated for type II SNe. We also have freedom when assigning a given model and a given latitude. Having selected two models and a specific model-latitude assignment to produce a 2D ejecta, the same setup is used for all epochs. No further adjustment is made to the setup during a sequence of polarization calculations. In the present paper, we focused on one hybrid model composed of model 1D-X1 (also named x1p5 in Hillier & Dessart 2019) and of model 1D-X2b, which deviates from model 1D-X1 by having enhanced 56 Ni mixing at large velocities. For an inclination of ∼ 70 deg, we find that this hybrid model can reproduce the SN 2012aw polarization at all epochs except the first two, without any adjustment. The model predicts the polarization is maximum in the line free region between Hα and O i 7774 Å, and also reproduces the variation in polarization across line profiles. The model predicts negligible polarization in line blanketed regions, while the observations suggest a nonzero residual polarization. However, the wavelength and time variation of the polarization is reproduced by the model, reflecting the internal changes in the corresponding ejecta (expansion, recession of the photosphere, cooling, recombination, reduction in density etc). As is well known (e.g., Jeffery 1991a), we also find that the polarization signatures are degenerate so other configurations can also reproduce the observations. 2D-X1-X2b at 107d and 67.5deg 2D-X1-X2B at 107d and 67.5deg Fig. 15. Left: Comparison of the mass density, 56 Ni mass fraction, and electron density for the 1D models used to build the 2D ejecta 2D-X1-X2b and 2D-X1-X2B. Right: Same as Fig. 13, but now comparing the results for the hybrid model 2D-X1-X2b (bipolar ejecta with an enhanced 56 Ni mass fraction along the polar regions) and the hybrid model 2D-X1-X2B (bipolar ejecta with an enhanced explosion energy along the polar regions). In both models the angle β 1/2 is 22.5 deg. [See Section 6.3 for discussion.]