Electron Acceleration in Blazars: Application to the 3C 279 Flare on 2013 December 20

Throughout the data comparison process, we use the numerical model, including the self-consistent ED (Paper II). We describe the ED ${N}_{e}(\gamma )$ in the frame of the blob using a steady-state Fokker–Planck equation,

$\begin{eqnarray}\begin{array}{rcl}0 & = & \displaystyle \frac{{\partial }^{2}}{\partial {\gamma }^{2}}\left(\displaystyle \frac{1}{2}\displaystyle \frac{d{\sigma }^{2}}{{dt}}{N}_{e}\right)-\displaystyle \frac{\partial }{\partial \gamma }\left(\left\langle \displaystyle \frac{d\gamma }{{dt}}\right\rangle {N}_{e}\right)\\ & & -\,\displaystyle \frac{{N}_{e}}{{t}_{\mathrm{esc}}}+{\dot{N}}_{\mathrm{inj}}\delta (\gamma -{\gamma }_{\mathrm{inj}}),\end{array}\end{eqnarray} \tag{ 1 }$

where $\gamma \equiv E/({m}_{e}{c}^{2})$ is the electron Lorentz factor, m_e is the electron mass, and c is the speed of light in a vacuum. Since this is a steady-state equation, we have set $\partial {N}_{e}/\partial t=0$. The acceleration and synchrotron energy loss timescales are appreciably shorter than the observed rise time of the flare (see Section 3.5), making the steady-state calculation appropriate in this case. We note that we use the term "electrons" here and throughout this paper to refer to both electrons and positrons.

In Equation (1), the energy-dependent particle escape timescale, t_esc, is related to the dimensionless escape parameter τ via

$\begin{eqnarray}&&{t}_{\mathrm{esc}}(\gamma )=\displaystyle \frac{\tau }{{D}_{0}\gamma },\end{eqnarray} \tag{ 2 }$

where

$\begin{eqnarray}&&\tau \equiv \displaystyle \frac{{R}_{b}^{{\prime} 2}{{qBD}}_{0}}{{m}_{e}{c}^{3}},\end{eqnarray} \tag{ 3 }$

in the Bohm limit (Paper I), and q is the fundamental charge. The Lorentz factor of the injected electrons is ${\gamma }_{\mathrm{inj}}$ in the particle injection rate

$\begin{eqnarray}&&{\dot{N}}_{\mathrm{inj}}=\displaystyle \frac{{L}_{{\rm{e}},\mathrm{inj}}}{{m}_{e}{c}^{2}{\gamma }_{\mathrm{inj}}},\end{eqnarray} \tag{ 4 }$

where ${L}_{{\rm{e}},\mathrm{inj}}$ is the electron injection luminosity. Both ${\gamma }_{\mathrm{inj}}$ and ${L}_{{\rm{e}},\mathrm{inj}}$ are implemented as free parameters, but the former is held constant in the subsequent analysis at the lower numerical grid limit to simulate a thermal particle source. We solved Equation (1) using the full Compton energy loss rate numerically, as described in Paper II.

In Equation (1), the broadening coefficient

$\begin{eqnarray}&&\displaystyle \frac{1}{2}\displaystyle \frac{d{\sigma }^{2}}{{dt}}={D}_{0}\,{\gamma }^{2},\end{eqnarray} \tag{ 5 }$

is consistent with hard-sphere scattering, where ${D}_{0}\propto {{\rm{s}}}^{-1}$ (Park & Petrosian 1995) is a free parameter. The drift coefficient expresses the mean rate at which each process contributes

$\begin{eqnarray}&&\left\langle \displaystyle \frac{d\gamma }{{dt}}\right\rangle ={D}_{0}\left[4\gamma +a\gamma -{b}_{\mathrm{syn}}{\gamma }^{2}-{\gamma }^{2}\displaystyle \sum _{j=1}^{J}{b}_{{\rm{C}}}^{(j)}H(\gamma {\epsilon }_{\mathrm{ph}}^{(j)})\right],\end{eqnarray} \tag{ 6 }$

where second-order Fermi acceleration occurs at a rate of

$\begin{eqnarray}&&{\dot{\gamma }}_{\mathrm{sto}}=4{D}_{0}\gamma ,\end{eqnarray} \tag{ 7 }$

and

$\begin{eqnarray}&&{\dot{\gamma }}_{\mathrm{ad}+\mathrm{sh}}\equiv {{aD}}_{0}\gamma \equiv {A}_{\mathrm{ad}+\mathrm{sh}}\gamma ,\end{eqnarray} \tag{ 8 }$

where a is a dimensionless free parameter that includes first-order Fermi acceleration and adiabatic cooling.⁵ The coefficients A_sh and A_ad represent the first-order Fermi acceleration rate and the adiabatic loss rate, respectively. The rate of synchrotron cooling is

$\begin{eqnarray}&&| {\dot{\gamma }}_{\mathrm{syn}}| \equiv {D}_{0}{b}_{\mathrm{syn}}{\gamma }^{2}=\displaystyle \frac{{\sigma }_{{\rm{T}}}{B}^{2}}{6\pi {m}_{e}c}{\gamma }^{2},\end{eqnarray} \tag{ 9 }$

where ${\sigma }_{{\rm{T}}}=6.65\times {10}^{-25}$ cm² is the Thomson cross section and B is the strength of the tangled, homogeneous magnetic field. Similarly, the Compton cooling rate for each individual component is

$\begin{eqnarray}\begin{array}{rcl}| {\dot{\gamma }}_{{\rm{C}}}| & \equiv & {D}_{0}{\gamma }^{2}{b}_{{\rm{C}}}^{(j)}H(\gamma {\epsilon }_{\mathrm{ph}}^{(j)})\equiv {\gamma }^{2}{B}_{{\rm{C}}}^{(j)}H(\gamma {\epsilon }_{\mathrm{ph}}^{(j)})\\ & = & {\gamma }^{2}\displaystyle \frac{4{\sigma }_{{\rm{T}}}{{\rm{\Gamma }}}^{2}{u}_{\mathrm{ph}}^{(j)}}{3{m}_{e}c}H(\gamma {\epsilon }_{\mathrm{ph}}^{(j)}).\end{array}\end{eqnarray} \tag{ 10 }$

Here ${b}_{{\rm{C}}}^{(j)}$ is a dimensionless constant related to Compton cooling for the different external radiation fields j, with energy densities ${u}_{\mathrm{ph}}^{(j)}$. The function H(y) is a complicated expression related to mitigation of energy losses with the full Compton cross section (Böttcher et al. 1997). We include a dust torus and 26 broad lines, for a total of J = 27 EC components (Finke 2016).

