text:
# Discovery of the optical counterpart and early optical observations of GRB9907121footnote 11footnote 1Based on observations collected at SAAO, Sutherland; ESO, Paranal and La Silla (ESO Programs 63.O-0618 and 63-O-0567); and AAT, Australia.

## 1. Introduction

The past two years have witnessed tremendous progress in our understanding of Gamma-Ray Burst (GRB) phenomena. Thanks to the ability of Beppo-Sax to provide arcminute sized error boxes for the bursts, the first optical counterpart was detected in February 1997 for GRB970228 (van Paradijs et al. 1997; Costa et al. 1997). The HST observations of this GRB provided the first clear indication that the GRB is associated with an external galaxy and is unrelated to its nuclear activity (Sahu et al. 1997a). Shortly thereafter, the redshift measurement for GRB 970508 proved beyond doubt that GRBs are extragalactic in nature (Metzger et al. 1997; Djorgovski et al. 1997). More than a dozen GRB optical counterparts have been detected since then, with redshifts as high as 3.42 in the case of GRB971214 (Kulkarni et al. 1998).

The extensive observations of GRB afterglows in X-ray, optical and radio wavelengths have been shown to be consistent with cosmological fire-ball models (e.g. Paczyński & Rhoads 1993; Mészáros & Rees 1997; Wijers et al. 1997). However, the exact cause of the GRBs has remained elusive, and progress in determining the actual mechanism causing the bursts has been slow. The two leading models for GRBs involve the collapse of a massive star (e.g. Woosley 1997; Paczyński 1998), and the merging of a neutron star with either another neutron star or a black hole (e.g. Eichler et al. 1989; Narayan, Paczyński, & Piran 1992; Sahu et al. 1997b). The “isotropic equivalent energy” of some GRBs is as high as 1054 ergs, which exceeds the energy equivalent of the total mass involved in the latter model, and hence is thought to favor the former. However, the beaming factor of the emission of the GRB in different wavebands can be high (in both models), so the energetics alone may not be conclusive in favoring one of the models over the other. Since the lifetime of a very massive star is of the order of a million years or shorter, the GRBs are expected to be within or close to star-forming regions in the massive-star collapse model. On the other hand, since the kick-velocity of a newly formed neutron star is of the order of 200 $`\text{km\hspace{0.17em}s}^1`$, and since the neutron star, on average, is already about 108 years old at the time of the burst, GRBs are generally expected to be far from star-forming regions in the merging neutron-star model (Bloom, Sigurdsson & Pols 1999). Thus the location of the GRB with respect to the star-forming regions in the host galaxy could be a distinguishing feature which can help in settling the question of the cause of the GRBs. However, if the presence of dense interstellar material is prerequisite for the onset of the optical afterglow, then OTs would be observed only in star-forming regions regardless of whether GRBs result from neutron-star mergers or from the collapse of massive stars. In any case, if some OTs are found far away from star-forming regions, neutron-star merger model would be favored. If all OTs are found in star-forming regions, the situation is less clear, and depending on the model, one may need to investigate other aspects to find distinguishing features. For example, one consequence of the massive-star collapse model is that the GRB should be accompanied by an underlying supernova, whose brightness variation is distinct from the power-law decay behavior of the OT and may be detectable.

There is a clear bimodality in the observed burst durations: long bursts have timescales of about 10 to 200 seconds, and short bursts have timescales of 0.1 to 1 second (Kouveliotou et al. 1995). All the afterglows that have been discovered so far belong to the subgroup with long bursts (Fishman 1999). All the OTs for which the host galaxies have been observed – although this sample is small – are found to be in galaxies whose spectral indicators suggest star-forming activity (e.g. Fruchter et al. 1999). So there is an interesting possibility that the long-duration GRBs are caused by the collapse of rapidly-rotating massive stars, and hence they occur in star-forming regions. The short-duration bursts may be of a different origin, and may be caused by the coalescence of two neutron stars. The detection of optical counterparts of short duration bursts, and their observations would be important to understand this scenario better.

GRB990712, being at $`z=0.434`$ (Galama et al. 1999a; Hjorth et al. 2000), is one of the closest GRBs observed so far and hence provides a good opportunity to determine the possible presence of an underlying supernova, the location of the GRB within the galaxy, and the luminosity of the galaxy in different wavebands. In this paper we present the optical observations of GRB990712 leading to the discovery of its optical counterpart, and subsequent optical imaging observations taken in different filters during the first $``$35 days after the burst.

## 2. Observations

GRB990712 was simultaneously detected on 1999 July 12.69655 with the Gamma-Ray Burst Monitor and the WFC unit 2 onboard the BeppoSAX satellite. The burst lasted for about 30 seconds in both the $`\gamma `$-ray (40 - 700 keV) and X-ray (2-26 keV) energy ranges and had a double-peaked structure. While its intensity in $`\gamma `$-rays is moderate, it exhibited the strongest X-ray prompt emission observed to date (Heise et al. 1999).

The initial Beppo-SAX position, with an error circle of 5 arcmin radius, was revised to one with an error radius of only 2 arcmin. Our first observations were made at this latter position with the SAAO 1m telescope at Sutherland, South Africa on July 12.87 UT, about 4.16 hours after the burst during time generally dedicated to PLANET microlensing observations (Albrow et al. 1998). The optical image was taken through an R filter, with an integration time of 900 sec. The detector was a $`1024\times 1024`$ SITe CCD with an image scale of 0.309 arcsec per pixel and a total field of view of $`5.3\times 5.3`$ square arcmin. Comparison of the image with the Digitized Sky Survey (DSS) image for the same field showed the presence of a new source well above the sky background, that was absent in the DSS image. This new source had a brightness of R = 19.4 $`\pm `$ 0.1 mag, which is about 2 magnitudes brighter than the limiting magnitude of the DSS image, indicating that the new source was the optical counterpart of the GRB (Bakos et al. 1999a,b). Fig. 1 shows the discovery image, along with the DSS image taken of the same region. Also marked in the figure are the reference stars used for photometric calibration of the OT. A spectrum of this source obtained a few hours later revealed emission and absorption lines which were used to derive redshift of 0.434 for this source (Galama et al. 1999a). Subsequent observations by various groups showed the decaying nature of this source further confirming the identification of the OT with the GRB.

We continued the observations with a multi-wavelength optical follow-up campaign using the telescopes at SAAO, ESO and AAO. A log of all the observations used in this analysis along with the derived magnitudes of the OT are listed in Table 1. The measurements include the acquisition images taken for the spectroscopic observations with the ESO VLT shortly after the discovery of the OT (Vreeswijk et al. 1999, in prep.).

### 2.1. Astrometry and Photometry of the Optical Counterpart

An astrometric solution of the field was carried out using 10 reference stars in the image which are marked in Fig. 1. Their coordinates, as taken from the USNO 2.0 catalog which uses the Hipparcos frame of reference, are listed in Table 2. The pixel centroids of the reference stars were determined using two-dimensional Gaussian fits. These centroids were combined with the USNO coordinates to determine an astrometric solution of the field. The resulting position of the optical counterpart is RA(2000) = 22h 31m 53.$`{}_{}{}^{s}061\pm 0^s.011`$, Dec(2000)= $`73^d24^{}28.^{\prime \prime }58\pm \mathrm{\hspace{0.17em}0}^{\prime \prime }.05`$.

Since sky conditions deteriorated due to clouds just after the discovery observations of the OT, no standard star observations could be taken on that night. The images were calibrated through observations of the standard stars F203 and F209 (Menzies et al. 1989) taken on 13th July at SAAO. All the photometric measurements were carried out using the IRAF aperture photometry task PHOT. The photometry of a few secondary standards in the field (Fig. 1) was carried out using the standard star observations, and the photometry of all GRB observations was then performed via these secondary standards. Since only two standard stars were observed, the extinction correction due to differential airmass between the GRB and the standard star observations were not applied in the preliminary analysis (Bakos et al. 1999a,b). Instead, the photometric measurements taken with the adjacent 50cm telescope were used to get the extinction coefficients, which resulted in a correction of –0.075 mag in V, –0.06 mag in R, and –0.02 mag in I. (This is consistent with the discrepancy pointed out by Kemp & Halpern 1999). Our final magnitudes of the secondary photometric standards are listed in Table 2, and Table 1 lists the derived magnitudes of the OT in different bands at different epochs including one extra measurement reported by Kemp and Halpern (1999).

## 3. Analysis

### 3.1. The Light Curve in Different Wavebands

All the photometric measurements listed in Table 1 are also shown in Fig. 2. The observed light curves in different bands clearly show a continuously-changing slope suggesting a additional component to the power-law decay of the OT as first noted by Hjorth et al. (1999a). In an approach similar to the one followed by Hjorth et al. (1999b), we have fitted this additional component in two ways: (i) with an underlying galaxy of constant brightness, and (ii) an underlying galaxy plus a supernova similar to GRB980425, appropriately scaled to the redshift of GRB990712.

First, a power law decline of the form $`f(t)t^\alpha `$ for the OT and a constant contribution from the background galaxy was used to fit the observed light curves simultaneously in different bands, with the same slope for all the bands. The resultant slope is $`0.97`$, and the total $`\chi ^2`$ is 47 for 39 d.o.f. (Relaxing the condition of the same slope in different bands does not alter the fit parameters significantly). The total $`\chi ^2`$ is, however, fairly flat (between 44 and 52) for any value of $`\alpha `$ between 0.95 and 1.02 beyond which the total $`\chi ^2`$ increases rapidly. Fig. 2 shows the best fit to the light curves. The derived magnitude of the underlying galaxy is not very sensitive to the exact value of $`\alpha `$ and are found to be $`\mathrm{V}=22.3\pm 0.05`$, $`\mathrm{R}=21.75\pm 0.05`$ and $`\mathrm{I}=21.35\pm 0.05`$ where the uncertainties represent the variation of the derived magnitudes as a result of changing $`\alpha `$ between 0.95 and 1.02.

One possible consequence of the massive-star model for the GRB is that one should observe an underlying supernova in the lightcurve of an OT. The physical connection of a GRB with a supernova (SN), first suggested by the discovery of a peculiar type Ic SN in the error box of GRB 980425 (Galama et al. 1998; Iwamoto et al. 1998), has been strengthened further by recent observations of other GRBs. Castro-Tirado and Gorosabel (1999) suggested that the light curve of GRB980326 resembled that of a SN, and indeed Bloom et al. (1999) showed that the late time light curve of GRB 980326 can be explained by an underlying SN1998bw type SN at a redshift of around unity. For the afterglow of GRB 970228, Reichart (1999) and Galama et al. (1999b) find that a power-law decay plus SN1998bw light curve redshifted to the distance of the burst, $`z=0.695`$ (Djorgovski et al. 1999), fits the observed light curve very well. In light of these new findings, we have also fitted the observed light curve of GRB 990712 assuming the presence of an underlying supernova.

In order to determine the contribution of the underlying supernova, we first calculate the expected V, R, I magnitudes of SN1998bw, placing it at the redshift of GRB 990712, $`z=0.434`$. This includes wavelength shifting and time profile stretching (both by a factor of 1+z), and rescaling the magnitudes for a distance corresponding to $`z=0.434`$, assuming $`\mathrm{\Omega }_0`$ = 0.2. To account for the wavelength shift correctly, we interpolate the redshifted UBVRI broadband flux spectrum of the SN for each bin in time with a spline fit, obtain the V, R, and I fluxes at their effective wavelengths, and convert these back to obtain magnitudes in the observer’s frame (Fukugita et al. 1995). The resultant SN light curves in different bands are shown at the bottom of Fig. 3. In order to fit the GRB light curve, the SN flux is subtracted from the observed GRB flux, and the residual fluxes are assumed to be due to the OT and the underlying galaxy. The procedure outlined earlier is then used to fit the light curves. The $`\chi ^2`$ in this case is 52 for 42 d.o.f., which is clearly higher than the model without the supernova. The resulting slope of the decay is $`0.96`$ which is quite similar but the magnitudes of the underlying galaxy are fainter: $`\mathrm{V}=22.7\pm 0.05`$, $`\mathrm{R}=22.25\pm 0.05`$ and $`\mathrm{I}=22.15\pm 0.05`$. The brightness of the OT of the GRB at a given time is also slightly lower, but the combined brightness of the OT and the underlying supernova (which is the quantity that can be measured if the OT can be resolved from the galaxy) is higher than the brightness of the OT in the absence of an underlying supernova.

