text:
# Modelling 20 years of synchrotron flaring in the jet of 3C 273

## 1 Introduction

The theory of synchrotron emission in relativistic jets was developed around 1980 mainly by Blandford & Königl (BK79 (1979)), Marscher (M80 (1980)) and Königl (K81 (1981)). In 1985, both Marscher & Gear (MG85 (1985), hereafter MG85) and Hughes et al. (HAA85 (1985)) proposed a model for the emission of a shock wave propagating down a simple relativistic jet. Both groups assumed that the jet is confined to a cone of constant opening angle and that the jet flow is adiabatic. The computer code of Hughes et al. (1989a ) was able to describe the low-frequency flux density and polarization variability of BL Lacertae (Hughes et al. 1989b ). The model of MG85 has however the advantage to include the effects of synchrotron and inverse-Compton energy losses of the electrons, which cannot be neglected at higher frequencies. It provides a simple explanation for the distinct components in the jet observed using very long baseline interferometry (VLBI). Multi-wavelength total flux measurements were found to be difficult to use to constrain the shock model of MG85, because the emission of all distinct features in the jet often overlap to form a nearly flat total spectrum, as illustrated by Marscher (M88 (1988)) for the quasar NRAO 140.

It is only since 1995 that the very well sampled total flux millimetre and radio light-curves of a few sources allowed to study the spectral evolution of individual synchrotron outbursts. By subtracting a quiescent spectrum assumed to be constant, Litchfield et al. (LSR95 (1995)) and Stevens et al. (SLR95 (1995), SLR96 (1996), SRG98 (1998)) could follow the early evolution of single synchrotron outbursts in 3C 279, PKS 0420$``$014, 3C 345 and 3C~273, respectively. These studies gave additional support to the MG85 model, but failed to constrain it strongly, mainly because they could only follow the evolution of an outburst until the onset of the next one.

To overcome this problem, Türler et al. (1999b , hereafter Paper I) proposed a new approach to derive the observed properties of synchrotron outbursts, which consists in a complete decomposition of long-term multi-wavelength light-curves into a series of self-similar flaring events. Two different approaches are presented in Paper I to describe the evolution of these events with both time and frequency. The so called “light-curve approach” is model-independent and describes empirically the shape of the light-curves of individual outbursts at different frequencies. A second approach based on three-stage shock models, like those of MG85 or Valtaoja et al. (VTU92 (1992)), describes directly the evolution of the flaring synchrotron spectrum. Finally, the “hybrid approach” of Türler (T00 (2000)) models strictly the shape of the synchrotron spectrum, but leaves the evolution of the spectral turnover as free as possible. The results of this third approach, based on self-similar outbursts as in Paper I, were found to be in very good agreement with the spectral evolution expected by the shock model of MG85.

After this last test, the next step, which is presented here, is to adapt the shock model of MG85 in order to describe physically both the average evolution of the outbursts and their individual specificity. This generalized shock model, described in Sect. 3, takes into account the effects of an accelerating or decelerating synchrotron source in a curved, non-conical and non-adiabatic jet. This model is confronted to the observations by fitting seventeen distinct outbursts simultaneously to thirteen submillimetre-to-radio light-curves of 3C 273. The proper parameterization to achieve this fit is described in Sect. 4 and the results of this light-curve decomposition are given in Sect. 5. In Sect. 6, we discuss the implications of our results, both on the global properties of the inner jet and on the peculiarities of individual outbursts. The main results of this study are summarized in Sect. 7.

Throughout this paper, we use the convention $`S_\nu \nu ^{+\alpha }`$ for the spectral index $`\alpha `$ and we use “$`\mathrm{lg}`$”, rather than “$`\mathrm{log}`$”, for the decimal logarithm “$`\mathrm{log}_{10}`$”, because of a lack of space in tables and long equations.

## 2 Observational material

The light-curves fitted here are part of the multi-wavelength database of 3C 273 presented by Türler et al. (1999a ). The 12 light-curves from 5 GHz to 0.35 mm are as described in Paper I, except that we now consider the observations up to 1999, including the most recent measurements from the Metsähovi Radio Observatory in Finland and from the University of Michigan Radio Astronomy Observatory (UMRAO).

We extend to lower frequencies the analysis of Paper I by adding a new light-curve at 2.7 GHz. This light-curve is constituted of observations from the Green Bank Interferometer (GBI) and from the 100 m telescope at Effelsberg in Germany (Reich et al. RRP98 (1998)). The fluxes obtained with the GBI at 2.7 GHz and 2.25 GHz include both 3C 273B (the inner jet) and 3C 273A (the hot spot at the far end of the jet), but these two components combine partially out of phase (Ghigo F.D., private communication). As a consequence, the GBI measurements are only part of the total flux of 3C 273 and have to be multiplied by a scaling factor.

To scale the GBI measurements at 2.25 GHz we use the contemporaneous Effelsberg single dish observations at 2.7 GHz. The nearly flat spectral index of 3C 273 at this frequency (Türler et al. 1999a ) allows us to neglect the small difference in frequency. We obtain that the GBI fluxes at 2.25 GHz have to be multiplied by 1.66 to fit the Effelsberg observations. A similar calibration of the earlier GBI measurements at 2.7 GHz is not possible due to the lack of contemporaneous single dish observations. We therefore let the scaling factor of these fluxes as a free parameter of our fit.111We obtain a best-fit value of $`1.12`$ for this factor with the jet model presented in Sect. 5. To smooth out the dips in the GBI light-curves (cf. Türler et al. 1999a ), we average all GBI observations into bins of 10 days. Finally, we do not consider the GBI measurements before 1980, to avoid a flux increase in the light-curve due to an outburst that started before 1979. We end up with a total of 5234 observational points to constrain the shock model.

## 3 Physical shock model

In the original shock model of MG85, the shock front was assumed to propagate with constant speed in a straight conical and adiabatic jet. Partial generalizations of this model were already derived by Marscher (M90 (1990)) and Marscher et al. (MGT92 (1992)) in the case of a bending jet, and by Stevens et al. (SLR96 (1996)) in the two cases of a straight non-adiabatic jet and a curved adiabatic jet. Here we further generalize the shock model of MG85 to account for the effects of an accelerating or decelerating shock front in a curved, non-conical and non-adiabatic jet. In Sect. 3.1, we describe the typical three-stage evolution of all outbursts, whereas in Sect. 3.2 we show how the physical conditions at the onset of the shock can influence the evolution of individual outbursts.

### 3.1 Typical three-stage evolution

Following MG85, we consider the synchrotron radiation emitted just behind a shock front in a cylindric portion of a jet having a radius $`R`$ across the jet and a length $`x`$ along the jet axis, as illustrated in Fig. 1. Within this volume, we assume that the magnetic field $`B`$ is uniform in strength and nearly random in direction and that the relativistic electrons have a power law energy distribution of the form $`N(E)=KE^s`$, with $`N(E)dE`$ being the number density of the electrons. Measured in the rest frame of the quasar, these synchrotron emitting electrons have a relativistic bulk velocity $`\beta =v/c`$, with $`c`$ being the speed of light, and a corresponding Lorentz factor $`\mathrm{\Gamma }=(1\beta ^2)^{1/2}`$. As a consequence, the synchrotron emission observed with an angle $`\theta `$ to the jet axis is Doppler boosted with a bulk Doppler factor given by $`𝒟=\mathrm{\Gamma }^1(1\beta \mathrm{cos}\theta )^1`$. We further assume that the jet opening half-angle $`\varphi `$ is smaller than $`\theta `$ at any time after the onset of the shock and that $`\theta `$ itself is small enough to verify $`\mathrm{sin}(\theta +\varphi )<1/\mathrm{\Gamma }`$, so that the line-of-sight depth of the emitting region is directly proportional to its thickness $`x`$ (Marscher et al. MGT92 (1992)).

The observed optically thin flux density $`S_\nu `$ and turnover frequency $`\nu _\mathrm{m}`$ of the self-absorbed synchrotron spectrum are then given by

