The physics of the optical light curve in supernovae

We present a new formula which models the rate of decline of supernovae (SN) as given by the light curve in various bands. The physical basis is the conversion of the flux of kinetic energy into radiation. The two main components of the model are a power law dependence for the radius--time relation and a decreasing density with increasing distance from the central point. The new formula is applied to SN 1993J, SN 2005cf, SN 1999ac, and SN 2001el in different bands.


Introduction
The light curve (LC) of supernovae (SN) at a given wavelength λ denotes the luminosity-time relation. The astronomers work in terms of apparent/absolute magnitude and therefore the LC in SN is usually presented as a magnitude versus time relation. We have two great astronomical classifications for the LC: type I SN and type II SN. The type I has a fast decrease in magnitude followed by a nearly linear increase. In luminosity terms, the SN has a fast increase followed by a nearly exponential decay. The type II has a fast decrease in magnitude followed by oscillations, type IIb, or a plateau, type IIp; a decay follows the plateau. In this complex morphology, we will always specify the type of SN under consideration. The luminosity is usually modeled by the formula where L and L λ,0 are the luminosity at time t and at t = 0 respectively, and τ is the typical lifetime, see Bowers and Deeming (1984). As an example, the radioactive isotope 56 Ni has τ = 8.767 days. On introducing the apparent magnitude m λ , the previous formula becomes where k ′ λ is a constant. The absolute magnitude M λ scales in the same way: where k ′′ λ is another constant. The observational fact that, as an example, in IC 4182 the LC has a half-life of 56 days, requires the production of 56 Co, see van Hise (1974). The previous formula is an empirical relation which is based solely on observations rather than theory. The theory for SNII LCs was first developed by Grasberg et al. (1971) and later analytically and numerically explored by Falk and Arnett (1973); Arnett (1980); Arnett and Fu (1989). A model for the luminosity in Hα of supernovae as a function of time can be found in Figure 7 of Chevalier and Fransson (1994). The LCs of type Ia SN have been explained (including the secondary maximum) by a time-dependent multigroup radiative transfer calculation, see Kasen (2006). A model for type II supernovae explosions has been built including progenitor mass, explosion energy, and radioactive nucleosynthesis, see Kasen and Woosley (2009). The model atmosphere code PHOENIX was used to calculate type Ia supernovae, see Jack et al. (2011). The previous works leave a series of questions unanswered or merely partially answered.
• Given the observational fact that the radius-time relation in young SNRs follows a power law, is it possible to find a theoretical law of motion which fits the observations?
• Can a model of an expansion in the framework of the thin layer approximation produce the observed radius-time relation?
• Can we express the flux of kinetic energy in the framework of an approximate law of motion and a medium characterized by a decreasing density?
• Can we parametrize the conversion of the flux of kinetic energy into total observed luminosity?
• Can we parametrize the fraction of conversion of the total luminosity into the optical bands?
In order to answer these questions, in Section 2.2 we analyze the existing equations of motion for SN 1993J as well a new adjustable equation. Section 3 reviews the basic formulas of synchrotron emission and reports the conversion of flux of kinetic energy into an observed band. Section 4 reports the application of the new formulas to different SNs in various bands.

The equation of motion
This section reviews three existing parameters for SNRs, the power law model and a new solution in the framework of the thin layer approximation.