It is important to note, however, that there is no evidence that the SNe possibly underlying the GRB afterglows have the same brightness or the same decay behavior. Our assumption that the underlying supernova is identical to SN1998bw is dictated by two reasons. First, this is the only SN associated with a GRB whose light curve has been monitored extensively. Second, the number of data points in our light curves does not allow us to vary the characteristics of the underlying SN light curve. So the fact that the $`\chi ^2`$ is higher in this case can be misleading since the true underlying supernova may be different from that of SN1998bw, making the resultant $`\chi ^2`$ higher for our simplified model. This chi-square analysis indicates that that either GRB990712 is not associated with an underlying supernova or the underlying supernova is much fainter than the supernova associated with GRB980425. However, the qualitative conclusion that the galaxy is expected to be fainter and the OT plus the SN is expected to be brighter in the presence of an underlying supernova is unlikely to change. Hence late time HST observations, in which the OT is well resolved so that the brightness of the OT and the galaxy can be estimated separately, or late time ground-based observations when the brightness of the OT is negligible, would greatly help in determining the presence of an underlying supernova.

In both the above scenarios (i.e. with and without an underlying supernova), the power-law index of the decay $`\alpha `$ is about $`1`$, making it one of the slowest decline rates of all the OTs observed so far. Since $`\alpha `$1 would lead to a divergence of the total energy integrated over time, the slope must steepen at later times. This has been observed for GRB990510 which declined with $`\alpha =0.82`$ at early times, later steepening with time (Harrison et al. 1999). Power-law decays are thought to arise from electrons shocked by the relativistic expansion of the debris into the ambient medium. In such a case, the information on the change of slope can be used to derive information on the cooling rate of the electrons in the post-shock region (see, e.g. Wijers et al. 1997; Sari, Piran, & Narayan 1998; Livio & Waxman 1999).

### 3.2. Spectral Energy Distribution of the OT

If the OT emission is due to the synchrotron radiation from the swept-up electrons in the post-shock region, then its energy distribution is expected to be of the general form $`f_\nu t^\alpha \nu ^\beta `$. In such a case, one expects a relationship between the power law index of the energy distribution of the electrons $`p`$, the spectral slope $`\beta `$, and the decay constant $`\alpha `$ (Sari et al. 1998). One must distinguish two cases: (i) both the peak frequency $`\nu _\mathrm{m}`$ and the cooling frequency $`\nu _\mathrm{c}`$ are below the optical/IR waveband, in which case $`p=(4\alpha +2)/3=1.97\pm 0.04`$ and $`\beta =p/2=0.99\pm 0.02`$, and (ii) $`\nu _\mathrm{m}`$ has passed the optical/IR waveband, but $`\nu _\mathrm{c}`$ has not yet done so, in which case $`p=(4\alpha +3)/3=2.31\pm 0.04`$ and $`\beta =(p1)/2=0.66\pm 0.02`$.

We can now compare the theoretically expected spectral slope ($`\beta `$) with the observed one. Although perhaps best determined by the spectroscopic measurements. our light curve analysis provides magnitudes of the OT in different bands, which can be directly used to derive $`\beta `$. We should note here that the brightness of the underlying galaxy and the possible contamination by an underlying supernova makes the determination of the magnitudes of the OT less reliable at later times, and our value of $`\alpha `$ is most likely dominated by the early part of the light curve when the OT was bright. Furthermore, $`\alpha `$ was kept fixed in $`all`$ the bands for the entire period of our observations, which implies that we cannot detect any spectral evolution of the OT. The $`synphot/calcphot`$ task in IRAF was used to determine the V$``$R and R$``$I values for a range of $`\beta `$ values, which were then compared with the colors derived from the light curve analysis to determine the spectral slopes. The V$``$R color implies a spectral slope of $``$1.8 $`\pm `$ 0.5, and the R$``$I color implies a spectral slope of $``$0.7 $`\pm `$ 0.1. Thus the slope derived from the V and R bands is closer to the case (i) mentioned above, and the slope derived from the R and I bands is closer to the case (ii). This is roughly consistent with the theoretical expectations and implies that, in the early part of the light curve, $`\nu _m`$ has passed the V and R and I bands, but $`\nu _c`$ has not passed the I band.

### 3.3. The Background Galaxy

The inferred magnitudes for an underlying galaxy suggest that the galaxy is relatively bright compared to the OT, and it may be possible to see the host galaxy directly in the images. The NTT images taken about 8 days after the burst had the best seeing ($``$ 1 arcsec) and the R images indeed show some hint of a slight extension. To further investigate this, all the images taken in B,V,R and I bands were co-added, giving rise to a single deep image of a total integration time of 40 min. From this combined image, we constructed a model point-spread-function (PSF) from a few isolated bright stars in the field using IRAF/DAOPHOT. This model PSF was used to subtract the point source contribution from the OT. The PSF-subtracted image shows some residual distribution elongated in the east-ward direction of the OT over a total of about 2 arcsec, which is probably the contribution of the underlying galaxy. Thus HST observations would show the detailed structure of the galaxy (Fruchter et al. in prep).

The derived magnitudes of the underlying galaxy were used to calculate its luminosity. The GRB is at a Galactic latitude and longitude of –40 and 315 degrees, respectively. The extinction models of Burstein & Heiles (1982) and Schlegel et al. (1998) yield E(B–V) of 0.015 and 0.033, respectively, for the position of the GRB. This extinction is small and similar to SN1998bw, which is neglected in our analysis. To derive the luminosity of the host galaxy, we use a redshift of z=0.434 which, for H0 = 70 $`\text{km\hspace{0.17em}s}^1`$and $`\mathrm{\Omega }`$=0.2, corresponds to a luminosity distance of 2160 Mpc, and a distance modulus of 41.67 mag. Applying the K-correction for z=0.434, the inferred luminosity of the galaxy is of the order of $`L_{}`$ (depending on the presence of an underlying supernova), where $`L_{}`$ corresponds to the luminosity of a typical galaxy at the ‘knee’ of the observed galaxy luminosity function (see, e.g., Lilly et al. 1999). If we exclude SN1998bw (which is not only ‘nearby’ ($`z=0.008`$), but for which the energy released in $`\gamma `$-rays is much smaller than any other GRB), then GRB990712 is the closest ‘cosmological’ GRB, and the apparent magnitude of its host galaxy is the brightest among all the GRB host galaxies observed so far. Compared to the host galaxy of GRB990713, which is the next brightest, the host galaxy of GRB990712 is about 0.25 magnitudes brighter in R if there is an underlying supernova, and more than a magnitude brighter if there is no underlying supernova.

We dedicate this paper to the memory of our colleague and friend Jan van Paradijs who passed away due to an illness on November 2, 1999, during the final stages of the preparation of this manuscript. Jan received several awards for his fundamental contributions in this field, including the 1998 Rossi prize of the American Astronomical Society. We will sorely miss his contagious enthusiasm, his brilliant insights, and his benevolent personality. This project is partly supported by the Chilean Fondecyt Project No. 3990024 and the Danish Natural Science Research Council (SNF). M. Dominik is supported by a Marie Curie Fellowship (ERBFMBICT972457).


format: text

text:
# NONTHERMAL EMISSION FROM A SUPERNOVA REMNANT IN A MOLECULAR CLOUD

## 1. INTRODUCTION

Supernovae and their remnants have long been identified as likely sites of galactic cosmic ray acceleration. Evidence for electron acceleration first came from observations of radio synchrotron radiation. The observation of nonthermal X-ray emission from SN 1006 by ASCA has provided evidence for electron acceleration up to $`100`$ TeV (Koyama et al. 1995; Reynolds 1998). The proton cosmic ray component is more difficult to observe, but it has long been recognized that pion decays from collisions with interstellar gas could give an observable flux of $`0.11`$ GeV photons. Initial estimates of $`\gamma `$-ray emission from supernova remnants (SNRs), concentrating on the component resulting from pion decays, were made by Chevalier (1977), Blandford & Cowie (1982), and Drury, Aharonian, & Völk (1994). When $`\gamma `$-ray emission was apparently detected from SNRs by CGRO, the $`\gamma `$-ray spectrum could not be fitted by a pure pion decay spectrum and some other component was needed (Esposito et al. 1996; Sturner et al. 1997; Gaisser, Protheroe, & Stanev 1998). In addition to the pion decays, the relevant processes are bremsstrahlung emission of relativistic electrons and inverse Compton emission. Gaisser et al. (1998) modeled these processes in detail in order to fit the observed $`\gamma `$-ray spectra of the remnants IC 443 and $`\gamma `$ Cygni. They assumed acceleration to a power law spectrum in the shock front and determined the spectral index, electron to proton ratio, and the upper energy cutoff.

The evolution of the relativistic component in a supernova remnant has been modeled by a number of groups. Sturner et al. (1997) assumed particle acceleration in the shock front with an $`E^2`$ energy spectrum and followed the evolution of the particle spectrum. Their SNR model was homogeneous: a factor 4 density jump at the shock front, constant density in the interior, and a constant magnetic field in the interior. Their shock acceleration model had a diffusion length parameter that was important for the maximum energy reached by the particles. They modeled the nonthermal emission from IC 443 as synchrotron emission in the radio and bremsstrahlung in $`\gamma `$-rays. de Jager & Mastichiadis (1997) dealt with the same processes in W44, as well as inverse Compton emission. They noted that the observed radio spectrum is flatter than would be expected from shock acceleration of newly injected particles and suggested that the particles originated from a pulsar in the supernova remnant. Recently, Ostrowski (1999) showed that second-order electron acceleration by the turbulent medium just after the shock could flatten the spectra to account for the radio observations of IC 443. Baring et al. (1999) presented calculations of the broad-band emission from nonlinear shock models of shell-type SNRs. They used Sedov adiabatic shock dynamics in a homogeneous medium and a Monte Carlo simulation of the particle acceleration taking into account the nonlinear shock structure. The set of models considered by Baring et al. (1999) covers the range of shock speeds 490 $`v_S`$ 4000 km s-1 and ambient medium number densities 10$`{}_{}{}^{3}n`$ 1 cm-3.

Massive stars that are the likely progenitors of core collapse supernovae are expected to be spatially correlated with molecular clouds. The remnants that are likely to be $`\gamma `$-ray sources in CGRO observations (Esposito et al. 1996) also show evidence for interaction with molecular gas. Chevalier (1999) recently studied the evolution of supernova remnants in molecular clouds and concluded that many aspects of the multiwavelength observations could be understood in a model where the remnants evolve in the interclump medium of a molecular cloud, which has density of $`525`$ H atoms cm-3, and become radiative at radii $`6`$ pc. The compression in the radiative shell is limited by the magnetic field. Molecular emission occurs when the radiative shell collides with molecular clumps. IC 443 is the remnant with the best evidence for high energy emission and also shows clear evidence for interaction of the shock with a molecular cloud (e.g., Burton et al. 1990; van Dishoeck et al. 1993; Cesarsky et al. 1999).

Most of the current models of particle acceleration in SNRs have dealt with adiabatic remnants in a homogeneous medium (see, however, the model of Boulares & Cox 1988 for the Cygnus Loop and the model of Jones & Kang 1993). In contrast, our aims here are to examine in detail the nonthermal emission of radiative supernova remnants in molecular clouds using a kinetic model of electron injection, acceleration and propagation, including the inhomogeneous structure deduced from multiwavelength SNR observations. We use global MHD (magnetohydrodynamic) models of supernova remnants interacting with molecular clouds as described by Chevalier (1999), including the theory of radiative shock structure by Shull & McKee (1979) and the kinetic model of electron injection and acceleration by shocks from Bykov & Uvarov (1999). The model described here is relevant to the mixed-morphology SNRs (Rho & Petre 1998), if they are interacting with molecular clouds. These remnants comprise a sizable fraction of the Galactic SNR population.