$`S_\nu `$ $``$ $`R^2xKB^{(s+1)/2}𝒟^{(s+3)/2}\nu ^{(s1)/2}\text{and}`$ (1)
$`\nu _\mathrm{m}`$ $``$ $`\left(xKB^{(s+2)/2}𝒟^{(s+2)/2}\right)^{2/(s+4)},`$ (2)
where $`K`$ and $`B`$ are measured in the comoving frame of the emitting plasma and $`S_\nu `$, $`\nu `$, $`\nu _\mathrm{m}`$, $`R`$, $`x`$ and $`𝒟`$ are measured in the observer’s frame.

The thickness $`x`$ of the emitting region behind the shock front is the crucial parameter that changes the outburst’s evolution from one stage to the other. If radiative losses are the dominant cooling process of the electrons accelerated at the shock front, $`x`$ is approximatively given by $`2v_{\mathrm{rel}}t_{\mathrm{cool}}`$, where $`v_{\mathrm{rel}}`$ is the excess velocity of the shock front relative to the emitting plasma and $`t_{\mathrm{cool}}`$ is the typical cooling time of the electrons.222The centre of the emitting region is given by the typical distance $`v_{\mathrm{rel}}t_{\mathrm{cool}}`$ from the shock front reached by the electrons in the upstream direction before loosing substantially their energy. The full width $`x`$ of the emitting region is thus about twice that distance.  During the first stage of the shock evolution Compton scattering is the dominant cooling process of the electrons and the Compton limited thickness $`x_1`$ of the emitting region is given by $`x_1u_{\mathrm{ph}}^1B^{1/2}𝒟^{1/2}\nu ^{1/2}`$, where the energy density of the synchrotron photons $`u_{\mathrm{ph}}`$ can be expressed as $`u_{\mathrm{ph}}K(B^{3s+7}R^{s+5})^{1/8}`$ (Eq. (13) of MG85).333This calculation assumes that only first order Compton scattering is important. The effect of higher order Compton losses would be to steepen even more the initial rise of the spectrum (Marscher et al. MGT92 (1992)).  By combining these two equations, we obtain444$`x_1`$ and $`x_2`$ in Eqs (3) and (4) are actually expressed in the comoving frame of the emitting plasma. The transformation of $`x_i`$ to the observer’s frame would be obtained by dividing these expressions by the Lorentz factor $`\mathrm{\Gamma }`$ to correct for length contraction and by multiplying them by the Doppler factor $`𝒟`$ to correct for front-to-back time delay (Marscher A.P., private communication). This tricky transformation is due to the non-spherical shape of the emitting slab, which appears rotated in the observer’s frame (cf. Marscher M87 (1987); Marscher et al. MGT92 (1992)). However, since the viewing angle $`\theta `$ is assumed to be smaller than about $`1/\mathrm{\Gamma }`$, $`𝒟`$ is proportional to $`\mathrm{\Gamma }`$ (e.g. Marscher M80 (1980)), so that these two corrections cancel out when using proportionalities.

$$x_1R^{(s+5)/8}K^1B^{3(s+1)/8}𝒟^{1/2}\nu ^{1/2}.$$
(3)
The shock evolution enters the synchrotron stage as soon as the photon energy density $`u_{\mathrm{ph}}`$ is equal to the magnetic field energy density $`u_B=B^2/(8\pi )`$. During this second stage, synchrotron radiation is the dominant energy loss process of the electrons. Due to the similarity between Compton and synchrotron expressions of the cooling time $`t_{\mathrm{cool}}`$ of the electrons, the synchrotron limited thickness $`x_2`$ of the emitting region is simply obtained by replacing $`u_{\mathrm{ph}}`$ in the expression of $`x_1`$ by $`u_B`$, so that $`x_2`$ becomes

$$x_2B^{3/2}𝒟^{1/2}\nu ^{1/2}.$$
(4)
During the final stage of the shock evolution synchrotron losses become less important than adiabatic expansion losses and the shock is assumed to remain self-similar so that $`x_3R`$.

By substituting $`x`$ in Eqs (1) and (2) by the corresponding expression of $`x_i`$ $`(i=1,2,3)`$ for each of the three stages $`i`$ of the shock evolution we obtain:

$`S_{\nu ,1}`$ $``$ $`R^{(11s)/8}B^{(s+1)/8}𝒟^{(s+4)/2}\nu ^{s/2}`$ (5)
$`S_{\nu ,2}`$ $``$ $`R^2KB^{(s2)/2}𝒟^{(s+4)/2}\nu ^{s/2}`$ (6)
$`S_{\nu ,3}`$ $``$ $`R^3KB^{(s+1)/2}𝒟^{(s+3)/2}\nu ^{(s1)/2}`$ (7)
$`\nu _{\mathrm{m},1}`$ $``$ $`R^{1/4}B^{1/4}𝒟^{(s+3)/(s+5)}`$ (8)
$`\nu _{\mathrm{m},2}`$ $``$ $`[K^2B^{s1}𝒟^{s+3}]^{1/(s+5)}`$ (9)
$`\nu _{\mathrm{m},3}`$ $``$ $`[RKB^{(s+2)/2}𝒟^{(s+2)/2}]^{2/(s+4)}.`$ (10)
Finally, by using these six last equations, it is straightforward to obtain the expressions for $`S_{\mathrm{m},i}S_{\nu ,i}(\nu _{\mathrm{m},i})`$ as

$`S_{\mathrm{m},1}`$ $``$ $`R^{11/8}B^{1/8}𝒟^{(3s+10)/(s+5)}`$ (11)
$`S_{\mathrm{m},2}`$ $``$ $`R^2[K^5B^{2s5}𝒟^{3s+10}]^{1/(s+5)}`$ (12)
$`S_{\mathrm{m},3}`$ $``$ $`[R^{2s+13}K^5B^{2s+3}𝒟^{3s+7}]^{1/(s+4)}.`$ (13)
We now assume that each of the quantities $`K`$, $`B`$ and $`𝒟`$ evolves as a power-law with the radius $`R`$ of the jet and we parameterize this evolution as

$$KR^kBR^b𝒟R^d.$$
(14)
If we replace these relations into the expressions for $`\nu _{\mathrm{m},i}`$ and $`S_{\mathrm{m},i}`$ (Eqs (8) to (13)), the turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ of the self-absorbed synchrotron spectrum will also evolve as a power-law with $`R`$ during each stage $`i`$ of the shock evolution as

$$\nu _{\mathrm{m},i}R^{n_i}\text{and}S_{\mathrm{m},i}R^{f_i}S_{\mathrm{m},i}\nu _{\mathrm{m},i}^{f_i/n_i},$$
(15)
where the exponents $`n_i`$ and $`f_i`$ are given by

$`n_1`$ $`=`$ $`(b+1)/4d(s+3)/(s+5)`$ (16)
$`n_2`$ $`=`$ $`[\mathrm{\hspace{0.17em}2}k+b(s1)+d(s+3)]/(s+5)`$ (17)
$`n_3`$ $`=`$ $`[\mathrm{\hspace{0.17em}2}(k1)+(b+d)(s+2)]/(s+4)`$ (18)
$`f_1`$ $`=`$ $`(11b)/8d(3s+10)/(s+5)`$ (19)
$`f_2`$ $`=`$ $`2[\mathrm{\hspace{0.17em}5}k+b(2s5)+d(3s+10)]/(s+5)`$ (20)
$`f_3`$ $`=`$ $`[\mathrm{\hspace{0.17em}2}s+135kb(2s+3)d(3s+7)]/(s+4).`$ (21)
Marscher (M87 (1987)) shows that even a geometrically thin source in a relativistic jet appears inhomogeneous when observed with a small angle $`\theta `$ between the jet axis and the line-of-sight. The optical depth of the source will therefore depend on the frequency, which has the effect of broadening the self-absorption turnover and leads to a lower optically thick spectral index $`\alpha _{\mathrm{thick}}`$ than the typical value of $`+5/2`$, which holds for a homogeneous synchrotron source. To estimate the value of the spectral index below the turnover frequency $`\nu _\mathrm{m}`$, we can cut the source into many self-similar cylindric portions of the jet having a length $`lx`$ proportional to their radius $`R`$ and chosen small enough for their synchrotron emission to be homogeneous. The emitted spectrum of each section will have a self-absorption turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ depending on $`R`$ according to (15) with $`i=3`$, because in this case $`lR`$ replaces $`x`$ in Eqs (1) and (2). The inhomogeneous source behind the shock front is therefore expected to have an optically thick spectral index of $`\alpha _{\mathrm{thick}}=f_3/n_3`$ due to the superimposition of the homogeneous spectra of the individual sections. However, the finite size of the emitting region should limit the frequency range over which this flatter spectral index pertains. We therefore expect a spectral break at a frequency $`\nu _\mathrm{h}`$ ($`<\nu _\mathrm{m}`$), at which the spectral index $`\alpha _{\mathrm{thick}}=f_3/n_3`$ returns to its homogeneous value of $`+5/2`$.