In the Thomson limit, $y\ll 1$, $H(y)\to 1$, the steady-state Fokker–Planck equation (Equation (1)) has the analytic solution (see Paper I for details)

$\begin{eqnarray}\begin{array}{rcl}{N}_{e}(\gamma ) & = & \displaystyle \frac{{\dot{N}}_{\mathrm{inj}}}{{{bD}}_{0}}\displaystyle \frac{{\rm{\Gamma }}(\mu -\kappa +1/2)}{{\rm{\Gamma }}(1+2\mu )}{e}^{({\gamma }_{\mathrm{inj}}-\gamma )b/2}{\gamma }_{\mathrm{inj}}^{-2-a/2}{\gamma }^{a/2}\\ & & \times \,\left\{\begin{array}{ll}{{ \mathcal M }}_{\kappa ,\mu }(b\gamma ){{ \mathcal W }}_{\kappa ,\mu }(b{\gamma }_{\mathrm{inj}}), & \ \ \gamma \leqslant {\gamma }_{\mathrm{inj}}\\ {{ \mathcal M }}_{\kappa ,\mu }(b{\gamma }_{\mathrm{inj}}){{ \mathcal W }}_{\kappa ,\mu }(b\gamma ), & \ \ {\gamma }_{\mathrm{inj}}\leqslant \gamma \end{array}\right.,\end{array}\end{eqnarray} \tag{ 11 }$

where $b\equiv {b}_{\mathrm{syn}}+{\sum }_{j=1}^{J}{b}_{{\rm{C}}}^{(j)}$, ${{ \mathcal M }}_{\kappa ,\mu }$ and ${{ \mathcal W }}_{\kappa ,\mu }$ are Whittaker functions (Slater 1960), with coefficients

$\begin{eqnarray}&&\kappa =2-\displaystyle \frac{1}{b\tau }+\displaystyle \frac{a}{2},\quad \mathrm{and}\ \quad \mu =\displaystyle \frac{a+3}{2}.\end{eqnarray} \tag{ 12 }$

We found the analytic Thomson regime solution, Equation (11), previously (Paper I). Here we provide a useful approximation,

$\begin{eqnarray}&&{N}_{e}^{\mathrm{app}}(\gamma )\approx \displaystyle \frac{{\dot{N}}_{\mathrm{inj}}}{{D}_{0}}{e}^{-b\gamma }\left[\displaystyle \frac{\tau {b}^{a+4}{\gamma }^{a+2}}{{\rm{\Gamma }}(a+4)}+\displaystyle \frac{1}{3+a}\displaystyle \frac{1}{\gamma }\right],\end{eqnarray} \tag{ 13 }$

in the ${\gamma }_{\mathrm{inj}}\lt \gamma$ regime (where most of our analysis takes place), and where $b\gamma \ll 1$ and (bτ)⁻¹ ≪ 1 are the simplifying assumptions, which are valid in all blazar analyses we have examined, both here and in Paper II. The precision of these assumptions is addressed further in Appendix A, however the simplified solution is relatively accurate for most applications, except where (bτ)⁻¹ ≳ 1.

The constraint $b\gamma \ll 1$ indicates that the rate of particle energy lost to synchrotron (Equation (9)) and Thomson processes (Equation (10), for $H(\epsilon \gamma )=1$)), where ($b={b}_{\mathrm{syn}}+{b}_{{\rm{C}}}^{(j)}$) must be smaller than the rate at which energy is gained by the second-order Fermi process (Equation (7)). Thus, the emitting blob must be in the acceleration region, which can be demonstrated by the energy budget for the flare presented in this work (Section 3.4). Similarly, the requirement (bτ)⁻¹ ≪ 1 is analogous to the expectation that the energy losses due to the synchrotron and Thomson processes outpace energy lost due to particle escape (i.e., the fast-cooling regime), which is apparent in the present application by examining the energy budget (Section 3.4).

The full derivation and complementing solution for $\gamma \gt {\gamma }_{\mathrm{inj}}$ are given in Appendix A. The main shape of the ED is governed by two components: one is driven by a balance between first- and second-order Fermi acceleration ($\propto {\gamma }^{a+2};$ orange curve in Figure 1) and the other is due to second-order Fermi acceleration ($\propto {\gamma }^{-1};$ green curve in Figure 1; see Appendix B). The ratio of first- to second-order Fermi acceleration acceleration (parameterized by $a={A}_{\mathrm{ad}+\mathrm{sh}}/{D}_{0}$) impacts the shared term. When a < −2, the lower energy power law is negative, and when a > −2, the $\propto {\gamma }^{a+2}$ term is increasing with increasing energy. In practice, the shape of the the term does not vary significantly for most analyses, but is important in simulating the flare examined here. The exponential cutoff is governed by the emission mechanisms, which are constrained to the Thomson limit for the analytic solution. Appendix A discusses the interpretation of the acceleration mechanisms in the simplified solution in further detail.

Zoom In Zoom Out Reset image size

The analytic solution and others derived from it, are useful for physical interpretation of the shape of the ED with regard to acceleration mechanisms. However, all of the data interpretation here is performed with the full numerical model using the full Compton cross section and energy losses.

We include thermal emission from the accretion disk (Shakura & Sunyaev 1973) and dust torus, in addition to emission from the jet from synchrotron, SSC, and EC of dust torus and BLR photons. The source 3C 279 has redshift z = 0.536, giving it a luminosity distance ${d}_{L}=9.6\times {10}^{27}\ \mathrm{cm}$ in a cosmology where (h, Ω_m, Ω_Λ) = (0.7, 0.3, 0.7).

The νF_ν disk flux is approximated as

$\begin{eqnarray}&&{f}_{{\epsilon }_{\mathrm{obs}}}^{\mathrm{disk}}=\displaystyle \frac{1.12}{4\pi {d}_{L}^{2}}\ {\left(\displaystyle \frac{\epsilon }{{\epsilon }_{\max }}\right)}^{4/3}{e}^{-\epsilon /{\epsilon }_{\max }},\end{eqnarray} \tag{ 14 }$

(Dermer et al. 2014), where $\epsilon ={\epsilon }_{\mathrm{obs}}(1+z)$ and ${m}_{e}{c}^{2}{\epsilon }_{\max }=10\,\mathrm{eV}$. The νF_ν dust torus flux is approximated as a blackbody,

$\begin{eqnarray}&&{f}_{{\epsilon }_{\mathrm{obs}}}^{\mathrm{dust}}=\displaystyle \frac{15{L}^{\mathrm{dust}}}{4{\pi }^{5}{d}_{L}^{2}}\displaystyle \frac{{\left(\epsilon /{\rm{\Theta }}\right)}^{4}}{\exp (\epsilon /{\rm{\Theta }})-1},\end{eqnarray} \tag{ 15 }$

