text:
# Optical and infrared observations of the supernova SN 1999el1footnote 11footnote 1Based on observations collected with the AZT-24 telescope (Campo Imperatore, Italy), the 1.8-m telescope of the Osservatorio Astronomico di Padova (Asiago, Italy), the TNT telescope operated by the Osservatorio Astronomico di Collurania-Teramo (Teramo, Italy), the Italian National Telescope TNG (La Palma, Canary Islands, Spain), the 1-m telescope of the Observatorio Astronómico Nacional (Tonantzintla, Mexico), the 2.1-m telescope of the Guillermo Haro Observatory (Cananea, Mexico).

## 1 Introduction

A significant fraction of the total number of Type II supernovae (SN II) show distinctive features which led Schlegel (1990) to introduce the class of SN IIn. These objects are characterized by a blue continuum near maximum and a slow temporal evolution of the spectrum; the slowly declining luminosity in the optical bands and the narrow $`\mathrm{H}\alpha `$ line profile superimposed on a broad base indicate the presence of diffuse lower velocity matter around the site of the explosion. Because of their relatively high luminosity and slowly declining light curve, the SN IIn are easier to discover than faint, narrow peak light curve objects. Their intrinsic frequency has been estimated by Cappellaro et al. (1997) to be 5% of all SN II. There is, however, a marked inhomogeneity among the events (Cappellaro & Turatto 2001), namely concerning the decline rate of the light curves, a matter which merits understanding. It has not yet been clarified to what extent individual differences among members of this type of SNe can be totally ascribed to the physical properties of pre-existing circumstellar matter (CSM).

In fact, recent attempts to explain the long term properties of visual light curves (LCs) with light echoes from circumstellar matter by Roscherr & Schaefer (2000) have not been fully successful. There are observational features which disagree with echo model predictions, forcing these authors to conclude that any pre-existing dust around the site of the explosion must be almost completely destroyed at the moment of the outburst. However, this contrasts with the observed late infrared emission from SN 1998S, which has been inferred to result from dust cooling (Fassia et al. 2000). Based on their observational data, Fassia et al. argue that dust was already present within the precursor wind of SN 1998S, and hence at least a fraction of it must have survived the explosion. This suggests that improved observations are needed in order to understand the environment of SN sites from this diagnostic tool. So far only a few of the models developed to test the effects of light echoes have taken into account the near infrared (NIR) light curves (see Emmering & Chevalier 1988).

On the other hand, the slow decline of light curves for Type IIn SNe appears to be well accounted for by simple models of ejecta-wind interactions (Terlevich 1994; Plewa 1995). At the same time, the peculiar shape of the emission lines has been also attributed to ejecta-wind interactions involving a two-component wind, characterized by dense clumps (or a dense equatorial wind) associated with a less dense spherically symmetric distribution (Chugai & Danziger 1994).

Therefore, several open questions remain, which can be formulated as follows:

1. How important is the issue of the inhomogeneity of this class of events and what causes this inhomogeneity?
2. What is the evolutionary scenario causing some massive precursors of SNe to explode within a pre-existing envelope of diffuse circumstellar matter and what physical mechanisms affect the observational properties?
3. Is dust present in the CSM before the explosion?
4. To what extent do echo processes and ejecta-wind interactions determine the observational properties of such SNe?

In this paper we present optical and infrared data for the Type IIn SN 1999el. The present work provides new photometric measurements and improves the observational constraints on this class of SNe. Observations and data reduction are presented in section 2. The results are presented and discussed in sections 3 and 4, and summarized in section 5.

## 2 Observations and data reduction

SN 1999el was discovered on October 20.45 UT (julian day 2451472) at mag $`15.4`$ (unfiltered) by the BAO SN Survey (Cao et al. 1999, IAUC 7288) $`23.4\mathrm{"}`$ east and $`9\mathrm{"}`$ south of the galactic nucleus of NGC 6951, an inclined ($`i=42\mathrm{°}`$) barred spiral galaxy hosting a type 2 Seyfert AGN, located at $`\alpha (1950)=20:36:36.59`$ and $`\delta (1950)=+65:55:46.00`$ (see Fig. 1). For a list of galactic properties see, e. g., Kohno, Kawabe, & Vila-Vilaró (1999) and references therein. The distance to NGC 6951 was recently discussed by Vinkó et al. (2001). The estimated distance moduli range from 31.85 mag (23.4 Mpc) to 32.59 mag (33 Mpc), where the lowest values are based on the Tully-Fisher relation and the highest ones on the multi-color light curve shape (MLCS) method applied to the SN 2000E (a Type Ia SN). We adopt our own determination ($`32.1\pm 0.4`$ mag) based on our observations of SN 2000E (see Valentini et al. 2002, in preparation).

Cao et al. (1999) report no trace of the SN on images taken on October 19.2 by KAIT (to a limiting magnitude of 18) and on October 19.44 by BAO, so we can assume October 20.00 as the zero point for the light curves, since it must be very close (within a few hours) to the shock breakout time. According to Cao et al. (1999), an image taken on October 21.44 indicates a magnitude of $`15.0`$ for the supernova, whereas on October 22.2 it appeared at a magnitude of $`15.2`$ (no band is indicated). It can be seen from Fig. 1 that SN 1999el is projected between two nearby foreground stars (c2 and c3) with a separation of no more than a couple of arcsec.

### 2.1 Optical data: acquisition, reduction and photometry

Optical photometry of SN 1999el has been obtained on images collected with several telescopes at different sites in the northern hemisphere. These were reduced using the standard techniques for de-biasing and flat-fielding.

One set of data comes from the 0.72-m Teramo-Normale Telescope (TNT), at the Teramo Observatory (Italy). The TNT is an f/14 Ritchey-Chretien reflector equipped with a Tektronic Tk512 CB1-1 front-illuminated $`512\times 512`$ pixels CCD. The pixel size of 27 $`\mu `$m yields a scale of $`0.46`$ arcsec/pixel. During each observational night we obtained 4 to 12 frames per filter in the B, R and I bands.

Other photometry was performed on images obtained with the 1-m telescope operated by the Observatorio Astronómico Nacional in Tonantzintla (Mexico), equipped in the f/15 Cassegrain port with a 1024$`\times `$1024 Thomson THX 31156 CCD and a Cousins $`R`$ filter. The field of view covered by this configuration is $`4.2^{}\times 4.2^{}`$. The chip was binned $`2\times 2`$ and read in standard mode (3.47 $`e^{}`$/pixel, 2.7$`e^{}`$/ADU). The effective on-chip integration of the SN field was split into 3 to 6 on-chip integrations of 300 s, offset from each-other by $`20^{\prime \prime }`$.

For imaging, also the 1.8-m telescope of the Osservatorio Astronomico di Padova sited in Asiago, Cima Ekar (Italy), equipped with AFOSC, was used. The instrument is a focal reducer type spectrograph/camera with a 1K $`\times `$ 1K Site Thinned CCD ($`24\mu `$m). The scale of $`0.473\mathrm{}`$/pixel gives a field of $`8\times 8`$ arcmin2.

Some photometric data were acquired using the 2.1-m telescope of the Guillermo Haro Observatory in Cananea (Mexico), equipped with the spectrophotometer LFOSC at the f/12 Cassegrain focus. Johnson $`V`$ and Cousins $`R`$ filters were used. The chip is a $`385\times 578`$ EEV P8693, with a field of view of $`6^{}\times 10^{}`$. This was read in the standard readout mode: 8$`e^{}`$/pixel, 4.7$`e^{}`$/ADU. The effective on-chip integration of the SN field was split into 2 to 3 integrations of 60 s in the $`R`$-band and of 120 s in the $`V`$-band.

Late time $`UBVRI`$ photometry was obtained with the Optical Imager Galileo (OIG) mounted on the 3.58-m Italian National Telescope (TNG, Roque de los Muchachos, La Palma, Canary Islands) and equipped with a mosaic of two thinned, back-illuminated EEV42-80 CCDs with $`2048\times 4096`$ pixels each (pixel size $`13.5`$ $`\mu `$m; pixel scale in $`2\times 2`$ binned mode $`0.144\mathrm{}`$/pixel). These nights were photometric and were used to calibrate the local sequence around the SN (1 to 8, see Fig. 1) which in turn was used for calibrating non-photometric nights. The instrumental color equations were obtained through observations of stars in the standard fields of Landolt (1992). Three of the local standards (1, 3 and 5) are a subset of those selected by Vinkó et al. (2001) to carry out differential photometry ($`BVRI`$) on SN2000E (their F1, B1, B3, respectively). The magnitudes adopted in our and their study for those sources agree within $`0.05`$ mag, with two exceptions: we found our star 1 to be $`0.13`$ mag brighter in $`R`$ and our star 5 to be $`0.13`$ mag brighter in $`I`$. The agreement can be considered good and confirms the validity of both our calibration and that of those authors. However, our measurements were obtained using a larger telescope (the TNG) under very good seeing and photometric conditions. This should have afforded us a better sky subtraction which, particularly at longer wavelengths, may be easily affected by nebular emission from the host galaxy.

During the SN 1999el follow-up observations, another SN event occurred within NGC 6951 (i. e., SN 2000E), so the CCD at the TNT was operated to image the 2 objects on the same frame. The exposure times were chosen so as to achieve the best signal-to-noise ratio for both supernovae, but avoiding saturation of the brightest calibration star in the field.

Since both supernovae are well embedded in the galaxy background, we obtained instrumental magnitudes by point-spread function (PSF) fitting using the DAOPHOT package in IRAF or a set of procedures developed for SN photometry in the same environment. For the TNT and Tonantzintla observations, the final instrumental magnitudes were obtained as an average of all valid values derived for a given night.

Errors have been obtained with artificial star experiments, i.e. placing stars of the same magnitude as the supernova in parts of the galaxy arm characterized by similar complexity as the SN site, then computing the r.m.s. of the derived magnitudes. For TNT images however, the uncertainty in the photometry has been estimated by evaluating the PSF fit residuals and the r.m.s. of the $`n`$ successive measurements obtained in a single night for the SN after calibration through the local standards. It is more difficult to assess the error due to the convolution of the PSFs of the SN and the nearby sources when their brightness becomes comparable and the seeing is poor, which is the case in the TNT photometry at late epochs. Hence, we omit measurements which appear uncertain because of the inferior quality of the images. Table 1 lists the results of optical photometry for SN 1999el.

### 2.2 Infrared data: acquisition, reduction and photometry

Most near-infrared (NIR) observations of SN 1999el were obtained at the AZT-24 1.1-m telescope in Campo Imperatore (Italy) with SWIRCAM, during a period spanning October 25, 1999 to May 26, 2000. SWIRCAM (D’Alessio et al. 2000) incorporates a $`256\times 256`$ HgCdTe NICMOS-3(-class) detector which, at the focus of AZT-24, yields a scale of $`1.04\mathrm{}`$/pixel, resulting in a field of view of $`4\times 4`$ arcmin2. The observations were performed through standard $`J`$ (1.25 $`\mu `$m), $`H`$ (1.65 $`\mu `$m) and $`K`$ (2.20 $`\mu `$m) broad-band filters. Any on-source image was obtained as a median composition of five frames with a 15 arcsec dithering and integration times of 60 s each at the $`J`$ and $`H`$ bands, and 120 s each at the $`K`$ band. Off-source frames were taken with a 10 arcmin offset with the same observational procedure and exposure times, but a larger (50 arcsec) dithering. Sky images obtained as a median of the original dithered off-source frames were subtracted from the single on-source images before composition.