A high-frequency spectral break is also expected in the optically thin part of the spectrum due to a change in the electron energy distribution induced by synchrotron and/or Compton losses. In the case of continuous injection or re-acceleration of electrons suffering radiative losses, the optically thin spectral index $`\alpha _{\mathrm{thin}}`$ is expected to steepen by a value of $`1/2`$ above a frequency $`\nu _\mathrm{b}`$ (Kardashev K62 (1962)). The break frequency $`\nu _\mathrm{b}`$ of the spectrum in the observer’s frame is related to the break energy $`E_\mathrm{b}`$ of the electron energy distribution as $`\nu _\mathrm{b}B𝒟E_\mathrm{b}^2`$ (e.g. Marscher M80 (1980)). For our analysis at submillimetre-to-radio frequencies, the evolution of $`\nu _\mathrm{b}`$ with the jet radius $`R`$ is only relevant during the final stage of the shock evolution. For adiabatic expansion in two dimensions, the energy $`E`$ of the electrons decreases with $`R`$ as $`ER^{2/3}`$ (e.g. Gear G88 (1988)). By using this relation and (14), we find that the evolution of the break frequency $`\nu _\mathrm{b}`$ during the adiabatic stage is given by

$$\nu _{\mathrm{b},3}R^{n_\mathrm{b}}\text{with}n_\mathrm{b}=(4/3+b+d).$$
(22)
Until now, we related the spectral turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ and the high-frequency spectral break $`n_\mathrm{b}`$ to the radius $`R`$ of the emitting region. But since we are interested in the temporal evolution of these quantities, we have to relate $`R`$ to the observed time $`t`$ after the onset of the outburst. According to the geometry of Fig. 1 and the basic principles of superluminal motion (e.g. Pearson & Zensus PZ87 (1987)), we obtain

$$t=\frac{(1+z)\mathrm{sin}\theta }{\beta _{\mathrm{app}}c}(LL_{\mathrm{on}})=\frac{(1+z)}{\beta c\mathrm{\Gamma }𝒟}(LL_{\mathrm{on}}),$$
(23)
where $`\beta _{\mathrm{app}}`$ is the apparent transverse velocity of the source in units of $`c`$ and $`L`$ measures the distance along the jet axis in the rest frame of the quasar. To continue working only with proportionalities, we are forced to assume that $`L_{\mathrm{on}}`$ is small with respect to $`L`$ during the major part of the shock evolution, so that $`LL_{\mathrm{on}}`$ tends rapidly towards $`L`$ during, or just after, the Compton stage. Under this assumption and by considering that $`\beta \mathrm{\Gamma }𝒟`$ is proportional to $`𝒟^2`$ for $`\mathrm{\Gamma }1`$ and $`\theta \stackrel{<}{}1/\mathrm{\Gamma }`$ (e.g. Marscher M80 (1980)), we have $`t𝒟^2L`$. Finally, by parameterizing the opening radius $`R`$ of the jet with the distance $`L`$ as $`RL^r`$ and by remembering that $`𝒟R^d`$ ((14)), we obtain

$$t𝒟^2R^{1/r}R^\rho \text{with}\rho =(2rd+1)/r.$$
(24)
With this relation, the three-stage evolution of (15) can now be expressed as a function of time $`t`$ rather than radius $`R`$, as

$$\begin{array}{cc}\nu _{\mathrm{m},\mathrm{i}}t^{\beta _i}\text{with}\beta _in_i/\rho \hfill & \\ S_{\mathrm{m},\mathrm{i}}t^{\gamma _i}\text{with}\gamma _if_i/\rho \hfill & \end{array}\}\frac{\gamma _i}{\beta _i}=\frac{f_i}{n_i},$$
(25)
where $`\beta _i`$ and $`\gamma _i`$ are defined as in Paper I and depend now on the five exponents $`s`$, $`r`$, $`k`$, $`b`$ and $`d`$ via the Eqs (16) to (21) for $`n_i`$ and $`f_i`$ $`(i=1,2,3)`$.

### 3.2 Specificity of individual outbursts

We described in the previous section how the spectral turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ of the emitted synchrotron spectrum behaves according to the physical properties of the jet and during each of the three stages of the shock evolution. What we describe here is how the transitions from one stage to the other are influenced by the physical quantities at the onset of the shock. The first transition from the Compton to the synchrotron stage is characterized by the condition $`u_{\mathrm{ph}}=u_B`$, which corresponds to $`x_1=x_2`$ and can be expressed by the proportionality

$$R_{1|2}^{(s+5)}K_{1|2}^8B_{1|2}^{3(s3)},$$
(26)
where $`K_{1|2}K(t_{1|2})`$, $`R_{1|2}R(t_{1|2})`$ and $`B_{1|2}B(t_{1|2})`$, with the subscript $`1|2`$ referring to the transition from the first to the second stage of the shock evolution at a time $`t_{1|2}`$ after the onset of the shock.

A similar expression can be derived for the second transition from the synchrotron to the adiabatic stage by imposing that $`x_2=x_3`$, which corresponds to $`R_{2|3}B_{2|3}^{3/2}𝒟_{2|3}^{1/2}\nu _{2|3}^{1/2}`$, where $`\nu _{2|3}\nu _\mathrm{m}(t_{2|3})`$, with the subscript $`2|3`$ referring to the transition from the second to the third stage of the shock evolution. By replacing $`\nu _{2|3}`$ in the condition above by the proportionality of (9) for $`\nu _{\mathrm{m},2}(t_{2|3})`$ – or equivalently by that of (10) for $`\nu _{\mathrm{m},3}(t_{2|3})`$ – we obtain

$$R_{2|3}^{(s+5)}K_{2|3}B_{2|3}^{2s+7}𝒟_{2|3}^1.$$
(27)
We can now further parameterize the proportionalities of (14), as

$$KK_{\mathrm{on}}R^kBB_{\mathrm{on}}R^b𝒟𝒟_{\mathrm{on}}R^d,$$
(28)
where the subscript “$`\mathrm{on}`$” refers to the onset of the shock. By including these relations expressed at the particular radii $`R=R_{1|2}`$ and $`R=R_{2|3}`$ into (26) and (27), respectively, we obtain

$$R_{1|2}^{\zeta _{1|2}}K_{\mathrm{on}}^8B_{\mathrm{on}}^{3(s3)}\text{and}R_{2|3}^{\zeta _{2|3}}K_{\mathrm{on}}B_{\mathrm{on}}^{2s+7}𝒟_{\mathrm{on}}^1,$$
(29)
where the exponents $`\zeta _{1|2}`$ and $`\zeta _{2|3}`$ are given by

$`\zeta _{1|2}`$ $``$ $`8k+3b(s3)(s+5)`$ (30)
$`\zeta _{2|3}`$ $``$ $`k+b(2s+7)d(s+5).`$ (31)
The expressions of (29) relate the radii $`R_{1|2}`$ and $`R_{2|3}`$ of the jet at which the two transitions $`1|2`$ and $`2|3`$ occur to the values of the physical quantities $`K`$, $`B`$ and $`𝒟`$ at, or just after, the onset of the shock. Any change of these quantities at the onset of the shock will therefore have the effect of displacing the position in the jet (via $`LR^{1/r}`$) at which the two transitions occur. A stronger magnetic field $`B_{\mathrm{on}}`$, for instance, will have the effect to prolong the synchrotron stage, by making it start slightly further upstream in the jet and finish much further downstream for a typical value of $`s2`$ and comparable values of $`\zeta _{1|2}`$ and $`\zeta _{2|3}`$.