The plan of our paper is as follows. In § 2, we treat the energy spectrum of energetic electrons as they pass through the regions associated with a radiative shock wave. We include the shock acceleration of particles from the thermal pool. In § 3, we present calculations of the nonthermal emission from the relativistic particles, including emission from shocks in clumps as well as emission from the radiative shock wave in the interclump region. Our model is specifically applied to the well-observed remnant IC 443, although we expect it to more generally apply to other supernova remnants in molecular clouds. We also examine the ionization by energetic particles and the energy involved in the relativistic particles. A discussion of results and future propects is in § 4. A detailed discussion of expected range of densities, shock velocities and magnetic fields determined from multiwavelength observations of IC 443, W44 and 3C391 was given in the paper by Chevalier (1999, and references therein).

## 2. EVOLUTION OF NONTHERMAL ELECTRONS

The structure of the flow expected in a radiative shock can be divided into the following regions (see, e.g., Draine & McKee 1993):
(I) preshock gas ionized and heated by UV (ultraviolet) emission and fast particles from the shock;
(II) shock transition region;
(III) postshock cooling layer;
(IV) cold shell of swept-up gas;
(V) hot low density gas interior to the shock.
Radiative shocks are subject to thermal instabilities (Chevalier & Imamura 1982; Bertschinger 1986) as well as dynamical instabilities (Vishniac 1983; Blondin et al. 1998). Although the effect of inhomogeneous magnetic fields on the stability of radiative shock flow is yet to be studied, MHD turbulence is expected in such a system. The transport and acceleration of nonthermal particles in the violent SNR environment are governed by the MHD flow and depend on the ionization structure and on the spectrum of MHD turbulence. Shocks with velocity $`v_{S7}`$ 1.1 generate sufficient ultraviolet UV radiation for full pre-ionization of hydrogen and helium (He+) in the preshock region, where the gas temperature is typically $`10^4`$ K (Shull & McKee 1979). Here, $`v_{S7}`$ is the shock velocity measured in units of $`100\mathrm{km}\mathrm{s}^1`$, and $`n_1`$ is the ion number density measured in units of 10 cm-3. Shocks with $`v_{S7}`$ 1.5 propagating in a medium with $`n_1`$ 1 photoionize gas ahead of the viscous jump to a depth $`N_i>\mathrm{\hspace{0.33em}3}\times 10^{17}`$ cm-2. The shock transition in such a case is collisionless and supercritical unless the preshock magnetic field is above $``$ 4$`\times `$10$`{}_{}{}^{5}n_{1}^{0.5}`$ G.

The geometry of the magnetic field is important for the efficiency of electron injection and acceleration. It has been argued (e.g., Levinson 1996; McClements et al. 1997) that quasi-perpendicular shocks are an efficient source of freshly injected electrons if their speed $`v_s\left(m_e/m_p\right)^{0.5}c`$; the injection efficiency drops sharply below that threshold. On the other hand, Bykov & Uvarov (1999) have shown that quasi-parallel MHD shocks of moderate Alfven Mach number $`M_a\left(m_e/m_p\right)^{0.5}`$ are suitable sites of electron injection and acceleration. The SNR shock waves in molecular clouds are expected to be in this regime of Mach number. Thus we shall discuss below the quasi-parallel portion of the shock (magnetic field inclinations $`\pi /4`$) propagating into a molecular cloud with $`v_{S7}1.5`$.

Nonresonant interactions of the electrons with fluctuations, generated by kinetic instabilities of the ions in the transition region of a quasi-parallel supercritical shock, play the main role in heating and preacceleration of the electrons. The electron heating and preacceleration (injection) occur in the collisionless shock front transition region on a scale $`\mathrm{\Delta }`$ of several hundred inertial lengths of the ions, $`l_i`$ (Bykov & Uvarov 1999). Here $`l_i=c/\omega _{pi}9\times 10^6n_1^{0.5}`$ cm and $`\omega _{pi}`$ is the ion plasma frequency. For a SNR shock in a molecular cloud, the scale $`\mathrm{\Delta }`$ 109 cm is much shorter then the scales of all the other regions.

To accelerate electrons injected in the shock transition region to relativistic energy $`E(p)`$, MHD turbulence should fill the acceleration region, which has a scale $`l_ak_m(p)/v_s`$. Here, $`k_m=\mathrm{max}\{k_i\}`$, where $`k_i(p)`$ is the diffusion coefficient of an electron of momentum $`p`$ and $`i=`$(I), (III) refers to the regions defined above. To provide the conditions for efficient transformation of the MHD flow power to the accelerated electrons, the MHD turbulent fluctuation spectrum in the vicinity of the shock should extend to scales $`lr_e(p)`$ (resonant scattering), where $`r_e(p)`$ is the gyroradius of an electron with momentum $`p`$. The MHD fluctuations responsible for the electron scattering are collisionless for particles with energies up to a GeV because the ion-neutral collision length is above 1014 cm for typical parameters of the radiative shock described below (a keV electron has $`r_e`$ 10$`{}_{}{}^{8}[B/10^6\mathrm{G}]_{}^{1}`$ cm).

To calculate the spectra of nonthermal electrons in the regions i = I $``$ IV, we use a kinetic equation for the nearly-isotropic distribution function $`N_i(z,p,t)`$:

$`{\displaystyle \frac{}{t}}N_i+u_i(z){\displaystyle \frac{}{z}}N_i{\displaystyle \frac{p}{3}}{\displaystyle \frac{}{p}}N_i\left({\displaystyle \frac{}{z}}u_i\right)=`$
$`k_i(p){\displaystyle \frac{^2N_i(z,p)}{z^2}}+{\displaystyle \frac{1}{p^2}}{\displaystyle \frac{}{p}}p^2D_i(p){\displaystyle \frac{N_i}{p}}+{\displaystyle \frac{1}{p^2}}{\displaystyle \frac{}{p}}[p^2L_i(p)N_i].`$ (1)
The Fokker-Planck type equation (1) takes into account diffusion and advection (with bulk velocity $`u_i[z]`$) of an electron in phase space due to particle interaction with MHD waves and the large scale MHD flow (Bykov & Toptygin 1993). Here, $`L_i(p)`$ is the momentum loss rate of an electron due to Coulomb collisions in a partially ionized plasma (e.g., Ginzburg 1979). The momentum diffusion coefficient $`D(p)`$ is responsible for second order Fermi acceleration, and $`k_i(p)`$ is the fast particle spatial diffusion coefficient. For low energy electrons (i.e., $`EE_C`$), Coulomb and ionization losses are important in the regions (I $``$ V), except for the narrow shock transition region (II) where acceleration is fast enough to overcome losses and nonthermal electron injection occurs. The characteristic energy $`E_C`$ depends on the plasma density and ionization state, magnetic field, and fast particle diffusion coefficients $`k_i(p)`$ and $`D_I(p)`$. The corresponding momentum $`p_{ci}`$ can be estimated as the point where the last two terms in equation (1) are equal. We used the energy $`E_C`$ calculated simultaneously with the electron distribution function as a convenient parameter to distinguish between different possible models of MHD turbulence in the postshock cooling layer of the radiative shock structure (see below in this section).

The diffusion coefficients $`k_I(p)`$ and $`D_I(p)`$ depend on the spectrum of the collisionless MHD turbulence, which is poorly understood. A plausible approximation for the diffusion coefficients $`k_i(p)`$ was assumed here. For all of the regions we used the following parameterization:

$$k_i(p)=k_{i0}\{\begin{array}{ccc}1,\hfill & p_Tpp_{}\hfill & \\ vp^a/v_{}p_{}^a,\hfill & p_{}pp_{}\hfill & \\ cp^2/(p_{}^{2a}p_{}^av_{})\hfill & p_{}pp_m.\hfill & \end{array}$$
(2)
Here, $`p_T`$ is the momentum of the upstream thermal electrons, $`p_m`$ is the momentum corresponding to the upper cut-off energy $`E_m`$ (see equation ), and $`p_{}`$ and $`p_{}`$ are defined in the next paragraph. The standard relation for the momentum diffusion coefficient $`D_i(p)`$ = $`p^2w_i^2/9k_i(p)`$ (e.g., Berezinsky et al. 1990) was used. In the low energy regime $`pp_{}`$, the particle transport is dominated by large scale turbulent advection (Bykov & Toptygin 1993). The large scale turbulence is due to stochastic MHD plasma motions on scales $`\mathrm{\Lambda }`$ longer than the particle mean free path due to resonant scatterings. The sources of MHD turbulent motions are the shock wave instabilities mentioned at the beginning of this section. The large scale vortex rms-velocity $`w_i`$ typically is a fraction of bulk speed $`u_i`$. The spatial diffusion is energy independent below $`p_{}`$ and $`k_{i0}w_i\mathrm{\Lambda }`$. For $`p_{}pp_{}`$, electrons are scattered by resonant MHD waves presumably generated by the streaming instability of shock accelerated particles (e.g., Blandford & Eichler 1987; Jones & Ellison 1991). In this case $`w_i`$ is close to the Alfven velocity $`v_a`$ and $`k_i(p)`$ is the electron diffusion coefficient due to resonant wave-particle interactions. If the power spectrum of magnetic field fluctuations of Alfvenic turbulence for the resonant fluctuations is approximated as $`dB_k^2/dkk^\theta `$ then the power-law index $`a`$ (equation ) is related to the index $`\theta `$ as $`a=2\theta `$. An important particular case corresponding to $`a=1`$ can be written as $`k_i(p)=(1/3)\eta _ivr_e(p)`$. The parameter $`\eta _i`$ determines the scattering “strength” and strong scattering (Bohm limit) corresponds to $`\eta _i`$ 1. The diffusion model with $`a=1`$ was successfully used for modeling the observed anomalous cosmic ray proton fluxes in the interplanetary medium (Ellison et al. 1999) and in the Monte Carlo simulations of nonlinear shocks in the shell type SNRs (e.g., Baring et al. 1999). The momentum $`p_{}`$ in our model can be estimated from the equation $`(1/3)\eta _ivr_e(p_{})=w_i\mathrm{\Lambda }`$.

MHD waves in the partially ionized plasma of the radiative shock flow are subject to ion-neutral damping (e.g., Völk et al. 1981). We define the momentum $`p_{}`$ from the condition $`r_e(p_{})=\lambda _{in}`$, where $`\lambda _{in}`$ is determined from the relation $`\lambda _{in}=2\pi v_a/\mathrm{\Gamma }_{in}`$ and $`\mathrm{\Gamma }_{in}`$ is the Alfven wave damping rate due to ion-neutral collisions (Kulsrud & Cesarsky 1971). For $`pp_{}`$, electrons are scattered by nonresonant small scale waves.

The limits imposed on the diffusive shock acceleration of particles due to ion-neutral MHD wave damping have been discussed by Draine & McKee (1993) and Drury et al. (1996). They considered the spectra of Alfven waves driven by the instabilities of shock accelerated ions in a partially ionized plasma. Drury et al. (1996) obtained the upper cut-off energy of a proton due to ion-neutral wave damping

$$E_mv_{S7}^3T_4^{0.4}n_n^1n_i^{0.5}𝒫_1^{CR}\mathrm{GeV},$$
(3)
where $`n_n`$ is the neutral particle density, $`n_i`$ is the ion density (both are measured in $`\mathrm{cm}^3`$), and $`𝒫_1^{CR}`$ is the total particle pressure normalized to 10% of the shock ram pressure. The ionization structure of the preshock region is inhomogeneous so the ion number density, temperature and other numbers in equation (3) must be taken at an upstream distance $`k_I(E_m)/v_S`$ from the collisionless shock transition. Application to our model yields $`E_m`$ 10 GeV if $`v_{S7}`$ ($`1.11.5`$) and $`k_I(E_m)10^{24}`$ cm2 s-1 because of the high ionization in the preshock region (I) and the postshock cooling layer (III) provided by the UV radiation. One can justify the same value of $`E_m`$ for the electrons because the synchrotron and inverse Compton losses are relatively unimportant for $`E_m`$ 10 GeV.

Region I. In the preshock region of the shock with $`v_{S7}>`$ 1.1, we assume a highly ionized plasma of temperature $`T_4`$ 1 up to depths $``$ a few times 1017 cm-2 upstream. For the IC 443 preshock region, we adopted a density $`25`$ cm-3 and $`v_{S7}`$ 1.5 (Fesen & Kirshner 1980). The magnetic field $`B_{}10^5`$ G (Chevalier 1999), while $`B_{}5\times 10^6`$ G (the inclination is assumed to be $`\pi /6`$). The Alfven velocity here is about 5 km s-1. We assume a conservative diffusion model in the shock upstream region where $`k_I(p)`$ has $`a`$ = 1 and the moderate scattering “strength” 10 $`\eta _I`$ 100 (Ellison et al. 1999).

Region II. In the shock transition region, injection and heating of the electrons occur due to nonresonant interactions with strong MHD fluctuations generated by the ions. Following the model developed by Bykov & Uvarov (1999) we take

$$D_{II}(p)p^2\overline{C}\left(\frac{\delta B}{B_0}\right)^2\left(\frac{v_a}{v}\right)^2\left(\frac{v}{l_i}\right),$$
(4)
where $`\overline{C}1`$ and $`\delta B>B_0`$. The Alfvenic Mach number of the shock, $`M_a>`$ 10, is moderate for our situation. Nevertheless, the effect of the accelerated particle pressure on the shock structure could be important (see, e.g., Jones & Ellison 1991; Baring et al. 1999). An exact treatment of the effect is not feasible in our model because it would require a kinetic description of the ion injection. We accounted for the effect by correction of the total shock compression ratio $`\delta _t`$ and by introducing an extended modified shock structure with a precursor and a subshock with compression ratio $`\delta _s<\delta _t`$ (e.g., Jones & Ellison 1991). In our case with Mach numbers $`M>`$ 15, the estimated extra compression was about 10% and thus $`\delta _t`$ 4.4. A subshock compression ratio $`\delta _s`$ 3.3 is consistent with the assumption of substantial heating of the precursor gas due to wave dissipation (Berezhko & Ellison 1999). These numbers are in agreement with those from a hybrid simulation of nonlinear shock structure (Giacalone et al. 1997). We found that the resulting high energy emission spectra are sensitive to the particular choice of precursor flow structure in the MeV regime. The keV and GeV emission is not sensitive to the structure if the compression ratios $`\delta _t`$ and $`\delta _s`$ are fixed.

Region III. The postshock cooling layer has about a column density $`3\times `$1017 cm-2 of highly ionized plasma with an initial density about $`\delta _t`$ times that of the preshock one. The Alfven velocity in the highly ionized portion of the region (III) is similar to that at the preshock region, $``$ 5 km s-1. Large scale vorticity with an amplitude $`<\mathrm{\hspace{0.33em}20}\mathrm{km}\mathrm{s}^1`$ on scales $`\mathrm{\Lambda }<\mathrm{\hspace{0.33em}10}^{14}`$ cm may occur here because of shell instabilities of the radiative shock with $`v_{S7}`$ 1.5 (see, e.g., section 5 in the review by Draine & McKee 1993). This would dominate the low energy electron propagation, providing a spatial diffusion coefficient of $`10^{20}`$ cm2 s-1 for the radiative shock parameters described above.

Coulomb losses are important here for electrons with energies below $`E_C`$ 20 keV if the large scale turbulent velocity has a substantial longitudinal component with $`w_{III}20\mathrm{km}\mathrm{s}^1`$. $`E_C>`$ 8 MeV is expected in the absence of a longitudinal component of large scale turbulence in the postshock cooling layer because the Coulomb losses are overcome by resonant interaction with Alfven waves. Since the Alfven velocity is relatively low ($`5\mathrm{km}\mathrm{s}^1`$), the value of $`E_C`$ is much higher in that case. An accurate description of MHD turbulence in the postshock cooling layer is not available now. Thus we considered both limiting cases described above and show the uncertainty introduced by the lack of data concerning the large scale turbulence properties in Fig. 1. To avoid overestimating the nonthermal emission, we used a conservative minimum value of $`E_C`$ = 120 keV for the radiative shock structure described above, although one could expect $`E_C`$ 20 keV in the most favorable case.

Region IV. The density and temperature in the radiative shell depend on the preshock density, magnetic field and abundances (Shull & McKee 1979; Chevalier 1999). From the model of Chevalier (1999) for the above preshock parameters, we adopted $`n`$ 275 cm-3 ($`N_H1.6\times `$ 1020 cm-2) and $`T`$ 102 K in the dense radiative shell. The mass of the dense radiative shell is $`1000`$ $`M_{}`$ and its magnetic field is $`6\times `$ 10-5 G. The ionization in the dense radiative shell is supported at the level of about a percent by shock accelerated particles penetrating into the shell (see section 3.5 for details). The phase velocity of the high frequency collisionless Alfven waves which are responsible for electron resonant scattering is about 50 km s-1; it is determined by the number density of the ions. The Coulomb loss rate of the electrons is then dominated by the neutral component.

We numerically calculated the electron distribution function in regions (I-IV). The results were then used in calculations of the nonthermal emission.

## 3. NONTHERMAL EMISSION FROM IC 443

The SNR IC 443 (G189.1+3.0) is a very good candidate for testing our model because the interaction of the SNR with a molecular cloud was established from radio and infrared line observations (e.g., DeNoyer 1979; Mufson et al. 1986; Burton et al. 1990; van Dishoeck et al. 1993; Richter et al. 1995; Claussen et al. 1997). The high energy fluxes are somewhat better determined than in other cases. Although we concentrate on IC 443, we briefly mention the case of W44, where observations of OH masers and CO clumps again provide convincing evidence of molecular cloud interaction (Claussen et al 1997; Frail & Mitchell 1998).

To model the nonthermal emission from the shell of IC 443 we integrated the local emissivities over the radiative shock structure (regions I $``$IV). As a parcel of gas containing accelerated nonthermal electrons evolves through the radiative shock structure including the dense shell, the local emissivity from bremsstrahlung and inverse Compton was calculated using standard theory (see Appendix A for details). The electron - ion bremsstrahlung emission differential cross sections (Bethe-Heitler formulae) were taken from Akhiezer & Berestetsky (1957) taking into account the Elwert factor, which is important for modeling the keV emission. The electron - electron bremsstrahlung contribution to the $`\gamma `$-ray emissivity was calculated using cross sections derived by Haug (1975). For the inverse Compton emissivity of IC 443, we used the same description of the background photon field as Gaisser et al. (1998). The hard emission production cross sections used in our calculations are similar to those used recently by other authors (Asvarov et al. 1990; Sturner et al. 1997; Gaisser et al. 1998; Baring et al. 1999). We also estimated the contribution of pion decays to the $`\gamma `$-ray emission using a model for the proton component and the cross sections from Dermer (1986). The resulting emission in the EGRET regime was roughly 7 times below that from electron bremsstrahlung, which is consistent with the results of Sturner et al. (1997). We do not discuss the pion component further.

The synchrotron emission of the electrons calculated using the scheme described by Ginzburg (1979) can be compared to the radio observations and to the upper limits at higher frequencies. We also included free-free absorption of low frequency radio waves (dotted line in Fig. 1). The synchrotron losses of a relativistic electron were included following Ginzburg (1979); bremsstrahlung energy losses are relatively unimportant in our case.

### 3.1. Energetic Nonthermal Emission from the Shell

The emission was integrated over the regions of the radiative shock in order to model the spectrum of the entire shell. Region IV is the dominant contributor to the synchrotron and bremsstrahlung emission. The results in Fig. 1 can be compared to the ASCA and EGRET observations of SNRs, as well as to the significant upper limits that have been set at TeV energies by Whipple observations (Buckley et al. 1998) and from the preliminary analysis of CGRO OSSE observations presented by Sturner et al. (1997). We tried different parameter sets to fit the observations.

In Fig. 1, $`\nu F_\nu `$ fits to the observations of IC 443 are presented for the model described in the previous section. The preshock density is taken to be $`n=25\mathrm{cm}^3`$ and the assumed distance to IC 443 is 1.5 kpc (e.g., Fesen & Kirshner 1980). Interstellar photoelectric absorption of the X-ray spectrum is accounted for with a line-of-sight column density $`N_H=2\times 10^{21}`$ cm-2 using the cross sections from Morrison & McCammon (1983).

For $`v_{S7}`$ = 1.5 in the diffusion model described by equation (2), we used $`k_{III0}=1.1\times 10^{20}`$ cm2 s-1, $`a`$ = 1.0, $`E(p_{}`$) = 1 MeV, and $`E(p_{}`$) = 20 GeV, which are compatible with theoretical estimates of resonant wave-particle interactions. The corresponding value $`\eta `$ 30. While in principle $`\eta `$ can be determined from observations of heliospheric shocks, single spacecraft observations do not provide direct values and are subject to ambiguous interpretation. Values of $`\eta 30`$ are consistent with modeling of the highly oblique solar wind termination shock (Ellison et al. 1999), but smaller values ($`\eta 310`$) have been inferred for interplanetary traveling shocks (Baring et al. 1997). On the other hand, very weak scattering ($`\eta 100`$) has been inferred from observations near corotating interaction regions (e.g., Fisk, Schwadron, & Gloeckler 1997) and in the solar wind when the interplanetary magnetic field is nearly radial (e.g., Möbius et al. 1998). We use $`\eta `$ 30 as a conservative value. If strong MHD turbulence is present in the upstream shock region, the $`\eta `$ values could be closer to 1. This might allow an electron acceleration model similar to that described earlier, but for slower shock velocities ($`<\mathrm{\hspace{0.33em}100}\mathrm{km}\mathrm{s}^1`$).

The diffusion propagation model is not unique. We obtained a similar fit to that presented in Fig. 1 for a diffusion model with $`k_{III0}`$ = 1020 cm2 s-1, $`a`$ = 0.7, $`E(p_{}`$) = 1 MeV, $`E(p_{}`$) = 20 GeV. As mentioned in §2, there is an important uncertainty concerning the lack of a quantitative model for the generation of large scale MHD turbulence in the postshock cooling layer. The uncertainty affects the hard X-ray spectral calculations. In view of this, we calculated the expected bremsstrahlung emission for two limiting cases covering the range of uncertainties. Curve 1 in Fig. 1 corresponds to fully developed Alfvenic large scale MHD turbulence in the postshock cooling layer ($`E_C`$ = 120 keV). Curve 2 corresponds to the case of a lack of large scale MHD turbulence ($`E_C`$ = 2 GeV). While the GeV bremsstrahlung emission (as well as synchrotron and inverse Compton emission) is similar in the 2 cases, the hard X-ray emission is sensitive to the uncertainty. We also obtained a similar fit to that presented in Fig. 1 for a model with $`v_{S7}`$ = 1.1 and a diffusion model described by $`k_{III0}=5\times 10^{19}`$ cm2 s-1, $`a`$ = 1.0, $`E(p_{}`$) = 1 MeV, and $`E(p_{}`$) = 10 GeV. The model with $`v_{S7}`$ 1.1 requires a somewhat higher power conversion efficiency than that with $`v_{S7}`$ 1.5 (see §3.6), but this value is preferred for the IC 443 radiative shock velocity (Chevalier 1999).