Flat-field frames were acquired at twilight using the differential flat technique that relies on the natural variation of the sky background. In this way, the frames with the highest mean counts are combined together (through median filtering) in one single image, as are those with the lowest mean counts, and the resultant lower signal image is subtracted from the resultant higher signal one so as to remove any biases and dark current. This is used to flat-field each image after sky-subtraction.

Deep $`JHK_\mathrm{s}`$ images of SN 1999el were obtained at the 3.58-m TNG telescope (La Palma) with ARNICA during the night of May 26, 2000. ARNICA is a $`256\times 256`$ pixels NIR camera (Lisi et al. 1996) and was matched to the TNG with a $`0.35`$ arcsec/pixel scale (yielding a $`90\mathrm{"}\times 90\mathrm{"}`$ field of view). Sequences of 4 dithered on-source and 4 off-source images were taken moving the telescope to a sky position after each on-source integration. Total on-source integration times amounted to 12 minutes at $`K_\mathrm{s}`$ and 8 minutes at $`J`$ and $`H`$. Five images per band of a group of standard stars (AS 40–5, AS 40–0 and AS 40–1; Hunt et al. 1998) were also obtained so as to locate the brightest of them (AS 40–5) at the center of the array and of each detector quadrant in turn. Flat-field frames in the $`K_\mathrm{s}`$ band were taken at sunset. Only the $`K_\mathrm{s}`$ frames were sky-subtracted after flat-fielding; off-source images (or the dithered ones for standard stars) were combined together and used for flat-fielding in the $`J`$ and $`H`$ bands. All images in each band were finally corrected for bad-pixels, registered and averaged together.

Late time NIR images were acquired with NICS at the TNG on December 14, 2000 (JD+893) during the commissioning of the camera. NICS is based on a $`1024\times 1024`$ Hawaii detector (Baffa et al. 2001) and we employed a $`0.25`$ arcsec/pixel scale resulting in a field of view of $`4.2\times 4.2`$ arcmin2. Dithered on-source and off-source frames were taken in the $`J_s`$ and $`K^{}`$ bands with total (on-source) integration times of $`1200`$ s and $`420`$ s, respectively. Data were reduced as described above, but no flat-field frames were available for $`K^{}`$ images.

The photometric measurements on SWIRCAM images were performed by PSF fitting using the DAOPHOT package in MIDAS in order to derive relative photometry between the supernova and a set of comparison stars properly selected in the field among the ones with the highest S/N ratio (1 to 9 in Fig. 1), of which those calibrated with the TNG data (see below) are a subset. For these, variability could be ruled out on the basis of an accurate image-by-image check of the fluctuations in their relative magnitudes. PSF fitting photometry was also needed owing to the partial convolution of the SN with the sources c1, c2 and c3 of Fig. 1.

The presence of the parent galaxy background has been taken into account by evaluating its influence on the derived photometry of the SN. In particular cases, i. e., when the S/N ratio was low, PSF fitting photometry was carried out using the ROMAFOT package in MIDAS, by means of which the background gradient due to the diffuse emission from the galaxy was accounted for by adopting a tilted plane.

Both aperture and PSF fitting photometry were performed on the ARNICA images for field stars using the package DAOPHOT in IRAF (Stetson 1987). However, we obtained more satisfactory PSF fits to SN 1999el using the package ROMAFOT in MIDAS. The group of 3 stars close, in projection, to the supernova, all resolved on the TNG images (see Fig. 1), is composed of, following the nomenclature used in the figure, c2 ($`K_\mathrm{s}=14.97\pm 0.03`$), $`2`$ arcsec east of the SN, c3 ($`K_\mathrm{s}=16.75\pm 0.32`$), $`2`$ arcsec west of the SN, and c1 ($`K_\mathrm{s}19`$). The observed standard stars were used to calibrate both the SN and the field stars 2 and 9 (see Fig. 1) which served as references for the photometry on SWIRCAM frames. $`K_\mathrm{s}`$ magnitudes were converted into $`K`$ magnitudes by extrapolating the continuum flux within the spectral region ($`0.1`$ $`\mu `$m wide) of $`K`$ not covered by the $`K_\mathrm{s}`$ filter, with corrections amounting to $`0.02`$ mag for the standard stars and to $`0.06`$ mag for SN 1999el (the uncorrected sources being fainter). Possible CO molecular line emission falling out of $`K_\mathrm{s}`$ but within $`K`$, and sometimes occurring in the case of Type IIn SNe, was therefore not detected and could not account for any fraction of the measured flux. The photometry of the standard stars for the 5 positions on the array indicates that residual errors of $`0.05`$, $`0.04`$ and $`0.02`$ mag in the $`J`$, $`H`$ and $`K_\mathrm{s}`$ bands, respectively, are still present owing to flat-field inaccuracies.

Differential photometry on the NICS images was performed by PSF fitting using ROMAFOT in MIDAS. By inter-comparing instrumental magnitudes obtained both for the field (and local standard) stars and for the close-by sources (c1-c3), we checked that even though not flat-fielded, $`K^{}`$ values where consistent at a $`0.10.2`$ mag level with ARNICA photometry throughout the whole final frame. We conservatively assumed uncertainties $`0.5`$ mag both for $`J_s`$ and for $`K^{}`$ and did not correct for the slight differences with the $`J`$ and $`K`$ bands. Table 2 lists the results of NIR photometry for SN 1999el.

### 2.3 Spectroscopy

Two spectra of the SN have been obtained a few hours apart in Asiago (November 7, 1999) and Cananea (November 8), i.e. about 18 days after the discovery. In Fig. 2 we show the latter obtained with a 5.5 Å/pixel grism, which for a $`3^{\prime \prime }`$ slit-width, yields a resolution of $`16`$ Å in the 4000–7200Å range. The total 4400 s exposure was split into 3 integrations to minimize the effect of cosmic rays. Wavelength calibrations were obtained using comparison spectra of HeArNe lamps, while the instrumental signature has been removed using a spectrophotometric standard star.

## 3 Results and analysis

### 3.1 Optical spectrum

The spectrum (Fig. 2) confirms the early finding by Filippenko (1999) that the SN shows lines with narrow components. Superimposed on the continuum there are strong, relatively narrow Balmer emission lines of H together with wider weak wings. The broad wings could be fitted with Gaussian profiles having FWHM of about 2500–3000 km s-1, while the narrow ones are barely resolved. Also narrow P-Cyg absorption components with minima displaced by about 600–900 km s-1 relative to the emission-line cores are present, indicating a slowly expanding shell of gas above the photosphere. Fe II lines have similar structures. Weak NaI D lines are visible.

In Fig.2 we compare the spectrum of SN 1999el with those of other H–dominated SNe at similar epochs showing signs of interaction with the CSM. The SN 1999el spectrum is very different from the spectra of SN 1988Z, the best studied Type IIn SN (Stathakis & Sadler 1991; Turatto et al. 1993; Chugai & Danziger 1994; Aretxaga et al. 1999). The line widths of SN 1988Z are much broader, implying expansion velocities of about 15000 km s-1. On the contrary the spectrum of SN 1999el closely resembles those of SNe 1999eb and 1995G (Pastorello et al. 2002, in preparation) at comparable age. They have similar continua, strong emission lines with similar profiles and narrow absorption components. It is the very narrow blue shifted absorption component particularly in the wing of the broad H$`\alpha `$ from the expanding envelope that is characteristic of a class of SNe, such as SN 1996L (Benetti et al. 1999), SN 1994aj (Benetti et al. 1998) and SN 1984E (Dopita et al. 1984). This absorption is commonly thought to arise from the CSM in the form of a wind.

### 3.2 Light curves

Figure 3 shows the LCs at optical and NIR bands; we draw attention to some particular features. Firstly, we see that a maximum of emission in the $`K`$ band has not been observed. However, it seems that LC maxima in the other NIR and optical bands have been recorded. There are indications of small delays in the time of maxima occurrence which increase with a decrease of wavelength. During the first 90 days, the decline rate is quite fast reminiscent of that of a Type II-L SN rather than that of a typical IIn SN, as already noted. After 90 days the LCs exhibit a much slower decrease and, in particular, that of the $`K`$ band flattens most. Since the $`K`$ point at JD+690 has been derived from a $`K_\mathrm{s}`$ measurement and any possible line contribution from CO emission is therefore excluded, the late behavior of the LC in $`K`$ is not due to CO line emission. Note that this point has been derived from TNG observations, so the SN and its nearby sources are well resolved. The NIR LCs appear to remain flat up to JD+893 ($`420`$ days from the shock outbreak), even if in the $`J`$ band the decrease seems slightly faster than in the $`K`$ band.

Figure 4 shows the time evolution of a series of color indices. We note the blueing in $`JH`$ and $`JK`$ during the first $`45`$ days (including $`10`$ days before the occurrence of maximum in the $`J`$ band). Afterwards, starting approximately around JD+550, a reddening in $`JK`$ and $`JH`$ colors begins to be evident. The increase in errors which is significant after JD+550 (affecting all colors) is due to the SN magnitude becoming comparable with that of the brightest of the nearby field stars, the two sources being difficult to deconvolve on the images from the smaller telescopes, even though we already noted that the $`K`$ point obtained at late times from TNG images (JD+690) is well established. This reddening effect is also evident as a flattening of the $`K`$ LC whereas $`J`$ and $`H`$ are still slightly steeper in their decline (see Fig. 3) and can be attributed to an infrared excess arising at later phases. SN 1999el exhibits a more standard behavior in $`BR`$, which increases after maximum and then flattens after $`70`$ days. This may be understood in terms of a cooling photosphere. In fact, the same trend seems to characterize all optical colors, with a change of $`0.5`$ mag in a period spanning $`200`$ days.

The behavior of the NIR emission can be followed on a color-color diagram ($`HK`$ vs. $`JH`$; see Fig. 5). Comparing the SN colors with the main sequence locus (solid line), the blackbody locus (dashed line) and the reddening law (arrow), a NIR excess is evident even at early epochs and this cannot be simply explained as reddened blackbody or stellar photospheric emission. After day 5, the colors evolve roughly along a track which appears almost parallel to the reddening law, evolving as if AV decreases with time. However, on day 219, a strong $`K`$ excess is present. The most likely cause for the NIR excess is continuum emission from hot dust pre-existing within the CSM of SN 1999el. Following Fassia et al. (2000), we can test this scenario in the following way. Given the (reddened) flux in the $`K`$ band for day 5, $`3.39\times 10^{12}`$ erg cm-2 s-1 $`\mu `$m-1, and the distance assumed in Sect. 3.4, we can determine the radius of a blackbody surface emitting the same luminosity in $`K`$ at 1500 K, roughly the dust sublimation temperature. It can be shown that reddening does not change the result by more than a 10–20%. This radius, $`1.5\times 10^{16}`$ cm, is comparable to the light travel distance from the explosion time to day 5, $`1.3\times 10^{16}`$ cm. Thus the two independently determined distances are consistent with the presence of dust at a radius where it has not sublimated and radiates at IR wavelengths as a result of heating by UV photons from the original outburst.

### 3.3 Reddening

In order to give an estimate of the absolute luminosity of SN 1999el it is necessary to determine its reddening. As can be seen in Fig. 2, the interstellar lines of NaI cannot be resolved and therefore add little useful information. Hence, we can only assume the galactic value as a lower limit. From the maps of Schlegel, Finkbeiner, & Davis (1998) we obtain AB $`=1.57`$ mag. This map appears to systematically overestimate the extinction where $`E(BV)0.15`$ mag, as pointed out by Arce & Goodman (1999). Since this occurs towards SN 1999el, as well, we expect a lower galactic contribution to the extinction. We dereddened optical and NIR magnitudes using the interstellar extinction law given by Cardelli, Clayton, & Mathis (1989) and assuming an AV ranging between $`1.18`$ (corresponding to AB $`=1.57`$ mag) and $`2.5`$ mag. These were then converted to fluxes according to Bessel (1979) and, in the NIR bands, the calibration given for the UKIRT standard. In Table 3 both apparent and absolute (dereddened) magnitudes at maxima are listed.