where again $\epsilon ={\epsilon }_{\mathrm{obs}}(1+z)$, and also Θ = k_BT_dust/(m_ec²), T_dust is the dust temperature, and k_B is the Boltzmann constant.

The jet blob νF_ν flux is computed using the ED solution to the electron Fokker–Planck equation (Equation (1)), ${N}_{e}^{{\prime} }({\gamma }^{{\prime} })$. We now add primes to indicate the distribution is in the frame co-moving with the jet blob.

The synchrotron flux

$\begin{eqnarray}&&{f}_{{\epsilon }_{\mathrm{obs}}}^{\mathrm{syn}}=\displaystyle \frac{\sqrt{3}\epsilon ^{\prime} {\delta }_{{\rm{D}}}^{4}{e}^{3}B}{4\pi {{hd}}_{{\rm{L}}}^{2}}{\int }_{1}^{\infty }d{\gamma }^{{\prime} }{N}_{e}^{{\prime} }({\gamma }^{{\prime} })\ R(x),\end{eqnarray} \tag{ 16 }$

where

$\begin{eqnarray}&&x=\displaystyle \frac{4\pi \epsilon ^{\prime} {m}_{e}^{2}{c}^{3}}{3{eBh}{\gamma }^{{\prime} 2}},\end{eqnarray} \tag{ 17 }$

and R(x) is defined by Crusius & Schlickeiser (1986). Synchrotron self-absorption is also included. The SSC flux

$\begin{eqnarray}\begin{array}{rcl}{f}_{{\epsilon }_{s}}^{\mathrm{SSC}} & = & \displaystyle \frac{9}{16}\displaystyle \frac{{\left(1+z\right)}^{2}{\sigma }_{{\rm{T}}}{\epsilon }_{s}^{{\prime} 2}}{\pi {\delta }_{{\rm{D}}}^{2}{c}^{2}{t}_{v}^{2}}\displaystyle {\int }_{0}^{\infty }\ d{\epsilon }_{* }^{{\prime} }\ \displaystyle \frac{{f}_{{\epsilon }_{* }}^{\mathrm{syn}}}{{\epsilon }_{* }^{{\prime} 3}}\\ & & \,\times \displaystyle {\int }_{{\gamma }_{1}^{{\prime} }}^{\infty }\ d{\gamma }^{{\prime} }\displaystyle \frac{{N}_{e}^{{\prime} }({\gamma }^{{\prime} })}{{\gamma }^{{\prime} 2}}{F}_{{\rm{C}}}\left(4{\gamma }^{{\prime} }{\epsilon }_{* }^{{\prime} },\displaystyle \frac{\epsilon }{{\gamma }^{{\prime} }}\right),\end{array}\end{eqnarray} \tag{ 18 }$

(e.g., Finke et al. 2008), where ${\epsilon }_{s}^{{\prime} }={\epsilon }_{s}(1+z)/{\delta }_{{\rm{D}}}$, ${\epsilon }_{* }^{{\prime} }\,={\epsilon }_{* }(1+z)/{\delta }_{{\rm{D}}}$, and

$\begin{eqnarray}&&{\gamma }_{1}^{{\prime} }=\displaystyle \frac{1}{2}{\epsilon }_{s}^{{\prime} }\left(1+\sqrt{1+\displaystyle \frac{1}{{\epsilon }^{{\prime} }{\epsilon }_{s}^{{\prime} }}}\right).\end{eqnarray} \tag{ 19 }$

The function F_C(p, q) was originally derived by Jones (1968), but had a mistake that was corrected by Blumenthal & Gould (1970). The EC flux (e.g., Georganopoulos et al. 2001; Dermer et al. 2009)

$\begin{eqnarray}\begin{array}{rcl}{f}_{{\epsilon }_{s}}^{\mathrm{EC}} & = & \displaystyle \frac{3}{4}\displaystyle \frac{c{\sigma }_{{\rm{T}}}{\epsilon }_{s}^{2}}{4\pi {d}_{L}^{2}}\displaystyle \frac{{u}_{* }}{{\epsilon }_{* }^{2}}{\delta }_{{\rm{D}}}^{3}\\ & & \times \,\displaystyle {\int }_{{\gamma }_{1}}^{{\gamma }_{\max }}d\gamma \displaystyle \frac{{N}_{e}^{{\prime} }(\gamma /{\delta }_{{\rm{D}}})}{{\gamma }^{2}}{F}_{{\rm{C}}}\left(4\gamma {\epsilon }_{* },\displaystyle \frac{{\epsilon }_{s}}{\gamma }\right),\end{array}\end{eqnarray} \tag{ 20 }$

where

$\begin{eqnarray}&&{\gamma }_{1}=\displaystyle \frac{1}{2}{\epsilon }_{s}\left(1+\sqrt{1+\displaystyle \frac{1}{\epsilon {\epsilon }_{s}}}\right)\ ,\end{eqnarray} \tag{ 21 }$

and u_* and ${\epsilon }_{* }$ are the energy density and dimensionless photon energy of the external radiation field, respectively. For the dust torus photons,

$\begin{eqnarray}&&{u}_{\ast }={u}_{{\rm{d}}{\rm{u}}{\rm{s}}{\rm{t}}}=2.2\times {10}^{-5}\left(\displaystyle \frac{{\xi }_{{\rm{d}}{\rm{u}}{\rm{s}}{\rm{t}}}}{0.1}\right){\left(\displaystyle \frac{{T}_{{\rm{d}}{\rm{u}}{\rm{s}}{\rm{t}}}}{1000{\rm{K}}}\right)}^{5.2}\,{\rm{e}}{\rm{r}}{\rm{g}}\,{{\rm{c}}{\rm{m}}}^{-3}\,,\end{eqnarray} \tag{ 22 }$

and

$\begin{eqnarray}&&{\epsilon }_{* }={\epsilon }_{\mathrm{dust}}=5\times {10}^{-7}\left(\displaystyle \frac{{T}_{\mathrm{dust}}}{1000\ {\rm{K}}}\right),\end{eqnarray} \tag{ 23 }$

consistent with Nenkova et al. (2008). Here ξ_dust is a free parameter indicating the fraction of disk photons that are reprocessed by the dust torus. For the BLR photons,

$\begin{eqnarray}&&{u}_{* }={u}_{\mathrm{line}}=\displaystyle \frac{{u}_{\mathrm{line},0}}{1+{\left({r}_{\mathrm{blob}}/{r}_{\mathrm{line}}\right)}^{\beta }},\end{eqnarray} \tag{ 24 }$