The model of electron injection from the cloud thermal pool described above is based on the ionization structure of a radiative shock, which is very sensitive to the shock velocity if $`v_{S7}`$ 1 (Shull & McKee 1979; Hollenbach & McKee 1989). In this respect, it is instructive to consider a model in which high energy electrons are not related to freshly injected particles, but are reaccelerated cosmic ray electrons (Blandford & Cowie 1982; Chevalier 1999; Cox et al. 1999). To model that case we suppose that the far upstream electron flux is just that observed for Galactic cosmic ray electrons near Earth, but extrapolated back to the energy $`E_{crm}`$ with the same slope as observed above a GeV. The calculations of $`\nu F_\nu `$ for this case are shown in Fig. 2. We consider the radiative shock structure with the parameters described in §2 (the same as were used for Fig. 1). The Galactic cosmic ray electron flux in the far upstream region was taken to be 3$`\times 10^2(E/\mathrm{GeV})^3`$ e± s-1 cm-2 GeV-1 sr-1 (e.g., Berezinsky et al. 1990) for $`EE_{crm}`$. We assumed a flat cosmic ray electron spectrum below $`E_{crm}`$, corresponding to flattening due to Coulomb losses in the interstellar medium. The maximum energy of cosmic ray electrons reaccelerated by the radiative shock is $`E_m`$ 10 GeV (see equation ), as in the case shown in Fig. 1. The energy $`E_{crm}`$ was considered as a free parameter here.

We do not distinguish between electrons and positrons because there is no difference between them in the MHD shock acceleration process considered above and most of the emission from high energy leptons involved in our model remains the same for both kind of particles. Because of interactions with the radiative shock structure, accelerated cosmic ray electrons could provide a good fit for both synchrotron radio and EGRET emission if $`E_{crm}`$ 10 MeV for the radiative shell parameters described above, but fall short for $`E_{crm}`$ 100 MeV.

The relatively low energies, $`E_{crm}`$ 10 MeV, required for the IC443 shell are a potential problem for this model. The radio synchrotron radiation from the Galaxy indicates that the cosmic ray electron spectrum flattens to $`N(E)E^2`$ below energies $`1`$ GeV (Webber 1983). If we use such a spectrum in our model, the nonthermal emission is significantly below that observed. However, there are uncertainties about betatron-type acceleration in the shock compression, about the magnetic field strength, and about the volume of the emitting region. Duin & van der Laan (1975) did successfully model the IC 443 radio flux with a cosmic ray compression model. We found that models with a higher magnetic field and/or larger volume could reproduce the radio flux, but that the $`\gamma `$-ray fluxes observed by EGRET could not be achieved. Cox et al. (1999) have recently concluded that cosmic ray compression by radiative shocks could provide a good fit to radio and gamma-ray emission from W44 if the cosmic ray electron population in the vicinity of the SNR is higher than that of measured in solar neighborhood. They used a somewhat higher magnetic field ($`B`$ 200 $`\mu `$G) in the radiative shell of W44 than the 60 $`\mu `$G field we assumed for the IC 443 shell.

The acceleration time is a constraining parameter in this model because of the relatively low shock speed $`v_{S7}1.5`$ in the dense medium. The test particle shock acceleration time $`t_a(E)`$ can be estimated from the equation

$$t_a=\frac{3}{v_Iv_{II}}_{p_0}^p\left(\frac{k_I(p)}{v_I}+\frac{k_{II}(p)}{v_{II}}\right)\frac{dp}{p}$$
(5)
where the velocities $`v_i`$ are measured in the shock rest frame (e.g., Axford 1981). For the radiative shock parameters described above, the acceleration time for $`E>`$ 10 GeV is about 3,000 years, for $`\eta `$ 30 (note that $`t_a\eta `$), consistent with the age of IC 443 estimated from the hydrodynamical model by Chevalier (1999). On the basis of Ginga observations of IC 443, Wang et al. (1992) suggested an age of about one thousand years and an association with the supernova of AD 837. The reasons for the low age are related to a thermal interpretation of the observed X-ray emission above 10 keV. We shall argue below for the nonthermal origin of the hard emission in a radiative shock model. A distinctive feature of the radiative shock SNR shell model with direct injection of the electrons from the thermal pool is a hard spectrum of nonthermal emission extending from keV to GeV energies (Fig. 1). There is also a hard emission component due to inverse Compton emission of relativistic electrons having a scale comparable to the scale of the whole radio image of the remnant.

### 3.2. Radio Emission from the Shell

Observations of the radio emission from IC 443 were performed with different instruments over the last 35 years (see, e.g., Erickson & Mahoney 1985; Green 1986; Claussen et al. 1997 and references therein). The 151 and 1419 MHz radio maps (Green 1986) reveal that the most intense flux comes from the shell of the remnant. We suggest that the eastern part of the shell is a radiative shock which has a high magnetic field and electron flux, dominating the synchrotron radiation from the remnant.

We calculated the synchrotron radiation from the remnant using the local emissivity given by equation (A17). The integrated spectrum from the remnant is dominated by the emission from the radiative shell (having magnetic field $`B6\times 10^5`$ G) if the magnetic field in the central parts of the remnant (region V) is moderate, $`B4\times 10^6`$ G. For modeling the shell emission, we took into account the interstellar free-free absorption at low radio frequencies, assuming a column density to IC 443 of $`N_H2\times 10^{21}`$ cm-2 and an average interstellar medium temperature $`T10^3`$ K. This column density is similar to that adopted for the recent analysis of ASCA data by Keohane et al. (1997) and is about four times less than the value used earlier by Petre et al. (1988) and Wang et al. (1992). We assumed that the ionization fraction is $`15`$%. The free-free absorption inside the cool radiative shell is not very important above 50 MHz because the emission measure is about 1 cm-6 pc. The solid and dot-dashed lines on Fig. 1 show the total synchrotron radiation from IC 443 (dominated by shell emission in our model) expected at the Earth compared to the observational data (Erickson & Mahoney 1985). Both of our models are illustrated: direct injection of thermal electrons (Fig. 1) and reacceleration of preexisting cosmic ray electrons by a radiative shock (Fig. 2). A radiative shock with total compression ratio $`\delta _t`$ = 4.4, including a shock precursor and subshock ($`\delta _s`$ = 3.3), was considered. An important feature of the radio observations of IC 443 and W44 is that their spectral indices are relatively flat, implying electron energy spectra flatter than those produced by standard strong shock acceleration in the test particle limit. In our model, the calculated radio spectrum fits that observed from the IC 443 shell because of the large shock compression ratio $`\delta _t4.4`$ and to a lesser extent because of the effect of second order Fermi acceleration discussed recently by Ostrowski (1999).

### 3.3. Nonthermal Emission from a Shocked Clump

Molecular clouds show complex structure with an ensemble of dense massive clumps embedded in the interclump medium (Blitz 1993). We considered the possibility of radiative shock interaction with dense molecular clumps, resulting in high energy emission. In IC 443 and W44, those regions do not appear to be significant sites of GHz radio continuum emission, but the clumps may contain most of the cloud mass and are possible sites of hard X-rays and MeV $`\gamma `$-ray emission. The interaction of a radiative shell with a molecular clump generates a dense slab bounded by two shocks (Chevalier 1999). The forward shock propagating into the clump is strong (see also the simulations of Jones & Kang 1993). Since the clump size is typically smaller or comparable to that of the radiative shock layer, it is embedded in a pool of nonthermal particles created by the radiative shock according to the scenario described above. Accelerated particles as well as ionizing radiation from the hot interior regions of the SNR should provide ionization of the external portions of a dense clump at the level $`>\mathrm{\hspace{0.33em}1}`$%. This is sufficient for a shock propagating into the dense clump with $`n`$ 104 $`\mathrm{cm}^3`$ and $`B`$ 0.1 mG (e.g., Claussen et al. 1997) to be collisionless. The percent ionization level provides Alfven Mach numbers $`M_a>`$ 1 and weak damping of collisionless Alfven waves up to wavelengths $``$ 1010 cm for a flow with velocity $`>`$$``$ 25 $`\mathrm{km}\mathrm{s}^1`$ in the clump. Then the injection and acceleration processes considered above would be efficient and the ionization in the preshock region would be maintained self-consistently at least at a percent level even for a shock propagating inside the dense clump.

The gas ionization by UV radiation in the radiative shock precursor is sensitive to the shock speed if $`v_{S7}`$ 1.1 (Shull & McKee 1979; Hollenbach & McKee 1989). This implies that a shock having $`v_{S7}`$ 1 propagating in a dense clump would be able to accelerate electrons only below $`E_m<`$ 1 GeV (see equation ). We shall estimate below also the gas ionization rate due to accelerated nonthermal particles. The lifetime of a pre-existing electron of MeV energy against ionization losses in the interclump medium ($`n`$ 10 cm-3) is about $`10^4`$ years. Unless a nearby source of fresh accelerated electrons is present, Galactic cosmic ray MeV electrons are unlikely to be present in the cloud. Direct injection of electrons from the thermal pool of the weakly ionized precursor would be a dominant source of shock accelerated electrons up to the cut-off energy $`E_m`$. Because of a lack of accelerated GeV electrons in that case ($`v_{S7}`$ 1) the source would not appear as a bright radio and hard EGRET $`\gamma `$-ray emission feature. Nevertheless, these shocks might manifest themselves as hard X-ray and MeV $`\gamma `$-ray emission sources with $`\nu F_\nu \nu ^a`$ $`(a>\mathrm{\hspace{0.33em}0.5})`$ up to a few MeV and with a soft spectrum at higher energies. A comparison of Fig. 3 with Fig. 1 shows that the hard X-ray flux at 7 keV and above could be dominated by the emission from the shocked clumps, which is consistent with ASCA observations of IC 443 presented by Keohane et al. (1997). The hard X-ray nonthermal emission flux may be variable on a timescale of years. A relatively high MeV emissivity is a distinctive feature of shocked clumps, but they are difficult to observe. They must be nearby (within a few hundreds of parsecs) to be detected with COMPTEL/OSSE CGRO.

For a forward shock of velocity $`v_{S7}`$ 1 propagating into a dense clump of number density $`n`$ 103 cm-3 with upstream ionization level of a few percent (at $`k_I[E_m]/v_S`$ distance), the upper cut-off energy $`E_m<`$ 0.5 GeV. This is about 40 times less than that for a radiative shock in the interclump medium. The spatial diffusion $`k_{i0}`$ below $`p_{}`$ scales as $`w_i\mathrm{\Lambda }`$. Since the scale of a postshock radiative cooling layer $`n^1`$ and $`\mathrm{\Lambda }`$ is a fraction of that, we have $`k_{III0}3\times `$ 1018 cm2 s-1 in the clump case.

The effective compression of the nonthermal electrons as well as the magnetic field by the system of two shocks bounding the slab could be high. The magnetic field in the clump is expected to be high. Claussen et al. (1997) derived line-of-sight magnetic field strengths of $``$0.2 mG, remarkably uniform on the scale of several parsecs. Relativistic electrons accelerated by a shock propagating into a clump should produce radio emission with a flat spectrum up to some hundreds of MHz, which is lower than the maximum radio frequency from the shell because of the lower $`E_m`$ in the clumps. The emission measure of a clump interacting with a radiative shock could be EM $`<`$$``$ 1,000 cm-6 pc, which is much higher than that of the shell. Thus, the internal free-free absorption could be important here providing a flatter radio spectrum in the 100s of MHz regime. This is possibly a reason for the observed spatial variations of the radio spectral index in IC 443 (e.g., Green 1986).