### 3.4 Blackbody fits

We fitted a blackbody law to the dereddened fluxes in order to get estimates of the photospheric temperature. Since not all bands are covered on the same days, we constructed sets of $`UBVRIJHK`$ fluxes complementing the available values with ones obtained through interpolation or using the results of spline fits to all data in a given wavelength. Then, we adopted a one-parameter procedure for $`\chi ^2`$ minimisation normalizing, for each selected day, all fluxes to that in the $`R`$-band in order to remove the dependence on radius and distance. $`R`$ was chosen because it is the better sampled band (see Table 1) and is likely to be much less affected by dust emission than the NIR bands. Two of the fits are shown in Fig. 6 as examples of earlier and later epoch emission.

The resulting photospheric temperatures are listed in Table 4 and shown in Fig. 7a for AV $`=1.18`$ and $`2.5`$ mag. The temperatures for AV $`=2.50`$ mag are similar to those derived in the case of SN 1998S (Fassia et al. 2000); an extinction as large as AV $`=3.00`$ mag drives the photospheric temperature at early phases towards values which appear too high. Hence, AV cannot much exceed $`2.5`$ mag and AV $`=1.18`$ mag can be safely assumed as a lower limit.

We performed several tests in order to evaluate the quality of our fits. Most important, we repeated our procedure normalizing all fluxes to that in the $`V`$ band. The differences between these newly obtained values and those listed in Table 4 increase with AV, ranging from $`<6`$% for AV $`=1.18`$ mag, to $`<20`$% for AV $`=2.50`$ mag. The most sensitive variations are observed on JD+495 (temperatures from $`V`$-normalized data hotter by 23% for AV $`=2.50`$ mag than those from $`R`$-normalized data) and on JD+506 (with a similar trend). Here, two kinds of problems exist. First, when $`\lambda T>>hc/k`$ a $`\lambda _0`$-normalized blackbody flux tends to become equal to $`(\lambda _0/\lambda )^4`$ and no longer depends on temperature. This occurs for $`T>>30000`$ K in the $`V`$-band and $`T>>40000`$ K in the $`U`$-band. Then, our fits are not reliable for very high temperatures. Second, differences in the results obtained by varying the normalization band may be related to the general increase in temperature found when excluding NIR fluxes from the fits (amounting, e. g., to up to 60% for AV $`=2.50`$ mag), although the trend shown in Fig. 7a remains unchanged. This reflects physical differences in the regions which emit the bulk of radiation at different wavelengths (particularly, in the optical and in the near-infrared). As an example, for SN 1998S Anupama, Sivarani, & Pandey (2001) when fitting blackbody curves to optical spectra find photospheric temperatures up to 30 % hotter than those obtained by Fassia et al. (2000) from photometric data including NIR fluxes. We also repeated our fitting procedure excluding fluxes in the $`U`$ band, where we have very few measurements. This yields somewhat higher temperatures, but only at early epochs (20 % on JD+486 and 60% on JD+490 for AV $`=2.50`$ mag); thereafter, differences amount to only a few percent. A major source of uncertainty is line emission in the photometric bands. Particularly at later stages ($`>100`$ days from LCs maxima), lines affecting the photometry in the $`I`$ and $`R`$ bands may develop (see Danziger et al. 1988) causing temperature underestimates.

The errors listed in Table 4 (calculated as the square root of the quadratic sum of the temperature residuals in the various bands) are largely dominated by those involving the NIR fluxes. We have already discussed the indication of a NIR excess from $`JHK`$ colors. This is further demonstrated by plotting the blackbody temperatures obtained scaling the dereddened fluxes in the NIR to the ratio of observed to blackbody fluxes in the $`R`$-band (which is equal to the square of the ratio of photospheric radius to distance). As shown in Fig. 7b, higher NIR temperatures are derived at both early and late epochs. This may be reconciled if the NIR fluxes are not photospheric in origin but come from larger (cooler) regions starting soon after the explosion and continuing. Interestingly, an excess in the $`J`$ band is exhibited until JD+560 and is quite evident during the period (JD+)520–550 (see also Fig. 3); it corresponds to the noted blueing at early days. The $`K`$-band excess at late epochs is very strong, so much so that the flux has been left out when evaluating the error for JD+604, JD+690 and JD+698. Fluctuations in the derived photospheric temperatures are appreciable in Fig. 7a and may be related to the undersampling of LCs in some bands (e. g., $`U`$ and $`V`$). Unaccounted for line emission may also contribute. The apparent warming during the period (JD+)550–600 coincides with the fading of the SN towards brightness values similar to those of the brightest and closest star and it may indicate that the two sources could not be accurately deconvolved at this epoch.

In conclusion, the photospheric temperatures derived from blackbody fits between (JD+)486–690 are mostly affected by uncertainties related to the unknown extinction, line emission within the photometric bands and the degree of NIR excess. We believe that assuming AV $`=1.18`$ mag and AV $`=2.50`$ mag should result in lower and upper limits whose errors are otherwise not larger than those intrinsic to other similar studies (such as Fassia et al. 2000). Most significantly, the temporal evolution appears well defined.

### 3.5 Physical parameters of the SN photosphere

The photospheric radius was determined from the ratio of blackbody to observed fluxes in the $`R`$-band. We slightly revised the distance to the host galaxy, assuming $`26.3\pm 4.8`$ Mpc, from the distance modulus ($`32.1\pm 0.4`$ mag) we derived from our own observations of SN 2000E. Once the photospheric radius and temperature are known, the absolute bolometric luminosity can be estimated through the Stefan-Boltzmann law. A photospheric velocity was also derived as the ratio of radius to elapsed time from epoch 0. These data are listed in Table 5 in the cases AV $`=1.18`$ mag and AV $`=2.5`$ mag and shown in Fig 8. Quoted errors include the uncertainties in photospheric temperatures and measured fluxes, but not on the distance, since this would affect all data similarly. Note that increasing the photospheric temperature, e. g., by assuming larger extinctions, produces smaller radii and velocities; this confirms that AV cannot much exceed the galactic value and rules out the occurrence of very high photospheric temperatures.

Qualitatively, the photospheric radius, velocity and the bolometric luminosity exhibit the same behavior revealed by Fassia et al. (2000) for SN 1998S. The photospheric radius increases during the first 30 days, stays roughly constant until $`80`$ days from epoch 0 and then declines. The exact amount of the variations depends on the assumed extinction. The decrease in the photospheric velocity during the period 30–80 days is probably caused by the recession of the hydrogen recombination front through the expanding ejecta. However, an in-depth comparison is prevented by the uncertainty in the extinction of SN 1999el. Again, similarly to what was found for SN 1998S by Fassia et al. (2000), at around $`100`$ days after outburst the light curve (see Fig. 8c) slows to the point of approximating the radioactive decay of 56Co. As shown in Fig. 8c, $`0.07`$ M of 56Ni (based on the relation given by Branch (1992) or direct comparison with SN 1987A) could power the luminosity at later epochs for $`A_\mathrm{V}=2.5`$ mag and can be assumed as an upper limit to the mass of 56Ni produced by the explosion.

## 4 Discussion

The early evolution (within $`100`$ days) of the light curves of SN 1999el indicates that this is one of the most rapidly declining Type IIn SN among those studied so far. It is useful to compare its features with those of SN 1998S, the other recent and extensively observed Type IIn SN (Fassia et al. 2000). In the $`B`$ band, SN 1999el exhibits a luminosity decrease even faster than that of SN 1998S and is more similar to a Type II linear SN than to previously studied Type IIn SNe, such as SN 1998Z or SN 1997ab. The light curves at various wavelengths of SN 1999el, with their steep decline, are closer to those of SNe 1994aj (Benetti et al. 1999) and 1996L (Benetti et al. 1998) than to that of SN 1995G (Pastorello et al. 2002), whose spectrum is on the contrary very similar to that of SN 1999el, as noted in Sect. 3.1. Therefore, the indication of interaction with the CSM arising from the spectrum is not supported by the rapidly declining light curve. Consequently, type IIn SNe seem to embrace a wide range of characteristics. Therefore the differences between SN 1999el and SN 1998S are worth noting.

It appears that the SN 1998S LCs attain a maximum in the NIR some days after those in the optical bands (however, the precise times of occurrence for the optical maxima are uncertain since measurements appear to have begun when the optical light curve had already started declining). Conversely, in the case of SN 1999el the maximum emission in the NIR $`H`$ and $`K`$ bands occurs before that in the optical bands, whose maxima almost coincide with that in the NIR $`J`$ band. The exact time of maxima at $`H`$ and $`K`$ bands is uncertain but they seem to have occurred some time before the beginning of our observations. The differences may be appreciated by looking at the color curves shown in Fig. 4. As previously noted, the evolution of $`JH`$ and $`JK`$ for SN 1999el exhibits an early appearance of a NIR excess (around day 5) and a NIR blueing effect in the first $`45`$ days. This seems observationally well established, since here the NIR LCs sample the SN when it was much brighter than the nearby sources c2 and c3 (see Fig. 1), and thus are easily measurable even though sometimes only partially resolved on SWIRCAM images. Although maxima at NIR wavelengths do appear to have been missed for SN 1998S, nevertheless there is no evidence for an early blueing effect in this case. Furthermore, a NIR excess from SN 1998S is not apparent at early epochs.

One possible origin of the quoted differences, in particular of the early NIR excess may be due to an IR echo from pre-existing dust in the vicinity of SN 1999el not evaporated at the moment of the shock outburst. This would suggest that dust clumps near SN 1999el are distributed differently than those near SN 1998S which may be devoid of nearby dust entirely. For SN 1998S this is supported by Leonard et al. (2000) whose spectroscopic observations suggest that the progenitor of SN 1998S underwent a major mass-loss episode that ended some 60 years before the explosion, followed by a weaker scarcely detectable mass-loss activity in the last seven years before the explosion. If a mass-loss wind associated with the SN 1999el progenitor had not ceased at the time of explosion, the density of the CSM around it would be much larger than in the case of SN 1998S.

In the above scenario, the interaction of the UV flash with the CSM around SN 1999el would have led to heating, destruction and sublimation of dust. The extent of the region where grain sublimation and destruction can occur may be estimated as (see, e. g., Waxman & Draine 2000)