The values of $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$, will also influence the times $`t_{1|2}`$ and $`t_{2|3}`$ after the onset of the shock at which the two transitions are observed, as well as the frequencies $`\nu _{1|2}`$ and $`\nu _{2|3}`$ and flux densities $`S_{1|2}`$ and $`S_{2|3}`$ of the self-absorption turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ at these times. By including the proportionalities of (28) into (24) for the observed time $`t`$ and into Eqs (9) and (12) for the evolution during the synchrotron stage of turnover frequency $`\nu _{\mathrm{m},2}`$ and the corresponding flux density $`S_{\mathrm{m},2}`$, we obtain for the first transition

$`t_{1|2}`$ $``$ $`𝒟_{\mathrm{on}}^2R_{1|2}^\rho `$ (32)
$`\nu _{1|2}`$ $``$ $`\left(𝒟_{\mathrm{on}}^{s+3}K_{\mathrm{on}}^2B_{\mathrm{on}}^{s1}\right)^{1/(s+5)}R_{1|2}^{n_2}`$ (33)
$`S_{1|2}`$ $``$ $`\left(𝒟_{\mathrm{on}}^{3s+10}K_{\mathrm{on}}^5B_{\mathrm{on}}^{2s5}\right)^{1/(s+5)}R_{1|2}^{f_2}.`$ (34)
Exactly the same expressions apply also to the transition $`2|3`$, because both transitions are related to the synchrotron stage.555We would have obtained other proportionalities by using the expressions of the Compton stage for the first transition or those of the adiabatic stage for the second transition. However, the end result of Eqs (36) and (37) would have been the same, because $`\nu _{1|2}\nu _{\mathrm{m},1}(t_{1|2})=\nu _{\mathrm{m},2}(t_{1|2})`$, and similarly for $`S_{1|2}`$, $`\nu _{2|3}`$ and $`S_{2|3}`$.

If we now substitute the expressions of (29) for $`R_{1|2}`$ and $`R_{2|3}`$, we obtain the proportionalities which define the place in the $`(S,\nu ,t)`$-space (cf. Fig. 5) where the transitions from one stage to the other occur according to the values of $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$. These proportionalities can be expressed as logarithmic shifts $`\mathrm{\Delta }\mathrm{lg}P=\mathrm{lg}P\mathrm{lg}P`$ from an average value $`\mathrm{lg}P`$ of the parameter $`P`$, which stands either for $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$ or for $`S`$, $`\nu `$ and $`t`$. The set of equations relating these shifts for the first transition can be written as

$`\mathrm{\Delta }\mathrm{lg}t_{1|2}`$ $`=`$ $`U_{t_{1|2}}\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}+V_{t_{1|2}}\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}+W_{t_{1|2}}\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}`$ (35)
$`\mathrm{\Delta }\mathrm{lg}\nu _{1|2}`$ $`=`$ $`U_{\nu _{1|2}}\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}+V_{\nu _{1|2}}\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}+W_{\nu _{1|2}}\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}`$ (36)
$`\mathrm{\Delta }\mathrm{lg}S_{1|2}`$ $`=`$ $`U_{S_{1|2}}\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}+V_{S_{1|2}}\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}+W_{S_{1|2}}\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}.`$ (37)
Similar expressions can be written for the second transition by replacing $`1|2`$ in Eqs (35) to (37) by $`2|3`$. The $`2\times 9`$ expressions of the parameters $`U`$, $`V`$ and $`W`$ obtained for both transitions are displayed in Table 1. There is a great symmetry among these parameters, which reflects the fact that the two equations

$`\mathrm{\Delta }\mathrm{lg}\nu _{2|3}\mathrm{\Delta }\mathrm{lg}\nu _{1|2}`$ $`=`$ $`\beta _2(\mathrm{\Delta }\mathrm{lg}t_{2|3}\mathrm{\Delta }\mathrm{lg}t_{1|2})\text{and}`$ (38)
$`\mathrm{\Delta }\mathrm{lg}S_{2|3}\mathrm{\Delta }\mathrm{lg}S_{1|2}`$ $`=`$ $`\gamma _2(\mathrm{\Delta }\mathrm{lg}t_{2|3}\mathrm{\Delta }\mathrm{lg}t_{1|2})`$ (39)
must be verified, in accordance with (25).

Finally, we can also derive the influence of $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$ on the high-frequency spectral break $`\nu _{\mathrm{b},2|3}\nu _\mathrm{b}(t_{2|3})`$, which is given by $`\nu _{\mathrm{b},2|3}𝒟_{\mathrm{on}}B_{\mathrm{on}}R_{2|3}^{n_\mathrm{b}}`$ (cf. (22)) and corresponds to the logarithmic shift

$$\mathrm{\Delta }\mathrm{lg}\nu _{\mathrm{b},2|3}=U_{\nu _\mathrm{b}}\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}+V_{\nu _\mathrm{b}}\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}+W_{\nu _\mathrm{b}}\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}},$$
(40)
with the expressions of Table 1 for $`U_{\nu _\mathrm{b}}`$, $`V_{\nu _\mathrm{b}}`$ and $`W_{\nu _\mathrm{b}}`$.

## 4 Parameterization

In the previous section, we described in detail the physical basis of the light-curve decomposition that we present in Sect. 5. Here, we are more concerned by the practical aspect of the decomposition and particularly by the proper definition of the 85 parameters used by the fit.

### 4.1 The synchrotron spectrum and its average evolution

According to Sect. 3.1, the evolution of a shock wave in a jet follows three distinct stages to which correspond different evolutions of the self-absorption turnover of the synchrotron spectrum emitted by the plasma behind the shock front. The evolution of the turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ during each of the three stages is fully determined through (25) by only five indices, which are the three indices $`k`$, $`b`$ and $`d`$ defined by (14), the index $`s`$ of the electron energy distribution $`N(E)`$ and the index $`r`$ of the relation $`RL^r`$ defining how fast the jet opens with distance $`L`$. In the model presented in Sect. 5, we impose the value of $`d`$ to be zero and leave the four other indices as free parameters of the fit. They are constrained by the slopes $`\beta _i`$ and $`\gamma _i`$ $`(i=1,2,3)`$ of the outburst’s evolution in the $`(\mathrm{lg}S,\mathrm{lg}\nu ,\mathrm{lg}t)`$-space (cf. Fig. 5) and are the most interesting parameters of the light-curve decomposition, because they describe physical properties of the jet. The two points $`(t_{1|2},\nu _{1|2},S_{1|2})`$ and $`(t_{2|3},\nu _{2|3},S_{2|3})`$ in this space at which the transitions $`1|2`$ and $`2|3`$ from one stage to the other occur are defined by only four of these six quantities, because the slopes $`\beta _2`$ and $`\gamma _2`$ can be used to determine $`\nu _{2|3}`$ and $`S_{2|3}`$ from the four other values.

Until now, we have used eight parameters, namely $`s`$, $`r`$, $`k`$, $`b`$, $`\mathrm{lg}t_{1|2}`$, $`\mathrm{lg}t_{2|3}`$, $`\mathrm{lg}\nu _{1|2}`$ and $`\mathrm{lg}S_{1|2}`$, to fully characterize the path in the $`(S,\nu ,t)`$-space followed in average by the maximum $`(\nu _\mathrm{m},S_\mathrm{m})`$ of the self-absorbed synchrotron spectrum. The shape of this spectrum is defined by the general expression (cf. Paper I)