We present in Fig. 3 the calculated $`\nu F_\nu `$ spectrum of nonthermal emission from a radiative shock - molecular clump interaction region. The forward shock velocity was high – $`100\mathrm{km}\mathrm{s}^1`$, the initial number density in the clump 1,000 cm-3, and the shocked clump radius 0.5 pc. The shock compression ratio was the same as that in the radiative shell model: $`\delta _s`$ 3 and $`\delta _t`$ 4.4. The ionization of gas in the shock upstream was about 5% (at the distances $`k_I[E_m]/v_S`$). To account for some uncertainty in the gas ionization, we present the radio spectra for two possible values of the emission measure EM of the shocked clump: EM = 110 cm-6 pc (solid line) and EM = 900 cm-6 pc (dashed line). The diffusion coefficient was $`k_{III0}3\times `$ 1018 cm2 s-1, $`E(p_{}`$) = 1 MeV. We illustrate the results for two possible maximal energies of the electrons (depending on the time of shock - clump interaction). Curve 1 in Fig. 3 corresponds to $`E_m=`$ 0.5 GeV and $`E_m=`$ 0.05 GeV for curve 2. The synchrotron radio spectrum for curve 2 is below the scale of the plot. The radio spectrum corresponding to curve 1 is very hard because of internal absorption in the shocked clump medium.

We presented above the results for nonthermal emission of a strong radiative shock of $`v_{S7}1`$ propagating into a molecular clump. For the case of radiative shock interaction with the molecular clumps in IC 443, lower shock velocities, $`v_{S7}`$ 0.3, are more realistic (e.g., Cesarsky et al. 1999; Chevalier 1999). The preshock ionization may be dominated mostly by nonthermal particles in such a shock because UV radiation is inefficient. A self-consistent model of MHD collisionless shock propagation into a dense clump requires simultaneous simulations of the MHD turbulence spectral properties and the nonthermal particle kinetics. Instead, we suppose that a collisionless shock transition exists for $`v_{S7}`$ 0.3 and $`n=10^4\mathrm{cm}^3`$ and assume the same diffusion model as described in §2. We calculated the upstream gas ionization, finding that the ionized fraction $`x>`$ 0.9 may hold up to depths about 1015 cm-2 due to shock accelerated electrons. The Mach number of the shock with $`v_{S7}`$ 0.3 is typically below 5 because of heating of upstream gas up to 104 K by accelerated particles and the compression ratio is lower than 4. The maximum energy of accelerated electrons in that case is typically below an MeV. We illustrate the calculated $`\nu F_\nu `$ spectrum of nonthermal emission from a 30 $`\mathrm{km}\mathrm{s}^1`$ velocity shock ($`\delta _t`$ = 3 ) interacting with a molecular clump of density 10$`{}_{}{}^{4}\mathrm{cm}_{}^{3}`$ by curve 3 in Fig. 3; this model corresponds to $`k_{III0}7\times `$ 1016 cm2 s-1 and $`E(p_{}`$) = 10 keV. Such collisionless shocks, if they exist, would provide a softer spectrum of nonthermal radiation than that for high velocity shocks with $`v_{S7}>`$ 1. A more comprehensive study is required to model the nonthermal emission from the low velocity shocks in detail. A search for hard X-ray emission correlated with molecular emission predicted by the simplified model described above is possible with high resolution instruments like Chandra and XMM.

### 3.4. Hard X-ray Emission