Some existing solutions
The Sedov-Taylor solution is where E is the energy injected into the process and t is time, see Sedov (1944); Taylor (1950a,b); Sedov (1959);McCray (1987). Our astrophysical units are: time, (t 1 ), which is expressed in years; E 51 , the energy in 10 51 erg; n 0 , the number density expressed in particles cm −3 (density ρ = n 0 m, where m = 1.4m H ). In these units, Equation (4) becomes The Sedov-Taylor solution scales as R(t) ≈ t 2/5 = t 0.4 . A second solution is connected with momentum conservation in the presence of a constant density medium, see Dyson, J. E. and Williams, D. A. (1997); Padmanabhan (2001); Zaninetti (2009). The astrophysical radius in pc as a function of time is where t 1 and t 0 are the time in years, R 0 is the radius in pc when t 1 = t 0 , and v 1 is the velocity in km s −1 when t 1 = t 0 . The thin layer solution in the presence of a constant density medium scales as t 0.25 . A relativistic solution of the thin layer approximation can be found in Zaninetti (2010). A sophisticated approach as given by Chevalier (1982a) and Chevalier (1982b) analyzes self-similar solutions with varying inverse power law exponents for the density profile of the advancing matter, R −n , and ambient medium, R −s . The previous assumptions give a law of motion R ∝ t n−3 n−s when n > 5.

The equation of motion as a power law
The equation of the expansion of an SNR can be modeled by a power law of the type where R is the radius of the expansion, t is the time, R 0 is the radius at t = t 0 , and α is an exponent which can be found from a numerical analysis. In order to find the unknown parameters, we analyzed the data of supernova SN 1993J , classified as type IIb, which began to be visible in M81 in 1993, see Ripero et al. (1993), and presented a circular symmetry for 4000 days, see Marcaide et al. (2009). Its distance is 3.63 Mpc (the same as M81), see Freedman et al. (1994). The velocity is As an example, Figure 1 reports the fit of SN 1993J . We have chosen this SN because: Figure 1: The theoretical radius as given by the power law fit represented by Equation (7) with α =0.828, R 0 =0.0087 pc and t 0 = 0.498 yr. The astronomical data of SN 1993J are represented by vertical error bars and are extracted from Table 1 in Marcaide et al. (2009).
• it presents a nearly spherical expansion, • the temporary radius of expansion has been measured for ≈ 10 yr in the radio band, seeMarcaide et al. (2009).
The observed radius-time relation of SN 1993J allows us to calibrate our model and the application of the least squares method through the FORTRAN subroutine LFIT from Press et al. (1992) allows finding α = 0.828 . Therefore the radius is growing more slowly than a free expansion with constant velocity, R ∝ t, but more quickly than the Sedov-Taylor solution, R ∝ t 0.4 , see Equation (4).