$$S_\nu =S_\mathrm{m}\left(\frac{\nu }{\nu _\mathrm{m}}\right)^{\alpha _{\mathrm{thick}}}\frac{1\mathrm{exp}(\tau _\mathrm{m}(\nu /\nu _\mathrm{m})^{\alpha _{\mathrm{thin}}\alpha _{\mathrm{thick}}})}{1\mathrm{exp}(\tau _\mathrm{m})},$$
(41)
where $`\tau _\mathrm{m}\frac{3}{2}\left(\sqrt{1\frac{8\alpha _{\mathrm{thin}}}{3\alpha _{\mathrm{thick}}}}1\right)`$ is a good approximation of the optical depth $`\tau _\nu `$ at the turnover frequency $`\nu _\mathrm{m}`$ and $`S_\mathrm{m}`$ is the real maximum of the spectrum.666Note that $`S_\mathrm{m}`$ as defined in Sect. 3.1 is not the real maximum, but the extrapolation $`S_\mathrm{m}^{\mathrm{thin}}S_\nu ^{\mathrm{thin}}(\nu _\mathrm{m})`$ down to $`\nu _\mathrm{m}`$ of the optically thin spectrum. The real maximum $`S_\mathrm{m}`$ is related to $`S_\mathrm{m}^{\mathrm{thin}}`$ by $`S_\mathrm{m}=S_\mathrm{m}^{\mathrm{thin}}(1\mathrm{exp}(\tau _\mathrm{m}))/\tau _\mathrm{m}`$, but this distinction is not important for this work, which is based on proportionalities.  We choose to impose the value of the optically thick spectral index $`\alpha _{\mathrm{thick}}`$ to $`f_3/n_3`$ according to the discussion in Sect. 3.1. The optically thin spectral index $`\alpha _{\mathrm{thin}}`$ is also determined by the model and has the value $`s/2`$ during the first two stages before flattening to $`(s1)/2`$ at the transition $`2|3`$ to the adiabatic stage, according to Eqs (5) to (7). This flattening by $`\mathrm{\Delta }\alpha _{\mathrm{thin}}=+1/2`$ is assumed to end at the time $`t_{2|3}`$ of the transition, but since it cannot be instantaneous we are forced to use one parameter to define the time $`t_\mathrm{f}`$ when it starts. Because this spectral change is expected to begin slowly before accelerating until the transition $`2|3`$, we choose to describe it with a logarithmic expression of the time $`t`$, as

$$\alpha _{\mathrm{thin}}(t)=\frac{s}{2}+\frac{1}{2}\frac{\mathrm{lg}(t/t_\mathrm{f})}{\mathrm{lg}(t_{2|3}/t_\mathrm{f})}\text{for}t_\mathrm{f}tt_{2|3}.$$
(42)
According to the considerations of Sect. 3.1, we further allow the spectrum defined by (41) to have a high- and a low-frequency spectral break as shown in Fig. 2. Since the exact shape of the breaks is difficult to assess (cf. Marscher M77 (1977); Band & Grindlay BG85 (1985)), we choose to keep them sharp, rather than smooth them arbitrarily, because sharp breaks have the advantage to define precisely the frequency at which they occur. The spectrum therefore simplifies to $`S_\mathrm{h}(\nu /\nu _\mathrm{h})^{2.5}`$ at the lowest frequencies $`(\nu <\nu _\mathrm{h})`$ and to $`S_\mathrm{b}(\nu /\nu _\mathrm{b})^{s/2}`$ at the highest frequencies $`(\nu >\nu _\mathrm{b})`$, where $`S_\mathrm{h}S_\nu (\nu _\mathrm{h})`$ and $`S_\mathrm{b}S_\nu (\nu _\mathrm{b})`$ are calculated by using (41). The frequency ratio $`\nu _\mathrm{h}/\nu _\mathrm{m}`$ of the low-frequency break to the spectral turnover is a free parameter of the fit and is assumed to remain constant throughout the outburst’s evolution. This ratio might actually slightly increase with time due to the increase of the thickness $`x`$ of the emitting region behind the shock front, but, for simplicity and because the effect might change from one stage to the other, we do not take this possible increase into account.

Since we do not extend the light-curve decomposition up to infrared frequencies, we do not need to consider a high-frequency spectral break during the two first stages of the outburst’s evolution, which would have the effect to further steepen the optically thin spectral index $`\alpha _{\mathrm{thin}}`$ to $`(s+1)/2`$. For frequencies above the break frequency $`\nu _\mathrm{b}`$, $`\alpha _{\mathrm{thin}}`$ remains therefore always at its value of $`s/2`$, even when it is flattening at lower frequencies. The evolution with time of the break frequency $`\nu _\mathrm{b}`$ is given by $`\nu _\mathrm{b}=\nu _{\mathrm{b},2|3}(t/t_{2|3})^{n_\mathrm{b}/\rho }`$, where $`\nu _{\mathrm{b},2|3}\nu _\mathrm{b}(t_{2|3})`$ is a free parameter of the fit and the ratio $`n_\mathrm{b}/\rho `$ is fixed by the values of $`r`$, $`b`$ and $`d`$ through Eqs (22) and (24).

### 4.2 Characteristics of different outbursts

The $`11`$ parameters defined above suffice to describe fully the spectral evolution of the average outburst shown in Fig. 5. The behaviour of the spectral turnover during each of the three stages of the evolution will be exactly the same for all individual outbursts. What can change from one outburst to the other is the position in the $`(\mathrm{lg}S,\mathrm{lg}\nu ,\mathrm{lg}t)`$-space where the two stage transitions occur. In Paper I, the two transition points were not allowed to move separately and therefore all outbursts were self-similar in the sense that they all had exactly the same evolution pattern in the $`(\mathrm{lg}S,\mathrm{lg}\nu ,\mathrm{lg}t)`$-space. Here, the distance between the two transition points can change and hence make the synchrotron stage relatively longer or shorter.

Although this new effect makes the light-curve modelling more complex, we do not introduce a further free parameter per outburst to describe it. Instead, we now give a physical origin to the three logarithmic shifts $`\mathrm{\Delta }\mathrm{lg}S`$, $`\mathrm{\Delta }\mathrm{lg}\nu `$ and $`\mathrm{\Delta }\mathrm{lg}t`$ in flux density $`S`$, frequency $`\nu `$ and time $`t`$ of Paper I, by replacing them by shifts of the values at the onset of the shock of the factor $`K`$ in the expression of the electron energy distribution $`N(E)`$, the magnetic field $`B`$ and the Doppler factor $`𝒟`$. The change of the three physical quantities $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$ results in different shifts in flux, frequency and time for the two transitions $`1|2`$ and $`2|3`$ according to Eqs (35) to (37). The shifts used in Paper I appear therefore to be quite unrealistic, because any such displacement in the $`(\mathrm{lg}S,\mathrm{lg}\nu ,\mathrm{lg}t)`$-space will also result in making the synchrotron stage longer or shorter.

Actually, the physical link between the shift of an outburst and its shape is able to solve some problems we had in Paper I (cf. Sect. 5), but does not improve the overall quality of the light-curve decomposition. To obtain a satisfactory fit, we have to add two outbursts, one in 1988 and one in 1991. Three other outbursts were added to extend the modelling of the light-curves to the most recent observations included here (cf. Sect. 2). We end up with seventeen outbursts, each characterized by four parameters, namely the date of its onset $`T_{\mathrm{on}}`$ and the three logarithmic changes $`\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}`$, $`\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}`$ and $`\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}`$ of the physical quantities $`K`$, $`B`$ and $`𝒟`$ at the onset of the shock. We therefore need a total of $`17\times 4=68`$ parameters to model the specificity of the seventeen individual outbursts.

### 4.3 Contributions from the outer jet and previous outbursts

The outer jet and especially its terminal hot spot called 3C 273A contributes significantly to the observed flux density at low radio frequencies, but this contribution can be considered as constant during the 20 years of this work. The spectrum of the hot spot 3C 273A is shown by Conway et al. (CGP93 (1993)) and corresponds roughly to a power-law with a spectral index $`\alpha `$ of $`0.85`$ above $`1`$ GHz, as described in Türler et al. (1999a ). As in Paper I, the fit of the light-curves by the series of outbursts is done above this constant power-law contribution.

Another important contribution at low frequencies comes from the superimposed decays of the outbursts peaking before 1979. In Paper I, we modelled this contribution by an exponential flux decay, which was appropriate because the decline of the individual outbursts was also modelled by an exponential decay in the light-curve approach, so that the characteristic $`\mathrm{e}`$-folding time was well constrained at all frequencies. It would be unsatisfactory to take here the same exponential decays, because the declining phase of a three-stage evolution implies that the final light-curve decay of each outburst is a power-law and not an exponential, as illustrated by the linear ending of the logarithmic light-curves of Fig. 5a.