where β ≈ 7.7 (Finke 2016), r_blob is the distance of the emitting blob from the BH (a free parameter). The line radii r_line and initial energy densities ${u}_{\mathrm{line},0}$ for all broad lines used are given by the Appendix of Finke (2016) relative to the Hβ line based on the composite Sloan Digital Sky Survey (SDSS) quasar spectrum of Vanden Berk et al. (2001). The parameters r_line and ${u}_{\mathrm{line},0}$ are determined from the disk luminosity using relations found from reverberation mapping, as described by Finke (2016) and in Paper II.

Since the analysis of the 2013 December 20 flare predicts very high-energy $\gamma$-rays, which could be attenuated, we include $\gamma \gamma$-absorption from the dust torus and BLR photons following Finke (2016). Absorption attenuates the emerging jet emission by a factor $\exp [-{\tau }_{\gamma \gamma }({\epsilon }_{1})]$, where ${\tau }_{\gamma \gamma }({\epsilon }_{1})$ is the absorption optical depth and ₁ is the dimensionless energy of the higher energy photon produced in the jet. We have computed the Doppler factor where $\gamma \gamma$ absorption with internal synchrotron photons becomes important (e.g., Dondi & Ghisellini 1995; Finke et al. 2008) and found that the minimum Doppler factor is much lower than the value used in our models here (Table 1). Therefore, we hereafter neglect internal synchrotron photoabsorption.

Table 1.  Free Model Parameters

Parameter (Unit)	Model A	Model B1	Model B2	Model B3
tvar (s)	1.5 × 104	3.5 × 104	2.3 × 103	9.0 × 102
B (G)	1.3	0.07	0.21	0.3
δD	30	18	49	70
rblob (cm)	1.9 × 1017	1.4 × 1017	1.6 × 1017	1.3 × 1017
ξdusta	0.1	0.1	0.1	0.1
Tdusta (K)	1470	1470	1470	1470
Ldiska (erg s−1)	1.0 × 1045	1.0 × 1045	1.0 × 1045	1.0 × 1045
D0 (s−1)	7.0 × 10−6	2.5 × 10−6	9.0 × 10−6	1.5 × 10−5
a	−4.1	+1.0	−2.0	−0.5
γinjb	1.01	1.01	1.01	2.01b
Linj (erg s−1)	8.8 × 1028	1.0 × 1032	8.5 × 1030	5.4 × 1030

Notes.

^aParameter held constant during analysis, although implemented as free. ^bThe positive slope in the ED at $\gamma \lesssim 10$ confounds the numerical normalization scheme, but the results are essentially equivalent to a Model B3 with γ_inj = 1.01 due to low particle injection number.

Download table as:  ASCIITypeset image

Emission from the BLR comes from a relatively narrow region at sub-parsec scales from the BH and different lines are produced at different radii (e.g., Peterson & Wandel 1999; Kollatschny 2003; Peterson et al. 2014), which can consist of concentric, infinitesimally thin, spherical shells for each emission line or concentric, infinitesimally thin rings Similarly, the dust torus can be modeled as a ring with an infinitesimally thin annulus or a more extended flattened disk with defined inner and outer radii. After testing each geometry, we find that in all cases $\gamma \gamma$-absorption is unimportant to model the energy range studied here, although attenuation by dust torus photons can have some effect at ≳800 GeV. The following analysis utilizes the concentric shell BLR and ring dust torus geometries, which are consistent with the emission calculations.

We first analyzed Epoch A with our model. The νF_ν SED data for Epoch A is shown in Figure 2 with our model result, and our model parameters are in Table 1. Figure 2 demonstrates that the $\gamma$-ray data is described by a 27 component EC model, including the dust torus and a stratified BLR. Scattering of dust torus photons is the largest single contributor to the production of $\gamma$-ray emission, although the sum of all BLR photon scattering is similar. The X-ray data are predominantly reproduced by the SSC process although EC/dust contributes heavily at harder energies. The IR-UV data is reproduced primarily by synchrotron radiation, which includes self-absorption at lower energies. Note that the radio emission is produced outside of the modeled region, so that the data are upper limits on the model.

Zoom In Zoom Out Reset image size

Figure 3 shows spectral data for Epoch B, which was observed during a 12 hr period on 2013 December 20, immediately following Epoch A. There are no X-ray data available during the 12 hr of peak flare for Epoch B, on 2013 December 20, and as shown in the analysis of Epoch A (Figure 2), the X-rays constrain primarily the SSC component of the spectrum. Thus, the variability timescale t_var and Doppler beaming factor ${\delta }_{{\rm{D}}}$ are not well constrained for this epoch. However, since the ED informs the SED and the parameters for the SED emission are all included in the ED as loss parameters, none of the jet components in the spectral model are fully independent from one another. We explore the parameter space with three models, each making different predictions for X-ray emission. Our model parameters for each one are in Table 1. Model B1 is the only one to reproduce the "ankle" in the $\gamma$-ray spectrum (i.e., the change in spectral index around 300 MeV), but requires a high SSC flux not previously observed for this source, leading to a very hard X-ray spectrum. For all of the Epoch B models, the $\gamma$-ray emission is dominated by scattering of BLR photons (unlike our model for Epoch A, where it was dominated by the scattering of dust photons). Model B2 is an intermediate possibility and predicts a more moderate X-ray spectrum. Model B3 has both the lowest flux and the lowest frequency peak for the SSC. Model B3 also has the smallest timescale for particle acceleration (Section 3.5), which is important since the flux-doubling timescale for the flare was quite short (∼ 2 hr). Model B3 has a high Doppler factor, indicating a very high blob velocity. All of the models for Epoch B have different ${\delta }_{{\rm{D}}}$ (and therefore ${\rm{\Gamma }}$) from the Epoch A model. However, this is not a problem, as the emission in Epochs A and B are likely produced by different blobs, which can be moving with different ${\rm{\Gamma }}$, since within a few days for the observer, the distance of a given emitting blob from the BH r_blob will change significantly. In Figure 3 the results of the Epoch A analysis are also shown for reference, demonstrating both the difference in the observed $\gamma$-ray spectra, as well as the possible changes to the X-rays. All of the B models have harder X-ray spectra for the flare than was observed during the quiescent period, which is qualitatively consistent with the optical and $\gamma$-ray hardness changes. Hayashida et al. (2015) estimate the Compton dominance A_C ≥ 300 for Epoch B, and thus classify the flare as a rare, extreme Compton flare. We predict Compton dominance values up to A_C = 1800 for Model B3 (Table 2).

Zoom In Zoom Out Reset image size

Table 2.  Calculated Parameters