An adjustable equation of motion
We assume that around the SNR the density of the interstellar medium (ISM) has the following two piecewise dependencies This assumption allows us to set up the initial conditions, otherwise the dependence ρ(R) ∝ ( R0 R ) d will have a pole at R = 0. At the moment of writing there is not a clear determination of the gradients around the SN and therefore d can be considered a free parameter.
In this framework the SN is moving and the density, which is at rest, decreases as an inverse power law with an exponent d which can be fixed from the observed temporal evolution of the radius. The mass swept, M 0 , in the interval The mass swept, M , in the interval 0 ≤ r ≤ R with r ≥ R 0 is Momentum conservation in the thin layer approximation requires that where v is the velocity at t and v 0 is the velocity at t = t 0 . The previous expression as a function of the radius is In this differential equation of first order in R, the variables can be separated and an integration term-by-term gives the following nonlinear equation An approximate solution of F N L (r) can be obtained assuming that Up to now, the physical units have not been specified, pc for length and yr for time are perhaps acceptable choices. With these units, the initial velocity v 0 is expressed in pc yr −1 and should be converted into km s −1 ; this means that v 0 = 1.02 10 −6 v 1 where v 1 is the initial velocity expressed in km s −1 .
The astrophysical version of the above equation in pc is where t 1 and t 0 are times in years, R 0 is the radius in pc at t 1 = t 0 and v 1 is the velocity at t 1 = t 0 in km s −1 . The approximate solution (15) has the following limit as t → ∞ where . ( On imposing where α is an observable parameter defined in Section 2.2. This means that the unknown parameter d can be deduced from the observed parameter α. More details on this model, as well as a relativistic version, can be found in Zaninetti (2011), where conversely the LC is not treated.

The energy cascade
This section contains the basic formula for the synchrotron emission, the transformation of the mechanical flux of energy into the observed luminosity, and the conversion of the predicted flux at a given wavelength to the apparent magnitude.

Synchrotron emission
In SNR we detect non-thermal emission with intensity (where ν is the frequency, λ the wavelength and β the power law index).
As an example in the case of SN 1993J after the transition from optically thick to optically medium, β becomes ≈ −0.6 after 2500 days, see Figure 8 in Martí-Vidal et al. (2011). The conversion of the flux of kinetic energy into synchrotron luminosity can be obtained by the following physical processes • Turbulent evolution in the advancing shock in the framework of both the Kolmogorov and Kraichnan spectrum, see Fan et al. (2010).
• Particle acceleration in a turbulent environment using a Monte Carlo approach for the diffusion and acceleration of the particles, coupled to a magnetohydrodynamics code in the SNR environment, see Schure et al. (2010).
• A model for the evolution of the magnetic field in the advancing layer, see Reynolds (2011).
• Diffusion of the relativistic electrons from the position of the advancing layer, see Section 6 in Zaninetti (2011).
The lifetime, τ syn , for synchrotron losses is where H is the magnetic field in Gauss and ν is the frequency of observation in Hertz, see Lang (1999). The outlined cascade of physical processes can work if the following inequalities are verified where t cas is the time scale of turbulence formation and t a is the time scale of electron acceleration. In the following we will assume that the synchrotron emission is the main source of luminosity. Two radioactive decays will be considered in Section 4.

Non-thermal and thermal emission
The synchrotron emission in SNRs is detected from 10 8 Hz of radio-astronomy to 10 19 Hz of gamma astronomy which means 11 decades in frequency. At the same time, some particular effects, such as absorption, transition from optically thick to optically thin medium, line emission, and the energy decay of radioactive isotopes ( 56 Ni, 56 Co) can produce a change in the concavity of the flux versus frequency relation, see the discussion about Cassiopea A in Section 3.3 of Eriksen et al. (2009). A comparison between non-thermal and thermal emission (luminosity and surface brightness distribution) can be found in Petruk and Beshlei (2007), where it is possible to find some observational tests which allow the estimation of the parameters characterizing the cosmic ray injection on supernova remnant shocks. At the same time, a technique to isolate the synchrotron emission from the thermal emission is widely used, as an example see X-limb of SN1006 Katsuda et al. (2010).

The temporal evolution
The density of kinetic energy, K, is where ρ is the density and V the velocity. In presence of an area A and when the velocity is perpendicular to that area, the flux of kinetic energy L m is which in SI is measured in J s −1 and in CGS in erg s −1 see formula (A28) in De Young (2002). In our case, A = 4πR 2 , which means where R is the instantaneous radius of the SNR and ρ is the density in the advancing layer. The source of synchrotron luminosity is assumed here to be proportional to the flux of kinetic energy. The density in the advancing layer is assumed to scale as R −d , which means that This last assumption is connected with the adjustable equation of motion which is derived in a decreasing density environment, see Section 2.3. On adopting this point of view, d is an unknown parameter which allows matching theory and observation. The temporal and velocity evolution are given by the power law dependencies of Equations (7) and (8) and therefore The synchrotron luminosity L λ and the observed flux S λ at a given wavelength λ are assumed to be proportional to the mechanical luminosity and therefore where S 0 is the flux when t = t 0 at a given wavelength λ obs . The apparent magnitude at a given color c, where c can be U, B, V, R or I, is m c = k c − 2.5 log 10 Se λ I λ dλ , where Se λ is the sensitivity function in the region specified by the wavelength λ, k c is a constant, and I λ is the energy flux reaching the earth. We now define a sensitivity function for a pseudo-monochromatic color system where δ denotes the Dirac delta function, see Bowers and Deeming (1984). In this system the apparent magnitude is m c = k c − 2.5 log 10 I λ .
On assuming that the intensity of emission and the flux of kinetic energy as given by (29) are directly proportional, we obtain where k c is a constant: and k b is a constant.
In the previous equations we have three unknowns: α, k c and d. In the case of SN 1993J the value of α is deduced from the data of the expansion. On fixing two times in the observed LC, t = t 0 and t = t 1 , we have two corresponding magnitudes m 0 and m 1 . The resulting nonlinear system of two equations in two unknowns can therefore be solved.
The (C 1 − C 2 ) color can be expressed as where k 12 is a constant and I λ is the energy flux reaching the earth. In a pseudo-monochromatic color system According to the previous equation, the color of a SN should be constant with time. As an example in the case of SN 2005cf (type Ia) (B-V) became stable after the first 120 days Pastorello et al. (2007). The constancy of the color has been obtained with the assumption that the spectral index is constant with time. The spectral index in the radio varies considerably but becomes constant, β ≈ −0.7, after ≈ 700 days, see Figure 8 in  Table 3 in Zhang et al. (2004). In other words, the constancy of the color can be applied after ≈ 700 days. More precisely, the observed luminosity at time t can be expressed introducing the initial mechanical luminosity, L m0 , defined as where the index 0 stands for the first measurement. The astrophysical version of the above equation is where n 1 is the initial number density expressed in units of particles cm −3 , R 1 is the initial radius expressed in units of pc, and v 10000 is the initial velocity expressed in units of 10,000 km s −1 . The spectral luminosity, L ν , at a given frequency ν is where S 0 is the flux observed at the frequency ν 0 and D is the distance. The total observed luminosity, L tot , is where ν min and ν max are the minimum and maximum frequencies observed. The total observed luminosity can be expressed as where ǫ is a constant of conversion from the mechanical luminosity to the total observed luminosity in the synchrotron emission. The fraction of the total luminosity deposited in a color f c is where ν c,min and ν c,max are the minimum and maximum frequency of a color. Table 3.3 presents some values of f c for the most important optical bands. At the time of writing, the number density in the advancing layer is unknown and we can therefore define ǫn 1 as the constant which allows adjusting theory and observations. About ǫ it should be said that by definition ǫ < 1. The rapid rise in intensity in a SN can be modeled by Equation (47) for the radiative transfer when a time dependent transition from optically thick to optically thin medium is considered. The solution of the radiative transfer equation for the specific intensity per unit frequency, I ν , at the end of an astrophysical object, is where τ ν is the optical depth, G ν is the source function, and I ν (0) the intensity beyond the astrophysical object, see equation (1.30) in Rybicki and Lightman (1991). On considering only the intensity of the object (I ν (0) = 0) the previous formula becomes where τ ν =1 represents the value at which the intensity is 63% of the source function. The temporal transition from optically thick to optically thin medium before the maximum can be modeled by imposing τ ν = t ta where t a is a typical time. This is an 'ad hoc' function that allows of modeling the transition before and after the maximum and the consequent change of concavity of the LC as function of time. The time t a can vary from the few seconds of a Gamma Ray Burst (GRB) to the few days of the optical bands.
A logarithmic form of Equation (45) introduces the apparent magnitude m ν where k ν is a constant.
We are now ready to introduce the two phase model which can be characterized by the following two piecewise dependencies 4 Applications Figure 2 reports the decay of the R magnitude of SN 1993J , which is type IIb, as well our theoretical curve. The theoretical temporal evolution of the Hα luminosity of SN 1993J as well the data are reported in Figure 3. The Hα luminosity which is derived from the i−band is also fitted by the model of Chevalier and Fransson (1994).
A second example is SN 2005cf (type Ia) which has been analyzed in Pastorello et al. (2007); Figure 4 and Figure 5 report the decay of the V and B magnitude of SN 2005cf as well our theoretical curve.
The (B-V) color evolution of SN 2005cf is reported in Figure 6. Only the second phase is reported. A third example is the sample of 44 type Ia supernovae which have been observed in the UBVRI bands, see Jha et al. (2006). We selected SN 1999ac in the U band and Figure 7 reports the LC as well our fit. A fourth example is SN 2001el (type Ia), which has been analyzed in Krisciunas et al. (2003). Figure 8 reports the decay of the V magnitude of SN 2001el . A comparison should be made with Figure 4 of Kasen (2006) which uses the decay of 56 Ni.
It is also interesting to plot the decay of the LC of SN 2001el , see Krisciunas et al. (2003), as given by two nuclear decay which, according to Equation (3), are straight lines, see Figure 9.    Table 6 of Pastorello et al. (200).  Table 6 of Pastorello et al. (2007).    Table 3 of Krisciunas et al. (2003) and the adopted distance modulus is µ = 31.4 according to Kasen et al. (2006). Figure 9: The V LC of SN 2001el (empty stars) in absolute magnitude, theoretical curve as given by Equation (3) when the radioactive decay of the isotope 56 Ni (τ = 8.757 d or T 1/2 = 6.07 d) was considered (full line), and theoretical curve of the radioactive decay of the isotope 56 Co (τ = 111.47 d or T 1/2 = 77.27 d) was considered (dashed line).

Conclusions
The SN's are classified as spherical SN, as an example SN 1993J , and as aspherical SN, see as an example Racusin et al. (2009) for SN 1987A . The theory here developed treats the spherical SN using classical dimensional arguments. The conversion of the flux of kinetic energy into luminosity after the maximum in the LC explains the curve of SNs in a direct form, see Equations (28) and (38) as well as in a logarithmic version, see Equation (33). The overall LC before and after the maximum can be built by introducing two different physical regimes, see Equation (47). The initial rise in intensity in the V-band is characterized by a typical time scale of t a ≈ 5 days and the decrease can be theoretically fitted for t ≈ 3500 days. This large range in time is also the great advantage of our model: the existing nuclear models cover ≈ 100 days, see Figure 2 in Leibundgut and Suntzeff (2003). The standard approach of formula (2), which predicts a linear increase in the apparent magnitude with time, does not correspond to the observations because the observed and theoretical magnitudes scale as m = a + b ln(t) where a and b are two constants. As an example, Figure (9) reports two commonly accepted sources which are the radioactive isotopes 56 Co, see Georgii et al. (2000); Pluschke et al. (2001); Georgii et al. (2002), and 56 Ni, see Truran et al. (2012); Dessart et al. (2012): the radioactive fit is acceptable only for the first few days. The application of the new formulas to three SNs in different bands gives acceptable results. As an example, Figure 2 reports the LC in the R-band for SN 1993J and Figure 3 reports the LC for the Hα of SN 1993J . An example of the two phase model as given by Equation (47) is reported in Figure (4) for SN 2005cf in the V-band. A careful analysis of the previous figures shows that the theoretical and observed curves present different concavities in the transition from small to large times. Similar results can be obtained assuming that all γ-rays produced by the decay of 56 Ni and 56 Co are converted into optical emission, see Figure 2 in Leibundgut and Suntzeff (2003). The observational fact that the initial velocity can be ≈ 30000 km s −1 requires a relativistic treatment that is necessary for future progress. The analysis here performed treats the SN as a single object and therefore is not connected with various types of recent cosmologies, see Astier and Pain (2012); Chavanis (2013); ElNabulsi (2013).
We conclude with a list of not yet solved problems: • The observational fact that the initial velocity can be ≈ 30000 km s −1 requires a relativistic treatment of the flux of kinetic energy that is left for future research; • The connection between the cosmic ray production and the γ-rays in SNR, see Dermer and Powale (2013), requires an analysis of the temporal behavior of the magnetic field.