The observed flux decay $`S_\nu ^{}(t)S_\nu (\nu ^{},t)`$ with time $`t`$ at a given frequency $`\nu ^{}`$ can be written as $`S_\nu ^{}(t)=S_\mathrm{m}(t)(\nu ^{}/\nu _\mathrm{m}(t))^{\alpha _{\mathrm{thin}}}`$, where the equality holds strictly if $`S_\mathrm{m}`$ is actually the extrapolation $`S_\mathrm{m}^{\mathrm{thin}}`$ of the optically thin spectrum (cf. footnote 6). During the final stage of the outburst’s evolution, the time dependence of this quantity can be expressed as $`S_\mathrm{m}(t)=S_\mathrm{m}^{}(t/t_\mathrm{m}^{})^{\gamma _3}`$ (cf. (25)), where $`t_\mathrm{m}^{}t_\mathrm{m}(\nu ^{})`$ is the time after the onset of the outburst when the spectral turnover $`(\nu _\mathrm{m},S_\mathrm{m})`$ passes at the frequency $`\nu ^{}`$, i.e. the time $`t`$ for which $`\nu _\mathrm{m}(t)=\nu ^{}`$. The corresponding equation for the turnover frequency is $`\nu _\mathrm{m}(t)=\nu _\mathrm{m}^{}(t/t_\mathrm{m}^{})^{\beta _3}`$. With these two equations, the light-curve decay at a given frequency $`\nu ^{}`$ during the adiabatic stage is given by

$$S_\nu ^{}(t)\left(\frac{t}{t_\mathrm{m}^{}}\right)^{\gamma _3\beta _3\alpha _{\mathrm{thin}}}\text{with}t_\mathrm{m}^{}=t_{2|3}\left(\frac{\nu ^{}}{\nu _{2|3}}\right)^{1/\beta _3}.$$
(43)
The quantity $`t_\mathrm{m}^{}`$ acts as a characteristic timescale of the power-law decay at a frequency $`\nu ^{}`$ and its frequency dependence describes how this timescale increases towards lower frequencies. In (43), $`\alpha _{\mathrm{thin}}`$ can be replaced by $`(s1)/2`$, if we do not consider a possible effect of the optically thin spectral break, which arises only at later times, once the break frequency $`\nu _\mathrm{b}`$ has reached the frequency $`\nu ^{}`$.

On this basis, we can model the contribution at a given frequency $`\nu `$ of the superimposed decays of the outbursts peaking before 1979 by the decay of an hypothetical outburst having an amplitude $`A_0(\nu )`$ at the date $`T_0=T_{\mathrm{on}}+t_\mathrm{m}(\nu )=1979.0`$ and then decaying as

$$S_\nu (T)=A_0(\nu )\left(1+\frac{T1979}{\mu t_\mathrm{m}(\nu )}\right)^{\gamma _3+\beta _3(s1)/2},$$
(44)
where we have used (43) with $`t`$ replaced by $`TT_{\mathrm{on}}=t_\mathrm{m}(\nu )+T1979`$ and we have introduced a factor $`\mu `$ ($`>1`$), which is frequency independent. There is no physical justification for introducing this factor $`\mu `$, but it was found to be a simple way to take into account the fact that the superimposition of several outbursts with different onset dates does decay with a longer characteristic timescale as a single outburst. For the whole range of frequencies covered by the thirteen light-curves considered in this work, $`A_0(\nu )`$ in (44) is the millimetre-to-radio spectrum of 3C 273 in 1979.0. This spectrum is modelled by a cubic spline parameterized at the four frequencies defined by $`\mathrm{lg}(\nu /\text{GHz})=0.3`$, $`0.5`$, $`1.0`$ and $`1.5`$ and is extrapolated to higher frequencies. We choose this low frequency range, because the influence of the spectrum $`A_0(\nu )`$ is the greatest in the radio domain and becomes negligible at submillimetre wavelengths (cf. Fig. 6). The value and the frequency dependence of the timescale $`t_\mathrm{m}(\nu )`$ is completely determined by the right part of (43) with $`\nu ^{}`$ replaced by $`\nu `$, so that the contribution from the outbursts prior to 1979 is only modelled by five free parameters, namely the values of $`A_0(\nu )`$ at the four above-mentioned frequencies and the frequency independent parameter $`\mu `$.

A final parameter of the fit is the GBI scaling factor mentioned in Sect. 2 and therefore the total number of free parameters used to adjust the 5234 observational points in the thirteen light-curves is $`11+68+6=85`$. We fit all light curves simultaneously, but usually with not more than ten of the $`85`$ parameters at a time. The four indices $`s`$, $`r`$, $`k`$ and $`b`$ are adjusted together, followed by the seven other parameters defining the spectral evolution of the average outburst. The specificity of individual outbursts is adjusted by a series of fits concerning only two or three outbursts at a time. Finally, we fit together the six remaining parameters, which define the initial flux decay and the GBI scaling factor. This procedure is repeated many times until the values of the parameters converge to their best-fit values. In total, several hundreds iterative fits are needed to obtain the light-curve decomposition of Fig. 3.

## 5 Results

We present here the results obtained by fitting the submillimetre-to-radio light-curves of 3C 273 according to the jet model described in Sect. 3 and with the set of parameters defined in Sect. 4. We choose to focus in this section on the results obtained by assuming that the synchrotron emitting plasma behind the shock front moves with constant speed and direction, so that the Doppler factor $`𝒟`$ remains constant during the outburst’s evolution. The index $`d`$ characterizing the decrease of $`𝒟`$ with jet radius $`R`$ (cf. (14)) was therefore fixed at a value of zero. All following figures and tables correspond to this particular jet model. In Sect. 6.1, we will however also discuss the results of an alternative jet model, in which the Doppler factor is free to vary, but we impose the jet to be conical.

The best-fit decomposition we could achieve with the jet model having a constant Doppler factor $`𝒟`$ is illustrated in Fig. 3. We note that the light-curves which are the best reproduced by the model are the radio light-curves, and especially those at 15 GHz and at 8 GHz. In the 22 GHz light-curve, the rapid variations from 1981.0 to 1983.5 are not well described by the model, and at 37 GHz, there is mainly a problem with the prominent peak of 1991.5. At 90 GHz, it is striking to note that the model reproduces quite well the features of the observed light-curve, but cannot reach a sufficiently high flux level. The same problem is also apparent at a wavelength of 2 mm and it is only at 1.1 mm that the model light-curve is again on average at the same flux level than the observations. Finally, at the highest frequencies, the model light-curve does not decrease enough between the outbursts to adjust the lowest flux observations.

This general discrepancy at millimetre and submillimetre frequencies is probably responsible for the relatively high reduced $`\chi ^2`$ value of $`\chi _{\mathrm{red}}^2\chi ^2/\text{d.o.f.}=20.5`$ obtained for the overall fit as compared to the values quoted in Paper I. Apart from this global problem, which will be discussed in Sect. 6.1, we note that the detailed structure of the light-curves is now better reproduced by the set of outbursts than it was in Paper I. In particular, the fast decay around 1985 is now more satisfactorily accounted for by the model. The rapid variations of the flare peaking near 1996.0 are also better described by the new parameterization. If we compare the high-frequency light-curves of Fig. 4 with those shown in Fig. 7b of Paper I, we see that the improved description of this outburst starting in 1995.8 is achieved with a faster decline just after the peak of the high-frequency light-curves. This faster decline results from the introduction of the new parameter $`t_\mathrm{f}`$, which defines the time when the optically thin spectral index begins to flatten and from the fact that this time is delayed with respect to the time $`t_{1|2}`$ of the first stage transition. This flattening, occurring now just at the end of the synchrotron stage, is in very good agreement with the shock model of MG85. The above-mentioned better description of the decay around 1985 is not due to the same effect. A satisfactory fit of 3C 273’s behaviour at this epoch could indeed only be achieved once we took into account the possibility that the stage transition points $`1|2`$ and $`2|3`$ could be shifted separately as it is the case here (cf. Sect. 4.2). The appropriate shape of the last three outbursts peaking before 1985 was obtained by reducing very much their intermediate synchrotron stage, as shown by the values in the last column of Table 3.

