Energy Conservation in the Thin Layer Approximation: II. The Asymmetric Classic Case for Supernovae Remnant ()
1. Introduction
The asymmetries observed in supernovae (SNs) and supernova remnants (SNRs) have been recently investigated by different approaches: on analysing a type Ia SN that exploded around CE 1900 [1] reported that it is asymmetric in the radio region but shows a bipolar morphology in the X region. The details of the deceleration of that SN in different directions were analysed by [2]. The symmetry of many SNRs in 6-cm and 20-cm Very Large Array images was measured by [3]. A three dimensional smoothed particle hydrodynamic simulation was used by [4] in order to explain the observed asymmetries in SNs. [5] implemented a three-dimensional MHD numerical code in which the initial mass distribution of the SN is asymmetric. The presence of an asymmetric circumstellar medium (CSM) responsible for the lack of narrow lines within the first two days of the explosion for the Type II-P supernova SN 2017gmr was suggested by [6]. The differences in the Si II line widths in type Ia supernovae with an asymmetric CSM were explained by [7]. An asymmetric explosion as seen from different viewing angles was suggested by [8] in order to explain some anomalies in the photospheric velocities. The model of energy conservation in the thin layer approximation for the expansion of SNRs was developed in [9] in the spherical case. Here we analyse the asymmetrical or non-spherical case. In order to accomplish this, Section 2 derives the basic differential equations which regulate the equation of motion for SNRs for four types of medium. Section 3 applies those analytical and numerical results to SN 1987A and SN 1006 and Section 4 analyses two models of the image formation for the two SNRs analysed.
2. Energy Conservation
A point in Cartesian coordinates is characterized by x, y and z, and the position of the origin is the centre of the explosive phenomena. The same point in spherical coordinates is characterized by the radial distance
, the polar angle
, and the azimuthal angle
. In the following, the profiles of the density considered are independent of the azimuthal angle and are functions of the distance z from the
plane. In spherical coordinates, the goal is to derive the instantaneous radius of expansion, r, as a function of the polar angle. The following approaches will model the CSM around the point of explosion of the SN which is characterized by a plane,
, which produces an up and down symmetry in the polar angle,
, and by a polar axis,
, which produces a right and left symmetry for the azimuthal angle
. The basic symmetries are now outlined
(1a)
(1b)
In other words, the numerical simulations will be developed in the first quadrant (
,
and then applied to the other three quadrants.
The conservation of kinetic energy in spherical coordinates along the solid angle
in the framework of the thin layer approximation is
(2)
where
and
are the swept masses at
and r, while
and v are the velocities of the thin layer at
and r. The above conservation law when written as a differential equation is
(3)
The velocity as a function of the radius is
(4)
In the following,
is the radius after which the density starts to decrease, r is the radius of expansion in spherical coordinates,
is the density at
,
is the Cartesian coordinate along the Z axis and
is the polar angle in spherical coordinates. The main astrophysical assumption adopted here is that the density of the CSM decreases as a function of the distance from the centre, due to the previous stellar winds. At the time of writing, there are no astronomical observations which outline such an effect, and therefore we test different profiles for the density and evaluate whether they are compatible with the observed sections of SNRs.
2.1. An Inverse Square Profile of the Density
We assume that the medium around the SN scales with the axial piecewise dependence
(5)
The mass
swept in the interval
is
The total mass
swept in the interval [0, r] is
The positive solution of Equation (2) gives the velocity as a function of the radius:
(6)
The differential equation which models the energy conservation is
(7)
The analytical solution of the above differential equation is
(8)
The above inverse square profile for density satisfies the basic symmetries as outlined in Equations (1a) and (1b).
2.2. Back Reaction for the Inverse Square Profile
The radiative losses per unit length are assumed to scale as
(9)
where
is a constant and
is the density in the thin advancing layer which is
. Inserting into the above equation the velocity to first order as given by Equation (6), the radiative losses
are
(10)
The sum of the radiative losses between
and r is given by the integral
(11)
The conservation of energy in the presence of the back reaction due to the radiative losses is
(12)
The analytical solution for the velocity to second order,
, is
(13)
where
(14)
The inclusion of the back reaction allows the evaluation of the SRS’s maximum length,
, which can be derived by setting this velocity to zero:
(15)
Figure 1 shows the SRS’s maximum length as a function of the constant of conversion
and the polar angle
. The above figure shows that the SNR cannot reach an infinite length at an infinite time but has a finite lifetime and length.
![]()
Figure 1. A 2D map of
in pc. The parameters are
and
.
2.3. A Power Law Profile of the Density
The medium is assumed to scale as
(16)
where
is a positive real number. The total mass
swept in the interval
is
The velocity
as a function of the radius is
(17)
The differential equation which governs the motion is
(18)
The above differential equation does not have an analytical solution; a Taylor expansion of order 3 covers a limited range in time
(19)
The above power law profile for the density satisfies the basic symmetries as outlined in Equations (1a) and (1b).
2.4. An Exponential Profile of the Density
The medium is assumed to scale as
(20)
The total mass
swept in the interval [0, r] is
(21)
The velocity
as a function of the radius is
(22)
The differential equation which governs the motion is
(23)
In the absence of an analytical solution, we present the Taylor expansion of order 3 for the trajectory.
(24)
The above exponential profile for the density satisfies the basic symmetries as outlined in Equations (1a) and (1b).
2.5. A Toroidal Profile of the Density
A torus is swept out by revolving a small circle of radius
about an axis lying in the same plane as the circle but outside it, say at a distance of
, see Figure 2. The above toroidal profile of the density, once the axis of revolution is identified with the polar axis, satisfies the basic symmetries as outlined in Equations (1a) and (1b). In Cartesian coordinates, the torus satisfies the equation
(25)
We now consider its intersection with the plane
and as a consequence the equation of the torus is now
(26)
We now evaluate the intersection between the torus and a straight line of equation which crosses the centre, (
),
(27)
where
is the polar angle in the spherical coordinates for the circle which represents the torus in the first quadrant. The line touches the torus at a single point at the critical value of the polar angle,
:
(28)
The above angle allows us to define the following three zones.
1)
: the line does not intersect the torus.
2)
: the line is tangent to the torus at one degenerate single point, with Cartesian coordinates (
), see Figure 3.
3)
: the line intersects the torus at two points, with Cartesian coordinates (
) and (
), see Figure 4.
The Cartesian coordinates at
are
(29a)
(29b)
The Cartesian coordinates of the two intersections when
are
(30a)
(30b)
![]()
Figure 3. The straight line which intersects the centre (red) and is tangent to the torus (blue) in the first quadrant.
![]()
Figure 4. The intersections between a straight line which intersects the centre (red), the torus (blue) to the tangent point (asterisk).
(30c)
(30d)
We now assume that the density of matter is
outside the torus and
inside the torus. The mass swept in the third zone requires a careful analysis. When
, the mass swept along a line before the intersection with the torus,
, is
(31)
where r is the momentary radius of expansion in spherical coordinates.
The mass swept when r is inside the torus in the third zone,
, is
(32)
The mass swept when r is outside the torus in the third zone is
(33)
The equation of motion is solved through the Euler method when the following recursive equations are solved
(34a)
(34b)
where
,
,
are the temporary radius, velocity, and total mass, respectively,
is the time step and n is the index. Due to the fact that the velocity is continuously updated, see Equation (34b), the method turns out to be stable.
3. Astrophysical Applications
The complex structure of SN 1987A, see the Hubble Space Telescope (ST) image in Figure 5, can be classified as a torus only, a torus plus two lobes, and a torus plus 4 lobes, see [10] [11]; we therefore speak of a strong asymmetry. The region connected with the radius of the advancing torus is here identified with our
equatorial region, in spherical coordinates,
. A useful resource for calibration
for calibration is the geometric section of SN 1987A which is given as a sketch in Fig. 5 of [12]. This geometric section was digitized and rotated in the x-z plane by −40˚, see Figure 6; it will be the astronomical section of reference in order to test the simulations.
SN 1006 started to be visible in 1006 AD and currently has a radius of 12.2 pc, see [13] [14]. The X-shape is shown in Figure 7 and the γ-ray shape (100 GeV) is shown in Figure 8.
More precisely, on referring to the above X-image, it can be observed that the radius is greatest in the north-east direction, see also the radio map of SN 1006 at 1370 MHz by [15], and the X-map in the 0.4 - 5.0 keV band of Fig. 1 in [13]. The following observed radii can be extracted:
in the polar direction and
in the equatorial direction. A geometric section of the above X-map was digitized and rotated in the x-z plane by −45˚, see Figure 9; it will be the test for the simulations based on the previous data we can speak of a weak asymmetry for SN 1006. The reliability of the models is evaluated in terms of the percentage reliability,
,
![]()
Figure 5. An ST image of SN 1987A in the year 1997. Credit is given to the Hubble Space Telescope.
![]()
Figure 6. Geometric section of SN 1987A in the x-z plane adapted by the author from Fig. of [12].
![]()
Figure 7. A Chandra X-ray (Red, Green, Blue) image of SN 1006 in the year 2013. Credit is given to the Chandra X-ray Observatory.
![]()
Figure 8. An high energy stereoscopic system (HESS) γ-ray image of SN 1006. Credit is given to the Max-Planck-Institut für Kernphysik.
![]()
Figure 9. Geometric section of SN 1006 in the x-z plane adapted by the author from our Figure 7.
(35)
where
is the theoretical radius of the SNR,
is the observed radius of the SNR, and the index j varies from 1 to the number of available observations.
3.1. Results for the Inverse Square Profile
In the case of an inverse square profile of the density we have an analytical solution as given by Equation (8). Figure 10 displays a cut of SN 1987A in the x-z plane.
A rotation around the z-axis of the previous geometric section allows building a 3D surface, see Figure 11.
Figure 10 displays a cut of SN 1987A in the x-z plane.
Figure 12 displays a cut of SN 1006 in the x-z plane and Figure 13 a 3D display.
The above results can be easily reproduced because we have an analytical expression for the trajectory.
3.2. Results for a Power Profile
In the case of a power law profile of the density we obtained a numerical solution for the differential Equation (18). Figure 14 presents the results for SN 1987A in the x-z plane and Figure 15 those for SN 1006.
The above results shows that a variable power law profile can model the observed sections of the two SNRS here simulated.
3.3. Results for an Exponential Profile
In the case of an exponential profile of the density we obtained a numerical solution for the differential Equation (23). Figure 16 presents the results for SN 1987A in the x-z plane and Figure 17 those for SN 1006.
![]()
Figure 10. Geometric section of SN 1987A in the x-z plane with an inverse square profile (green points) and the observed profile (red stars). The parameters
,
,
,
and
give
.
![]()
Figure 11. 3D surface of SN 1987A with parameters as in Figure 10 in the framework of an inverse square profile. The three Euler angles are
,
and
.
![]()
Figure 12. Geometric section of SN 1006 in the x-z plane with an inverse square profile (green points) and the observed profile (red stars). The parameters
,
,
,
and
give
.
![]()
Figure 13. 3D surface of SN 1006 with parameters as in Figure 12, inverse square profile. The three Euler angles are
,
and
.
![]()
Figure 14. Geometric section of SN 1987A in the x-z plane with a power law profile (green points) and the observed profile (red stars). The parameters
,
,
,
,
and
give
.
![]()
Figure 15. Geometric section of SN 1006 in the x-z plane with a power law profile (green points) and the observed profile (red stars). The parameters
,
,
,
,
and
give
.
![]()
Figure 16. Geometric section of SN 1987A in the x-z plane with an exponential profile (green points) and the observed profile (red stars). The parameters
,
,
,
and
give
.
![]()
Figure 17. Geometric section of SN 1006 in the x-z plane with an exponential profile (green points) and the observed profile (red stars). The parameters
,
,
,
and
give
.
The above results shows that an exponential profile for the density is comparable with the observed sections of the two SNRS here simulated.
3.4. Results for a Toroidal Profile
In the case of a toroidal profile of the density, we obtained a numerical solution for the two recursive Equations (34a) and (34b). Figure 18 presents the results for SN 1987A in the x-z plane and Figure 19 for SN 1006.
A toroidal region around the point of the explosion with density greater than the surrounding medium produces theoretical sections which are comparable with the observed sections of the two SNRS here simulated.
4. Theory of the Image
The astronomical images are given by cut or 2D maps for the intensity of emission. The shape of the intensity of emission is a combination of different processes which are as follows.
1) An assumption on the transfer equation: we adopt an optically thin medium.
2) Type of emission: we choose non thermal emission.
3) Geometry of the radiative zone: the zone of emission resides in the thin advancing shell, which in our case is not spherical.
More details on the above assumptions can be found in Section 7 of [16].
![]()
Figure 18. Geometric section of SN 1987A in the x-z plane for a toroidal profile (green points) and the observed profile (red stars). The parameters
,
,
,
,
,
,
,
and
give
.
![]()
Figure 19. Geometric section of SN 1006 in the x-z plane for a toroidal profile (green points) and the observed profile (red stars). The parameters
,
,
,
,
,
,
,
and
give
.
4.1. How to Build an Image
The first model assumes that the intensity of the image is proportional to the length along the line of sight within the emitting region. The derivation of this length can be complicated, due to the complex morphology of the analysed object and to the infinite points of view of the observer. Figure 20 presents a sketch of the emitting layer, some lines of sight, and the position of the observer. This sketch allows building a geometric model for the theoretical image formed between the two layers. The second model for the intensity of the observed radiation assumes a synchrotron source for luminosity proportional to the flux of kinetic energy,
,
(36)
where A is the considered area, V the velocity and
the density, see Formula (A28) in [17]. In our case
, where
is the considered solid angle along the chosen direction. This means
(37)
where R is the instantaneous radius of the SNR and
is the density in the advancing layer in which the synchrotron emission takes place. The observed luminosity along a given direction can be expressed as
(38)
![]()
Figure 20. First quadrant for the layer between the upper (red curve) and the lower (green curve) section of SN 1987A in the x-z plane with an inverse square profile when
. Parameters as in Figure 10. The observer is at infinity of the x-axis and two lines of sight, s1 and s2, are marked.
where
is a constant of conversion from the mechanical luminosity to the observed luminosity in the considered band. The numerical algorithm which allows us to build a complex image is now outlined.
· An empty (value = 0) memory grid
which contains 4003 pixels is considered.
· We first generate an internal 3D surface by rotating the ideal image of 180˚ around the polar direction and a second external surface at a fixed distance
from the first surface. As an example, we fixed
, where R is the momentary radius of expansion. The points on the memory grid which lie between the internal and external surfaces are stored in memory on
with a variable integer number according to Formula (43).
· Each point of
has spatial coordinates
which can be repre-sented by the following
matrix, A,
(39)
The orientation of the object is characterized by the Euler angles
and therefore by a total
rotation matrix, E, see [18]. The matrix point is represented by the following
matrix, B,
(40)
· The intensity map is obtained by summing the points of the rotated images along a particular direction.
· The effect of the insertion of a threshold intensity,
, given by the observational techniques, is now analysed. In both models the threshold intensity can be parametrized by
, the maximum value of intensity characterizing the map.
More details on the theory of the image can be found in [16].
4.2. Results for the First Model
The radiation of the SN is supposed to originate from between a lower boundary of radius
, given as an example by
of Equation (8), and an upper boundary with radius given by the equation
(41)
The above layer is visualized in Figure 21 for SN 1987A when an inverse square profile is adopted.
The image of SN 1987A is formed between the two layers. The data of the two layers are stored in memory as Cartesian coordinates
(42a)
(42b)
(42c)
(42d)
where n is the considered index. In order to build the first model we assume that z varies between
and
where N represents the number of points. The image is assumed to be proportional to the distance between the x coordinate on the upper layer and the x coordinate on the lower layer as a function of the position z on the polar axis, see Figure 22. In principle, astronomers should produce these observational cuts in order to compare the observations with the theory here developed.
The same algorithm allows building a map of the theoretical intensity for the rotated image of SN 1987A, see Figure 23. The observed image of SN 1987A should be digitized in order to make a comparison with the theory.
4.3. Results for the Second Model
This second model allows simulating particular effects, such as the triple ring system of SN 1987A, see Figure 24, where three zones or holes with theoretical intensity under the threshold value are visible. The enigma of the three holes is solved. Characteristic features, such as the “jet appearance” visible in some maps for SN 1006, see the X-ray map 3 and γ-ray map 3, are theoretically modeled in Figure 25. The enigma of the jet, which is a cone, in a nearly spherical expansion is solved.
![]()
Figure 21. The layer between the upper (green points) and the lower (red points) section of SN 1987A in the x-z plane with an inverse square profile when
. Parameters as in Figure 10.
![]()
Figure 22. Cut of the mathematical intensity I, as a function of the polar z-axis for SN 1987A with an inverse square profile. Parameters as in Figure 10 and
.
![]()
Figure 23. Map of the theoretical intensity I, for SN 1987A with an inverse square profile made by 400 × 400 pixels. Parameters as in Figure 10.
![]()
Figure 24. Model map of SN 1987A rotated in accordance with the observations, for an inverse square medium with parameters as in Figure 10. The three Euler angles characterizing the orientation are
,
and
. In this map
.
![]()
Figure 25. Model map of SN 1006 rotated in accordance with the γ-ray observations with HESS, for a toroidal medium. Physical parameters as in Figure 19. The three Euler angles characterizing the orientation of the observer are
,
and
. In this map
.
![]()
Table 1. Synoptic percentage reliability for the best model of SNRs with different profiles for the density.
5. Conclusions
Models
We have adopted four models in the framework of the conservation of energy for the thin layer approximation and applied them to two SNRs; the results for the percentage reliability, see Formula (35) are shown in Table 1. The best fit for the asymmetric section of SN 1987A is represented by the exponential model for decreasing density in the polar direction and by the toroidal profile of the density for SN 1006.
Image theory
Two models for the intensity of the image in an SNR have been presented, which both require an analytical or numerical function for the advancing radius in terms of the polar angle. The appearance of the triple ring in SN 1987A, see Figure 24, and the jet-feature in SN 1006, see Figure 25, have been simulated.
Back reaction
The derivation of an approximate form for the losses, see Equation (10), allows deriving a finite rather than infinite length for an SNR, see Equation (15).
Acknowledgements
At the time of writing, an animated version of Figure 11 is visible at http://personalpages.to.infn.it/~zaninett/image/sn1987a_animation.gif.