Parameter (Unit)	Model A	Model B1	Model B2	Model B3
${R}_{\mathrm{Ly}\alpha }$ (cm)	2.7 × 1016	2.7 × 1016	2.7 × 1016	2.7 × 1016
RHβ (cm)	1.0 × 1017	1.0 × 1017	1.0 × 1017	1.0 × 1017
${R}_{b}^{{\prime} }$ (cm)	8.5 × 1015	1.2 × 1016	2.4 × 1015	1.2 × 1015
${\phi }_{j,\min }$ ($^\circ$)	1.3	2.4	0.4	0.5
PB (erg s−1)	8.2 × 1044	1.6 × 1042	4.5 × 1042	5.0 × 1042
Pe (erg s−1)	8.6 × 1045	2.2 × 1046	9.0 × 1045	4.3 × 1045
ζea	1.0 × 101	1.4 × 104	2.0 × 103	8.5 × 102
μab	3.8	8.8	3.6	1.7
uext (erg cm−3)	3.5 × 10−4	9.6 × 10−4	5.8 × 10−4	1.3 × 10−3
udust (erg cm−3)	1.6 × 10−4	1.6 × 10−4	1.6 × 10−4	1.6 × 10−4
uBLR (erg cm−3)	1.9 × 10−4	8.0 × 10−4	4.2 × 10−4	1.1 × 10−3
AC	4.7	1600	800	1800
Ljet (erg s−1)	5.0 × 1044	5.4 × 1046	4.5 × 1045	3.1 × 1045
σmax	4.3	2.9	2.2	1.8

Notes.

^aζ_e = u_e/u_B = P_e/P_B is the equipartition parameter, where ζ_e = 1 indicates equipartition. ^b ${\mu }_{a}\equiv ({P}_{{B}}+{P}_{e})/{P}_{a}$ is the ratio of jet to accretion power. Magnetically arrested accretion explains values of μ_a ≲ a few.

Download table as:  ASCIITypeset image

Throughout the simulations, we hold constant several parameters that the model can in principle vary (see Table 1). The parameters L_disk, ξ_dust, and T_dust are not expected to vary significantly on timescales of days. The Lorentz factor of particles injected into the base of the blob ${\gamma }_{\mathrm{inj}}$ is generally taken to be near unity as these particles are expected to originate from the thermal population in the accretion disk. However the numerical machinery can produce incomplete solutions for positive ED slopes near ${\gamma }_{\mathrm{inj}}$ when ${\gamma }_{\mathrm{inj}}/{\gamma }_{\max }$ is very small (and for negative ED slopes near ${\gamma }_{\mathrm{inj}}$ as ${\gamma }_{\mathrm{inj}}/{\gamma }_{\max }$ goes to 1). Thus, for Model B3, we used ${\gamma }_{\mathrm{inj}}=2$. However, since the number of particles injected into the blob is low compared to the total number of particles in the ED at the injection energy, the ED solution is effectively the same.

Models A, B2, and B3 have stronger second-order Fermi acceleration (as indicated by a larger D₀), while in Model B1 a first-order Fermi process dominates over adiabatic cooling (${A}_{\mathrm{ad}+\mathrm{sh}}\gt 0$). Notably, all of the models are able to represent the available broadband multiwavelength data, suggesting that it is possible for Fermi acceleration processes to meet the acceleration requirements for the observed emission.

In Figure 3, the data suggest that the peak luminosity of the synchrotron emission decreases during the flare, which might indicate a decrease in the magnetic field B for a similar particle distribution. In Model B1, this effect is exaggerated, with the particle distribution being much more energetic than in Model A, and thus the magnetic field strength is especially low (see Table 1).

The dust and BLR energy densities (u_dust and u_BLR, respectively) represent the energy in the incident photon fields in the AGN rest frame, and that are available for EC scattering. The combined BLR energy density is elevated by a factor of a few in Model B1 from the Model A values, which are within the previously observed range. The elevated energy density of EC photon sources should be expected of an orphaned $\gamma$-ray flare for a similar ED, since the Compton dominance A_C is high (see Table 2). However, in the case of Model B1, the requisite external energy densities are lower because more of the scattering energy is provided by the elevated energy in the ED. Thus Model B1 has more moderate energy requirements from the environment outside the jet than Models B2 and B3.

Models B2 and B3 employ more moderate estimates of the SSC flux, and the magnitude and spectral index in the X-rays are closer to those observed during Epoch A (Figure 3). In both Models B2 and B3, the magnetic field B ≲ 0.3 G (Table 1), is consistent with the analysis in Hayashida et al. (2015). These higher (than Model B1) values for the magnetic field strength produce similar synchrotron simulations because in Models B2 and B3, there is less energy in the jet electrons (Figure 4) and a lower ${\delta }_{{\rm{D}}}$. Additionally, the lower energy ED requires higher ${\delta }_{{\rm{D}}}$ and higher energy densities of the external radiation fields due to dust and the BLR to maintain the Compton dominance in the EC portion of the simulation.

Zoom In Zoom Out Reset image size

Since the ED for each SED is produced independently from the others, and is an integral part of the simulation, it is instructive to look at the shapes produced for each model. Figure 4 has ED curves for each SED model in Figure 3, with the same color scheme. Model A (gray) appears as a relatively simple power law with an exponential cutoff, which is due to the balance between first-order Fermi (including adiabatic expansion) and second-order Fermi accelerations, where the ratio between the two a = −4.1 < a_critical = −2, thus $N(\gamma )$ is decreasing with increasing $\gamma$. The second-order Fermi acceleration dominated portion of the solution is subdominant in Epoch A (acceleration dependencies in the ED are derived in Appendix B).

In each of the Epoch B model EDs, both terms from Equation (13) are apparent, but they are arranged differently in Figure 4 than in Figure 1. The second-order Fermi term, which was cut off around ${\gamma }_{\max }\sim {10}^{5}$ in Figure 1 is cut off at ${\gamma }_{\max }\sim {10}^{7}$ for Model B1 and ${\gamma }_{\max }\sim {10}^{6.5}$ for Models B2 and B3 (Figure 4). It is interesting to note that for Model B1, the second-order Fermi component dominates at 10² ≲ γ ≲10⁵. All of the cutoffs occur at much higher Lorentz factors than Model A, indicating that acceleration is providing more energy to the particles. So, regardless of the particular simulation, the flare requires more high-energy particles than does the preceding quiescent period. The feature in the EDs at $\gamma \approx {10}^{3.5}$ is due to the first-order Fermi/adiabatic/second-order Fermi term in Equation (13). Model B2 has a = −2, which produces a slope of 0 in that term. Models B1 and B3 have a > −2, and produce positive slopes at $\gamma \lesssim {10}^{3.5}$ in the balanced term. Positive slopes indicate that acceleration outpaces emission in that energy range. However, the combined term acceleration is damped by emission mechanisms at much lower energies than the second-order Fermi dominated term for all three Epoch B models.