$$R_0=3.7\times 10^{19}(Q_{\mathrm{abs}}L_{49}a_{01}^1)^{\frac{1}{2}}\mathrm{cm}$$
(1)
where $`Q_{\mathrm{abs}}`$ is the grain absorption efficiency at UV wavelengths, $`L_{49}`$ the maximum luminosity in the range 1–7.5 eV in units of $`10^{49}`$ erg s-1 (optical-UV flash) and $`a_{01}`$ the mean size of grains in units of $`0.1`$ $`\mu `$m. Assuming $`Q_{\mathrm{abs}}1`$ and $`a_{01}1`$, $`R_010^{17}`$ cm is obtained for a maximum luminosity of $`10^{44}`$ erg s-1 (the greatest value indicated in Table 5), a distance which is two orders of magnitude larger than the largest photospheric radius attained by the supernova (see Table 5), but also somewhat greater than the estimate based on the occurrence of the NIR excess very early (see Sect. 3.2). Since the grain distruction and sublimation due to the UV flash imply also the unlocking of Fe ions (if they exist in the dust), we have been tempted to ascribe the early 45 days blueing to \[FeII\] line emission. The \[FeII\] lines which may affect $`J`$ and $`H`$ bands result from the transitions involving metastable levels $`a^6`$D9/2$`a^4`$D7/2 ($`\lambda =1.26`$ $`\mu `$m, A$`=5.6\times 10^3`$ s-1) and $`a^4`$F9/2$`a^6`$D7/2 ($`\lambda =1.64`$ $`\mu `$m, A$`=1.9\times 10^3`$ s-1), collisionally excited by free electrons (see Mouri, Kawara, & Taniguchi 2000, and references therein). Very large line integrated fluxes (up to 0.5–1 times the continuum flux integrated on the $`K`$ band; see Fig. 5) would be required to produce the observed changes in colors, indicating that \[FeII\] lines alone cannot account for it. In addition, since relative luminosity of the main \[FeII\] lines in the J and H windows is constant, the contribution of this mechanism to the color evolution is definitively ruled out.

If the early IR excess of SN 1999el is due to the presence of a CSM generated by a continuous wind of the progenitor star, it might be reasonable to expect its signature in the optical light curves. This seems not to be the case. Roscherr and Schaefer (2000) analysed theoretically the modifications of optical LCs supernovae produced by echo processes. They assume spherical symmetry for the dust and find that critical parameters are $`r_{min}`$, $`r_{max}`$ and $`\tau `$, where $`r_{min}`$ and $`r_{max}`$ are the inner and outer radii of the dust shell and $`\tau `$ the optical depth in the $`R`$ band for a purely radial photon trajectory, a quantity directly related with the total mass of dust. Roscherr and Schaefer (2000) find that, in the case of thick shells, $`\tau `$ is the more sensitive parameter. Changes in $`\tau `$ effect changes in the luminosity of light maxima in various bands whose amounts depend on the wavelength range. Even the time at which the light maximum is attained may be shifted (by up to several days; Roscherr, private communication) according to its wavelength. Unfortunately, these computations do not address the NIR bands. The decline rate is also strongly affected by changes in $`\tau `$ at levels once again dependent on wavelength. In this respect a relevant parameter is the $`\beta _{100,B}`$ index which corresponds to the luminosity drop in 100 days in the $`B`$ band. Therefore the $`\beta _{100,B}`$ index can provide information on the presence of a CSM around the SN. Observational data provide a value of $`\beta _{100,B}`$ greater than 4 mag for both SN 1999el and 1998S. A comparison with the theoretical results quoted above suggests the surprising result that the expected amount of dust surrounding both SN 1999el and SN 1998S should be vanishingly small. (Or that almost all the pre-existing dust is evaporated at the moment of the UV flash).

For SN 1999el this result is in contrast with the suggestion that the early IR excess is due to an IR echo caused by near-by dust.

We suggest that the failure of the $`\beta _{100,B}`$ index to display the presence of CSM around SN 1999el may be ascribed to the assumption made to derive theoretically this index, namely the assumption of spherically distributed dust.

The theoretical results of Emmering and Chevalier (1988) on the IR echo light originated by an asymmetric CSM distribution around a SN progenitor clearly illustrate how many parameters influence the observed results. Indeed they find that IR light curves can substantially change by varying the geometrical properties of a given non-spherical (but still axi-symmetric) CSM distribution and/or its orientation relative to the line of sight. Even if not explicitly accounted for in their work it can be nevertheless deduced that additional parameters, such as a non-radial density profile of the dust, play an additional pivotal role.

We intend in the future to provide a tentative model of the dust distribution around SN 1999el . Finally we note that both for SN 1998S (Fassia et al. 2000) and SN 1999el the observations do not rule out a scenario of dust formation after the explosion. In the case of SN 1999el an occurrence like this would imply that the signature of dust appears both at early and also at late stages. In particular, for SN 1999el assuming AV $`=1.18`$ mag, the dereddened $`K`$ magnitude on JD+691 (see Table 2) yields an observed flux of $`1.2\times 10^{13}`$ erg cm-2 s-1 $`\mu `$m-1. As in Sect. 3.2, we can determine the radius of a blackbody-emitting surface for a temperature of 1500 K, which is roughly the sublimation temperature for dust grains. This amounts to $`3\times 10^{15}`$ cm, of the order of the largest photospheric radius attained (see Table 5), a distance that might be covered by dust formed in the ejecta moving with a velocity $`1600`$ km s-1.

Therefore the scenario we propose embraces the evolution of the progenitor star in that cold, high-luminosity, mass loss from evolved massive stars may proceed in episodes with asymmetric distributions. It is however conceivable that single mass-loss episodes are not completely stochastic in duration, periodicity and intensity. Over long temporal intervals, possibly covering several cycles of activity, it is possible to define mean wind properties, such as, for example, the dependence of the mean mass-loss rate on the main structural parameters such as temperature and luminosity. If the SN explosion occurs during a long pause in mass-loss we observe a normal Type II SN, otherwise the typical features of the Type IIn SNe will appear with varying degrees of visibility.

## 5 Conclusions

Optical and NIR light curves have been presented for the Type IIn SN 1999el. These data, in comparison with those of SN 1998S, suggest that the progenitor star of SN 1999el was still undergoing a probably asymmetric mass loss episode at the moment of the explosion. Dust was already present in the expanding circumstellar envelope. Dust lying in the inner $`10^{16}`$ cm evaporated at the outburst.

Later, on, the evolution of NIR colors seems to be affected by the echo light from a pre-existing asymmetrically expanding envelope. In the absence of mass loss SN 1999el would have been presumably observed as a bona-fide linear Type II SN. The way in which mass loss occurs in the last evolutionary phases of massive stars seems to be responsible for various observed departures of Type II SNe from the standard II-L and II-P behavior.

This paper is partially based on observations made with the Italian Telescopio Nazionale Galileo (TNG) operated on the island of La Palma by the Centro Galileo Galilei of the CNAA (Consorzio Nazionale per l’Astronomia e l’Astrofisica) at the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofisica de Canarias. The infrared data collected with the AZT-24 are products of the Supernova Watch-dogging InfraRed Telescope (SWIRT), a joint project of the Astronomical Observatories of Collurania-Teramo (Italy), Pulkovo (Russia) and Rome (Italy). We thank F. Ghinassi and F. Mannucci for acquiring the latest NIR images with NICS. We warmly thank Bruce Roscherr for advice and suggestions about the echo light issue.


format: text

text:
# Abstract

## Abstract

Magnetization reversal in magnetic nanostructures is investigated numerically over time-scales ranging from fast switching processes on a picosecond scale to thermally activated reversal on a microsecond time-scale. A simulation of the stochastic Landau-Lifshitz equation of motion is used as well as a time quantified Monte Carlo method for the simulation of classical spin systems modeling magnetic Co nanoparticles. For field pulses larger than the Stoner-Wohlfarth limit spin precession effects govern the reversal behavior of the particle while for lower fields a magnetization reversal is only possible when it is assisted by thermal fluctuations.

## 1 Introduction

The miniaturization of magnetic structures plays an important role for fundamental research as well as for technical applications. This leads to an incremental interest in the understanding of the behavior of small magnetic particles and structures down to the nanometer scale. But with decreasing size of the magnetic system thermal activation becomes relevant. Hence, much effort is focused now on the understanding of magnetization dynamics at finite temperatures. .

In the following, we briefly describe numerical techniques for the study of magnetization dynamics in nanostructures, modeled as classical spin systems, where a finite temperature is taken into account. First, we will focus on the underlying model and the two basic methods, Langevin dynamics  and Monte Carlo simulation . Then we study the magnetization reversal in Co nanoparticles. Starting with the deterministic spin dynamics on short time scales which plays a crucial role in high-speed data storage , we go on to the probabilistic long-time behavior where a thermally assisted reversal can occur even for magnetic fields below the coercive field .

## 2 Classical Spin Model

The micromagnetic properties of a nanoparticle can be described using a model of classical magnetic moments which are localized on a given lattice. Such a spin model can be motivated following different lines: on the one hand it is the classical limit of a quantum mechanical, localized spin model — the Heisenberg model . On the other hand, a classical spin model can also be interpreted as the discretized version of a micromagnetic continuum model , where the charge distribution for a single cell of the discretized lattice is approximated by a point dipole .

The interpretation as an atomic model restricts the use of computer simulations to the investigation of rather small systems of only a few million atoms - corresponding to particle sizes of only a few nanometers. On the other hand, within a continuum model, the space might be discretized on a much larger length scale, as compared to an atomic distance. However, in continuum theory usually a constant absolute value of the magnetization vector is assumed, an assumption which fails for higher temperatures since the space averaged magnetization breaks down when approaching the critical temperature. Hence, one expects correct thermal properties only in the limit of small cell sizes of the order of atomic distances.

In the following, let us consider a classical three dimensional Heisenberg Hamiltonian for localized spins,

$``$ $`=`$ $`J{\displaystyle \underset{ij}{}}𝐒_i𝐒_j\mu _s𝐁{\displaystyle \underset{i}{}}𝐒_id_z{\displaystyle \underset{i}{}}(S_i^z)^2`$ (1)
$`w{\displaystyle \underset{i<j}{}}{\displaystyle \frac{3(𝐒_i𝐞_{ij})(𝐞_{ij}𝐒_j)𝐒_i𝐒_j}{r_{ij}^3}},`$
where the $`𝐒_i=𝝁_i/\mu _s`$ are three dimensional magnetic moments of unit length representing atomic magnetic moments. The first sum is the ferromagnetic exchange of the moments with the coupling constant $`J`$. The second sum is the coupling of the magnetic moments to an external magnetic field $`B`$, the third sum represents a uniaxial anisotropy, here, for $`d_z>0`$ favoring the $`z`$ axis as easy axis of the system, and the last sum is the dipolar interaction where $`w=\mu _0\mu _s^2/(4\pi a^3)`$ describes the strength of the dipole-dipole interaction. The $`𝐞_{ij}`$ are unit vectors pointing from lattice site $`i`$ to $`j`$ and $`r_{ij}`$ is the distance between these lattice sites in units of $`a`$. The dipole-dipole interaction can be computed efficiently with the help of fast Fourier transformation methods . One should however note that in a Monte Carlo simulation with a single-spin flip algorithm the FFT method is an approximation the implementation of which was described in details before .

The transformation of the above introduced atomic parameters to the material parameters which are usually used in a continuum model is given by $`J=2aA_x`$ where $`A_x`$ is the exchange energy, $`d_z=Ka^3`$ where $`K`$ is the anisotropy energy density, and $`\mu _s=M_sa^3`$ where $`M_s`$ is the spontaneous magnetization.

## 3 Landau-Lifshitz Equation and Langevin Dynamics

In the short time limit spin precession is important which can be taken care of by studying the corresponding equation of motion. The basic numerical approach which includes thermal activation, is the direct numerical integration of the Langevin equation of the problem. In order to obtain thermal averages one has to calculate many of these trajectories starting with the same initial conditions taking an average over these trajectories for the quantities of interest. This method is referred to as the Langevin dynamics formalism .

The underlying equation of motion for a magnetic system is the Landau-Lifshitz-Gilbert (LLG) equation,