IC 443 was a target of X-ray observations with HEAO 1 (Petre et al. 1988), Ginga (Wang et al. 1992), ROSAT (Asaoka & Aschenbach 1994) and ASCA (Keohane et al. 1997). ASCA GIS observations discovered the localized character of the hard X-ray emission (Keohane et al. 1997). Most of the 2-10 keV GIS photons came from an isolated emitting feature and from the southeast elongated ridge of hard emission. The ridge is also coincident with the 95% confidence error circle of the EGRET CGRO source (Esposito et al. 1996), while the isolated feature is outside the circle (Keohane et al. 1997). The integrated 2-10 keV flux observed by ASCA GIS was (5$`\pm 1)\times `$ 10-11 erg cm-2 s-1, representing more than 90% of the total flux from the SNR (Keohane et al. 1997). The hard X-ray ridge, as defined by Keohane et al. (1997), is only a part of the more extended radio shell discussed in the previous section. The flux density at 7 keV from the isolated emitting feature was about 4$`\times `$ 10-5 photons cm-2 s-1 keV-1 and that from the ridge was 2$`\times `$ 10-5 photons cm-2 s-1 keV-1. HEAO A-2 measured a flux at 2–10 keV of (7$`\pm 1)\times `$ 10-11 erg cm-2 s-1 (Petre et al. 1988), while Ginga measured a 2–20 keV flux of about 9$`\times `$ 10-11 erg cm-2 s-1 (Wang et al. 1992). These numbers are consistent with the assumption that the integrated flux is dominated by the extended soft component coming from the central parts of the remnant.

The soft X-ray 0.2–3.1 keV surface brightness map of IC 443 from the Einstein Observatory (Petre et al. 1988) shows bright features in the northeastern part of the remnant. The presence of nearly uniform X-ray emission from the central part of the remnant is a characteristic feature of mixed-morphology SNRs (Rho & Petre 1998). It might be due to the effect of thermal conduction (Cox et al. 1999; Shelton et al. 1999), although Harrus et al. (1997) argued for an alternative radiative shock model for W44, another mixed-morphology SNR. That feature in our model corresponds to the emission from hot (T$`<\mathrm{\hspace{0.33em}10}^7`$ K), low density gas interior to the shock (region V).

In the radiative shock model with direct injection from the thermal pool and second order Fermi acceleration in postshock cooling layer, most of the hard X-ray emission ($`>\mathrm{\hspace{0.33em}7}`$ keV) comes from the shell - radiative shock structure with a rising $`\nu F_\nu `$ X-ray spectrum (curve 1 in Fig. 1 for the shell and curves 1 and 2 in Fig. 3 for a shocked molecular clump). In addition to the shell, there is a more extended source of hard X-rays from inverse Compton emission. It is the dominant source of hard X-rays in the scenario with cosmic ray electron reacceleration (Fig. 2) and in the case of a lack of second order Fermi acceleration in the radiative shell (see curve 2 in Fig. 1). In the scenario with reacceleration of cosmic ray electrons, it is impossible to obtain the 7 keV flux observed by ASCA from the IC 443 ridge. Thus, a resolved hard X-ray image of the IC 443 remnant may be very informative in determining the electron acceleration scenario.

The isolated emitting feature in the southern part of IC 443 detected by ASCA and brightest around 7 keV, has a flat spectrum low frequency radio continuum source (Green 1986) and is close to a region of excess H2 emission. High spatial resolution observations of this isolated emitting feature are required to decide between a low-luminosity pulsar nebula and shock interaction models (e.g., radiative shock interaction with a dense molecular clump, Fig. 3). Observations of the hard emission with the sensitivity and resolution available with BeppoSAX, Chandra XRO, XMM and GLAST would be required to constrain the extension and spectral properties of the hard X-ray emission from IC 443.

### 3.5. Gas Ionization by Nonthermal Electrons

Nonthermal particles accelerated by a MHD collisionless shock wave provide an efficient ionizing agent. We calculated above the nonthermal emission generated by accelerated energetic electrons; one can also calculate the ionization produced by the same electrons in the molecular cloud. This might connect the high energy observations with radio, IR and optical emission from shocked atomic and molecular gas. To estimate the gas ionization due to accelerated electrons in the radiative shock structure we used the electron spectra calculated with the kinetic model of electron acceleration and propagation described in §2 and electron ionization cross sections from Appendix B. In Fig. 4, we present the calculated ionization rate $`\zeta _e`$ and estimated ionized fraction $`x`$ in the shell of a radiative shock due to accelerated electrons for the case of $`v_{S7}=1.5`$ for both models of large-scale turbulence in the postshock cooling layer considered in §3.1.

There should also be contributions to $`x`$ from accelerated nucleons and from UV radiation which are not included in Fig. 4. We accounted for only radiative recombination, assuming a small molecular fraction in the radiative shell at column densities below 2 $`\times 10^{20}`$ cm-2.

The solid curves in both panels of Fig. 4 correspond to the case of fully developed Alfvenic turbulence in the postshock cooling layer, which provides efficient second order Fermi acceleration making Coulomb losses relatively unimportant above $`E_C`$ = 120 keV. This case corresponds to curve 1 in Fig. 1 for emission from the shell with substantial hard X-ray emission. The dotted lines in Fig. 4 correspond to the case of a lack of large scale MHD turbulence where the Coulomb losses overcome the second order Fermi acceleration up to $`E_C`$ = 2 GeV. The nonthermal emission expected in this case is illustrated by curve 2 in Fig. 1. The radiative shell ionization by nonthermal electrons is sensitive to the turbulent structure of the postshock cooling layer and there is a correlation between the hard X-ray spectrum and the ionization structure of the radiative shell. We shall consider the effect of nonthermal particles on the ionization structure of a shocked clump elsewhere.

### 3.6. Energy in Nonthermal Electrons

The energy of the electron component for the parameter set used to compute the spectra in Fig. 1 is substantial. The power required for electron acceleration and maintenance for curve 1 in Fig. 1, $`1.5\times 10^{37}`$ erg s-1, is high because of the strong Coulomb losses of keV electrons in the dense ionized medium, but it is consistent with the estimate of the total radiative losses from IC 443, dominated by infrared emission (Mufson et al. 1986; Burton et al. 1990). A similar power is required for both models 1 and 2 in Fig. 3, where radiative shock - molecular clump interaction is illustrated. For curve 2, the power requirements are less, $`6\times 10^{36}`$ erg s-1. The implied efficiency of power conversion from the MHD shock flow (which can be estimated as 3 $`\times 10^{38}`$ erg s-1) to the nonthermal electrons is about 5% in these cases. The model with $`v_{S7}1`$ would provide emission similar to that shown in Fig. 1, but requires about three times higher efficiency. This is higher than the electron acceleration efficiency estimated from GeV regime cosmic rays observed near the Earth, but evolved SNRs interacting with dense molecular gas probably cannot be considered as the main source of the observed Galactic cosmic rays because of the relatively low maximum energies of accelerated particles.

The pressure of the nonthermal electron component downstream of the shock (region IV) is $`𝒫_e`$ 7$`\times `$ 10-11 erg cm-3. The magnetic pressure in the dense shell is about 1.5 $`\times `$ 10-10 erg cm-3, which is higher than the thermal gas pressure. The uniform ($``$ 60 $`\mu `$G) and stochastic magnetic field components dominate the total pressure in the shell.

Let us consider the energy requirements for the scenario with electron injection from preexisting cosmic rays (Fig. 2). The nonthermal electron pressure in this case is $`𝒫_e`$ 5.2$`\times `$ 10-11 erg cm-3 (for $`E_{crm}`$ 8 MeV), similar to that for the model with injection from the thermal pool, because it is dominated by GeV particles. An advantage of this scenario is that much less power is required. It is about 4 $`\times `$ 1035 erg s-1, which is less than 10% of that for the scenario with injection from the thermal pool. The difference is because Coulomb losses are unimportant for the high energy electrons involved in that scenario. However, a potential problem is the relatively low $`E_{crm}`$ required. The lifetime of a 20 MeV electron against ionization losses in a neutral medium of number density $``$ 1 cm-3 is about 106 years. This implies that the source of MeV cosmic ray electrons should exist within about 100 pc of the molecular cloud if the diffusion coefficient is $`3\times `$ 1026 cm2 s-1 at these energies.

## 4. DISCUSSION AND FUTURE PROSPECTS

For individual supernova remnants, the ambient density is an important parameter that can be estimated from multiwavelength observations (Chevalier 1999). Another important parameter for modeling the nonthermal particles in a SNR in a molecular cloud is the collisionless MHD turbulence spectrum, particularly in the radiative shock cooling layer. A substantial level of MHD collisionless turbulence could overcome Coulomb losses in the dense plasma downstream from the radiative shock. An accurate model of MHD turbulence in the postshock cooling layer is not available now. Thus we considered both limiting cases described in the § 2 and show in Fig. 1 the uncertainty introduced by the lack of data concerning the MHD turbulence properties. Coulomb losses are important in the postshock cooling region only for electrons with energies below $`E_C`$ 20 keV if the large scale turbulent velocity has a substantial longitudinal component of $`20\mathrm{km}\mathrm{s}^1`$. A higher value ($`E_C>`$ 8 MeV) is expected in the absence of a longitudinal component of large scale turbulence in the postshock cooling layer because the Coulomb losses are overcome by resonant interaction with Alfven waves. Since the Alfven velocity is relatively low ($`5\mathrm{km}\mathrm{s}^1`$), the value of $`E_C`$ is much higher in that case. We used a conservative minimum value of $`E_C`$ = 120 keV for the radiative shock structure described above, although one could expect $`E_C`$ 20 keV in the most favorable case. The nonthermal emission of SNR in a molecular cloud (especially in the hard X-ray regime) and ionization structure of the radiative shock are sensitive to the MHD turbulence model. Spatially resolved observations of the nonthermal emission from SNRs may be able to constrain models of MHD turbulence.

All of the SNR candidates from the EGRET list: IC 443, $`\gamma `$ Cyg, W44, and Monoceros (Esposito et al. 1996) are old remnants interacting with molecular clouds. The $`\gamma `$-ray emission from electrons accelerated by the radiative shock calculated for the IC 443 parameters presented in Figs. 1 and 2 is in good agreement with that observed by EGRET (Esposito et al. 1996), as well as with the upper limits established by Whipple (Buckley et al. 1998). The high energy $`\gamma `$-rays ($`>\mathrm{\hspace{0.33em}50}`$ MeV) from the radiative shell should be spatially correlated with the radio emission. The EGRET telescope detected an extended excess (95% confidence circle $`42^{}`$ in radius) from the Monoceros SNR correlated with the 1.42 GHz radio emission (Esposito et al. 1996). This is in accordance with the model of a leptonic origin of GeV emission from extended SNRs interacting with clouds (de Jager & Mastichiadis 1997). GLAST (Gamma-ray Large Area Space Telescope) will be an excellent instrument for future $`\gamma `$-ray observations, because it will have the capability of spatially resolving the $`\gamma `$-ray emission from SNRs.

Spatially resolved observations of W44, IC 443, 3C391 and some other mixed-morphology SNRs from the list given by Rho & Petre (1998) with BeppoSAX, Chandra XRO, XMM would be very valuable tools to test our model. The hard X-ray detector (HXD) aboard the forthcoming ASTRO-E mission could be used for observations of the hard 10–700 keV continuum with a field of view of 0.8 FWHM at 60 keV and 2.8 at 500 keV from the extended hard X-ray emission from IC 443, W44 and the Monoceros SNR predicted by the radiative shock model (see Fig. 1 and Fig. 3). Spatially resolved spectra are needed to distinguish the shell emission from the shocked clump emission.

Molecular clumps interacting with moderately fast radiative shocks are also expected to be sources of hard X-rays and MeV $`\gamma `$-rays (up to 100s of MeV). Nonthermal continuum radio emission (100s of MHz) with a sharply rising $`\nu F_\nu `$ spectrum (Fig. 3) and a time dependent cut-off frequency is expected in the fast shock model. Due to potentially substantial internal free-free absorption, the spectrum of radio emission from a localized clump might constrain the ionized gas density in the clump. The MeV $`\gamma `$-ray spectrum of a localized clump may be resolved with the forthcoming $`\gamma `$-ray missions INTEGRAL (e.g., Schönfelder 1999; Winkler 1999) and GLAST. With an expected angular resolution about 12 FWHM, imager IBIS aboard INTEGRAL would allow detection of hard X- ray ( $`>`$$``$ 50 keV) emission from molecular clouds within $``$ 1 kpc. A comparison of hard emission spectra with the radio spectrum can provide valuable information about the density and magnetic field in a clump. A possible variability on a timescale of a few years for clump hard emission and radio emission (on a longer time scale) could further constrain the model.

For low velocity (below $`30\mathrm{km}\mathrm{s}^1`$) shocks interacting with a dense ($`10^4\mathrm{cm}^3`$), magnetized ($`B<`$ 0.5 mG) molecular clump, one may expect hard X-ray emission below the MeV regime correlated with the regions of molecular emission of shocked gas. Nonthermal radio continuum is not expected in the case of low velocity shocks because of a lack of accelerated relativistic electrons at sufficiently high energies.

Finally, we note that the high density of energetic particles in the vicinity of a shock wave in a molecular cloud can affect the ionization and thermal properties of the gas. Rich molecular spectra have been observed from IC 443 and 3C391 (Reach & Rho 1999). It will be interesting to see whether there is a signature of the presence of nonthermal particles that can be discerned from the molecular spectra.

A.M.B thanks Hans Bloemen for useful discussions. We are grateful to M.G. Baring for comments and for making available the corrections to the e-e bremsstrahlung formulae and to the referee for very constructive comments. The work of A.M.B and Yu.A.U was supported by the INTAS grant 96-0390 and that of R.A.C. by NASA grant NAG5-8232.

## Appendix A Appendix A: Radiative Processes

Given the electron distribution function, the high energy emission flux $`J(E_\gamma )`$ at a distance $`R`$ from the source can be calculated from the equation:

$$J(E_\gamma )=\frac{1}{4\pi R^2}\underset{V}{}𝑑V\frac{dn_\gamma (E_\gamma ,\stackrel{}{r})}{dt}$$
(A1)
where $`dn_\gamma (E_\gamma ,\stackrel{}{r})/dt`$ is the emissivity, which here includes bremsstrahlung, synchrotron radiation and the inverse Compton effect. We generally assumed a pure hydrogen composition in our calculations.

High energy particles penetrating through a partially ionized medium produce photons due to interactions with atoms, electrons and ions. We assumed a Maxwellian distribution for the ambient matter and in most cases took the target ambient particles to be at rest. Then, the emissivity can be calculated from

$$\frac{dn_\gamma (E_\gamma )}{dt}=4\pi \left(\underset{E_\gamma }{\overset{\mathrm{}}{}}𝑑E_eN_p\sigma _{ep}(E_\gamma ,E_e)J(E_e)+\underset{E_{min}}{\overset{\mathrm{}}{}}𝑑E_eN_e\sigma _{ee}(E_\gamma ,E_e)J(E_e)\right),$$
(A2)
where $`J(E_e)`$ is the electron flux measured in units s-1 cm-2 keV-1 sr-1 , $`\sigma _{ep}`$ and $`\sigma _{ee}`$ are differential cross sections of one photon emission due to electron-proton and electron-electron interactions integrated over angles, and $`E_{min}`$ is the minimum energy of incoming electron that could generate a photon with energy $`E_\gamma `$. For e-p interactions, $`E_{min}=E_\gamma `$ is a good approximation. The correction due to finite proton mass is less than the uncertainty in the cross section. In the case of e-e interactions, we used the following equation for $`E_{min}`$:

$$E_\gamma =E_{min}m_ec^2/\left[2m_ec^2+E_{min}\sqrt{E_{min}(E_{min}+2m_ec^2)}\right].$$

The inverse bremsstrahlung contribution is negligible unless there is a very high (far in excess of 100) ratio of protons to electrons.

We used the following cross section approximations in our calculations:

### A.1. e-e bremsstrahlung.

Haug (1975) has obtained a general formula for the e-e bremsstrahlung cross section, but it is more convenient to use simplified approximations. We used relativistic and non-relativistic asymptotic expressions, which were matched in the intermediate energy region (see Baring et al. 1999).

In the relativistic limit, we used approximate formulae given by Baier et al. (1967) with a correction factor from Baring et al. (1999):111Equations (A4) and (A6) are somewhat different from that given by Baring et al. (1999). M.G. Baring (private communication) made available the corrections to their original formulae.

$$\sigma _{ee}(E_e,E_\gamma )=(\sigma _1+\sigma _2)A(ϵ_\gamma ,\gamma _e),$$
(A3)
where $`\gamma _e=(E_e+m_ec^2)/m_ec^2`$, $`ϵ_\gamma =E_\gamma /m_ec^2`$,

$$\sigma _1(\gamma _e,ϵ_\gamma )=\frac{4r_o^2\alpha }{ϵ_\gamma }\left[1+\left(\frac{1}{3}\frac{ϵ_\gamma }{\gamma _e1}\right)\left(1\frac{ϵ_\gamma }{\gamma _e1}\right)\right]\left[\mathrm{ln}\left(2(\gamma _e1)\frac{\gamma _e1ϵ_\gamma }{ϵ_\gamma }\right)\frac{1}{2}\right],$$
(A4)
$$\sigma _2=\frac{r_o^2\alpha }{3ϵ_\gamma }\{\begin{array}{cc}[16(1ϵ_\gamma +ϵ_\gamma ^2)\mathrm{ln}\left(\frac{\gamma _e}{ϵ_\gamma }\right)\frac{1}{ϵ_\gamma ^2}+\frac{3}{ϵ_\gamma }4+4ϵ_\gamma 8ϵ_\gamma ^2\hfill & \\ 2(12ϵ_\gamma )\mathrm{ln}(12ϵ_\gamma )(\frac{1}{4ϵ_\gamma ^3}\frac{1}{2ϵ_\gamma ^2}+\frac{3}{ϵ_\gamma }2+4ϵ_\gamma )],\hfill & \hfill ϵ_\gamma \frac{1}{2},\\ \frac{2}{ϵ_\gamma }\left[\left(4\frac{1}{ϵ_\gamma }+\frac{1}{4ϵ_\gamma ^2}\right)\mathrm{ln}(2\gamma _e)2+\frac{2}{ϵ_\gamma }\frac{5}{8ϵ_\gamma ^2}\right],\hfill & \hfill ϵ_\gamma >\frac{1}{2},\end{array}$$
(A5)
$`r_o=e^2/m_ec^2`$ is the classical electron radius, $`\alpha =e^2/\mathrm{}c`$, and

$$A(ϵ_\gamma ,\gamma _e)=1\frac{10}{3}\frac{(\gamma _e1)^{\frac{1}{5}}}{\gamma _e+1}\left(\frac{ϵ_\gamma }{\gamma _e}\right)^{\frac{1}{3}}.$$
(A6)
According to Baring et al. (1999), these formulae agree with those from Haug (1975) with an accuracy $`10\%`$ for electrons of energy $`5`$ MeV .

In the non-relativistic limit, we used equations from Garibyan (1952) and Fedyushin (1952) with correction coefficients from Baring et al. (1999):

$$\sigma _{ee}=\frac{4r_o^2\alpha }{15ϵ_\gamma }F\left(\frac{4ϵ_\gamma }{\gamma _e^21}\right)$$
(A7)
where the photon energy is in the range $`0<ϵ_\gamma <\frac{1}{4}(\gamma _e^21)`$ and $`F(x)`$ ($`0<x<1`$) has the form:

$$\begin{array}{ccc}F(x)\hfill & =\hfill & B(\gamma _e)\left[17\frac{3x^2}{(2x)^2}C(\gamma _e,x)\right]\sqrt{1x}\hfill \\ & & +\left[12(2x)\frac{7x^2}{2x}\frac{3x^4}{(2x)^3}\right]\mathrm{ln}\left(\frac{1+\sqrt{1x}}{\sqrt{x}}\right).\hfill \end{array}$$
(A8)
The correction coefficients are:

$$B(\gamma _e)=1+\frac{1}{2}(\gamma _e^21);C(\gamma _e,x)=\frac{10x\gamma _e\beta _e(2+\gamma _e\beta _e)}{1+x^2(\gamma _e^21)}$$
(A9)
Baring et al. (1999) found less than a $`10\%`$ difference between equation (A7) and the expressions of Haug (1975) for electron energies $`E_e500`$ keV.

### A.2. e-p bremsstrahlung.

In the non-relativistic limit, we took the e-p bremsstrahlung cross section from Akhiezer & Berestetsky (1957):

$$\sigma _{ep}=\frac{16r_0^2Z^2\alpha m_e^2c^2}{3p_1^2}\frac{4\pi ^2\zeta _1\zeta _2}{(e^{2\pi \zeta _1}1)(1e^{2\pi \zeta _2})}\frac{1}{ϵ_\gamma }\mathrm{ln}\left(\frac{p_1+p_2}{p_1p_2}\right),$$
(A10)
which is valid if $`\zeta _11`$; here $`\zeta _i=\alpha Z/\beta _i`$.

For the relativistic limit, we used Born approximation formulae from Akhiezer & Berestetsky (1957):

$$\begin{array}{cc}\sigma _{ep}=\hfill & r_o^2Z^2\alpha \frac{p_2}{p_1}\{\frac{4}{3}2ϵ_1ϵ_2\frac{p_1^2+p_2^2}{p_1^2p_2^2}+m_e(\frac{\eta _1ϵ_2}{p_1^3}+\frac{\eta _2ϵ_1}{p_2^3}\frac{\eta _1\eta _2}{p_1p_2})+\hfill \\ & L[\frac{8}{3}\frac{ϵ_1ϵ_2}{p_1p_2}+\frac{ϵ_\gamma ^2}{p_1^3p_2^3}(ϵ_1^2ϵ_2^2+p_1^2p_2^2)+\frac{m_e^2ϵ_\gamma }{2p_1p_2}(\eta _1\frac{ϵ_1ϵ_2+p_1^2}{p_1^3}\eta _2\frac{ϵ_1ϵ_2+p_2^2}{p_2^3}+\frac{2ϵ_\gamma ϵ_1ϵ_2}{p_1^2p_2^2})]\},\hfill \end{array}$$
(A11)
where $`L=\mathrm{ln}\left(\frac{p_1^2+p_1p_2ϵ_1ϵ_\gamma }{p_1^2p_1p_2ϵ_1ϵ_\gamma }\right)=2\mathrm{ln}\left(\frac{p_1^2+p_1p_2m_e^2}{m_eϵ_\gamma }\right)`$, $`\eta _1=\mathrm{ln}\left(\frac{ϵ_1+p_1}{ϵ_1p_1}\right)=2\mathrm{ln}\left(\frac{ϵ_1+p_1}{m_e}\right)`$, and $`\eta _2=2\mathrm{ln}\left(\frac{ϵ_2+p_2}{m_e}\right)`$. In the velocity range $`Z\alpha (\beta _2^1\beta _1^1)`$, we included the Elwert factor from Koch & Motz (1959): $`f_E=\frac{\beta _1(1exp[2\pi Z\alpha /\beta _1)])}{\beta _2(1exp[2\pi Z\alpha /\beta _2)])}`$.

In equation (A11), $`c=1`$ is assumed. In the intermediate energy range, we matched equation (A10) and equation (A11).

### A.3. Synchrotron radiation.

We used standard formulae for synchrotron radiation in our calculations (e.g., Ginzburg 1979). The emissivity can be written as:

$$\begin{array}{ccc}E_\gamma \frac{dn_\gamma (E_\gamma )}{dt}\hfill & =\hfill & \frac{\sqrt{3}\alpha e}{m_ec}𝑑\chi 𝑑EN(E,\chi )H\mathrm{sin}^2(\chi )F\left(\frac{\nu }{\nu _c}\right),\hfill \\ F(x)\hfill & =\hfill & x\underset{x}{\overset{\mathrm{}}{}}K_{5/3}(\eta )𝑑\eta ,\hfill \end{array}$$
(A12)
where $`\chi `$ is the angle between the electron velocity and the magnetic field, $`N(E,\chi )`$ is the electron distribution function over energy and angle, $`\nu _c=\frac{3eH\mathrm{sin}(\chi )}{4\pi m_ec}\left(\frac{E}{m_ec^2}\right)^2`$, and $`K_{5/3}(x)`$ is the McDonald function. Integrating equation (12) over angle under the assumption of a chaotic magnetic field orientation, we obtain the equation

$$E_\gamma \frac{dn_\gamma (E_\gamma ,\stackrel{}{r})}{dt}=\frac{\sqrt{3}\alpha e}{2\pi m_ec}𝑑EN(E,\stackrel{}{r})H_{}F\left(\frac{\nu }{\nu _c}\right),$$
(A13)
where $`N(E,\stackrel{}{r})`$ is the electron distribution function as a function of energy and position and $`H_{}=\sqrt{2/3}H(\stackrel{}{r})`$ is the averaged magnetic field.

### A.4. Inverse Compton emission

High energy electrons colliding with photons result in high energy photon production (e.g., Jones 1968; Gaisser et al. 1998; Sturner et al. 1997; Baring et al. 1999). If the electron energy $`E_e\frac{m_e^2c^4}{4E_\gamma }\frac{610^{10}}{E_\gamma (\text{eV})}`$eV ($`E_\gamma `$ is the energy of the emitted photon), then the formulae in the Thompson limit are (e.g., Ginzburg 1979):

$$\frac{dn_\gamma (E_\gamma )}{dt}=\sqrt{3}\pi \sigma _Tm_ec^2𝑑E_{ph}\sqrt{\frac{E_\gamma }{E_{ph}}}N_{ph}(E_{ph},\stackrel{}{R})J_e(m_ec^2\sqrt{\frac{3E_\gamma }{4E_{ph}}},\stackrel{}{R}),$$
(A14)
where $`\sigma _T=\frac{8\pi }{3}\left(\frac{e^2}{m_ec^2}\right)^2`$ is the Thompson cross section, $`N_{ph}(E_{ph})`$ is the background photon energy spectrum, and $`J_e(E)`$ is the incoming electron flux.

If the electron energy $`E\frac{m_e^2c^4}{4E_\gamma }\frac{610^{10}}{E_\gamma (\text{eV})}`$eV, we must use the Klein-Nishina formulae (integrated over angles):

$$\frac{dn_\gamma (E_\gamma )}{dt}=4\pi 𝑑E_e𝑑E_{ph}N_{ph}(E_{ph})\sigma _{KN}(E_\gamma ,E_e,E_{ph})J(E_e),$$
(A15)
$$\sigma _{KN}(E_\gamma ,E_e,E_{ph})=\frac{2\pi r_o^2}{ϵ_{ph}\gamma _e^2}\left[2q\mathrm{ln}(q)+1+q2q^2+\frac{q^2(1q)\mathrm{\Gamma }^2}{2(1+q\mathrm{\Gamma })}\right],$$
(A16)
$$\text{ where }\mathrm{\Gamma }=4ϵ_{ph}\gamma _e,q=\frac{ϵ_\gamma }{(\gamma _eϵ_\gamma )\mathrm{\Gamma }},\text{ }0q1.$$

In the above, $`ϵ_\gamma =E_\gamma /m_ec^2`$, $`ϵ_{ph}=E_{ph}/m_ec^2`$, $`\gamma _e`$ is the electron Lorentz factor, and $`r_o=e^2/m_ec^2`$ is the classical electron radius. The interstellar background photon spectrum, $`N_{ph}(E_{ph},R)`$, was adopted from Mezger et al. (1982) and Mathis et al. (1983) taking into account infrared data from Saken et al. (1992) for the case of IC443 (see also Gaisser et al. 1998).

### A.5. Free-free absorption

Radio waves propagating in the partially ionized thermal plasma are subject to thermal free-free absorption. For $`h\nu kT`$ (the Rayleigh – Jeans regime), the free-free absorption coefficient $`\alpha _\nu ^{ff}`$ (cm-1) can be calculated from

$$\alpha _\nu ^{ff}=\frac{4e^6}{3m_ekc}\left(\frac{2\pi }{3km_e}\right)^{1/2}T^{3/2}Z^2n_en_i\nu ^2\overline{g_{ff}}.$$
(A17)
Using an approximate formula for the Gaunt factor $`\overline{g}_{ff}(\nu ,T)`$ from Rybicki & Lightman (1979) appropriate to our parameter range, we obtained a simplified expression for the optical depth $`\tau (\nu ,T)`$

$$\tau 0.022\nu _{100}^{2.1}T_{100}^{1.34}EM,$$
(A18)
where the radio wave frequency $`\nu _{100}`$ is measured in units of 100 MHz, the plasma temperature in units of 100 K is $`T_{100}`$, and the emission measure $`EM`$ is in units of cm-6 pc.

## Appendix B Appendix B: Electron Ionization

The primary ionization rate $`\zeta _e`$ due to electron impact was calculated from the equation:

$$\zeta _e=4\pi \underset{E_{min}}{\overset{\mathrm{}}{}}𝑑E\sigma _{eion}(E)J_e(E),$$
(B1)
where $`J_e(E)`$ is the electron flux measured in units s-1 cm-2 keV-1 sr-1 , $`\sigma _{eion}`$ is the differential cross section for hydrogen atom primary ionization by electron-atom interactions integrated over atomic electron states, and $`E_{min}`$ is the minimum energy of the ionizing incoming electron.

In the regime below 1 keV, we obtained the following fit to the ionization cross section $`\sigma _{eion}`$ using experimental data from Fite & Brackmann (1958):

$$\sigma _{eion}(E)=\pi a_0^2A_1\mathrm{exp}(bA_4b^3A_5b^5A_6)/E(\mathrm{eV})^{A_3},$$
(B2)
where $`b=A_2/E(\mathrm{eV})`$ and $`E`$(eV) is the incoming electron energy measured in eV. The fitting coefficients are $`A_1`$=76.11, $`A_2`$=14.34, $`A_3`$=0.89, $`A_4`$=3.82, $`A_5`$=–2.55, and $`A_6`$=4.5. $`a_0=\mathrm{}^2/m_ee^2`$ is the Bohr radius.

At high energies ($`>\mathrm{\hspace{0.33em}1}`$ keV), the Born approximation provides good accuracy (e.g., Mott & Massey 1965):

$$\sigma _{eion}=_0^{\kappa _{max}}_{K_{min}}^{K_{max}}I_{0\kappa }(K)𝑑K𝑑\kappa ,$$
(B3)
where

$`I_{0\kappa }(K)`$ $`=`$ $`{\displaystyle \frac{2\pi \mathrm{\hspace{0.17em}2}^{10}\kappa }{\text{}a_0^2k^2K}}{\displaystyle \frac{\mu ^6\left[\text{}K^2+(\mu ^2+\kappa ^2)/3\right]}{\left[\mu ^4+2\mu ^2\left(\text{}K^2+\kappa ^2\right)+\left(\text{}K^2\kappa ^2\right)^2\right]^3}}`$ (B4)
$`\times `$ $`\mathrm{exp}\left[{\displaystyle \frac{2\mu }{\kappa }}\mathrm{arctg}\left({\displaystyle \frac{2\mu \kappa }{K^2\kappa ^2+\mu ^2}}\right)\right]{\displaystyle \frac{1}{1\mathrm{exp}\left(2\pi \mu /\kappa \right)}},`$
$`\mu =Z/a_0`$, $`k(E)`$ is the incoming electron wave number, $`k_\kappa `$ is the outgoing electron wave number, and $`\kappa `$ is the emitted atomic electron wave number. In addition, $`K=|\stackrel{}{k}\stackrel{}{k}_\kappa |`$, $`K_{max}=k+k_\kappa `$, and $`K_{min}=kk_\kappa `$.

We found a convenient approximation providing a reasonably accurate fit (better than 10% above 1 keV and better than 5% above 3 keV) to equations (B3) and (B4):

$$\sigma _{eion}(E)=\pi a_0^2P_1E(eV)^{P_2},$$
(B5)
where $`P_1`$=97.72 and $`P_2`$=0.93.

The gas ionization fraction $`x=n_i/(n_i+n_n)`$ (we consider here a hydrogen plasma) can be expressed through the ionization rate $`\zeta `$ and recombination rate $`\alpha _r(T)`$ as:

$$x=\left(\text{}s^2/4+s\right)^{0.5}s/2$$
(B6)
where $`s=\zeta /(\alpha _rn)`$. We used the radiative recombination rates given by Spitzer (1978) and dissociative recombination rates from Schneider et al. (1994).