The approximation Equation (13) is based on the analytic derivation in the Thomson limit. In the numerical solution, we include the full Compton cross section, and it is possible to separate the loss mechanisms. Figure 5 shows the rate of energy gain and loss at each Lorentz factor for Model B1. The red acceleration curve in Figure 5 is given by the sum of acceleration rates for first- and second-order Fermi processes (Equations (8) and (7), respectively). Acceleration dominates at Lorentz factors $\gamma \lesssim {10}^{3.5}$, where it is intersected by the blue Compton curve, which corresponds to the first turnover in the ED (Figure 4). The Compton loss rate (Equation (10)) changes at higher energies due to the Klein–Nishina effect, which allows it to trace the acceleration rate through $\gamma \sim {10}^{5}$, causing a decline in the ED over the same range (Figure 4). The synchrotron loss rate (Equation (9)) is subdominant until $\gamma \gtrsim {10}^{6}$, at which point it provides the definitive exponential cutoff to the ED.

Zoom In Zoom Out Reset image size

Our model gives r_blob ∼ (5–7) × 10¹⁷ cm during previous quiescent and flare states of 3C 279 (Paper II). If the X-ray emission does not significantly change between Epochs A and B, then the emitting region size ${R}_{b}^{{\prime} }$ may decrease as acceleration and emission increase. In all our models the minimum jet opening angle has a reasonable value.

It is clear within our analysis that first- and second-order Fermi accelerations are sufficient to power the observed flare, even assuming that the source particles are from a thermal distribution. Thus, no other forms of acceleration are necessary to explain the behavior. However, we briefly explore reconnection because it is part of the existing conversation in the literature of this particular flare's acceleration. Magnetic reconnection can in principle describe a rapid flare with very short variability timescale and a hard particle spectrum, both of which are observed. It requires a magnetization parameter, $\sigma \gg 1$. This can be constrained by (Paper I)

$\begin{eqnarray}&&\sigma \lt {\sigma }_{\max }=\displaystyle \frac{3{D}_{0}{R}_{b}^{{\prime} }}{c},\end{eqnarray} \tag{ 25 }$

based on the maximum Larmor radius of an electron that fits inside the blob. The model calculations of σ_max can be found in Table 2. In all models, σ_max ∼ 1. Since particle reconnection generally requires σ ≫ 1 , a lack of particle acceleration by reconnection is consistent with our models. First- and second-order Fermi acceleration are sufficient to explain the emission during the epochs explored here.

Equipartition is used in a wide range of astrophysical analyses, and can be used as a simplifying assumption for otherwise poorly constrained parameters (e.g., Dermer et al. 2014). There is an uncertainty involved, since the power in protons (either accelerated or "cold") is poorly constrained from the modeling (e.g., Beck & Krause 2005). In the discussion in the rest of this section, we neglect the presence of protons in the jet; however, this uncertainty should be kept in mind. In our models for 3C 279 we find the power in the field P_B is not always equivalent to the power in the jet electrons P_e, although P_B ∼ P_e within an order of magnitude or two for previously examined epochs (Paper II). Our Model A analysis shows the quiescent jet was close to equipartition (Table 2), where the equipartition parameter ${\zeta }_{e}\equiv {P}_{e}/{P}_{B}=1$ represents an equipartition state. However, in all of the B models P_B ≪ P_e (Table 2), and the jet is strongly electron dominated. The larger the frequency integrated flux, the larger the equipartition parameter for Epoch B. This is consistent with the analysis by Hayashida et al. (2015) for the Epoch B flare; they also found the jet was matter dominated. A matter dominant jet might be indicative of a larger than usual influx of electrons into the emitting region. An analysis for Epoch B was attempted, in which P_B ≈ P_e, but the jet opening angle was unphysically large, and therefore it is not presented here.

Neglecting protons, the total jet power P_B+e = P_e + P_B. For a maximally rotating BH, one expects that the accretion power P_a = L_disk/0.4, giving ${P}_{{\rm{a}}}=2.4\times {10}^{45}\ {\mathrm{ergs}}^{-1}$ for our models. For all our models except Model B1, ${\mu }_{a}\equiv {P}_{B+e}/{P}_{a}$ is a factor ∼ a few (Table 2). This is consistent with extracting spin from a BH with a magnetically arrested accretion disk (Tchekhovskoy et al. 2011). However, for Model B1, μ_a = 8.8, which may be too large too extract from the Blandford & Znajek (1977) process.

The total jet luminosity due to particle radiation (e.g., Finke et al. 2008)

$\begin{eqnarray}&&{L}_{{\rm{j}}{\rm{e}}{\rm{t}}}=\displaystyle \frac{2\pi {d}_{L}^{2}}{{{\rm{\Gamma }}}^{2}}{\int }_{0}^{\infty }\displaystyle \frac{d\epsilon }{\epsilon }\,\left({f}_{{\epsilon }_{{\rm{o}}{\rm{b}}{\rm{s}}}}^{{\rm{s}}{\rm{y}}{\rm{n}}}+{f}_{{\epsilon }_{s}}^{{\rm{S}}{\rm{S}}{\rm{C}}}+\displaystyle \sum _{j=1}^{J}{f}_{{\epsilon }_{s}}^{{\rm{E}}{\rm{C}},(j)}\right),\end{eqnarray} \tag{ 26 }$

is reported in Table 2. The radiative efficiency L_jet/P_B+e is expected to be <1. This is indeed the case, for all models except for Model B1.

Based on the excessive μ_a and large radiative efficiency (L_jet/P_B+e > 1), we believe Model B1 is unphysical and cannot explain Epoch B.

The injection luminosity L_inj of Epoch A is similar to other injection luminosities we have found for 3C 279 in other epochs (Paper II). The injection luminosities for the B models are somewhat higher due primarily to the increased rate of particle injection. This supports the idea of an influx of particles from the accretion disk area instigating the flare.

One benefit of the Fokker–Planck equation (or transport equation) formalism is the conservation of energy. We can compute the rate of energy gains and losses in electrons (in the co-moving frame: that is, injected into the system (${P}_{\mathrm{inj}}^{{\prime} }$); escaping the system (${P}_{\mathrm{esc}}^{{\prime} }$); accelerated by the second-order Fermi acceleration process (${P}_{\mathrm{sto}}^{{\prime} }$); accelerated by first-order Fermi acceleration or lost by adiabatic expansion (${P}_{\mathrm{sh}+\mathrm{ad}}^{{\prime} }$); lost by synchrotron (${P}_{\mathrm{syn}}^{{\prime} }$) or EC radiation (${P}_{\mathrm{EC}}^{{\prime} }$)). The first-order Fermi acceleration and adiabatic losses are taken together, since in practice separating them introduces an unconstrainable free parameter in our formalism. (See Paper II for the details on how these rates are calculated.)