$$\frac{𝐒_i}{t}=\frac{\gamma }{(1+\alpha ^2)\mu _s}𝐒_i\times \left[𝐇_i(t)+\alpha \left(𝐒_i\times 𝐇_i(t)\right)\right],$$
(2)
with the gyromagnetic ratio $`\gamma =1.76\times 10^{11}(\mathrm{T}s)^1`$, the dimensionless damping constant $`\alpha `$, and the internal field $`𝐇_i(t)=𝜻_i(t)/𝐒_i`$. Langevin dynamics is introduced here in form of the noise $`𝜻_i(t)`$ which represents thermal fluctuations, with $`𝜻_i(t)=0`$ and $`\zeta _i^\eta (t)\zeta _j^\theta (t^{})=2\delta _{ij}\delta _{\eta \theta }\delta (tt^{})\alpha k_\mathrm{B}T\mu _s/\gamma `$ where $`i,j`$ denote once again lattice sites and $`\eta ,\theta `$ Cartesian components.

The LLG equation with Langevin dynamics is a stochastic differential equation with multiplicative noise. For this kind of differential equation a problem arises which is called the Itô-Stratonovich dilemma . As a consequence, different time discretization schemes may converge to different results with decreasing time step. As was pointed out in  the multiplicative noise in the Langevin equation above has to be treated by means of the Stratonovich interpretation. An appropriate discretization scheme leading to a Stratonovich interpretation is the Heun method  which is used in the following.

## 4 Monte Carlo Methods

In the long time limit only spin relaxation and thermal fluctuations are relevant which can be studied very conveniently using Monte Carlo methods with quantified time step . Within a Monte Carlo approach  trajectories in phase space are calculated following a master equation  for the time development of the probability distribution $`P_s(t)`$ in phase space,

$$\frac{\text{d}P_s}{\text{d}t}=\underset{s^{}}{}(P_s^{}w_{s^{}s}P_sw_{ss^{}}).$$
(3)
Here, $`s`$ and $`s^{}`$ denote different states of the system and $`w_{s^{}s}`$ is the transition rate for a change from a state $`s^{}`$ to a state $`s`$. These rates have to fulfill the condition

$$\frac{w_{ss^{}}}{w_{s^{}s}}=\mathrm{exp}\left[\frac{E(S)E(S^{})}{k_\mathrm{B}T}\right].$$
(4)
The master equation describes exclusively the coupling of the system to the heat bath . Hence, only the irreversible part of the dynamics of the system is considered including only the relaxation and the fluctuations. A Monte Carlo simulation does not include the energy conserving part of the equation of motion. Hence, no precession of magnetic moments will be found.

Monte Carlo approaches in general have no physical time associated with each step of the algorithm. However, recently a time quantified Monte Carlo method was proposed in  and later successfully applied to different model systems . Here, the interpretation of a Monte Carlo step as a realistic time interval $`\mathrm{\Delta }t`$ was achieved by a comparison of one step of the Monte Carlo process with a time interval of the LLG equation in the high damping limit. We will use this algorithm in the following. The trial step of this algorithm is a random movement of the magnetic moment within a cone with a given size $`r`$ with

$$r^2=\frac{20k_BT\alpha \gamma }{(1+\alpha ^2)\mu _s}\mathrm{\Delta }t.$$
(5)
In order to achieve this efficiently one constructs a random vector with constant probability distribution within a sphere of radius $`r`$. This random vector is added to the initial moment and subsequently the resulting vector is normalized .

Using this algorithm one Monte Carlo step represents a given time interval $`\mathrm{\Delta }t`$ of the LLG equation in the high damping limit as long as $`\mathrm{\Delta }t`$ is chosen appropriately (for details see ).

## 5 Precessional Reversal in Co Nanoparticles

In the following we consider Co nanoparticles, where the material parameters are $`A_x=1.310^{11}`$ J/m, $`K=6.810^5`$ J/m2 and $`M_s=1.4310^6`$ A/m. For simplicity we simulate a simple cubic lattice with atomic distance $`a=0.25`$ nm. Our simulation starts with a spin configuration where all magnetic moments point into the $`z`$ direction, aligned with the easy axis, and with the $`z`$ component of the external magnetic field antiparallel to the magnetization so that the system is in an unstable, or at least metastable state.

In sufficiently small particles the magnetic moments rotate coherently during the magnetization reversal. A quantitative description of coherent rotation in ellipsoidal single domain particles was developed by Stoner and Wohlfarth . Depending on the angle between the applied field $`𝐁`$ and the $`z`$ (easy) axis of the system, the coercive field $`B_c`$ varies following the so-called Stoner-Wohlfarth asteroid . Under an angle of $`45^{}`$ it is $`B_c=d_{\mathrm{eff}}V/\mu _\mathrm{s}`$ where $`V`$ is the volume of the particle and $`d_{\mathrm{eff}}`$ is an effective anisotropy constant. In our case it is $`B_c0.7`$ T.

First, we are interested in fast switching processes where the applied field is higher than the coercive field and the reversal is dominated by spin precession. We simulate the system with Langevin dynamics as described before. Figure 1 shows the time dependence of the magnetization of our Co particle in the low damping limit ($`\alpha =0.1`$) for a very small ellipsoidal shaped particle of length $`L=4`$nm and diameter $`L=2`$nm.

The $`z`$ and the $`x`$ components of the magnetization are shown as well as its absolute value. The magnetic field $`𝐁`$ is set under an angle of $`45^{}`$ to the $`z`$ axis within the $`yz`$ plane so that the response of the system to the external field sets in directly. The wavering magnetization of the system clearly follows from the precession of the spins.

Note, that the precession time of our system is not simply given by the precession time of a single spin in an external field ($`\tau _p=2\pi (1+\alpha ^2)/\gamma B32`$ ps in our case). Instead, the whole internal field is relevant for the spin precession, also the contributions from the dipolar field, the exchange and the anisotropy. This internal field is not constant in time and in general it is non-homogenous within the system. However, in the very small particle which we consider here, the internal field is sufficiently homogenous so that the magnetization moves coherently. Hence, the absolute value of the magnetization remains constant in time as is shown in the figure. Note also, that even after the new stable state is reached the magnetization still keeps on oscillating around its equilibrium value, driven by thermal fluctuations.

## 6 Thermally Activated Reversal

Let us now turn to the case $`B<B_c`$, where the reversal process can only occur when it is thermally activated. Since we are now interested in the long time and high damping limit ($`\alpha =4`$) where the behavior of the particle is governed by thermal fluctuations we can use Monte Carlo methods. Figure 2 shows the typical time dependence of the magnetization of the same Co particle as before. The field $`𝐁`$ is set antiparallel to the initial state here, so that the zero-temperature coercive field is given by $`B_c=2d_{\mathrm{eff}}V/\mu _\mathrm{s}`$ which in our case is $`B_c1.4`$ T.

As one can see, the $`z`$ component of the magnetization remains nearly constant for a time period which is rather long as compared to the previous simulation. Then, suddenly, the magnetization drops and its sign changes. As before the absolute value of the magnetization is constant in time and the reversal mechanism is mainly a coherent rotation. The value of the switching time is approximately 6.6 ns in our simulation. However, this thermally activated switching is not a deterministic process as it was the short time dynamics studied before, where the switching followed mainly from the deterministic part of the equation of motion. Instead, the thermal activation process here is a probabilistic event. The probability distribution $`P(t_s)`$ for switching events taking place after a time $`t_s`$ follows an exponential law ,

$$P(t_s)\mathrm{exp}(t_\mathrm{s}/\tau ),$$
(6)
in the limit of large time-scales. Here, $`\tau `$ is a characteristic time scale,

$$\tau =\tau _0\mathrm{exp}\left(\mathrm{\Delta }E/k_\mathrm{B}T\right),$$
(7)
where $`\tau _0`$ is a prefactor and $`\mathrm{\Delta }E`$ an energy barrier which both are related to a certain reversal mechanism (see e. g.  for analytically determined prefactors and energy barriers in different systems and  for numerical work on this subject). In general the prefactor may depend on the system parameters, the temperature, the applied magnetic field and the damping constant.

For the case of a Stoner-Wohlfarth particle with the applied field parallel to the easy axis the energy barrier has the form

$$\mathrm{\Delta }E=d_{\mathrm{eff}}V\left(1\frac{B}{B_c}\right)^2.$$
(8)
This energy barrier as well as the prefactor have been calculated by Brown  under the assumption that all magnetic moments are parallel, so that the system behaves like one single magnetic moment.

For a further analysis we extract the energy barrier which governs the reversal process from our numerical data. Figure 3 shows the temperature dependence of the characteristic time, i. e. the mean switching time, obtained from our simulations. The slope of the solid line corresponds to the theoretical value of the energy barrier for a reversal by coherent rotation obtained from Eq. 8. Obviously, it is in very good agreement with our numerical data for low enough temperatures.

## Acknowledgments

This work was supported by the Deutsche Forschungsgemeinschaft (SFB 491 and project NO290/1) and by the European Union (COST action P3, working group 4).


text:
# Computational Geometry Column 43

## Abstract

The concept of pointed pseudo-triangulations is defined and a few of its applications described.

A *pseudo-triangle* is a planar polygon with exactly three convex vertices. Each pair of convex vertices is connected by a reflex chain, which may be just one segment. Thus, a triangle is a pseudo-triangle. A *pseudo-triangulation* of a set $`S`$ of $`n`$ points in the plane is a partition of the convex hull of $`S`$ into pseudo-triangles using $`S`$ as a vertex set. Under the name of “geodesic triangulations,” these found use in ray shooting \[CEG+94\], and then for visibility algorithms by Pocchiola and Vegter \[PV95\], who named them (after a dual relationship to pseudoline arrangements) and studied many of their properties \[PV96\]. Recently the identification of *minimum pseudo-triangulations* by Streinu \[Str00\] has generated a flurry of new applications that we selectively sample here.

Minimum pseudo-triangulations have the fewest possible number of edges for a given set $`S`$ of points. The term *pointed pseudo-triangulation* is gaining prominence because of the following theorem111 Thm. 3.1(2) of \[Str00\] (where the property is called “acyclicity”); Lem. 1 in \[KKM+01\].  (see Fig. 1):

###### Theorem 1

A pseudo-triangulation of $`S`$ is minimum iff every point $`pS`$ is *pointed* in the sense that its incident edges span less than $`\pi `$ (i.e., they fall within a cone apexed at $`p`$ with aperture angle $`<\pi `$).

Pointed pseudo-triangulations have remarkable regularity in that each contains exactly $`n2`$ pseudo-triangles. In this they parallel triangulations of simple polygons and are unlike triangulations of point sets, which may have from $`n2`$ to $`3(n2)`$ triangles. Another indication of their well-behavedness is the recent result \[KKM+01\] that every set of $`n`$ points in general position has a pointed pseudo-triangulation with maximum vertex degree of $`5`$, a property not enjoyed by triangulations (of either polygons or point sets). This bounded degree may prove useful for designing data structures employing pseudo-triangulations, such as those used for collision detection \[KSS02\].

We mention two disparate applications of pointed pseudo-triangulations. Two years ago a long-standing open problem was settled by Connelly, Demaine, and Rote when they proved that polygonal linkages cannot lock in the plane: open polygonal chains can be straightened, and closed polygonal chains can be convexified \[CDR00\] \[O’R00\]. Their unlocking motions are *expansive*: no distance between any pair of vertices decreases. Their proof relies on proving the existence of a solution to differential equations describing this motion. Although not difficult to solve numerically in practice,222 http://db.uwaterloo.ca/~eddemain/linkage/animations/.  theirs does not constitute an algebraic proof. This was provided by Streinu \[Str00\]. She proved that a pointed pseudo-triangulation is rigid when viewed as a bar and joint framework, and that removal of a convex hull edge yields a mechanism with one degree of freedom (1-dof) that executes an expansive motion. Adding bars to a given collection of polygonal chains to produce a pointed pseudo-triangulation, and using a useful “flipping” property of pseudo-triangulations, leads ultimately to a piecewise-algebraic straightening/convexifying motion. These 1-dof expansive motions have a rich structure that is just now being elucidated \[RSS01\].