The best-fit values of the adjusted parameters and the corresponding values of other interesting quantities are displayed in Table 2. These values define the three-dimensional evolution of the average outburst in 3C 273, as shown in Fig. 5. We note that the decline of the turnover flux $`S_\mathrm{m}`$ with turnover frequency $`\nu _\mathrm{m}`$ is now flatter than obtained with the three-stage approach of Paper I. Instead of $`\gamma _3/\beta _3=1.14`$, we now find $`\gamma _3/\beta _3=0.88`$. This difference could be due to the introduction of the high-frequency break in the spectrum (cf. Fig. 2). Another difference is the steeper rise of $`S_\mathrm{m}`$ with $`\nu _\mathrm{m}`$ by nearly a factor 2 ($`\gamma _1/\beta _1`$ is now of $`1.83`$ instead of $`0.99`$). This might be due to the introduction of the low-frequency spectral break, but it must be noticed that this initial part of the evolution is the less constrained stage of the outburst’s evolution. Its behaviour could therefore be imposed by the best-fit values of $`s`$, $`r`$, $`k`$ and $`b`$, which are more strongly constrained by the two other stages of the evolution. Nevertheless, both the steeper rise and the flatter decline obtained here are in better agreement with the original outburst’s evolution proposed by MG85.

The value of $`2.05`$ that we obtain for the index $`s`$ of the electron energy distribution $`N(E)`$ corresponds to an optically thin spectral index of $`\alpha _{\mathrm{thin}}=(s1)/20.53`$ for the final adiabatic stage of the outburst’s evolution. These values for $`s`$ and $`\alpha _{\mathrm{thin}}`$ are slightly higher than those obtained in Paper I, but are still well below the values of $`s=2.4`$ and $`\alpha _{\mathrm{thin}}=0.7`$ measured by MG85 for the strong flare of 1983. The parameter $`K`$ of the electron energy distribution $`N(E)`$ is found to evolve with the width of the jet $`R`$ approximately as $`KR^{3.0}`$, which is just a bit steeper than the decrease $`KR^{2(s+2)/3}R^{2.7}`$ expected if the jet flow was adiabatic (e.g. Gear G88 (1988)). The decrease of the magnetic field $`B`$ with $`R`$ is found to be between the two extreme cases of $`BR^2`$ expected if $`B`$ was parallel to the jet axis and $`BR^1`$ expected if $`B`$ was transverse to the jet (e.g. Begelman et al. BBR84 (1984)). The obtained index of $`b=1.58`$ suggests that the parallel $`B_{}`$ and the perpendicular $`B_{}`$ components of the magnetic field are roughly equal. Finally, the best-fit value of $`r=0.82`$ obtained with this model having a constant Doppler factor $`𝒟`$ suggests that the inner jet of 3C 273 is not really conical, but tends to open slightly less with distance $`L`$ along the jet axis (cf. Fig. 1).

The individual characteristics of the outbursts are presented in Table 3. They are expressed by the logarithmic shifts $`\mathrm{\Delta }\mathrm{lg}P=\mathrm{lg}P\mathrm{lg}P`$ from an unknown average value $`\mathrm{lg}P`$ of the physical quantities $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$ at the onset of the shock. The dispersions $`\sigma `$ of the seventeen values in the columns of Table 3 are $`\sigma 0.40`$ for $`\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}`$; $`\sigma 0.29`$ for $`\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}`$; and $`\sigma 0.09`$ for $`\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}`$. We note that the Doppler factor $`𝒟`$ is the quantity with the smallest relative changes from one outburst to the other. To study the possibility that the logarithmic shifts of the three physical quantities $`K_{\mathrm{on}}`$, $`B_{\mathrm{on}}`$ and $`𝒟_{\mathrm{on}}`$ are correlated, we applied a Spearman rank-order test (e.g. Press et al. PTV92 (1992)) to these three data sets. We find out that there is a very significant anti-correlation between the values of $`\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}`$ and those of $`\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}`$. The probability that a stronger anti-correlation could occur by chance is less than $`\mathrm{6\hspace{0.17em}10}^5`$. A weaker trend, with a probability of non-correlation of less than 2 %, suggests that $`\mathrm{\Delta }\mathrm{lg}𝒟_{\mathrm{on}}`$ is positively correlated with $`\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}`$. On the other hand, no significant correlation is found between $`\mathrm{\Delta }\mathrm{lg}K_{\mathrm{on}}`$ and $`\mathrm{\Delta }\mathrm{lg}B_{\mathrm{on}}`$. These relations are discussed in Sect. 6.2.

## 6 Discussion

In Paper I, we discussed the relationship between the onset of an outburst and the ejection of a new VLBI component in the jet and we concluded that the identified outbursts do fairly well correspond to the VLBI knots. The model presented here does not change this conclusion and therefore we postpone a further discussion of this relationship until more recent VLBI measurements of 3C 273 are published. We focus, in Sect. 6.1, on the global properties of the inner jet, before discussing, in Sect. 6.2, the peculiarities of individual outbursts.

### 6.1 On the global properties of the inner jet

We have shown in Fig. 3 the decomposition of a series of light-curves. It is, of course, also possible to show the corresponding decomposition of the observed spectrum at different epochs. As an example, we give in figure Fig. 6 the spectral dissection of the two most different submillimetre-to-radio spectra observed in 3C 273 during the last 20 years. Until now, such a spectral dissection could only be achieved with VLBI observations and the best example of this is probably the decomposition by Marscher (M88 (1988)) of the radio spectrum of the relatively high redshift quasar NRAO 140. Actually, the great similarity between the dissection of NRAO 140, based on a single epoch VLBI image, and the dissection of Fig. 6, based on multi-wavelength variability observations, gives strong support to the idea that VLBI knots, observed in the jet, are physically linked to the outbursts seen in the light-curves.

The kind of spectral shape shown in Fig. 6a, with a relative maximum of the spectrum in the millimetre range, was observed in several blazars by Brown et al. (BRG89 (1989)) and interpreted as the evidence of two synchrotron components contributing to the spectrum. The two decompositions of Fig. 6, suggest that these two components are not of distinct nature, but that the flaring component progressively becomes part of the more slowly varying component as it evolves towards lower frequencies. It means that, apart from the contribution of the outer jet, the quiescent radio emission of 3C 273 can be entirely attributed to the superimposed decays of earlier outbursts.

We have assumed in this work that the emission of the inner jet is negligible outside of a small region containing the shocked plasma just behind the shock front (cf. Fig. 1). It means that we have implicitly made the hypothesis that the emission of the unresolved VLBI core is due to the emission of the one or two most recent outbursts observed during the earlier phases of their evolution. This is a rather extreme assumption, since it is usually believed that the base of the jet is also contributing to the emission, because of the higher electron and magnetic field density at the narrow end of the jet. The present work shows that the main features of submillimetre-to-radio emission of 3C 273 can be understood as being only due to shock waves. There are, however, some indications that there is still some place left for an underlying contribution of the inner jet. The presence of an inner jet contribution is suggested by the relatively poor fit at millimetre and submillimetre wavelengths and especially by the fact that it is not the shape of the variations which is the greatest problem, but much more the average flux level reached by the model, which is too low at some frequencies and too high at others (cf. Sect. 5). The idea behind this is that the model seems to make a compromise between the opposite requests of, on the one hand, the millimetre observations asking to flatten the optically thin spectral index to have more flux during the slow decline of the outbursts and, on the other hand, the submillimetre observations asking to steepen this index, in order to decrease the flux level at the end of the outbursts.

It is possible that this problem could be solved by adding an underlying nearly constant spectrum with a maximum around a frequency of 90 GHz, in association with a steeper optically thin spectral index for the individual outbursts. This would imply a higher value for the index $`s`$ of the electron energy distribution, in better agreement with the value of $`2.4`$ inferred from the measure of the millimetre-to-infrared slope of the flaring spectrum during the initial stage of the 1983 outburst (MG85). Such an underlying contribution from the inner jet to the millimetre and submillimetre emission of 3C 273 could possibly be related to the so called $``$ component identified in the blue-bump of 3C 273 by Paltani et al. (PCW98 (1998)). Indeed, this slowly varying component, which dominates the optical spectrum, can be interpreted as synchrotron emission and might reveal the blazar-like characteristics of 3C 273.