The component powers resulting from our simulations for the epochs examined here can be found in Table 3. The relative error percentages found by adding up the powers are low (within expected numerical errors), indicating the expected energy balance.

Table 3.  Power in the Physical Components

Variable (Unit)	Model A	Model B1	Model B2	Model B3
${P}_{\mathrm{esc}}^{{\prime} }$ (erg s−1)	−6.0 × 1031	−5.3 × 1035	−9.6 × 1033	−2.9 × 1033
${P}_{\mathrm{inj}}^{{\prime} }$ (erg s−1)	8.8 × 1029	1.0 × 1032	8.6 × 1030	1.3 × 1031
${P}_{\mathrm{sto}}^{{\prime} }$ (erg s−1)	3.2 × 1043	1.9 × 1044	7.2 × 1042	1.4 × 1042
${P}_{\mathrm{sh}+\mathrm{ad}}^{{\prime} }$ (erg s−1)	−2.4 × 1043	+4.7 × 1043	−3.6 × 1042	−1.8 × 1041
${P}_{\mathrm{syn}}^{{\prime} }$ (erg s−1)	−2.4 × 1041	−6.3 × 1041	−1.3 × 1040	−3.0 × 1039
${P}_{\mathrm{EC}}^{{\prime} }$ (erg s−1)	−8.6 × 1041	−2.4 × 1044	−3.6 × 1042	−1.3 × 1042
%σerr	2.9	0.03	0.03	0.04

Download table as:  ASCIITypeset image

Since in all cases ${P}_{\mathrm{sto}}^{{\prime} }\gt {P}_{\mathrm{sh}+\mathrm{ad}}^{{\prime} }$, the particles are primarily accelerated by the second-order Fermi process, rather than first-order Fermi acceleration. This is consistent with our results for other SEDs in Papers I and II. First-order Fermi acceleration dominates over adiabatic losses for only Model B1, since this is the only simulation where ${P}_{\mathrm{sh}+\mathrm{ad}}^{{\prime} }\gt 1$. This is consistent with this being the only model with ${A}_{\mathrm{ad}+\mathrm{sh}}\gt 0$ (Table 1). For the other models, adiabatic losses are an important energy loss mechanism. For Model A, adiabatic losses dominate over radiative losses by a factor of ∼100. For Models B2 and B3, they are the same order of magnitude as radiative losses. We note the Epoch B models have a much larger ${P}_{\mathrm{EC}}^{{\prime} }/{P}_{\mathrm{syn}}^{{\prime} }$ than the Epoch A simulation, as expected, since Epoch B has a much larger A_C than Epoch A.

The injection energy is not an important contribution to the energy budget. The escape power is always slightly larger in magnitude than the injection power because the escape occurs in the Bohm limit, meaning higher energy particles are preferentially lost from the electron population (Paper I).

An additional benefit of the particle transport method we employ is the ability to compare the timescales for each physical process, as expressed by the coefficients in Equations (7)–(10), which have units of [s⁻¹] (Paper I). Of particular interest during the flare is the acceleration timescale, because the flare duration is short (∼12 hr) and the flux-doubling timescale is rapid (∼2 hr). Hence the acceleration mechanism providing energy to the flare must act on commensurate timescales (e.g., Hayashida et al. 2015; Paliya et al. 2016).

The mean timescale for the second-order Fermi acceleration of electrons via hard-sphere scattering with magnetohydrodynamic (MHD) waves is computed in the frame of the observer using

$\begin{eqnarray}&&{t}_{\mathrm{sto}}=\displaystyle \frac{1+z}{{\delta }_{{\rm{D}}}}\displaystyle \frac{1}{4{D}_{0}},\end{eqnarray} \tag{ 27 }$

where the factor of 4 comes from the derivative separating the broadening and drift coefficient components of the second-order Fermi up-scattering. First-order Fermi acceleration is closely linked with adiabatic expansion in the transport model during analysis of the SED, since both processes have the same energy dependence (Paper I). However, we constrain the first-order Fermi acceleration timescale by assuming that the energy loss rate due to adiabatic expansion during the flare, ${A}_{\mathrm{ad}}^{\mathrm{flare}}$, is the same as during quiescence, ${A}_{\mathrm{ad}}^{\mathrm{quiescent}}$. Thus, the first-order Fermi acceleration timescale during the flare can be constrained using Model A as the limiting case, which yields

$\begin{eqnarray}&&{t}_{\mathrm{sh}}\lesssim \displaystyle \frac{1+z}{{\delta }_{{\rm{D}}}}\displaystyle \frac{1}{{A}_{\mathrm{ad}\ +\ \mathrm{sh}}^{\mathrm{flare}}-{A}_{\mathrm{ad}}^{\mathrm{quiescent}}},\end{eqnarray} \tag{ 28 }$

where ${A}_{\mathrm{ad}}^{\mathrm{quiescent}}\lesssim {A}_{\mathrm{sh}+\mathrm{ad}}^{\mathrm{quiescent}}$ (values in Table 1). The total acceleration timescale, t_acc, depends on both the first- and second-order Fermi timescales, t_sh and t_sto, respectively, via

$\begin{eqnarray}&&{t}_{\mathrm{acc}}=\displaystyle \frac{1}{{t}_{\mathrm{sh}}^{-1}+{t}_{\mathrm{sto}}^{-1}},\end{eqnarray} \tag{ 29 }$

where the value of t_sh is an upper limit given by Equation (28). The values provided in Table 4 for t_acc are therefore also upper limits. Thus, first- and second-order Fermi acceleration rates, as included in the model, are sufficiently rapid to produce the observed flux-doubling timescale for the 2013 flare.

Blazar jets are thought to contain standing shocks (e.g., Marscher et al. 2008). A superluminal blob passing through a standing shock will increase the amount of first-order Fermi acceleration (e.g., Marscher 2012), which is consistent with the Epoch B models, compared with our model for Epoch A (Table 1). MHD simulations demonstrate that first-order Fermi acceleration gives rise to higher levels of stochastic turbulence downstream of the shock (e.g., Inoue et al. 2011), which agrees with the increase in second-order Fermi acceleration during the flare analyses (Table 1). Both the first- and second-order Fermi acceleration contribute to the higher maximum electron Lorentz factor as well as the higher number of high-energy particles (see Figure 4). These are important to the higher frequency position of each emission component in the SED.