Second, pointed pseudo-triangulations were used to make an advance on an open problem posed by Urrutia \[Urr00\]: Can a simple polygon of $`n`$ vertices be illuminated by $`cn`$ interior $`\pi `$-floodlights placed at vertices, with $`c<1`$ (for sufficiently large $`n`$)? For several years the lower bound on $`c`$ was $`3/5`$, but no one could prove that $`c`$ was smaller than $`1`$. Starting from a pointed pseudo-triangulation of the polygon, Speckmann and Tóth have recently established \[ST01\] that $`(2/3)n`$ lights suffice, finally breaking the $`c=1`$ barrier.333 See http://cs.smith.edu/~orourke/TOPP/, Problem 23, for further developments.

We close with an open problem, posed in \[RRSS01\]. For a given set $`S`$ of points, let $`\mathrm{\#}T`$ be the number of triangulations of $`S`$, and $`\mathrm{\#}PPT`$ be the number of pointed pseudo-triangulations of $`S`$. They conjecture that $`\mathrm{\#}T\mathrm{\#}PPT`$ for all $`S`$, with equality only when the points of $`S`$ are in convex position. This has been established for all sets of at most $`9`$ points \[BKPS01\].


text:
# Scale-dependent rigidity of polymer-ornamented membranes

## Abstract

We study the fluctuation spectrum of fluid membranes carrying grafted polymers. Contrary to usual descriptions, we find that the modifications induced by the polymers cannot be reduced to the renormalization of the membrane bending rigidity. Instead we show that the ornamented membrane exhibits a scale-dependent elastic modulus that we evaluate. In ornamented lamellar stacks, we further show that this leads to a modification of the Caillé parameter characterizing the power-law singularities of the Bragg peaks.

The characterization of polymer-membrane interactions is a fundamental issue of colloidal science. In cosmetics, pharmaceutics or detergency, many formulations are suspensions of self-assembled surfactant bilayers with polymers added for performance, processing, conditioning or delivery vandepas94 . Likewise, lipid bilayers form the walls of living cells and liposomes, and host a great variety of macromolecules for coating protection, ion exchange and mechanical reinforcement. albertsbook . In many instances, the polymers are end-tethered to the soft interfaces. Grafting is easily achieved experimentally by using polymer chains that carry hydrophobic groups. Typical examples are provided by the so-called PEG-lipids, hydrophilic chains of polyethylene glycol covalently linked to a double-tail phospholipid molecule. Recent studies have shown that polymers grafted to bilayers can induce gelation warriner96  or other phase changes yang98 ; castro99  in liquid lamellar phases. They stabilize monodisperse vesicles joannic97  and modify the geometry of monolamellar ringsdorf91 ; cates91  and multilamellar frette99  cylindrical vesicles. They also lead to drastic changes in the structure and phase behaviour of ternary amphiphilic systems endo00 .

As first explained by Canham and Helfrich, fluid membranes are fluctuating objects seifert97 . In thermal equilibrium, the surface assumes all the possible shapes allowed by the geometry with the associated Boltzmann probability $`𝒫\mathrm{exp}\{_0/(k_BT)\}`$. The Hamiltonian $`_0`$ is a quadratic function of $`c_1`$ and $`c_2`$, the two principal curvatures at any given point of the surface,

$$_0=𝑑S\left\{\frac{\kappa }{2}(c_1+c_22c_0)^2+\overline{\kappa }c_1c_2\right\}.$$
(1)
The spontaneous curvature $`c_0`$ vanishes for symetric bilayers. The Gaussian rigidity $`\overline{\kappa }`$ plays an important role in the determination of the topology of the system, and the integral $`𝑑S\overline{\kappa }c_1c_2`$ does not depend on the particular shape of the membrane. The amplitude of the thermal fluctuations is controled by the bending rigidity $`\kappa `$, that ranges from a few $`k_BT`$ to a few tens of $`k_BT`$. Conversely, the analysis of the height correlations in a membrane system allows for the determination of the constitutive rigidity. Many techniques were specifically developped to extract $`\kappa `$ by comparing experiments and theoretical predictions for vesicles, lamellar stacks, bicontinuous phases and other geometries. These methods have so far been of limited application in polymer-membrane systems because of the lack of prediction for the fluctuation spectrum of bilayers in the presence of polymers. In this Letter, we make a first step to reach this gap by investigating the role of grafted chains on the fluctuation spectrum of fluid membranes. We first analyze the subtle interplay between local monomer concentration and surface curvature. We show that the chains induce additional height correlations that can be interpreted as the result of effective, scale-dependent elastic coefficients. For the specific case of an ornamented lamellar phase of membranes, we also determine the Caillé parameter that characterizes the power-law singularities around the Bragg peaks.

In a mean-field description, ornamented membranes are assumed to have new parameters $`\kappa _{\mathrm{𝑒𝑓𝑓}}=\kappa +\mathrm{\Delta }\kappa `$ and $`\overline{\kappa }_{\mathrm{𝑒𝑓𝑓}}=\overline{\kappa }+\mathrm{\Delta }\overline{\kappa }`$ different from those of the bare membrane. We consider $`𝒩`$ polymers grafted to each side of a membrane of area $`𝒮`$, the chains being free to diffuse along the membrane. The surface coverage $`\sigma =𝒩/𝒮`$ is low enough to avoid polymer/polymer interactions: $`\sigma <\sigma ^{}R_F^2`$, with $`R_F`$ the Flory radius of the chains. In this limit of small density, the excess surface energy is proportional to $`\sigma `$ and one anticipates corrections of the order of $`\mathrm{\Delta }\kappa \mathrm{\Delta }\overline{\kappa }k_BT(\sigma R_F^2)`$ in both good and theta solvent. The effect of solvent quality and polymer architecture resides in the numerical prefactors of the variations. For Gaussian chains tethered by one extremity, the calculations can be carried out explicitly and one gets $`\mathrm{\Delta }\kappa =\left(1+\frac{\pi }{2}\right)k_BT\sigma R_g^2`$ and $`\mathrm{\Delta }\overline{\kappa }=2k_BT\sigma R_g^2`$ hiergeist96 ; marquesfournier , with $`R_g=(Na^2/6)^{1/2}`$ the radius of gyration of a Gaussian chain of $`N`$ monomers of size $`a`$. At the overlap concentration $`\sigma ^{}`$, the variations are of the order of the thermal energy $`k_BT`$. However, for most practical situations, the chains can freely move along the membrane surface: the degrees of freedom associated with the end positions are annealed, as opposed to the quenched results of references hiergeist96 ; marquesfournier . A thorough study of the system reveals that averaging the free energy over the anchor position is not equivalent to averaging the partition function itself. As a consequence, the response to bending is lowered if one correctly accounts for the tendency of the chains to be at the top of the membrane ondulation, leading to bickelprep

$`\mathrm{\Delta }\kappa `$ $`=`$ $`k_BT\sigma R_g^2,`$ (2)
$`\mathrm{\Delta }\overline{\kappa }`$ $`=`$ $`2k_BT\sigma R_g^2.`$ (3)
The reduction of the effective modulus is due to the linear coupling between positions and spontaneous curvature and actually follows the general argument first recognized by Leibler leibler86 .

These mean-field results miss however a crucial feature, namely that the polymer-membrane interactions strongly depend on the length scale. We first consider the partition function $`𝒵_N`$ of a single Gaussian chain and construct a perturbative expansion in the limit of small surface deformations around the flat plane. In absence of overhangs and inlets, the surface profile can be described by a single-valued function $`h(𝐫)`$, where $`𝐫`$ spans the reference plane. We expand $`𝒵_N`$ in a series $`𝒵_N^{(0)}+𝒵_N^{(1)}+𝒵_N^{(2)}+\mathrm{}`$, where $`𝒵_N^{(i)}`$ is of order $`h^i`$ and $`𝒵_N^{(0)}(z)z/(\sqrt{\pi }R_g)`$ denotes the partition function of a chain tethered at a small distance $`z`$ of a flat, impenetrable wall. It has been shown recently that the first-order term of the expansion is related to the pressure field applied by the grafted chain to the surface breidenich00 ; bickel00 . At distance $`r=\sqrt{x^2+y^2}`$ from the grafting point, this entropic pressure decays with the scaling form $`p(r)k_BTr^3`$ ($`arR_g`$) in both good and theta solvents, and then vanishes sharply beyond $`R_g`$ bickel01 .

We now evaluate the contribution of the tethered polymers to the fluctuation spectrum of the membrane. If one introduces the Fourier Transform (FT) of the height profile $`\stackrel{~}{h}(𝐪)=𝑑𝐫e^{i\mathrm{𝐪𝐫}}h(𝐫)`$, then the energy $`_0`$ of the pure membrane reads

$$_0=\frac{1}{2}\frac{d𝐪}{(2\pi )^2}\stackrel{~}{h}(𝐪)\stackrel{~}{h}(𝐪)\kappa q^4.$$
(4)
Contrary to the previous mean-field results, we assert here that the presence of the chains entails non-local contributions to the elastic description of the membrane. This can be shown by extending the perturbative scheme up to second order in $`h`$, a calculation that we perform analytically for a Gaussian chain. In this case, the single-chain partition function obeys the Edwards equation doibook

$$\left(\frac{}{n}\frac{a^2}{6}^2\right)𝒵_n(𝐫,z)=0,$$
(5)
where $`(𝐫,z)`$ is the position of the constrained extremity. The linearity of eq. (5) implies that each term of the development obeys an Edwards equation. The solution of the successive orders are coupled through the (Dirichlet) boundary condition on the wall $`𝒵_N(𝐫,h(𝐫))=0`$. The solution of eq. (5) obeying the impenetrability condition is then determined recursively, and is given to each order by bartonbook

$$𝒵_N^{(i)}(𝐫,z)=\frac{a^2}{6}_0^N𝑑n𝑑𝐫^{}\frac{G_{Nn}^{(0)}}{z^{}}(𝐫,z;𝐫^{},0)𝒵_n^{(i)}(𝐫^{},0),$$
(6)
where the Green function $`G_N^{(0)}`$ of eq. (5) satisfies the Dirichlet condition on the horizontal plane bickel01 . This solution corresponds to a “loop expansion” of $`𝒵_N`$: the term of order $`i`$ enumerates the conformations that encounter the surface $`i`$ times, as depicted on figure 1. Finally, the statistical weight $`𝒵_N`$ is averaged for all possible positions of the anchoring point along the surface.

In general, the membrane carries a finite polymer concentration $`\sigma `$. We neglect for sake of simplicity the interactions between the chains, whose origin is twofold: in addition to the usual steric repulsions, the mean deformation resulting from the polymer pressure gives rise to attractive, deformation-sharing interactions bickel01 . However, the latter do not modify the fluctuation spectrum of the membrane. Our results are thus exact in the low-density limit $`\sigma R_g^21`$ where direct polymer interactions can be neglected. Assuming $`𝒩`$ chains on each side of the membrane, the total partition function of the system is given by the path integral

$`𝒵`$ $`=`$ $`{\displaystyle \frac{1}{(𝒩!)^2}}{\displaystyle 𝒟[h]e^{\beta _0}\left(\frac{1}{𝒮}𝑑𝐫𝒵_N\right)^{2𝒩}}`$ (7)
$`=`$ $`{\displaystyle \frac{1}{(𝒩!)^2}}{\displaystyle 𝒟[h]e^{\beta _{eff}}},`$
where the factor $`1/𝒩!`$ reflects the fact that particles are indistinguishable. The properties of the decorated surface are described by the new Hamiltonian $`_{eff}=_0+\mathrm{\Delta }`$, with