As mentioned in Sect. 5, the obtained value of $`1.6`$ for the index $`b`$ suggests that the parallel $`B_{}`$ and the perpendicular $`B_{}`$ components of the magnetic field $`B`$ are roughly equal. The simplest interpretation of this result is that the magnetic field $`B`$ is turbulent behind the shock front, leading to a re-isotropization of the field. Another possibility, but which requires much more fine tuning, is that $`B_{}`$ is dominant in the underlying jet and that the shock amplifies just enough the $`B_{}`$ component (see MG85) to end with $`B_{}B_{}`$. The obtained value of $`r0.8`$ suggests that the jet opens less with distance $`L`$ along the jet axis than if it was conical. Such a behaviour was recently observed in M 87 (Virgo A, 3C 274) on an image of its very inner jet (Junor et al. JBL99 (1999)). A non-conical jet can arise due to an accelerating jet flow (Marscher M80 (1980)), but this interpretation would be in contradiction with the constant Doppler factor $`𝒟`$ assumed in this model. Another possibility is that there is an external pressure which confines the jet, so that the jet does not expand freely in two dimensions. Such a collimation could be of magnetic origin, either in the case of a magnetically self-confined jet (cf. Begelman et al. BBR84 (1984)) or in the case of the two-flow model (Sol et al. SPA89 (1989)), which assumes that a highly relativistic electron-positron beam is confined within the magneto-hydrodynamic structure of a mildly relativistic electron-proton jet.

The discussion above concerns the results obtained with a shock model in which we assume a constant Doppler factor $`𝒟`$ by imposing a zero value for the index $`d`$. We obtain a comparable $`\chi ^2`$ of the fit using an alternative model in which we assume that the jet is conical and thus that the index $`r`$ equals 1. For this model, we obtained very similar values for most parameters given in Table 2, including the indices $`s`$ and $`b`$, which therefore seem to be strongly constrained by the light-curve decomposition. What changes with this conical jet model are the values of the indices $`k`$ and $`d`$. The value of $`k`$ decreases to $`2.72`$, which is now almost the value of $`k_{\mathrm{ad}}=2(s+2)/32.70`$ corresponding to an adiabatic jet flow and the best-fit value for $`d`$ is $`0.10`$, which suggests either that the emitting region is slightly decelerating while it travels down the jet or that the jet is bending away from the line-of-sight. Even if the shock front is not decelerating, a slight deceleration of the centre of the emitting region is a natural consequence of the continuous increase of its thickness $`x`$ during the whole shock evolution. Therefore, even for a straight jet, this positive value of $`d`$ might not be in contradiction with a value of $`k`$ typical for an adiabatic jet flow. On the other hand, it is well established that the path projected on the sky followed by the VLBI components in the parsec-scale jet of 3C 273 is curved (e.g. Abraham et al. ACZ96 (1996)). It is therefore also possible that the slight decrease of the Doppler factor $`𝒟`$ suggested by the conical jet model is related to this overall bending of the jet.

The fact that both the model with a constant Doppler factor $`𝒟`$ and the conical jet model give an equally good fit suggests that the light-curve decomposition cannot uniquely determine the values of the three parameters $`r`$, $`k`$ and $`d`$. Indeed, by leaving all these parameters free to vary, we could not achieve a better fit than obtained with one of them fixed. It means that besides the two models discussed above, jet models with other combinations of the values of $`r`$, $`k`$ and $`d`$ cannot be excluded by our light-curve decomposition.

### 6.2 On the peculiarities of individual outbursts

Concerning the specificity of individual outbursts, we noted in Paper I that short- and long-lived outbursts in 3C 273 are usually not peaking at the same frequency. Long-lived outbursts were found to peak at lower frequencies and we proposed that this relationship might be related to the distance down the jet at which the shock forms. The proper interpretation of the peculiarities of individual outbursts is now complicated by the fact that the outbursts are not self-similar anymore, because of the different behaviour of the two stage transitions $`1|2`$ and $`2|3`$ when the physical quantities $`K`$, $`B`$ and $`𝒟`$ change (cf. Sect. 3.2). For instance, relatively high values of $`K_{\mathrm{on}}`$ and low values of $`B_{\mathrm{on}}`$ result in a shorter synchrotron stage. This effect is particularly pronounced in the three successive outbursts between 1982 and 1984 (cf. Table 3). The apparent absence of a flat peaking stage in 3C 345 (Stevens et al. SLR96 (1996)) could therefore be due to a relatively low magnetic field $`B`$ as compared to the factor $`K`$ of the electron energy distribution.

We note in Sect. 5 that there is apparently no correlation from one outburst to the other between the two quantities $`K_{\mathrm{on}}`$ and $`B_{\mathrm{on}}`$ at the onset of the shock. This result suggests that differences from one outburst to the other are not primarily due to changes of the compression ratio $`\eta `$ of the shock, because in this case we would expect correlated $`K_{\mathrm{on}}`$ and $`B_{\mathrm{on}}`$ variations due to their similar dependence on $`\eta `$ (cf. MG85). Changes of $`K_{\mathrm{on}}`$ which are independent of $`B_{\mathrm{on}}`$ could be possible if the acceleration process of the electrons crossing the shock front was not adiabatic. In this case, the energy gain $`\xi `$ of the electrons would not simply be related to the compression ratio as $`\xi =\eta ^{1/3}`$ (MG85), but might have different behaviours from one outburst to the other. Since the normalization $`K_{\mathrm{on}}`$ of the electron energy distribution is proportional to $`\xi ^{s1}`$ (MG85), changes of the electron energy gain $`\xi `$ for different shocks in 3C 273 could lead to the observed changes of $`K_{\mathrm{on}}`$. This interpretation is supported by the significant anti-correlation found between $`K_{\mathrm{on}}`$ and the Doppler factor $`𝒟_{\mathrm{on}}`$, because a non-adiabatic shock wave converts bulk kinetic energy into internal energy (e.g. Begelman et al. BBR84 (1984)). Outbursts with a small value of $`𝒟`$ and a great value of $`K`$, would therefore be associated with more efficient shocks, for which a greater fraction of their bulk kinetic energy, as measured by $`𝒟`$, would be converted into internal energy by increasing the average energy gain $`\xi `$ of the electrons.

## 7 Summary and conclusion

This work presents a phenomenological model describing the evolution of synchrotron outbursts expected to be emitted by shock waves in a relativistic jet. This model, which is a generalization of the original shock model of MG85, is then confronted to the very well sampled long-term submillimetre-to-radio light-curves of 3C 273. Many iterative fits are needed to adjust the more than 5000 observations with a series of seventeen successive outbursts defined by 85 parameters of the model. This fitting procedure allows us not only to define the average outburst’s evolution in 3C 273, but also to study the peculiarities of each individual outburst during the last 20 years. The main results are the following:

* The quiescent low frequency emission can be understood as being entirely due to the superimposition of slowly decaying outbursts, which flared a few years before, and a constant contribution from the hot spot (3C 273A) of the outer jet.
* The values of the indices describing how the physical quantities characterizing the jet evolve with its opening are found to be in general agreement with the simple jet model considered by MG85.
* The best-fit values of these indices suggest however either that the jet opens slightly less than a conical jet if we assume that the emitting region has a constant bulk Doppler factor or, alternatively, that this Doppler factor decreases if we impose the jet to be conical.
* In both cases, the magnetic field in the emitting region behind the shock front seems to be rather turbulent and the jet flow is found to be nearly adiabatic, especially for the conical jet model.
* The peculiarities of individual outbursts can be understood as being due to changes at the onset of the shock of the magnetic field strength $`B`$, the normalization $`K`$ of the electron energy distribution and the Doppler factor $`𝒟`$.
* Shocks with a high value of $`K`$ have usually a lower value of $`𝒟`$. This anti-correlation might be related to the shock efficiency to convert bulk kinetic energy into internal energy of the plasma behind the shock front.

In the past, shock models were found to be difficult to test and constrain with total flux measurements. We show here that very strong observational constraints can be derived from long-term multi-wavelength monitoring campaigns. This new ability, together with interferometric imaging techniques and numerical simulations of shocks in relativistic jets (e.g. Gómez et al. GMM97 (1997)) should lead in a near future to important progresses in the understanding of the physics involved in relativistic jets.

###### Acknowledgements.

We thank H. Teräsranta for providing us previously unpublished measurements of 3C 273 at 22 and 37 GHz from the Metsähovi Radio Observatory. This research has made use of data from the University of Michigan Radio Astronomy Observatory which is supported by funds from the University of Michigan.