As particle energy increases, so does the Larmor radius (${r}_{{\rm{L}}}=\gamma {m}_{e}{c}^{2}{q}^{-1}{B}^{-1}$). A particle with a larger Larmor radius, will travel preferentially closer to the edge or sheath of the jet (e.g., Hillas 1984). If the magnetic field is radially dependent (stronger near the jet core), then the apparent magnetic field of the jet may be lower by a factor of a few when more particles spend more time near the edge of the emitting region in the jet. Massaro et al. (2004) discuss a log-parabolic particle distribution (which can be formed by second-order Fermi acceleration; Tramacere et al. 2011) as one in which the confinement efficiency of a collimating magnetic field decreases with increasing gyro-radius. This physical interpretation can explain the smaller magnetic field (Table 1) and the greater loss to escaping particles (Table 3), during the flare.

In Model B1, the first-order Fermi acceleration contributes more than in Models B2 and B3. This is consistent with the smaller bulk Lorentz factor (recall we assume ${\delta }_{{\rm{D}}}={\rm{\Gamma }}$, Table 1), since increased first-order Fermi acceleration indicates a given particle is scattered through the shock front more times. Thus, even larger Larmor radii may indicate that emitting particles occupy some of the jet sheath, where the magnetic field is significantly lower, which explains the 95% drop in field strength in the analysis of the flare (Table 1). While the interpretation is somewhat more extreme for Model B1, the second-order Fermi acceleration coefficient and variability timescale are more similar to Model A than are Models B2 and B3. However, Model B1 has an unphysical radiative efficiency (Table 2), suggesting that the simulation of extreme SSC emission employed therein is not an appropriate representation of the data.

Models B2 and B3 are also well described by a large influx of material injected into the base of the jet, which may strengthen the effects of particle acceleration at a shock while temporarily weakening the comparative power of the magnetic field, without initiating widespread reconnection. It would be particularly interesting to study high-energy polarization in blazar flares, especially as several new instruments are in the planning stages. X-ray and $\gamma$-ray polarimetry can more definitively separate leptonic from hadronic models, elucidate the shock versus magnetic reconnection debate, and provide a new avenue through which to explore the geometry of the acceleration/emission region (e.g., Dreyer & Böttcher 2019).

Besides reporting on the 2013 December epochs that we model here (and other epochs that we do not consider), Hayashida et al. (2015) also modeled these epochs with the BLAZAR code (Moderski et al. 2003). They used a double broken power-law ED to model Epoch A, and a broken power-law distribution to model Epoch B. Similar to our modeling, Hayashida et al. (2015) found that EC/dust dominated the $\gamma$-ray emission during Epoch A and EC-BLR dominated during Epoch B. For Epoch A, they found a larger magnetic field, lower ${\rm{\Gamma }}$, and larger r_blob. They modeled Epoch B with two sets of parameters. Their model parameters for this epoch are generally similar to ours for our Models B2 and B3, although or magnetic field for our Model B1 is quite lower compared to their models. Our models for Epoch B have similar r_blob, but larger ${\rm{\Gamma }}$ for our Models B2 and B3.

Asano & Hayashida (2015) model the 2013 December flare with a time-dependent model (Asano et al. 2014) that treats particle acceleration quite similar to ours. They model both the SED and the $\gamma$-ray light curve. In their model, scattering of a UV field (presumably representing a BLR) dominates the $\gamma$-ray emission and they have ${\rm{\Gamma }}=15$ and r_blob ∼ 6 × 10¹⁶ cm, so their emitting region is closer to the BH than in our models (to a degree consistent with our use of a stratified BLR, Paper II). Their model successfully reproduces the observed LAT γ-ray light curve, although the high ${\rm{\Gamma }}$ may be a problem for models of jet acceleration by magnetic dissipation.

The 2013 flare was further examined by Paliya et al. (2016) using time-dependent lepto-hadronic and two-zone leptonic models. Those authors examined a 3 day window that included the flare, and noted a 3 hr flux-doubling timescale. They adopted a smooth broken power law as the particle distribution, and found that the isolated flare could represent γ-ray emission from a smaller blob, while the rest of the spectrum was produced in a larger volume. Another possibility is that proton acceleration could power the enhanced $\gamma$-ray emission. The inclusion of significant non-flare data in the analysis of Paliya et al. (2016) raises concerns about identifying the physics of the flare specifically. Hence the shape of the ED and the nature of the associated particle acceleration mechanism(s) during the 2013 flare from 3C 279 are still open questions.

Our model from Paper II included a numerical steady-state solution to a particle transport equation that included first- and second-order Fermi acceleration and particle escape. Particles in this single homogeneous blob can lose energy to adiabatic expansion, synchrotron radiation, and inverse-Compton scattering of incident photon fields, including SSC, EC/disk, EC/dust, and EC/BLR, where the BLR is composed of 26 individual lines radiating at infinitesimally thin concentric shells. We added to this $\gamma \gamma$-absorption due to the accretion disk, dust torus, and stratified BLR consistent with the emission and Compton scattering geometry (Finke 2016). We applied our model to the SED of the extremely Compton-dominant flare of 3C 279 observed on 2013 December 20 (Hayashida et al. 2015). The preceding 3 days of quiescent data are similarly analyzed as a baseline for the flare simulations. We derive a simplified version of the Thomson regime approximation (Paper I), which assists in the interpretation of the ED. Our primary results are as follows:

• 1.  
  It is possible to simulate the acceleration in the flare SED with reasonable levels of only first- and second-order Fermi processes (with the latter process dominating in our models). Acceleration by reconnection is not needed, and the maximum magnetization parameter (σ) found from our model parameters is consistent with this.
• 2.  
  It is possible for the BLR to be the dominant EC component without significant $\gamma \gamma$ absorption from BLR photons.
• 3.  
  The quiescent period displays electron and field powers near equipartition, while the flare is strongly electron dominated.
• 4.  
  There is insufficient energy available for the X-rays to have undergone a flare comparable to the $\gamma$-ray flare simultaneously.
• 5.  
  The simplified ED analysis clarifies that first- and second-order Fermi acceleration can influence different components of the overall ED.
• 6.  
  Based on our modeling of Epoch B, the νF_ν flux at ≈10 MeV is likely $\lt {10}^{-9}\ {\mathrm{ergs}}^{-1}\ {\mathrm{cm}}^{-2}$ (i.e., Model B1 is highly unlikely). The emission in this energy range could be probed by a future $\gamma$-ray mission such as e-Astrogam or the All-sky Medium Energy Gamma-ray Observatory (AMEGO), which could provide further constraints on blazar SEDs.

This model has been informative in the analysis of this particularly unusual flare. We plan to use the same model to analyze other intriguing blazar behaviors. Additionally the simplified analytic ED can be applied to astrophysical jets more broadly, where both first- and second-order Fermi acceleration are expected to contribute.