$$\mathrm{\Delta }=\frac{1}{2}\frac{d𝐪}{(2\pi )^2}\stackrel{~}{h}(𝐪)\stackrel{~}{h}(𝐪)G(𝐪).$$
(8)
The algebra will be detailed elsewhere bickelprep , and we now discuss the results. Integrating out the degrees of freedom of the grafted chains generates the additional contribution

$$G(𝐪)=\frac{2k_BT\sigma }{R_g^2}\left\{e^{q^2R_g^2}1+q^2R_g^2\right\}.$$
(9)
The correlation function $`G(𝐪)`$ is displayed on figure 2. At length scales larger than the polymer size $`qR_g1`$, one can extract a correction to the bending rigidity that reads $`\mathrm{\Delta }\kappa =lim_{q0}[q^4G(𝐪)]=k_BT\sigma R_g^2`$. We retrieve in this limit the mean-field correction (2), as expected. Interestingly, the polymer-membrane interactions also induce a tension-like term at short distances: $`G(𝐪)2k_BT\sigma q^2`$ for $`qR_g1`$. The correlations between chain and membrane conformations are stronger inside the polymer mushroom: as a consequence, the fluctuations of the surface are smoothed out at lenght scales smaller than the polymer size. The coefficient of this effective surface tension corresponds to the two-dimensional perfect gas pressure (recall that both sides are ornamented with $`𝒩`$ chains). Note that our results contrast with the usual interfacial behaviour, where tension is dominant at large scales whereas bending effects appear only at shorter scales. We should emphasize that the results presented here do not depend on the nature of the interface: they can be used to study a wide class of surfaces, provided that grafting of macromolecules is feasible.

We now apply our predictions to the lyotropic phase $`L_\alpha `$. This phase consists in a stack of regularly spaced membranes, that form a quasi-crystalline structure in one direction while retaining their fluid properties in the perpendicular plane. In pure water, the smectic order is generally stabilized by long-range electrostatic interactions. At high salt concentration or for neutral surfactant molecules, the stability of the phase arises from the entropic repulsions between fluctuating membranes: we will refer to this system as to the Helfrich system. In a continuum elastic theory, the multilayer membrane system is described by the displacement field $`u(𝐫)`$ in the $`z`$ direction normal to the layers. The standard theory of a pure $`L_\alpha `$ phase is based on the smectic energy density pggprostbook

$$\frac{F}{V}=\frac{1}{2}B\left(\frac{u}{z}\right)^2+\frac{1}{2}K\left(_{}^2u\right)^2,$$
(10)
the subscript  referring to the coordinates along the layers. $`K`$ is the smectic curvature modulus and is directly related to the membrane bending rigidity through the relation $`K=\kappa /d`$, with $`d`$ the average layer spacing. The smectic compressibility modulus $`B`$ reflects the interlamellar interactions, and is given in a Helfrich system by $`B=9\pi ^2(k_BT)^2/(64\kappa d^3)`$.

The lyotropic phase is a good candidate to probe experimentally the fine modifications to the properties of the pure phase induced by the addition of macromolecules tsapis01 ; sens01 . X-ray or neutron scattering experiments allow to measure the exponent $`\eta `$ that characterizes the power-law singularities around the Bragg peak: $`I(q)|qq_0|^{1+\eta }`$, with $`q_0=2\pi /d`$ safinya86 . The “strength” of the smectic order varies conversely to $`\eta `$. For pure membrane sytems, the FT of the average heigth fluctuations is $`|u(𝐪)|^2=Vk_BT/(Bq_z^2+Kq_{}^4)`$, which leads to the classical prediction first derived by Caillé caille72

$$\eta =q_0^2\frac{k_BT}{8\pi \sqrt{KB}}.$$
(11)
Likewise, a decorated membrane has according to eq. (9) a scale-dependent bending rigidity $`\kappa (𝐪)=\kappa (1+\alpha g_D(q^2R_g^2))`$, where we define the dimensionless parameter $`\alpha =\frac{k_BT}{\kappa }\sigma R_g^2`$, and $`g_D(x)=\frac{2}{x^2}(e^x1+x)`$ is the Debye function debye . This is certainly the most fundamental step for forecasting the properties of membrane-polymer complexes, but one should keep in mind that some macromolecular interactions such as repulsion between different membranes of the lamellar stack cannot be reduced only to the modification of the single membrane parameters. Notwithstanding this point, the height fluctuations are now given by

$`|u(𝐪)|^2_{pm}={\displaystyle \frac{Vk_BT}{Bq_z^2+Kq_{}^4\left(1+\alpha g_D(q_{}^2R_g^2)\right)}}`$ (12)
$`|u(𝐪)|^2\alpha Kq_{}^4g_D(q_{}^2R_g^2){\displaystyle \frac{Vk_BT}{\left(Bq_z^2+Kq_{}^4\right)^2}}`$
at lowest order in $`\alpha `$. It ensues that the $`q_{}`$ corrections lead to the same power-law divergence of the scattering intensity with new exponent

$$\eta _{pm}=\eta \left(1\frac{k_BT}{2\kappa }\sigma R_g^2\right).$$
(13)
The low-$`q`$ part of $`G(𝐪)`$ (*i.e.*, the correction to $`\kappa `$) actually provides the main contribution to the integrals leading to (13). Providing that the interlamellar distance is not affected, the tethered chains stiffen the lamellar order, that leads to a sharpening of the Bragg peak. A decrease of $`\eta `$ with polymer concentration is indeed observed in numerous experiments — see for instance warriner96 ; castro99 . However, the measured variation of $`\eta `$ is much sharper than the predicted linear behaviour at low coverage. As explained above, the chains are also likely to induce non-trivial variations of the compressibility $`B`$. A simple explanation of these modifications has been proposed by Castro-Roman *et al.* castro99 : a reduction in $`\eta `$ would be due to the polymer layer that effectively increases the membrane thickness. Within this picture, however, it is still difficult to explain why the membrane thickness varies so sharply with polymer coverage. Our study suggests that this could be caused by a non-homogeneous chain distribution, arising from the coupling between position and curvature. This naturally leads to a higher density of polymers in curved regions, which also correspond to the contact zones between the fluctuating bilayers. The concentration at which $`\eta `$ is significantly modified should then be as lower as the coverage needed to roughly have one polymer by membrane contact point.

In conclusion, we have developed a local description of the fluctuation spectrum of fluid membranes ornamented with grafted polymers. Contrary to current mean-field predictions, we find that the effect of the grafted chains cannot be reduced to an effective increase of the membrane bending rigidity. Instead, we show that the elastic coefficient is a scale-dependent quantity. On scales much larger than the polymer size, the bending rigidity is indeed renormalized. This low-$`q`$ correction is the dominant contribution to the so-called Caillé exponent measured by X-ray experiments in Lα phases. On scales smaller than the polymer size, polymer-membrane interactions smooth the fluctuations, leading to a tension-like contribution to the correlation function. Our results were obtained for Gaussian chains, a good representation for polymers in $`\theta `$-solvents. It is not clear at this point whether or not this description is robust with respect to the introduction of monomer-monomer excluded volume interactions. In particular, the “tension” regime that develops at short length scales might exhibit a different, non-trivial form. Indeed, the Debye function $`g_D`$ is known to scale as $`g_Dq^{1/\nu }`$, with $`\nu =1/2`$ in theta solvent and $`\nu 3/5`$ in good solvent conditions. If the relation $`Gq^4g_D`$ still holds for chains with excluded volume, one should then expect $`G(𝐪)`$ to behave like $`G(𝐪)q^{7/3}`$ for large wavevectors $`𝐪`$. Work on this issue is currently under progress.

###### Acknowledgements.

We gratefully acknowledge E. Blokhuis, R. Bruinsma, A. Lau and P. Sens for inspiring discussions. This work was supported by the Chemistry Department of the CNRS, under AIP “Soutien aux Jeunes Equipes”.


text:
# 1 Introduction

## 1 Introduction

The recent publication of an experimental realization of certain artificial materials — that are effectively isotropic, homogeneous, and possess negative real permittivity and permeability in some frequency range (Shelby, Smith, & Schultz, 2001) — has created a big splash in the electromagnetics research community. Although the first samples of these materials are quite primitive (Lakhtakia, 2001), there is enough experimental as well as theoretical evidence (Dewar, 2001; McCall, Lakhtakia, & Weiglhofer, 2001) along with the promise of applications (Pendry, 2001), that concerted theoretical research in advance of its experimental counterpart is justified.

This communication focuses on two aspects of the planewave responses of planar slabs whose permittivity and permeability have negative real parts at the frequency of interest. The first aspect comprises the magnitudes and phases of the reflection and the transmission coefficients; and the second encompasses the longitudinal shifts experienced by beams on total reflection. An $`\mathrm{exp}(i\omega t)`$ time–dependence is implicit, with $`\omega `$ as the angular frequency.

## 2 Theoretical Preliminaries

Suppose the region $`0zL`$ is occupied by an isotropic, homogeneous, dielectric–magnetic medium with relative permittivity $`ϵ_a=ϵ_a^{}+iϵ_a^{\prime \prime }`$ and relative permeability $`\mu _a=\mu _a^{}+i\mu _a^{\prime \prime }`$. The medium being passive, the principle of energy conservation requires that $`ϵ_a^{\prime \prime }0`$ and $`\mu _a^{\prime \prime }0`$. The half–spaces $`z<0`$ and $`z>L`$ are vacuous. A plane wave is incident on the face $`z=0`$ of the slab; accordingly, a plane wave is reflected into the half–space $`z<0`$ and another plane wave is transmitted into the half–space $`z>L`$.

### 2.1 s–polarization

Let the incident plane wave be $`s`$–polarized. Then, the electric field phasor in the two half–spaces can be expressed as follows:

$$𝐄(z)=\underset{¯}{u}_y\mathrm{exp}(i\kappa x)\{\begin{array}{cc}\mathrm{exp}(ik_0\mathrm{cos}\theta _{inc}z)+r_s\mathrm{exp}(ik_0\mathrm{cos}\theta _{inc}z),\hfill & z0\hfill \\ t_s\mathrm{exp}\left[ik_0\mathrm{cos}\theta _{inc}(zL)\right],\hfill & zL\hfill \end{array}.$$
(1)
Here and hereafter, $`k_0=\omega \sqrt{ϵ_0\mu _0}`$ is the wavenumber in vacuum; $`ϵ_0`$ and $`\mu _0`$ are the permittivity and the permeability of vacuum, respectively; $`\eta _0=\sqrt{\mu _0/ϵ_0}`$ is the intrinsic impedance of vacuum; the quantity

$$\kappa =k_0\mathrm{sin}\theta _{inc}$$
(2)
involves $`\theta _{inc}`$ as the angle of incidence measured from the $`+z`$ axis; $`\underset{¯}{u}_x`$, $`\underset{¯}{u}_y`$ and $`\underset{¯}{u}_z`$ are cartesian unit vectors; while $`r_s`$ and $`t_s`$ are the reflection and the transmission coefficients, respectively.

Following standard practice, and after defining the 2$`\times `$2 matrixes

$$[P_s]=\left[\begin{array}{cc}0& \omega \mu _0\mu _a\\ \omega ϵ_0ϵ_a+\frac{\kappa ^2}{\omega \mu _0\mu _a}& 0\end{array}\right],[M_s]=\mathrm{exp}\left\{i[P_s]L\right\},$$
(3)
we obtain the following matrix equation for the boundary value problem described in the previous two paragraphs:

$$t_s\left[\begin{array}{c}1\\ \eta _0^1\mathrm{cos}\theta _{inc}\end{array}\right]=[M_s]\left[\begin{array}{c}r_s+1\\ (r_s1)\eta _0^1\mathrm{cos}\theta _{inc}\end{array}\right].$$
(4)
This equation’s solution yields $`r_s`$ and $`t_s`$, from which the remittances $`R_s=|r_s|^2`$ and $`T_s=|t_s|^2`$ can be calculated.

### 2.2 p–polarization

If the incident plane wave is $`p`$–polarized, the magnetic field phasor outside the slab $`0zL`$ can be written as

$$𝐇(z)=\underset{¯}{u}_y\eta _0^1\mathrm{exp}(i\kappa x)\{\begin{array}{cc}\mathrm{exp}(ik_0\mathrm{cos}\theta _{inc}z)+r_p\mathrm{exp}(ik_0\mathrm{cos}\theta _{inc}z),\hfill & z0\hfill \\ t_p\mathrm{exp}\left[ik_0\mathrm{cos}\theta _{inc}(zL)\right],\hfill & zL\hfill \end{array}.$$
(5)
The solution of the boundary value problem involves the calculation of the coefficients $`r_p`$ and $`t_p`$ from the matrix equation

$$t_p\left[\begin{array}{c}\mathrm{cos}\theta _{inc}\\ \eta _0^1\end{array}\right]=[M_p]\left[\begin{array}{c}(r_p1)\mathrm{cos}\theta _{inc}\\ (r_p+1)\eta _0^1\end{array}\right],$$
(6)
where

$$[M_p]=\mathrm{exp}\left\{i[P_p]L\right\},[P_p]=\left[\begin{array}{cc}0& \omega \mu _0\mu _a\frac{\kappa ^2}{\omega ϵ_0ϵ_a}\\ \omega ϵ_0ϵ_a& 0\end{array}\right].$$
(7)
The remittances $`R_p=|r_p|^2`$ and $`T_p=|t_p|^2`$ can be obtained thereafter.

## 3 Reflection and Transmission Analysis

Both matrix equations (4) and (6) can be stated together compactly as

$$t[f]+r[M][g]=[M][h],[M]=\mathrm{exp}\left\{i[P]L\right\},$$
(8)
where $`[f]`$, $`[g]`$ and $`[h]`$ are 2$`\times `$1 matrixes.

Now, $`ϵ_a^{}>0`$ and $`\mu _a^{}>0`$ if the region $`0zL`$ is filled with a conventional medium. Suppose the transformation

$$\left\{ϵ_a^{}ϵ_a^{},\mu _a^{}\mu _a^{}\right\}$$
(9)
occurs; i.e., the conventional medium is replaced by a counterpart the real parts of whose permittivity and permeability are exactly negative of those of the conventional medium, but both $`ϵ_a^{\prime \prime }`$ and $`\mu _a^{\prime \prime }`$ remain unchanged.

### 3.1 Magnitude Invariance and Phase Reversal

According to (3) and (7), the transformation

$$\left\{[P][P]^{}\right\}$$
(10)
must follow, with the asterisk denoting the complex conjugate. As shown in the Appendix, we get

$$\left\{[M][M]^{}\right\}$$
(11)
in consequence. The matrixes $`[f]`$, $`[g]`$ and $`[h]`$ in (8) being real–valued, (11) in turn mandates the transformation

$$\left\{rr^{},tt^{}\right\}$$
(12)
for the reflection and the transmission coefficients, and

$$\left\{RR,TT\right\}$$
(13)
for the remittances.

The symmetry inherent in the relationship (9)$``$(12) is of significance for experiments on planar, isotropic dielectric–magnetic slabs. When the signs of the real parts of both the permittivity and the permeability of a medium are altered, the planewave reflection and transmission coefficients remain unchanged in magnitude but suffer a phase reversal. The alteration cannot be detected if only the planewave remittances are measured.

### 3.2 Goos–Hänchen Shifts

A major consequence of the relationship (9)$``$(12) can be appreciated when the incident plane wave is replaced by a beam. By definition, a plane wave has infinite transverse extent; but a beam is of finite transverse extent. Newton had conjectured that if a beam were to impinge on a planar interface with an electromagnetically rarer medium and total reflection were to occur, the reflected beam would emerge with a longitudinal displacement (along the $`+x`$ axis). Such displacements are often called the Goos–Hänchen shifts after the two scientists who reported the first measurements in 1947 (Goos & Hänchen, 1947). Subsequently, these shifts have been estimated as well as measured for planar interfaces between many different pairs of materials (Lakhtakia, 2002; Lotsch, 1970a,b, 1971a,b), and several applications have emerged as well (Jackson, 1999; Chauvat et al., 2000; de Fornel, 2000).

Analytical treatment of the Goos–Hänchen shifts essentially requires the solution of the same problem as dealt with in Section 2 (Haibel, Nimtz, & Stahlhofen, 2001). The vacuous half–spaces are, however, filled with a medium possessing $`ϵ_b`$ and $`\mu _b`$ as its relative permittivity and relative permeability, both positive and real–valued. Equations (1)–(7) still apply, but with (i) $`k_0`$ replaced by $`k_0\sqrt{ϵ_b\mu _b}`$, and (ii) $`\eta _0`$ by $`\eta _0\sqrt{\mu _b/ϵ_b}`$. The relationship (9)$``$(12) still holds.

——————————————————————————————————–

——————————————————————————————————–

With the assumption of a Gaussian profile for the incident beam, and after denoting the reflection coefficient $`r=|r|\mathrm{exp}(i\phi )`$, the longitudinal shift is estimated as

$$d=\frac{\phi }{\kappa }.$$
(14)
This quantity is meaningful, in the present context, only when $`\theta _{inc}`$ exceeds the critical angle; i.e., when total reflection occurs. A depth of penetration equal to $`(d/2)\mathrm{cot}\theta _{inc}`$ is also sometimes calculated.

Two Goos–Hänchen shifts are possible in the present scenario: $`d_s`$ and $`d_p`$ for $`s`$– and $`p`$–polarized beams, respectively. Their characteristics are well–known (Jackson, 1999; de Fornel, 2000) and do not bear repetition here. It suffices to note here that

$$\left\{ϵ_a^{}ϵ_a^{},\mu _a^{}\mu _a^{}\right\}\left\{d_sd_s,d_pd_p\right\},$$
(15)
this relationship being an immediate consequence of the relationship (9)$``$(12) and the definition (14).

The result (15) can be understood with the help of Figure 1 as follows: When total reflection occurs and the electromagnetically rarer medium is of the conventional type (i.e., $`ϵ_a^{}>0`$ and $`\mu _a^{}>0`$), the reflected beam appears as if the incident beam has skipped forward. On the contrary, according to (15), the reflected beam must appear to skip backward, if the electromagnetically rarer medium is of the kind reported by Shelby, Smith, & Schultz (2001) (i.e., $`ϵ_a^{}<0`$ and $`\mu _a^{}<0`$). No violation of causality would occur in the latter instance, because some time must elapse during the traversal of the distance that is apparently skipped. This result is also consistent with the propagation velocity being oppositely directed to the energy velocity when $`ϵ_a^{}<0`$ and $`\mu _a^{}<0`$ (McCall, Lakhtakia, & Weiglhofer, 2001).

## 4 Conclusions

Examining the alteration in the planewave response of an isotropic, dielectric–magnetic, homogeneous, planar slab when the real parts of both the permittivity and the permeability are changed from positive to negative, we arrived at the following deductions:

* The planewave reflection and transmission coefficients suffer phase reversals.
* The planewave remittances remain unaffected.
* When conditions for total reflection prevail, the longitudinal shift experienced by a beams reverses its direction.

## Appendix

Let a matrix $`[P]`$ of size $`m\times m`$ be diagonalizable. Then it can be decomposed as

$$[P]=[U][G][U]^1,$$
(16)
where the matrix $`[G]=\mathrm{diag}(g_1,g_2,\mathrm{}g_m)`$ contains the consecutive eigenvalues of $`[P]`$, while the successive columns of $`[U]`$ are the corresponding eigenvectors of $`[P]`$. It follows that (Hochstadt, 1975)

$$\mathrm{exp}\left\{i[P]z\right\}=[U]\mathrm{exp}\left\{i[G]z\right\}[U]^1.$$
(17)
Let the matrix $`[Q]=[P]^{}`$ so that

$$[Q]=[U]^{}\left([G]^{}\right)\left([U]^1\right)^{}.$$
(18)
Consequently,

$`\mathrm{exp}\left\{i[Q]z\right\}`$ $`=`$ $`[U]^{}\mathrm{exp}\left(i[G]^{}z\right)\left([U]^1\right)^{}`$ (19)
$`=`$ $`\left([U]\mathrm{exp}\left\{i[G]z\right\}[U]^1\right)^{}`$
$`=`$ $`\left(\mathrm{exp}\left\{i[P]z\right\}\right)^{}.`$
Chauvat, D., O. Emile, F. Bretenaker, & A. Le Floch. 2000. Direct measurement of the Wigner delay associated with the Goos–Hänchen effect. Phys. Rev. Lett. 84:71–74.

de Fornel, F. 2000. Evanescent waves — From Newtonian optics to atomic optics. Berlin: Springer; pp. 12–18.

Dewar, G.A. 2001. Candidates for $`\mu <0`$, $`ϵ<0`$ nanostructures. Int. J. Mod. Phys. B 15:3258–3265.

Goos, F. & H. Hänchen. 1947. Ein neuer und fundamentaler Versuch zur Totalreflexion. Ann. Phys. Lpz. 1:333–346.

Haibel, A., G. Nimtz, & A.A. Stahlhofen. 2001. Frustrated total reflection: the double–prism revisited. Phys. Rev. E 63:047601.

Hochstadt, H. 1975 Differential equations — A modern approach. New York, NY: Dover Press; p. 57.

Jackson, J.D. 1999. Classical electrodynamics, 3rd ed. New York: Wiley; pp. 306–309.

Lakhtakia, A. 2001. An electromagnetic trinity from “negative permittivity” and “negative permeability”. Int. J. Infrared Millim. Waves 22:1731–1734.

Lakhtakia, A. 2002. Truncation of angular spread of Bragg zones by total reflection, and Goos–Hänchen shifts exhibited by chiral sculptured thin films. Arch. Elektr. Über. 56:000–000 (accepted for publication).

Lotsch, H.K.V. 1970a. Beam displacement at total reflection: The Goos–Hänchen effect, I. Optik 32: 116–137.

Lotsch, H.K.V. 1970b. Beam displacement at total reflection: The Goos–Hänchen effect, II. Optik 32: 189–204.

Lotsch, H.K.V. 1971a. Beam displacement at total reflection: The Goos–Hänchen effect, III. Optik 32: 299–319.

Lotsch, H.K.V. 1971b. Beam displacement at total reflection: The Goos–Hänchen effect, IV. Optik 32: 553–569.

McCall, M.W., A. Lakhtakia, & W.S. Weiglhofer. 2001. The negative index of refraction demystified. Report No. 2001/30. Department of Mathematics, University of Glasgow, UK.

Pendry, J. 2001. Electromagnetic materials enter the negative age. Physics World 14(9):47–51.

Shelby, R.A., D.R. Smith, & S. Schultz. 2001. Experimental verification of a negative index of refraction. Science 292:77–79.