text:
# The evolution of the high energy tail in the quiescent spectrum of the soft X–ray transient Aql X-1

## 1 Introduction

The large luminosity swing of transient X–ray binaries allows the sampling a variety of physical conditions that are inaccessible to accreting compact objects in persistent sources. The very low luminosity that characterises the quiescent state of neutron star soft X–ray transients, SXRTs, ($`L_X10^{32}10^{33}\mathrm{erg}\mathrm{s}^1`$) opens up the possibility of studying these old, fast spinning neutron stars in different and yet unexplored regimes such as accretion onto the neutron star magnetosphere (propeller), resumed millisecond radio pulsar activity, and/or low-level atmospheric emission from the cooling of the neutron star in between the accretion intervals of the outbursts (e.g. Campana et al. 1998a; Brown et al. 1998; Rutledge et al. 2002b).

In recent years the quiescent properties of a handful SXRTs have been studied in some detail. The main outcome of these investigations is that the quiescent X–ray spectra of SXRTs display a soft component plus a hard (power-law) component contributing a comparable flux in the 0.5–10 keV band (Campana 2001; Bildsten & Rutledge 2000; Wijnands 2001). The soft component has been frequently modelled with a black body model of 0.1–0.3 keV temperature and few km radius. Especially promising is the idea that the soft component of SXRTs may be produced from the cooling of the neutron star heated during the repeated outbursts (van Paradijs et al. 1987; Stella et al. 1994; Campana et al 1998a). The theory of deep crustal heating by pycnonuclear reactions compares well with the observations (Brown et al. 1998; Campana et al. 1998a; Rutledge et al. 1999; Colpi et al. 2001). In particular, Rutledge et al. (1999) fitted neutron star atmospheric models to the soft component of quiescent spectra of SXRTs and derived slightly smaller temperatures (0.1–0.3 keV) and larger radii (10–15 km, consistent with the neutron star radius) than those inferred from simple black body fits.

Observationally, the hard component is well described by a power law tail. In the quiescent spectrum of Aql X-1 and Cen X-4 observed by ASCA and BeppoSAX this component is statistically significant (Asai et al. 1996, 1998; Campana et al. 1998b, 2000) with photon index in the 1–2 range. The same power law is needed (even if not statistically significant) in the analysis of Chandra data in order to achieve an emitting radius of the cooling component consistent with the neutron star radius (otherwise the inferred radius would be smaller; Rutledge et al. 2001a, 2001b). The nature of this hard component is still uncertain. Models range from Comptonization to Advection/Convection Dominate Accretion Flow (ADAF/CDAF) to shock emission from the neutron star that resumed its radio pulsar activity in quiescence. The latter model envisages a situation similar to that of the eclipsing radio pulsar PSR B1259–63 or of the ‘black widow’ pulsar PSR B1957+20: a shock at the boundary between the relativistic MHD wind from the radio pulsar and the matter outflowing from the companion star (Tavani & Arons 1997; Tavani & Brookshaw 1991; Campana et al. 1998a). For the model to explain the observed luminosity in the hard power law component of quiescent SXRTs, some few percent of the pulsar spin-down luminosity must be converted into shock emission111Ongoing deep searches in the radio band have not yet revealed any steady or pulsed emission from quiescent SXRTs (Burgay et al. 2003); however free-free absorption due to matter in the binary system might be an important limiting factor in these searches (Stella et al. 1994; Burgay et al. 2003).. The shock emission model predicts synchrotron emission with power law photon indexes in the 1.5–2 range and extending over a wide range of frequencies. Indirect indications for the presence of this emission also from UV observations of Cen X-4 with HST revealing a flat spectrum (i.e. $`\mathrm{\Gamma }2`$) which matches well the extrapolated X–ray power law component (McClintock & Remillard 2000). Power law indexes outside the above range are indications of strong (inverse-Compton) cooling.

In a way similar to what is routinely done in the optical, a promising tool to probe the X–ray emitting regions is through a rough eclipse mapping technique (e.g. Horne 1985). Chandra observations of the eclipsing SXRT 4U 2129+47 were the first to exploit the potential of this technique by looking at the extension of the emitting regions through eclipses (Nowak, Heinz & Begelman 2002). During eclipses the soft component gets totally eclipsed whereas the hard component is too faint to be revealed. This implies an upper limit on the emission size of $`\stackrel{<}{}10\%`$ the orbital separation (Nowak et al. 2002). The inclination of Aql X-1 has been estimated, from ellipsoidal variations in the R and I band lightcurves, to be greater than 36 degrees (Welsh, Robinson & Young 2000).

In this paper we investigate in more detail the quiescent spectrum of one of the best studied SXRT sources: Aql X-1. We take advantage of four Chandra exposures (Rutledge et al. 2002a) and one 76 ks long (unpublished) BeppoSAX exposure. All these data were collected just after the November 2000 outburst. Based on the Chandra data Rutledge et al. (2002a) claimed that the soft component decreased by $`50\%`$ over three months, then increased by $`35\%`$ in one month, and then remained constant ($`<6\%`$ change) over the last month. The variability of these observations was ascribed to an intrinsic variability of the soft component hinting to accretion onto the neutron star surface. Here we discuss in more detail these observations together with the BeppoSAX long exposure, probing the shock emission model.

In Section 2 we deal with the data. In Section 3 we describe the spectral fitting and related results. Discussion and conclusions are reported in Section 4.

## 2 Data

### 2.1 Chandra

Aql X-1 was observed by Chandra after the November 2000 outburst on four occasions (see Table 1). Observations were carried out with the backside-illuminated ACIS-S detector (S3) at an off-axis position of $`4^{}`$ and limited read-out area of 1/8 (achieving a time resolution of 0.44 s) in order to limit problems connected to pile-up (see Rutledge et al. 2002a). As expressly required, Aql X-1 fell on the same physical pixels in order to avoid problems with the CCDs quantum efficiency (see Rutledge et al. 2002a). For the analysis we use CIAO 2.2.1 with CALDB 2.15 (these are later versions than what used by Rutledge et al. 2002a). In all observations, we extracted the source counts from an elliptical $`4.5^{\prime \prime }\times 3^{\prime \prime }`$ region centered on source with a position angle matching the source. Background photons were extracted from an annular region with inner and outer radii of $`10^{\prime \prime }`$ and $`20^{\prime \prime }`$, respectively. Data were extracted, using psextract, into pulse-invariant (PI) spectra. We grouped all the spectra to have (at least) 30 photons per channel. We also corrected all the ancillary (arf) files with the recently released corrarf tool to account for the continuous degradation in the ACIS CCD’s quantum efficiency.

### 2.2 BeppoSAX

We analysed data from the two imaging instruments on board the BeppoSAX satellite: the Low Energy Concentrator Spectrometer (LECS; 0.1–10 keV, Parmar et al. 1997) and the Medium Energy Concentrator Spectrometer (MECS; 1.6–10.5 keV, Boella et al. 1997). Non-imaging instruments provided only upper limits. Only two of the three MECS units were operating at the time of the observations. LECS data were collected only during satellite night-time resulting in shorter exposure times. For a summary of the observations see Table 1. The observation took place on 2001 April 14, observing Aql X-1 for a net exposure time of 76 ks with the MECS and 30 ks for the LECS.

Products were extracted using the FTOOLS package (v. 5.1). LECS and MECS events were extracted from a circle of $`4^{}`$ radius. The background was subtracted using spectra from blank sky files at the same detector coordinates (after checking that the background of the observation was comparable). We rebinned the LECS and MECS spectra in order to have 80 counts per spectral bin each.

## 3 Overall spectral analysis

### 3.1 Spectral model

Rutledge et al. (2002a) analysed the four Chandra spectra together. Spectra are different, e.g. a hard power-law tail has been detected only during two of the four observations. To account for these differences Rutledge et al. (2002a) considered a model made by an absorbed atmosphere model plus a power law and let vary single parameters (fixed all the other) trying to account for the variations. They found that differences in the spectra cannot be explained as entirely due to either a changing power law flux and/or index or to a variable column density. On the contrary, differences can be acceptably (in terms of $`\chi ^2`$ statistics) explained as entirely due to a (non monotonic) temperature variability of the thermal emission component. These variations cannot be explained within the deep crustal heating model and the authors suggested that these might originate from quiescent accretion onto the neutron star surface.

As discussed in the introduction, we want to test here the hypothesis that quiescent emission of SXRTs (and Aql X-1 in particular) are produced by a soft thermal component, likely arising from cooling of the neutron star, plus a power law hard tail, arising from shock emission due to an active millisecond pulsar. In theory, the soft component is steady on relatively short time whereas the hard component can likely vary depending on the geometry and density the outflowing matter. Hydrodynamical simulations (Brookshaw & Tavani 1993) as well as radio observations of millisecond radio pulsars (MSPs) in binary systems with a sizeable mass transfer show complex geometries and, more importantly, variations from one orbital cycle to another. The example of the recently discovered MSP PSR J1740–5340 (D’Amico et al. 2001; Ferrario et al. 2001) is enlighting. This is a 3.7 ms MSP orbiting every 32.5 hr a main sequence companion. The source is located in the globular cluster NGC 6397 at 2.5 kpc. The radio pulsar gets partially and totally eclipsed over a wide range of orbital phases. It emits X–rays as observed by Chandra (Grindlay et al. 2001) likely arising from shock emission. As testified by this source, the MSP is eclipsed for large part of the orbit and variations from orbit to orbit are seen.

This case motivates us to consider a spectral model for fitting the quiescent X–ray spectra of Aql X-1 made by a soft thermal component from the entire neutron star, a variable power law component (the strength of which depends on the interaction with the surrounding matter) and a variable column density due to variations intrinsic to the source (over a fixed interstellar column density).

### 3.2 Spectral analysis

Spectral analysis was carried out with the XSPEC software (v. 11.2.0). The spectral model we adopted consist of a (fixed) cooling spectrum (we use here the hydrogen atmosphere model by Gänsicke, Braje & Romani 2002, hyd\_spectra.mod in XSPEC) plus a variable power law. A variable column density was also adopted (tbabs in XSPEC, Wilms, Allen & McCray 2000). We fitted all the five spectra (4 Chandra and 1 BeppoSAX, LECS plus MECS) together. For the BeppoSAX data we used the public response matrices available in January 2000. During the fit a variable normalisation factor between LECS and MECS was included to account for the mismatch in the absolute flux calibration of the BeppoSAX instruments (Fiore, Guainazzi & Grandi 1999). A variable normalisation factor was also included to account for mismatch between BeppoSAX and Chandra222We find a ratio of the BeppoSAX over the Chandra normalization of $`1.9_{0.6}^{+1.4}`$, consistent with previous determinations (e.g. Piro et al. 2001).. The spectrum provides a statistically acceptable description of the entire data set with a reduced $`\chi ^2`$ of 1.00 (null hypothesis probability $`49.0\%`$, see Table 2 and Fig. 1). This indicates that, at least at a first sight, all the data concerning the quiescence following the November 2000 outburst of Aql X-1 are consistent with a model made by a cooling neutron star and a variable power law and absorption components.

A detailed theory of the shock emission mechanism has been developed by Tavani & Arons (1997) tailored on the young radio pulsar PSR B1259–63 orbiting a Be star. A more detailed discussion on SXRTs is within Campana et al. (1998a). The expected spectrum is a power law spectrum with photon index in the 1.5–2 range, with a positive correlation between the quantity of matter at the shock region (possibly traced by the column density) and the power law index. This behaviour has been observed in PSR B1259–63 which showed a photon index of $`2`$ at periastron and hardened towards apastron. However, the column density to PSR B1259-63 is larger than the one to Aql X-1 and we expect that variations in the column density are therefore highly suppressed. This correlation might provide us with a further check on the shock emission mechanism. In Fig. 2 we show the column density and power law indexes for the five observations. A correlation is indeed present suggesting that this mechanism might be at work. We note however that the column density and the power law index are correlated parameters and for this reason we derived the errors ($`68\%`$ c.l. in Fig. 2), for the two parameters together ($`\mathrm{\Delta }\chi ^2=2.30`$). Correlation is tight. A fit with a constant provides a $`\chi _{\mathrm{red}}^2=5.4`$ with a null hypothesis probability of $`2\times 10^4`$ whereas a linear fit gives $`\chi _{\mathrm{red}}^2=0.2`$ (with an F-test probability of $`99.8\%`$). We also carried out a weighted linear correlation test finding a correlation probability of $`r_\mathrm{w}=0.9`$ ($`95\%`$ probability). A further correlation can be tested between the power law flux and the column density. A weighted Pearson correlation test gives $`r=0.8`$ ($`90\%`$ probability).

A similar correlation has also been found in the quiescent emission of the transient black hole candidate V404 Cyg (Kong et al. 2002). They however suspected that the correlation is not intrinsic to the source, but it is an artifact of the fitting process. In their case the slope of the correlation is nearly (to within $`5\%`$) the same as the slope of the major axis of the parameter confidence contours. Moreover, they did not found correlation between flux and column density. This is contrary to our findings. In Fig. 2 we also included the confidence countors showing that this alignement effect is not present.

We consider then a cut-off power law model in addition to the neutron star atmosphere one. The overall fit is as good as the others ($`\chi _{\mathrm{red}}^2=1.05`$, null hypothesis probability of $`35\%`$). The cut-off energy is only loosely constrained and only an upper limit for each observation can be put. Drawing contour plots in the cut-off energy – column density plane we find an anticorrelation between the two quantities: the lower the energy of the cut-off the larger the column density. From a Pearson correlation test we obtain $`r=0.9`$, corresponding to a significance of $`96\%`$. This might indicate an (almost) constant luminosity that is shared among different quantities of matter.

To further test this idea we fit the data with a different model with a smaller emission at low energies. We consider a model consisting of a Comptonized component (COMPTT, Titarchuk 1994) plus an atmosphere component. We found a large degree of freedom with this model which forced us to freeze the input soft photon temperature to the temperature of the neutron star atmosphere (that is the same for all the observations) and take the same amount of optical depth for all the observations (that turns out to be 7.7) and leave the plasma temperature free to vary. The fit is as good as the one with the power law ($`\chi _{\mathrm{red}}^2=1.00`$). We derive an anti-correlation between the column density and the plasma temperature. A weighted Pearson correlation test gives $`r=0.9`$ ($`93\%`$ probability).

## 4 Discussion

Rutledge et al. (2002a) analysing Chandra data of the Aql X-1 quiescent phase after the November 2000 outburst found a variable flux and X–ray spectrum. They interpreted these variations in terms of variations of the neutron star effective temperature, which changed from $`130_5^{+3}`$ eV (C1), down to $`113_4^{+3}`$ eV (C2), and finally increased to $`118_4^{+9}`$ eV (C3, C4). Interestingly, during observation C4 they also found short-term variability (at $`32\%`$ rms) and a possible absorption feature near 0.5 keV (even if this feature can also be explained as due to a time-variable response in the ACIS detector; Rutledge et al. 2002a).

Short term variability is a powerful tool for the study of the emission mechanism(s) responsible for the SXRTs quiescent emission. A factor of 3 variability over timescales of days (Campana et al. 1997) and $`40\%`$ over 4.5 yr (Rutledge et al. 2001b) has been reported in Cen X-4. Several other neutron star systems have also been found to be variable in quiescence by factors of 3–5 (e.g. Rutledge et al. 2000) but data have been collected over several years and with different instruments. Chandra data on Aql X-1 are the first that show a clear luminosity variation and, more importantly, an increase during quiescence. No known mechanism associated with crustal heating can account for this variability (Rutledge et al. 2002a).

Here we approach the same Chandra data plus an unpublished long BeppoSAX observation of Aql X-1 in quiescence carried out in the same period to probe a different spectral model. Deep crustal heating (Brown et al. 1998; Rutledge et al. 1999; Colpi et al. 2001) has been proposed as a physically sound mechanism powering the soft component of the quiescent spectra of SXRTs. Several mechanisms have been proposed to explain the hard tail component often observed in quiescent SXRTs. One of them, physically motivated by the recent observations of the MSP PSR J1740–5340 (D’Amico et al. 2001; Grindlay et al. 2001), relies on the shock emission between the relativistic MSP wind and matter outflowing from the companion (Tavani & Arons 1997; Campana et al. 1998a). These two components are not exclusive. With this physical scenario, we fit the spectra fixing the soft component for all the observations and leave free to vary the hard component and the column density (that changes by a factor of $`5`$). This model is consistent with the entire dataset ($`\chi _{\mathrm{red}}^2=1.00`$, for 109 d.o.f. and with a null hypothesis probability of $`49.0\%`$). Fitting the same dataset with the best fit model by Rutledge et al. (2002a) with the addition of the BeppoSAX data, we obtain a slightly worse fit ($`\chi _{\mathrm{red}}^2=1.17`$, for 113 d.o.f. and with a null hypothesis probability of $`11.1\%`$). An equally good fit is provided by a neutron star atmosphere plus a Comptonization component (COMPTT). This model provides a high degree of freedom with a $`\chi _{\mathrm{red}}^2=0.99`$ fit (for 104 d.o.f. and with a null hypothesis probability of $`50.5\%`$) obtained fixing the soft Wien temperature to atmosphere temperature and the plasma temperature kept the same in all the observations. We conclude that the scenario proposed is (at least) equally well consistent with the data, meaning that a shock emission scenario can account for the spectral variability observed in Aql X-1. We also note that Rutledge et al. (2002) found $`32\%`$ (rms) variability in observation C4. In their case the power law component contributed only $`12\%`$ of the flux. From our fit the power law component contributes to $`38\%`$ of the total flux, so it can in principle account for all of the short-term variability.

Despite the low number of points, spectral parameters derived for the power law index show some correlation with the column density (interpreted as a measure of the variable mass around the system, over a fixed interstellar amount) as well as with the power law flux. This correlation might be expected in the shock emission scenario (Tavani & Arons 1997). What is now expected is the large value of the power law index in the last observations. This might then provide an indication of a different regime in the system, possibly underlying a larger inverse Compton cooling. The hard part of the spectrum is in fact consistent also with a thermal bremsstrahlung spectrum.

Further observations can shed light on this new interesting field, namely variability in the quiescent phase of SXRTs, which up to now has been often unconsidered.

We thank an anonymous referee and G. Ghisellini for useful comments.


format: text

text:
# Simultaneous surface acoustic wave and surface plasmon resonance measurements: electrodeposition and biological interactions monitoring

## I Introduction

Surface plasmon resonance (SPR) is a well accepted direct detection technique for monitoring biological processes liedberg ; gizeli1 . While ellipsometry is another well known method for analyzing thin film properties kasemo , its use in liquid medium for monitoring biochemical reactions is made difficult by the varying environment through which the probing light beam has to propagate. In the Kretschmann configuration, the laser generating the SPR evanescent wave is only propagating through the substrate, leading to a better control over the influence of the various buffer solutions used during a protein adsorption experiment.

The use of various acoustic wave devices for monitoring bound mass changes in liquid media is well known gizeli2 . Love mode devices, based on the propagation of a guided shear acoustic wave, present sensitivity improvements over the more usual quartz crystal microbalance gizeli3  as well as a compatibility with measurements in liquids.

We take advantage of the unique geometrical setup of the surface acoustic wave (SAW) device which leaves the backside of the quartz wafer free of electrodes to inject a laser in order to generate an evanescent surface plasmon on the gold coated sensing area. Such a setup enables simultaneous estimates of the bound mass and dielectric surface properties changes during electrochemical and biochemical reactions occurring on the sensing electrode bailey ; johan1 .

Estimates of water content in protein layers is an important topic in the development of biosensors, since the detection sensitivity is directly related to the number of active sites on the surface of the sensor to which the receptor molecules are bound. When physical measurements of this bound layer are performed, such as via SPR or acoustic measurements, a large water content will lead to overestimates of the potential detection limit of the biosensor due to an overestimate of number of molecules bound to the surface. From a more fundamental point of view, determining water level content should provide more accurate physical parameters of the protein layer itself such as density and optical index by allowing the subtraction of the influence of water once it has been identified.

## II Experimental data

We use a modified commercial SPR instrument (IBIS II, IBIS Technologies BV, The Netherlands) to detect the surface plasmon resonance angle after replacing the gold coated glass slide by a Love mode SAW device (Fig. 1). The excitation laser in this instrument has a wavelength of 670 nm. All reflected intensity vs. angle curves (6 o angle span recorded on a 200 pixels CCD array) were recorded and later fitted by a polynomial to extract with high accuracy the position of the dip. The acoustic wave device is made of a 500 $`\mu `$m thick ST-cut quartz wafer on which 200 nm thick sputtered $`Al`$ interdigitated electrodes are patterned. The surface is coated by a 1.13 $`\mu `$m PECVD silicon dioxide layer acting as a guiding layer, and the 4.9$`\times `$5.4 mm2 sensing area is coated with 10 nm $`Ti`$ and 50 nm $`Au`$. This area acts both as a working electrode for electrochemistry or a grounded electrode during biochemical experiments, as well as a supporting layer for the surface plasmon resonance generation. The influence of the SAW device substrate over the detection of the SPR is limited to interference patterns due to the birefringence of the quartz, and the optical index mismatch between the quartz and the deposited silicon dioxide layer. The former effect is reduced by orienting the optical axis of quartz so that it is in the plane defined by the normal to the sensing surface and the wavevector of the laser, thus minimizing the optical index difference between the ordinary and the extraordinary axis. The remaining interference patterns lead to fringes with a low enough contrast that the surface plasmon peak is easily identified and tracked during the experiments.

The phase of the SAW delay line is monitored at a fixed frequency using a network analyzer HP4396A at 123.200 MHz. The phase shift is converted to a frequency shift thanks to the linear phase to frequency relationship recorded in the Bode plot. The observed phase shift leads to a frequency shift (as would be observed in a phase locked loop configuration) which in turn can be converted to an adsorbed mass change through the sensitivity of the device.

We first electrodeposited copper on the surface in order to calibrate the mass sensitivity of the acoustic wave device in liquid medium (Fig. 2) me . At the same time, the SPR displayed resonance angle shifts due to the varying potential kolb . When the voltage applied by the potentiostat is above 0.2 V with respect to the $`Cu`$ pseudo reference electrode, the SPR angle slightly shifts due to the electroreflectance effect described by Kötz et al. kolb . Below 0.2 V, under potential deposition starts depositing a mono layer of $`Cu`$ atoms on the $`Au`$ surface as visible both in the phase shift of the SAW device and as a reversal of the trend of the angle shift of the SPR. Below 0 V, a rough $`Cu`$ layer is deposited, leading to a loss of the SPR peak and a large phase shift of the SAW device due to the large added mass. The rough, discontinuous film of copper clusters a few tens to hundreds of nanometers high was observed by in-situ AFM imaging (data not shown me ). Simulating the influence of sub-wavelength metallic structures on the surface plasmon resonance is a difficult task ebbesen ; JACS  which we have not attempted to solve. Simultaneously, the SAW device displays a frequency drop $`\mathrm{\Delta }f`$ due to the added mass $`\mathrm{\Delta }m`$ to the sensing electrode, leading with the combination of the current measurement from the potentiostat to a direct estimate of the mass sensitivity $`S=\frac{\mathrm{\Delta }f}{f_0}\frac{A}{\mathrm{\Delta }m}`$ where $`A=0.49\times 0.54`$ cm2 is the area of the sensing working electrode and $`f_0=123.200`$ MHz the base frequency of the SAW device. The electrodeposited mass is deduced from the current measurement $`I`$ with a sampling period $`\delta t`$: $`\mathrm{\Delta }m=\frac{_jI_j\times \delta t}{F}\times \frac{M}{n_e}`$ where $`F=96485`$ C is the Faraday constant, $`M`$ the molar weight of the electrodeposited metal and $`n_e`$ the number of electrons exchanged for each metallic atom reduced. The resulting sensitivity for the Love mode device under consideration here is 150$`\pm `$15 cm2/g.

Once the mass sensitivity of the SAW is calibrated, we adsorbed on the same device a crystalline monolayer of S-layer proteins (100 $`\mu `$g/ml of the protein SbpA of Bacillus sphaericus CCM 2177 in 0.5 mM tris, 10 mM $`CaCl_2`$, pH=9 buffer) review ; erika  while monitoring simultaneously with the SAW device the mass variation and with SPR the surface dielectric change. The adsorption and desorption steps, the latter being possible with the use of 2 % $`NaOCl`$, was repeated several times (Fig. 3). A property of this protein which makes it suitable for calibration of a new instrument is that it only forms a monolayer and will not stack to multiple layers even at high concentrations.

## III Discussion

The parameters required for analyzing the SAW data are the layer thickness $`d`$ and the layer density $`\rho _{layer}`$. The parameters required for analyzing the SPR data are the common layer thickness $`d`$ and the optical index $`n`$ of the adsorbed layer. We thus have three free parameters and only two measurements, namely the adsorbed mass per unit area and an SPR resonance angle shift which is a function of both $`n`$ and $`d`$. The two options are to add one more measurement parameter, such as an SPR dip position as a function of wavelength which would lead to a unique identification of $`d`$ and $`n`$ multil , or to reduce the number of variables by adding the assumption that the layer is made of a homogeneous mixture of a proportion $`x`$ of protein and $`1x`$ of water. The layer density $`\rho _{layer}`$ is then $`\rho _{layer}=x\rho _{protein}+(1x)\rho _{water}`$ where $`\rho _{protein}`$ is assumed to be equal in the $`1.2`$ g/cm3 to 1.4 g/cm3 range caruso ; kasemo ; hook2 ; marsch  and $`\rho _{water}=1`$ g/cm3 is the density of water. Similarly, we use then the optical index of the layer $`n_{layer}=xn_{protein}+(1x)n_{water}`$ with the optical index of the protein layer assumed to be in the $`n_{protein}=1.45`$ to $`1.465`$ range hook2 ; marsch ; liedberg  and that of water $`n_{water}=1.33`$. The remaining unknown parameters are then the layer thickness $`d`$ and the water content proportion $`1x`$. By simulating a stack of planar multilayers (glass $`n_{glass}=1.518`$, 1200 nm silicon dioxide $`n_{SiO_2}=1.45`$, 2 nm titanium $`n_{Ti}=2.76+i3.84`$ palik , 50 nm gold $`n_{Au}=0.14+i3.697`$ palik , proteins, water) for the angle shift as a function of water content and layer thickness, one obtains a set of pairs of values for these two parameters compatible with the observed angle shift (Fig. 4, bottom). Then, by calculating the mass per unit area $`\mathrm{\Delta }m=\rho d`$, a comparison with the experimentally observed adsorbed mass as seen from the phase shift of the SAW device leads to a unique set of parameters both compatible with the optical index change and the adsorbed mass (Fig. 4, top). In the case under consideration here, an observed angle shift of 380 to 400 mo and an adsorbed mass of 540 to 580 ng/cm2 is only compatible with $`x=0.75\pm 0.15`$ and $`d=4.7\pm 0.7`$ nm. The thickness result is compatible with atomic force microscope (AFM) measurements in liquid erika  while the mass per unit area is compatible with morphological data obtained by electronic microscopy review . While the main source of uncertainty is due to the wide possible values of the density and optical index reported in the literature, the experimental results display good reproducibility. The water content is lower than that observed by other authors for different kinds of proteins kasemo ; marsch  but compatible with the 30 to 78% value cited in Ref. stenberg . The dense packing of protein is explained by the regular arrangement of identical S-layer subunits in the p4 lattice weygand1 ; weygand2 .

## IV Conclusion

We have shown here how the combination of Love mode SAW device with SPR provides advantageous combinations of information on the bound mass as well as on the dielectric changes when monitoring protein adsorption. The resulting protein layer thickness and protein content percentage, respectively $`d=4.7\pm 0.7`$ nm and $`x=75\pm 15`$ %, is in agreement with independent AFM estimates.

## V Acknowledgments

The SAW transducers were fabricated within the framework of a Belgian PhD scholarship program (Fonds pour la formation à la Recherche dans l’Industrie et dans l’Agriculture – FRIA). We wish to thank R. Giust (LOPMD, Besançon, France) for kindly providing SPR simulation routines.


text:
# A Micro-Canonical Description of Hadron Production in Proton-Proton Collisions

## Abstract

A micro-canonical treatment is used to study particle production in pp collisions. First this micro-canonical treatment is compared to some canonical ones. Then proton, antiproton and pion 4$`\pi `$ multiplicities from proton-proton collisions at various center of mass energies are used to fix the micro-canonical parameters $`E`$ and $`V`$. The dependences of the micro-canonical parameters on the collision energy are parameterised for the further study of pp reactions with this micro-canonical treatment.

1Institute of Particle Physics, Central China Normal University, Wuhan, China

2Institut fuer Theoretische Physik, J. W. Goethe-Universitaet, Frankfurt am Main, Germany

3Laboratoire SUBATECH, Universite de Nantes - IN2P3/CNRS \- Ecole des Mines de Nantes, Nantes, France

In recent years, the use of thermal statistical models to describe particle ratios for heavy ion collisions has received great attention. Thermal models are even used to describe proton-proton and $`\mathrm{e}^+\mathrm{e}^{}`$ collisions. The corresponding volumes for these reactions are very small ($`V25fm^3`$). For small volumes, energy and flavour conservation become important. Therefore a micro-canonical treatment should be employed. With the increase of the volume, the micro-canonical results should converge towards the canonical one.

The following micro-canonical treatment  is employed. Following the general philosophy of statistical approaches to hadron production, we suppose that the result of a high energy collision can be considered as a distribution of “clusters”, “droplets”, or “fireballs”, which move relative to each other. Here, we are only interested in $`4\pi `$ particle yields and average transverse momenta, and the distribution of clusters may be identified with one single “equivalent cluster”, being characterised by its volume $`V`$ (the sum of individual proper volumes), its energy $`E`$ (the sum of all the cluster masses), and the net flavour content $`Q=(N_uN_{\overline{u}},N_dN_{\overline{d}},N_sN_{\overline{s}})`$. The basic assumption is that a cluster, characterised by $`V`$, $`E`$, and $`Q`$, decays “statistically” according to phase space. More precisely, the probability of a cluster to hadronize into a configuration $`K=\{h_1,\mathrm{},h_n\}`$ of hadrons $`h_i`$ is given by the micro-canonical partition function $`\mathrm{\Omega }(K)`$ of an ideal, relativistic gas of the $`n`$ hadrons $`h_i`$,

$$\mathrm{\Omega }(K)=\frac{V^n}{(2\pi \mathrm{})^{3n}}\underset{i=1}{\overset{n}{}}g_i\underset{\alpha 𝒮}{}\frac{1}{n_\alpha !}\underset{i=1}{\overset{n}{}}d^3p_i\delta (E\mathrm{\Sigma }\epsilon _i)\delta (\mathrm{\Sigma }\stackrel{}{p}_i)\delta _{Q,\mathrm{\Sigma }q_i},$$

with $`\epsilon _i=\sqrt{m_i^2+p_i^2}`$ being the energy, and $`p_i`$ being the 3-momentum of particle $`i`$. The term $`\delta _{Q,\mathrm{\Sigma }q_i}`$ ensures flavour conservation; $`q_i`$ is the flavour vector of hadron $`i`$. The symbol $`𝒮`$ represents the set of hadron species considered: we take $`𝒮`$ to contain the pseudoscalar and vector mesons $`(\pi ,K,\eta ,\eta ^{},\rho ,K^{},\omega ,\varphi )`$ and the lowest spin-$`\frac{1}{2}`$ and spin-$`\frac{3}{2}`$ baryons $`(N,\mathrm{\Lambda },\mathrm{\Sigma },\mathrm{\Delta },\mathrm{\Sigma }^{},\mathrm{\Xi }^{},\mathrm{\Omega })`$ and the corresponding antibaryons. $`n_\alpha `$ is the number of hadron species $`\alpha `$, and $`g_i`$ is the degeneracy of particle $`i`$.

For the present investigation, we limit the number of hadrons in $`𝒮`$ to 54. Strange particles are produced according to phase space, i.e. without applying any artificial suppression factor. The results are compared to that of the canonical calculations using the approach of Becattini et al. , but without strangeness suppression and the number of hadrons which can be produced is limited to the same 54 species allowed in our micro-canonical treatment. In principle it is possible to include more hadrons in the micro-canonical ensemble. This makes a detailed comparison more difficult, because the less known decay channels of the additional hadrons have to agree.

For micro-canonical treatment, there exists a critical volume $`V_c`$, above which the physics results, i.e. particle yield per volume and average transverse momentum, become independent of volume and coincide to the canonical ones. In the following we check the volume effect in particle yields, compare the micro-canonical results to canonical one and try to find the critical volume $`V_c`$. To do this calculation, we keep the same baryon density as input, i.e. the baryon number is two (a proton-proton collision) for the volume $`25\mathrm{fm}^3`$ which is usually used by canonical calculations (dashed lines), baryon number one for volume $`12.5\mathrm{fm}^3`$ (solid lines) and baryon number four for volume $`50\mathrm{fm}^3`$ (dotted lines), c.f. Figs 1, 2, 3. The results of a canonical calculation with the parameters fitted to describe the particle yields observed in pp collisions at 27.4 GeV  are marked as dot points. The parameters used in this canonical calculations are $`V=25.5\mathrm{fm}^3`$ and $`T=162\mathrm{MeV}`$. In contradistinction to the calculation with this parameter which is presented in ref.  here the strange particles are not suppressed by a $`\gamma _s`$ factor and only the 54 hadron species mentioned above are produced. Therefore this canonical calculation can directly be compared with our micro-canonical approach. The average energy of the clusters obtained in this calculation is 8.74 GeV, resulting in an average energy density of $`\epsilon =0.342\mathrm{GeV}/\mathrm{fm}^3`$. More detailed explanation of Figs 1, 2,3 can be found in .

Figs 1, 2, 3 tell us, the lighter is the partcile, the smaller is the critical volume $`V_c`$. 90% newly-produced hadrons from the collisions are pions while the critical volume $`V_c`$ is not bigger than $`12.5\mathrm{fm}^3`$. However, for the heavy particles, i.e. $`\mathrm{\Delta }^{++}`$ and $`\mathrm{\Omega }`$, the critical volume $`V_c`$ is bigger than $`50\mathrm{fm}^3`$. The canonical results agree well with the micro-canonical one at the same volume, $`25\mathrm{fm}^3`$. Therefore we raise the question if it is safe to study the heavy particles from small reaction systems with canonical models.

Canonical calculations get volume $`25\mathrm{fm}^3`$ from fitting particle yields. Here we also fit the $`4\pi `$ particle yields from pp collisions to fix the micro-canonical parameters $`E`$ and $`V.`$ The micro-canonical model does not have any strangeness suppression factor. It fails to describe non-strange and strange hadrons at the same time. So we only consider non-strange hadron production here. To fit the data, we minimize

$$\chi ^2=\frac{1}{\alpha }\underset{j=1}{\overset{\alpha }{}}\frac{[n_{\mathrm{exp},j}(\sqrt{s})n_{\mathrm{th},j}(E,V)]^2}{\sigma _j^2}$$

to determine the parameter *E* and *V* for at each collision energy $`\sqrt{s}`$, where $`n_{\mathrm{exp},j}`$, $`n_{\mathrm{th},j}`$ and $`\sigma _j`$ are respectively the experimental multiplicity, multiplicity from the micro-canonical treatment and experimental variance of particle $`j`$. Both the experimental and theoretical multiplicities here are after decays. We take the parameterization of energy dependence of $`\pi ^+`$, $`\pi ^{}`$, proton and antiproton instead of the data to fix the micro-canonical parameters $`E`$ and $`V`$ for a convenient reason.

The minimum $`\chi ^2`$ method produces the micro-canonical parameters $`E`$ and $`V`$ dependence on collision energy $`\sqrt{s}`$ as shown in Fig. 4. $`E`$ and $`V`$ increase monotonously with $`\sqrt{s}`$ . We parameterize the micro-canonical parameters $`E`$ (GeV) and $`V`$ ($`\mathrm{fm}^3`$) dependence on the collision energy $`\sqrt{s}`$ (GeV) for some further study of pp cillisions with this canonical treatment:

$`E`$ $`=`$ $`3.8+3.76\mathrm{ln}\sqrt{s}+6.4/\sqrt{s}`$ (1)
$`V`$ $`=`$ $`30.04+14.93\mathrm{ln}\sqrt{s}0.013\sqrt{s}`$ (2)
Note that $`V_{\mathrm{min}}=1\mathrm{fm}^3`$ has been used for Eq. 2 at very low energies. The Eq. 1 itself is the same as charged hadron excitation function in pp collisions. This agrees well with the results in heavy ion collsions where the freeze-out energy per particle $`E/N=1\mathrm{GeV}`$, taking into account that the increase of the multiplicity from decay has roughly the same fraction of neutral particles in the total hadron multiplicity.

Acknowledgment

We would like to thank F. Becattini for many discussions and for providing us with the calculation with a restricted set of hadrons and without strangeness suppression.


text:
# Lower bounds on the complexity of simulating quantum gates

## I Introduction

Understanding quantum dynamics is at the heart of quantum physics. Recent ideas from quantum computation have stimulated interest in studying the physical resources needed to implement quantum operations. In addition to a qualitative understanding of what resources are necessary, we would like to quantify the resource requirements for universal quantum computation and other information processing tasks. Ultimately, we would like to understand the minimal resources that are necessary and sufficient to implement particular quantum dynamics.

As a first step towards answering these questions, it has been shown that there is a sense in which all entangling dynamics are qualitatively equivalent. In particular, it has been shown that any $`n`$-qudit two-body Hamiltonian capable of creating entanglement between any pair of qudits is, in principle, universal for quantum computation, when assisted by arbitrary single-qudit unitaries Jones and Knill (1999); Leung et al. (2000); Dodd et al. (2002); Wocjan et al. (2002); Dür et al. (2001); Bennett et al. (2002); Nielsen et al. (2002); Vidal and Cirac (2002). Thus, any particular entangling two-qudit Hamiltonian can be used to simulate any other, provided local unitaries are available. This suggests that such dynamics are a fungible physical resource.

Having established the qualitative equivalence of all entangling dynamics, we would like to quantify their information processing power. In particular, it is interesting to consider the minimal time required to implement a unitary operation, $`U`$, on a two-qubit system, using a fixed Hamiltonian, $`H`$, and the ability to intersperse fast local unitary operations on the two qubits. This problem was studied by Khaneja, Brockett and Glaser Khaneja et al. (2001), who found a solution using the theory of Lie groups. Their results, although giving a solution in principle, are neither explicit about the form of the minimal time, nor do they explain how to construct all elements of the time-optimal simulation. Further work by Vidal, Hammerer, and Cirac Vidal et al. (2002), from a different point of view, resulted in an explicit formula for the minimal time, and gave a constructive procedure for minimizing that time (see also Hammerer et al. (2002), where an alternate proof is given by the same authors).

The purpose of the present paper is to give a simplified proof that the formula of Vidal, Hammerer and Cirac is, in fact, a lower bound on the simulation time. Note that the difficult part of Vidal et al. (2002); Hammerer et al. (2002) was proving the lower bound; finding a protocol to meet the lower bound was comparatively easy.

The main advantages of our proof are its simplicity and conceptual clarity, as compared to the ingenious, but rather complex, arguments in Khaneja et al. (2001); Vidal et al. (2002); Hammerer et al. (2002). This simplicity is achieved by making use of a powerful result from linear algebra, Thompson’s theorem. We expect that Thompson’s theorem might be useful for many other problems in quantum information theory. A second advantage of using Thompson’s theorem is that it does not rely on special properties of two-qubit unitary operators. Therefore, essentially the same arguments give a lower bound on the time required to implement an $`n`$-qubit unitary operation using a fixed $`n`$-qubit interaction Hamiltonian, and fast local unitary operations.

Our approach to the proof of the lower bound has its roots in the framework of *dynamic strength measures* for quantum operations Nielsen et al. (2003). The dynamic strength framework is an attempt to develop a quantitative theory of the power of dynamical operations for information processing. The idea is to associate with a quantum dynamical operation, such as a unitary operation $`U`$, a quantitative measure of its “strength.” In Nielsen et al. (2003) it was shown that such strength measures can be used to analyze the minimal time required for the implementation of a quantum operation. The present paper takes a similar approach, but instead of using a single real number to quantify dynamic strength, we use a vector-valued measure. This can also be compared to the analysis of optimal simulation of Hamiltonian dynamics using a set of several strength measures Childs et al. (2003).

Our paper is structured as follows. Section II reviews some background material on majorization, Thompson’s theorem, and the structure of the two-qubit unitary matrices. The main result, the lower bound on optimal simulation, is proved in Section III. We conclude in Section IV by presenting our generalization of the lower bound to $`n`$ qubits and suggesting some directions for future work. In addition, an appendix gives a procedure for calculating a canonical decomposition of two-qubit unitary gates.

## II Background

This section reviews the relevant background needed for our proof. Section II.1 reviews the basic notions of majorization, introduces Thompson’s theorem, and explains how to use Thompson’s theorem and majorization to relate properties of a product of unitary operators to properties of the individual unitaries. Section II.2 introduces the *canonical decomposition*, a useful representation theorem for two-qubit unitary operators, and Section II.3 presents an analogous decomposition for Hamiltonians.

### II.1 Majorization and Thompson’s theorem

Our analysis uses the theory of majorization together with Thompson’s theorem. More detailed introductions to majorization may be found in Nielsen and Vidal (2001), Chapter 2 and 3 of Bhatia (1997), and in Marshall and Olkin (1979); Alberti and Uhlmann (1982).

Suppose $`x=(x_1,\mathrm{},x_D)`$ and $`y=(y_1,\mathrm{},y_D)`$ are two $`D`$-dimensional real vectors. The relation *$`x`$ is majorized by $`y`$*, written $`xy`$, is intended to capture the intuitive notion that $`x`$ is *less ordered* (i.e., more disordered) than $`y`$. To make the formal definition we introduce the notation $``$ to denote the components of a vector rearranged into non-increasing order, so $`x^{}=(x_1^{},\mathrm{},x_D^{})`$, where $`x_1^{}x_2^{}\mathrm{}x_D^{}`$. Then $`x`$ is majorized by $`y`$, that is, $`xy`$, if

$$\underset{j=1}{\overset{k}{}}x_j^{}\underset{j=1}{\overset{k}{}}y_j^{}$$
(1)
for $`k=1,\mathrm{},D1`$, and the inequality holds with equality when $`k=D`$.

To connect majorization to Hamiltonian simulation, we use a result of Thompson relating a product of two unitary operators to the individual unitary operators. Recall that an arbitrary pair of unitary operators can be written in the form $`e^{iH}`$ and $`e^{iK}`$, for some Hermitian $`H`$ and $`K`$. Thompson’s theorem provides a representation for the product $`e^{iH}e^{iK}`$ in terms of $`H`$ and $`K`$:

###### Theorem 1 (Thompson Thompson (1986)).

Let $`H,K`$ be Hermitian matrices. Then there exist unitary matrices $`U,V`$ such that

$$e^{iH}e^{iK}=e^{i(UHU^{}+VKV^{})}.$$
(2)
The proof of Thompson’s theorem in Thompson (1986) depends on a result conjectured earlier by Horn Horn (1962). A proof of this conjecture had been announced and outlined by Lidskii Lidskii (1982) at the time of Thompson’s paper. However, remarks in Thompson (1986) suggest that Lidskii (1982) did not contain enough detail to be considered a fully rigorous proof. Fortunately, a proof of Horn’s conjecture has recently been fully completed and published. See, for example, Fulton (2000); Knutson (2002) for reviews and references.

Thompson’s theorem may be related to majorization using the following theorem of Ky Fan:

###### Theorem 2 (Ky Fan Fan (1949); Bhatia (1997)).

Let $`H,K`$ be Hermitian matrices. Then $`\lambda (H+K)\lambda (H)+\lambda (K)`$, where $`\lambda (A)`$ denotes the vector whose entries are the eigenvalues of the Hermitian matrix $`A`$, arranged into non-increasing order.

Combining the results of Ky Fan and Thompson, we have the following:

###### Corollary 3.

Let $`H,K`$ be Hermitian matrices. Then there exists a Hermitian matrix $`L`$ such that

$$e^{iH}e^{iK}=e^{iL};\lambda (L)\lambda (H)+\lambda (K).$$
(3)
We will not apply this corollary directly, but we have included it here because it captures the spirit of our later argument, combining the Thompson and Ky Fan theorems to relate the properties of a product of unitaries to the individual unitaries themselves. Corollary 3 can be regarded as a vector-valued analogue of the chaining property for dynamic strength measures used in Nielsen et al. (2003) to establish lower bounds on computational complexity.

### II.2 The canonical decomposition of a two-qubit gate

The *canonical decomposition* is a useful representation theorem characterizing the non-local properties of a two-qubit unitary operator. It was proved by Khaneja, Brockett, and Glaser Khaneja et al. (2001) using ideas from Lie theory. Kraus and Cirac Kraus and Cirac (2001) have given a constructive proof using elementary notions, while Zhang *et al.* Zhang et al. (2003) have discussed the decomposition in detail from the point of view of Lie theory. The decomposition states that any two-qubit unitary $`U`$ may be written in the form

$$U=(A_1B_1)e^{i(\theta _xXX+\theta _yYY+\theta _zZZ)}(A_2B_2),$$
(4)
where $`A_1`$, $`A_2`$, $`B_1`$, $`B_2`$ are single-qubit unitaries, and the three parameters, $`\theta _x`$, $`\theta _y`$, and $`\theta _z`$ characterize the non-local properties of $`U`$.111Prior to Khaneja et al. (2001), Makhlin Makhlin (2002) gave a proof that the non-local properties of $`U`$ are completely characterized by $`\theta _x`$, $`\theta _y`$ and $`\theta _z`$, but did not write down the canonical decomposition explicitly. Without loss of generality, we may choose the local unitaries to ensure that

$$\frac{\pi }{4}\theta _x\theta _y|\theta _z|,$$
(5)
and we refer to the set of parameters chosen in this way as the *canonical parameters* for $`U`$. We will see below that these parameters are unique. We define the canonical form of $`U`$ to be

$$U_c:=(A_1^{}B_1^{})U(A_2^{}B_2^{});$$
(6)
up to local unitaries, $`U_c`$ is equivalent to $`U`$. It will be convenient to assume through the remainder of this section that $`U`$ has unit determinant. This is equivalent to requiring that $`A_1,A_2,B_1,B_2`$ can all be chosen to have unit determinant.

The canonical parameters turn out to be crucial to results about simulation of two-qubit gates. If

$$U_c=e^{i(\theta _xXX+\theta _yYY+\theta _zZZ)}$$
(7)
is the canonical form of $`U`$, then we define the *non-local content*, $`\varphi (U)`$, of $`U`$ by $`\varphi (U):=\lambda (H_U)`$, where

$$H_U:=\theta _xXX+\theta _yYY+\theta _zZZ.$$
(8)
Explicitly, the components of $`\varphi (U)`$ are

$`\varphi _1`$ $`=`$ $`\theta _x+\theta _y\theta _z`$ (9)
$`\varphi _2`$ $`=`$ $`\theta _x\theta _y+\theta _z`$ (10)
$`\varphi _3`$ $`=`$ $`\theta _x+\theta _y+\theta _z`$ (11)
$`\varphi _4`$ $`=`$ $`\theta _x\theta _y\theta _z.`$ (12)
We now outline a simple procedure to determine the canonical parameters of a two-qubit unitary operator. Our explanation initially follows Hammerer et al. (2002) and Leifer et al. (2002). However, as explained below, there is an ambiguity in the procedure described in those papers, related to the fact that the logarithm function has many branches. Our procedure resolves this ambiguity.

To explain the procedure, we need to introduce a piece of notation, and explain a simple observation about single-qubit unitary matrices. The *spin flip* operation on an arbitrary two-qubit operator is defined as

$$\stackrel{~}{M}:=(YY)M^T(YY),$$
(13)
where $`Y`$ is the Pauli sigma $`y`$ matrix, and the transpose operation is taken with respect to the computational basis. Note that the spin flip operation may also be written as $`\stackrel{~}{M}=M^T`$, where the transpose is taken with respect to a different basis, the *magic basis* Hill and Wootters (1997),

$`{\displaystyle \frac{|00+|11}{\sqrt{2}}};`$ $`i{\displaystyle \frac{|00|11}{\sqrt{2}}};`$
$`i{\displaystyle \frac{|01+|10}{\sqrt{2}}};`$ $`{\displaystyle \frac{|01|10}{\sqrt{2}}}.`$ (14)
The observation about single-qubit unitary matrices that we need is the following. Let $`U`$ be any single-qubit unitary matrix with unit determinant. Then

$`UYU^T=Y,`$ (15)
where the transpose is taken in the computational basis. This simple identity is easily verified.

Now suppose $`U`$ is an arbitrary two-qubit unitary with unit determinant. By definition of the spin flip, and substituting the canonical decomposition, we have

$`U\stackrel{~}{U}`$ $`=`$ $`(A_1B_1)U_c(A_2B_2)(YY)`$ (16)
$`\times (A_2^TB_2^T)U_c(A_1^TB_1^T)(YY).`$
By the identity Eq. (15) we see that

$`U\stackrel{~}{U}=(A_1B_1)U_c(YY)U_c(A_1^TB_1^T)(YY).`$ (17)
Using the fact that $`YY`$ commutes with $`XX`$, $`YY`$, and $`ZZ`$, we see that $`YY`$ commutes with $`U_c`$, and thus

$`U\stackrel{~}{U}=(A_1B_1)U_c^2(YY)(A_1^TB_1^T)(YY).`$ (18)
Finally, applying Eq. (15) again gives

$`U\stackrel{~}{U}=(A_1B_1)U_c^2(A_1^{}B_1^{})`$ (19)
Eq. (19) suggests a procedure to determine the canonical parameters for $`U`$, based on the observation that

$$\lambda (U\stackrel{~}{U})=\lambda (U_c^2)=(e^{2i\varphi _1},e^{2i\varphi _2},e^{2i\varphi _3},e^{2i\varphi _4}),$$
(20)
where the $`\varphi _j`$ are related to the canonical parameters $`\theta _x,\theta _y`$ and $`\theta _z`$ by Eqs. (9)–(12). It is tempting to conclude that one can determine $`\theta _x,\theta _y,\theta _z`$ from the eigenvalues of $`U\stackrel{~}{U}`$, simply by taking logarithms and inverting the resulting linear equations. Indeed, such a conclusion is reached in Hammerer et al. (2002) and Leifer et al. (2002), using arguments similar to those just described. Unfortunately, determining the canonical parameters is not quite as simple as this, because $`ze^{iz}`$ is not a uniquely invertible function. In particular, $`e^{iz}=e^{i(z+2\pi m)}`$, where $`m`$ is any integer, so there is some ambiguity about which branch of the logarithm function to use in calculating the canonical parameters. In fact, we prove later that no one branch of the logarithm function can be used. However, these considerations do allow us to reach the following conclusion:

###### Lemma 4.

Let $`U`$ be a two-qubit unitary. Then there exists a Hermitian $`H`$ such that

$`U\stackrel{~}{U}=e^{2iH},\lambda (H)=\varphi (U).`$ (21)
Moreover, if $`H`$ is any Hermitian matrix such that $`\lambda (U\stackrel{~}{U})=\lambda (e^{2iH})`$ then it follows that $`\lambda (H)=\varphi (U)+\pi \stackrel{}{m}`$, where $`\stackrel{}{m}`$ is some vector of integers.

Although this lemma is sufficient to prove our later results, there is in fact a simple method for exactly calculating the canonical parameters. Because there are many applications of the canonical decomposition, we describe this method in the appendix. The method will not be needed elsewhere in the paper.

### II.3 The canonical form of a two-qubit Hamiltonian

Finally, we introduce one additional concept, the *canonical form* of a two-qubit Hamiltonian, $`H`$ Dür et al. (2001). Any two-qubit Hamiltonian $`H`$ can be expanded as

$$H=\underset{j,k=0}{\overset{3}{}}h_{jk}\sigma _j\sigma _k.$$
(22)
Then let

$$H^{}:=\frac{H+\stackrel{~}{H}}{2}=\underset{j,k0}{}h_{jk}\sigma _j\sigma _k.$$
(23)
That is, $`H^{}`$ is just the Hamiltonian that results when the local terms in $`H`$ are removed. It is not difficult to show that $`H`$ and $`H^{}`$ are interchangeable resources for simulation in the sense that, given fast local unitaries, evolution according to $`H`$ for a time $`t`$ can be simulated by evolution according to $`H_c`$ for a time $`t`$, and vice versa. Furthermore, by doing appropriate local unitaries, it can be shown Dür et al. (2001) that simulating $`H^{}`$ (and thus $`H`$) is equivalent to simulating the canonical form of $`H`$,

$`H_c=h_xXX+h_yYY+h_zZZ,`$ (24)
where $`h_xh_y|h_z|`$. Once again, $`H`$ and $`H_c`$ are interchangeable resources for simulation.

Note that the three parameters $`h_x,h_y,h_z`$ are completely characterized by the three degrees of freedom in $`\lambda (H_c)=\lambda (H+\stackrel{~}{H})/2`$, just as the three canonical parameters $`\theta _x,\theta _y,\theta _z`$ are completely characterized by the three degrees of freedom in $`\lambda (U_c^2)=\lambda (U\stackrel{~}{U})`$.

## III Simulation of two-qubit gates

We now return to the main purpose of the paper, proving results about the time to simulate a unitary gate using entangling Hamiltonians and fast local gates. We aim to prove the following result:

###### Theorem 5 (Vidal, Hammerer, Cirac Vidal et al. (2002); Hammerer et al. (2002), cf. Khaneja, Brockett and Glaser Khaneja et al. (2001)).

Let $`U`$ be a two-qubit unitary operator, and let $`H`$ be a two-qubit entangling Hamiltonian. Then the minimal time required to simulate $`U`$ using $`H`$ and fast local unitaries is the minimal value of $`t`$ such that there exists a vector of integers $`\stackrel{}{m}`$ satisfying

$`\varphi (U)+\pi \stackrel{}{m}{\displaystyle \frac{\lambda (H+\stackrel{~}{H})}{2}}t.`$ (25)
Note further that only two vectors of integers need to be checked, $`\stackrel{}{m}=(0,0,0,0)`$ and $`\stackrel{}{m}=(1,1,1,1)`$, since all the other possibilities give rise to weaker constraints on the minimal time, $`t`$ Vidal et al. (2002); Hammerer et al. (2002). The difficult part of the proof of Theorem 5 is the proof that Eq. (25) is a lower bound on the simulation time, $`t`$, and it is this part of the proof that we focus on simplifying. The proof that this lower bound may be achieved follows from standard results on majorization, and we refer the interested reader to Vidal et al. (2002); Hammerer et al. (2002) for details.

To prove that Eq. (25) constrains the minimal time for simulation, we begin by characterizing the canonical decomposition of a product of unitary matrices. Let $`\mathrm{\Lambda }(U):=\lambda (U\stackrel{~}{U})`$, and define the equivalence relation $`AB`$ for Hermitian matrices $`A`$ and $`B`$ iff $`\lambda (A)=\lambda (B)`$. Then we have:

###### Lemma 6.

Let $`U_j`$ be unitary matrices, and let $`H_j`$ be Hermitian matrices such that $`U_j\stackrel{~}{U}_j=e^{2iH_j}`$. Then there exist Hermitian matrices $`K_j`$ such that $`H_jK_j`$, and

$$\mathrm{\Lambda }(U_N\mathrm{}U_1)=\lambda (e^{2i(K_1+\mathrm{}+K_N)}).$$
(26)
###### Proof.

We induct on $`N`$. The result is trivial for $`N=1`$, so we need only consider the inductive step. Using the fact $`\lambda (AB)=\lambda (BA)`$, we have

$$\mathrm{\Lambda }(U_{N+1}\mathrm{}U_1)=\lambda (\stackrel{~}{U}_{N+1}U_{N+1}U_N\mathrm{}U_1\stackrel{~}{U}_1\mathrm{}\stackrel{~}{U}_N).$$
(27)
By the inductive hypothesis there exist Hermitian $`K_j^{}`$ such that $`H_jK_j^{}`$ and

$$\lambda (U_N\mathrm{}U_1\stackrel{~}{U}_1\mathrm{}\stackrel{~}{U}_N)=\lambda (e^{2i(K_1^{}+\mathrm{}+K_N^{})}).$$
(28)
Therefore, $`U_N\mathrm{}U_1\stackrel{~}{U}_1\mathrm{}\stackrel{~}{U}_N=e^{2i(K_1^{\prime \prime }+\mathrm{}+K_N^{\prime \prime })}`$, for some $`K_j^{\prime \prime }H_j`$. Observe also that

$$\stackrel{~}{U}_{N+1}U_{N+1}U_{N+1}\stackrel{~}{U}_{N+1}=e^{2iH_{N+1}},$$
(29)
and thus $`\stackrel{~}{U}_{N+1}U_{N+1}=e^{2iK_{N+1}^{\prime \prime }}`$ for some $`K_{N+1}^{\prime \prime }H_{N+1}`$. It follows by substitution that

$$\mathrm{\Lambda }(U_{N+1}\mathrm{}U_1)=\lambda (e^{2iK_{N+1}^{\prime \prime }}e^{2i(K_1^{\prime \prime }+\mathrm{}+K_N^{\prime \prime })}).$$
(30)
Applying Thompson’s theorem gives

$$\mathrm{\Lambda }(U_{N+1}\mathrm{}U_1)=\lambda (e^{2i(K_1+\mathrm{}+K_{N+1})})$$
(31)
for some $`K_jK_j^{\prime \prime }H_j`$, which completes the inductive step of the proof. ∎

Given this result, it is straightforward to complete the proof of Eq. (25).

###### Proof.

Write $`U`$ in the form

$`U=e^{iHt_1}V_1e^{iHt_2}V_2\mathrm{}V_{k1}e^{iHt_k},`$ (32)
where $`t_1,\mathrm{},t_k`$ are times of evolution, $`t=t_1+\mathrm{}+t_k`$ is the total time for simulation, and $`V_j`$ are local unitaries. Without loss of generality, we may assume $`H`$ is in canonical form. Applying Lemma 6, we obtain

$`\mathrm{\Lambda }(U)=\lambda (e^{2i(H_1t_1+\mathrm{}+H_kt_k)})`$ (33)
where $`H_jH`$ for each $`j`$. Here we have used the observation $`V_j\stackrel{~}{V}_j=1`$, so all the contributions from local unitaries vanish. It follows from Lemma 4 that

$`\varphi (U)+\pi \stackrel{}{m}=\lambda (H_1t_1+\mathrm{}+H_mt_m),`$ (34)
and using Ky Fan’s theorem gives

$`\varphi (U)+\pi \stackrel{}{m}\lambda (H)(t_1+\mathrm{}+t_m),`$ (35)
which is Eq. (25), as desired. ∎

## IV Discussion

In this paper, we have provided a simplified proof of a lower bound on the time required to simulate a two-qubit unitary gate using a given two-qubit interaction Hamiltonian and local unitaries. The bound follows easily from standard results on majorization together with Thompson’s theorem on products of unitary operators.

Although we have described canonical decompositions of two-qubit gates in some detail, we note that our proof does not actually require properties of the decomposition unique to two qubits. In fact, it is straightforward to prove an analogue of Eq. (25) for an $`n`$-qubit system. For an $`n`$-qubit operator $`M`$, suppose we define a generalized spin flip $`M\stackrel{~}{M}`$, where $`\stackrel{~}{M}`$ is the transpose operation in a basis such that, whenever $`M`$ is local, $`M`$ is orthogonal, i.e., $`M\stackrel{~}{M}=I`$. It is not difficult to construct examples of such bases, at least when $`n`$ is even. An example is the basis obtained by rotating the computational basis using the transformation $`(IiY^n)/\sqrt{2}`$, for $`n`$ even. This basis change gives $`\stackrel{~}{M}=Y^nM^TY^n`$, where the transpose is taken in the computational basis, and thus this operation generalizes the transpose in the magic basis. In this general setting the following lower bound on the time required to implement an $`n`$-qubit gate holds:

###### Corollary 7.

Let $`U`$ be an $`n`$-qubit unitary operator, and let $`H`$ be an $`n`$-qubit Hamiltonian. Then the time required to simulate $`U`$ using $`H`$ and fast local unitaries satisfies

$`{\displaystyle \frac{1}{2}}\mathrm{arg}\lambda (U\stackrel{~}{U})+\pi \stackrel{}{m}{\displaystyle \frac{\lambda (H+\stackrel{~}{H})}{2}}t.`$ (36)
for some vector of integers $`\stackrel{}{m}`$.

The proof follows simply by taking the arguments of both sides of Eq. (33) and applying Ky Fan’s theorem. All steps leading up to Eq. (33) remain valid for $`n`$-qubit systems using the above definition of the generalized spin flip.

Unfortunately, we have not found any interesting examples with $`n>2`$ for which Eq. (36) provides a nontrivial lower bound on the time required to implement some quantum gate. It would be interesting to construct cases where Eq. (36) (or some similar condition) does give a nontrivial constraint on multipartite gate simulation. One might imagine that such techniques could be used to prove circuit lower bounds on certain quantum computations, although it does not seem likely that such bounds would be especially strong, given the well-known difficulty of this problem.

###### Acknowledgements.

We thank Aram Harrow and Tobias Osborne for helpful discussions, and Andrew Doherty for an informative seminar on related problems. AMC received support from the Fannie and John Hertz Foundation, and thanks the University of Queensland node of the Centre for Quantum Computer Technology for its hospitality. AMC was also supported in part by the Cambridge–MIT Institute, by the Department of Energy under cooperative research agreement DE-FC02-94ER40818, and by the National Security Agency and Advanced Research and Development Activity under Army Research Office contract DAAD19-01-1-0656. Finally, we acknowledge the hospitality of the Caltech Institute for Quantum Information, where this work was completed. This work was supported in part by the National Science Foundation under grant EIA-0086038.

## Appendix: A method for computing the canonical parameters of a two-qubit unitary gate

In this appendix, we describe a method for computing the canonical parameters of a two-qubit unitary, based on the discussion in Section II.2. The key is to take logarithms in just the right way. From Eqs. (5) and (9)–(12), we see that

$$\frac{3\pi }{2}2\varphi _12\varphi _22\varphi _32\varphi _4\frac{3\pi }{2}.$$
(37)
It is not difficult to find examples where the first or last inequality is saturated, so no single fixed branch of the logarithm function can be used to determine the $`\varphi _j`$. One might hope instead that there exists a method for choosing a different branch for each particular $`U`$, so that the corresponding $`2\varphi _j`$ lie within that branch. However, even this is not possible in general. To understand this, note that

$$2\varphi _12\varphi _4=4(\theta _x+\theta _y).$$
(38)
In cases where $`\theta _x=\theta _y=\pi /4`$, we have $`2\varphi _12\varphi _4=2\pi `$, in which case the values $`2\varphi _j`$ do not lie in any one branch.

We now show how to compute the $`\varphi _j`$. The idea is that we can first take the argument of the eigenvalues in Eq. (20) over some fixed branch. Then we can systematically determine which of the resulting values have been shifted by $`2\pi `$ from the value $`2\varphi _j`$ (due to an incorrect branch) and correct these values accordingly.

Let $`S_j`$, $`j=1,\mathrm{},4`$ be defined as follows:

$$2S_j=\mathrm{arg}(e^{2i\varphi _j}).$$
(39)
That is, $`2S_j`$ are the arguments of the eigenvalues of $`U\stackrel{~}{U}`$, where we take the argument over the branch $`(\frac{\pi }{2},\frac{3\pi }{2}]`$, so that the $`S_j`$ are contained in the interval $`(\frac{\pi }{4},\frac{3\pi }{4}]`$. Considering the range of values that $`\varphi _j`$ may take, from Eq. (37), and the particular branch we are using, it is clear that:

$$S_j=\{\begin{array}{ccc}\varphi _j+\pi \hfill & & \text{if }\varphi _j\frac{\pi }{4}\hfill \\ \varphi _j\hfill & & \text{otherwise.}\hfill \end{array}$$
(40)
From Eqs. (9)–(12) we have

$$\varphi _1+\varphi _2+\varphi _3+\varphi _4=0.$$
(41)
Combining Eqs. (40) and (41), we see that

$$S_1+S_2+S_3+S_4=\pi n,$$
(42)
where $`n`$ is the number of $`\varphi _j`$ that are less than or equal to $`\frac{\pi }{4}`$. Possible values for $`n`$ are $`0,1,2`$ and $`3`$ (all four $`\varphi _j`$ cannot simultaneously be $`\frac{\pi }{4}`$, since that would contradict Eq. (41)). Since the $`\varphi _j`$ obey the ordering in Eq. (37), then the $`n`$ values of $`\varphi _j`$ that are less than or equal to $`\frac{\pi }{4}`$ are $`\varphi _4,\mathrm{},\varphi _{4n+1}`$, and the remaining $`4n`$ values greater than $`\frac{\pi }{4}`$ are $`\varphi _1,\mathrm{},\varphi _{4n}`$. Thus, using Eq. (40), we see that the set of values $`S_j`$ consist of $`n`$ “shifted” $`\varphi _j`$ values

$$\varphi _4+\pi ,\mathrm{},\varphi _{4n+1}+\pi ,$$
(43)
and $`4n`$ “non-shifted” values of $`\varphi _j`$

$$\varphi _1,\mathrm{},\varphi _{4n}.$$
(44)
Furthermore, all of the shifted values in (43) are no less than any of the non-shifted values in (44). This is shown by combining Eq. (5) with Eq. (38), giving $`\varphi _1\varphi _4\pi `$, which when combined with Eq. (37) implies that $`\varphi _j\varphi _k+\pi `$ for all $`j,k`$, as required. Therefore, the largest $`n`$ values of $`S_j`$ are guaranteed to be the values in (43). Thus subtracting $`\pi `$ from the largest $`n`$ values of $`S_j`$, gives us $`\varphi _4,\mathrm{},\varphi _{4n+1}`$, and the the remaining $`4n`$ values of $`S_j`$ give us $`\varphi _1,\mathrm{},\varphi _{4n}`$.

In summary, the nonlocal parameters $`\theta _x,\theta _y`$ and $`\theta _z`$ may be computed as follows. Find the arguments of the eigenvalues of $`U\stackrel{~}{U}`$ over the branch $`(\frac{\pi }{2},\frac{3\pi }{2}]`$. Call these values $`2S_j`$. Calculate $`n=(S_1+S_2+S_3+S_4)/\pi `$. Replace the $`n`$ largest values of $`S_j`$ by those values minus $`\pi `$. The resulting values, when placed in nonincreasing order, are equal to $`(\varphi _1,\varphi _2,\varphi _3,\varphi _4)`$. The parameters $`\theta _x,\theta _y`$ and $`\theta _z`$ are then found by inverting Eqs. (9)–(12).


text:
# Submillimetre observations of a sample of broad absorption line quasars

## 1. Introduction

A small fraction ($`0.2`$Hewett & Foltz 2003; Becker et al. 2000) of quasars display broad absorption features in their spectra, with widths of several thousand $`\mathrm{km}\mathrm{s}^1`$, blue-ward of the prominent broad emission lines. Seen only in permitted lines, these absorption troughs are attributed to bulk outflows from the quasar’s heart. Typically, the BALs are seen in high ionization lines, such as $`\mathrm{Ly}_\alpha `$, $`\mathrm{C}_{\mathrm{IV}}\lambda 1549`$ and $`\mathrm{N}_\mathrm{V}\lambda 1240`$, although more rarely some systems also display low ionization lines such as $`\mathrm{Mg}_{\mathrm{II}}\lambda 2800`$; these latter systems comprise $`10\%`$ of the BALQSO population and are referred to as LoBALQSOs, whereas HiBALQSOs display solely high ionization absorption features. The BAL phenomenon in quasars appears to be part of a much broader manifestation of absorption outflows from AGN; lower velocity and narrower UV absorption features are detected more readily in Seyfert spectra and similar spectral features from ‘associated absorbers’ are also detected in some (mainly steep-spectrum) radio-loud quasars (see Crenshaw, Kraemer & George 2002 for a review). Indeed, the recently released catalog of BALQSOs from the Sloan Digital Sky Survey (Reichard et al., 2003) reveals numerous quasars possessing weak, difficult to classify absorption features, hinting at the possible ubiquitous nature of intrinsic absorption outflows.

The BAL phenomenon, and in particular its relation to the overall quasar population, has been the subject of debate for a number of years. Two distinct competing interpretations for the occurrence of BALs amongst quasar spectra are the orientation and evolution hypotheses. According to the orientation interpretation (Weymann, et al., 1991; Schmidt & Hines, 1999), all quasars possess BAL outflows with a restricted covering factor, so that the frequency of detection simply translates to the rate at which our line-of-sight intercepts the outflow. According to the evolution hypothesis (Briggs, Turnshek, & Wolfe, 1984; Voit, Weymann, & Korista, 1993), the incidence rate of the BAL phenomenon is interpreted as the duration of a phase of a quasar’s natural life cycle.

Observational evidence supporting the orientation hypothesis comes largely from spectral comparisons of BAL and non-BAL quasars (Weymann, et al., 1991) and polarization studies (Hines & Wills, 1995; Goodrich & Miller, 1995; Cohen et al., 1995; Hutsemekers, Lamy, & Remy, 1998; Schmidt & Hines, 1999). Evidence in favour of the evolution hypothesis comes largely from IRAS-selected BALQSOs (e.g. Boroson & Meyers, 1992), particularly those showing evidence for recent mergers or close interactions, consistent with proposed evolutionary scenarios for quasars (e.g. Sanders et al., 1988). Further support to the evolution hypothesis has been provided indirectly by radio observations of BALQSOs, which are inconsistent with simple orientation schemes (Becker et al., 2000).

Both hypotheses, however, also rely on observations indicating that: 1. BALQSOs are substantially more reddened than non-BAL quasars (Hall, et al., 1997; Reichard et al., 2003); 2. their continua are virtually indistinguishable from those of non-BAL quasars once dust extinction is taken into account (e.g. Yamamoto & Vansevicius, 1999); and 3. LoBALQSOs are much redder than HiBALQSOs (Sprayberry & Foltz, 1992). In the orientation scheme, these observations support the idea that BAL outflows become increasingly more dusty as the viewing angle becomes more inclined, with LoBALs seen along sight-lines grazing the putative nuclear dusty torus which plays a key role in orientation-based unification schemes for AGN. In evolutionary scenarios, the merging of two gas-rich systems fuels both vigorous star formation and a central AGN, resulting in the strong dust emission associated with an Ultraluminous Infrared Galaxy (ULIRG) phase (Sanders et al., 1988; Voit, Weymann, & Korista, 1993; Sanders & Mirabel, 1996; Canalizo & Stockton, 2001). The overall evolutionary path in this scenario is gas-rich mergers $``$ ULIRG $``$ LoBALQSO $``$ HiBALQSO $``$ unobscured quasar.

Since evolutionary models implicitly predict that BALQSOs should exhibit enhanced dust emission with respect to non-BAL quasars, one method of distinguishing between the two competing hypotheses is to measure submm luminosities. To date, only a handful of BALQSOs have been observed at submm wavelengths (e.g. Hughes, Dunlop, & Rawlings, 1997; Lewis et al., 1998; Page et al., 2001). In this paper, we present submm photometry observations of seven BALQSOs in order to determine how their dust emission properties compare to that of the general quasar and AGN population.

## 2. Observations

Our target sample was drawn from the literature111 This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. and consisted of a range of BAL systems. Given weather constraints, seven of our sample was observed. These were PSS 1537+1227 222 Currently unpublished, the position of this quasar was deduced from images available in the press as 15 37 42, 12 27 44 (J2000), 0840+3633 (Becker et al., 1997), 1104-0004 (Brotherton et al., 1998), 1556+3517 (Najita et al., 2000), 1053-0058 (Brotherton et al., 1998) and LBQS 0059-2735 (Hazard et al., 1987; Weymann, et al., 1991), whose submm flux are presented in this paper The details of our observations of Hawaii 167 (Cowie, et al., 1994; Egami et al., 1996) are presented elsewhere (Lewis & Chapman, 2000).

This sample was observed in May 2000 with the Submillimetre Common-User Bolometer Array (SCUBA) on the James Clerk Maxwell Telescope 333The James Clerk Maxwell Telescope is operated by The Joint Astronomy Centre on behalf of the Particle Physics and Astronomy Research Council of the United Kingdom, the Netherlands Organisation for Scientific Research, and the National Research Council of Canada. on Mauna Kea, Hawaii. We used the PHOTOMETRY three bolometer chopping mode described in Chapman et al. (2000) and Scott et al. (2000) to keep the source in a bolometer throughout the observation. This mode has the additional advantage of allowing a check on the apparent detection of a source over three independent bolometers. While the 450$`\mu `$m and 850$`\mu `$m arrays are illuminated simultaneously, the bolometer alignment is not perfect, and we did not include the 450$`\mu `$m offbeams in our final flux estimate, except to check that the source had offbeam flux consistent with the detection in the primary bolometer.

The observations incorporate chopping (7.8125 Hz) and nodding (every 9 seconds), and the final flux density in an individual bolometer is a double-difference with $`N=\mathrm{\hspace{0.17em}3}`$ beams. The central beam has an efficiency of unity and the two off beams have

$$ϵ=0.5\mathrm{exp}\left(\frac{d^2}{2\sigma _\mathrm{b}^2}\right),$$
(1)
where $`d`$ is the angular distance of the off-beam centre from the source, and $`\sigma _\mathrm{b}`$ is the Gaussian half-width of the beam. For the secondary bolometer the beam efficiency is simply 0.5. However, distortion in the field results in our chosen third bolometer being slightly offset from the source position, resulting in a beam efficiency of 0.44. Our 2- and 3-beam measured fluxes agree within $`<1\sigma `$ of the primary beam measurement in all non-detection cases. For the two sources we claim as detections, H167 and LBQS0059, the detection significance increases after folding in the negative flux density from the two offbeam pixels.

The effective integration time on source varied from 1200 s to 4800 s for the seven objects in our sample. The secondary was chopped with a throw of 52 arcsec to keep the source on bolometer at all times. Pointing was checked before and after the observation on blazars and a sky-dip was performed to measure the atmospheric opacity directly. The rms pointing errors were 1.4 arcsec, while the average atmospheric zenith opacities at 450$`\mu `$m and 850$`\mu `$m were 1.7 and 0.21 respectively. The data were reduced using the Starlink package SURF (Scuba User Reduction Facility, Jenness & Lightfoot 1998) and our own reduction routines to implement the three bolometer chopping mode. Spikes were first carefully rejected from the data, followed by correction for atmospheric opacity and sky subtraction using the median of all the array pixels, except for obviously bad pixels and the source pixels. The data were then calibrated against standard planetary and compact H ii  region sources, observed during the same night.

## 3. Results

Table 1 summarizes the submm photometry and dust properties of our sample of BALQSOs. There is a strong redshift dependence upon the submm properties of these BAL systems, with only upper limits for all systems below $`z1.59`$, despite deep and uniform noise limits achieved for this sample. The inferred associated dust masses, FIR fluxes and star formation rates are substantially higher for the two highest-$`z`$ sources compared to the values calculated from the upper limits for the lower-$`z`$ sources. Whilst our sample is small, these results nevertheless support the notion that the high redshift universe is considerably more dusty than the local universe. Our results also corroborate other separate submm observations of three different BALQSOs, all with $`S_{850\mu \mathrm{m}}>5`$ mJy and all lying at $`z>1.6`$ (Hughes, Dunlop, & Rawlings, 1997; Lewis et al., 1998; Page et al., 2001).

How do the submm properties of our BALQSO sample compare to those of other AGN at high redshift? Priddey et al. (2003) obtained 850$`\mu \mathrm{m}`$ fluxes for a sample of 57 optically bright quasars in the range $`1.5<z<3.0`$, overlapping with the redshift range of our BALQSO sample. With slightly worse noise limits, their study detected nine targets at $`>3\sigma `$ significance, a similar detection rate to this sample. The BALQSOs presented here have 850$`\mu \mathrm{m}`$ properties very similar to the Priddey et al. (2003) bright quasar sample, although only a slight trend with redshift is found in their sample, and is consistent with negative K-correction effects.

At higher redshifts, however, the dependence of FIR luminosity on redshift becomes more evident in samples of quasars (e.g. McMahon et al., 1999; Isaak et al., 2002) and even more so for radio galaxies (e.g. Archibald, et al., 2001). In the Archibald, et al. (2001) study of 47 radio galaxies between $`1<z<5`$, 20 of the sources lie in the redshift range $`1.2<z<2.35`$, overlapping with our BALQSO sample, and 4 are detected at $`3\sigma `$, a level similar to both our BALQSO sample and the Priddey et al. (2003) bright quasar survey. The FIR fluxes of the detected radio galaxies are broadly similar to those in both quasar samples. Archibald, et al. (2001), however, also found a dramatic increase of FIR flux with redshift $`[L_{\mathrm{FIR}}(1+z)^3]`$, a much stronger evolution than that seen in the quasar samples. This trend is interpreted as increasingly younger stellar populations associated with the radio galaxies at earlier epochs. Finally, 850$`\mu \mathrm{m}`$ fluxes for a sample of X-ray absorbed AGN (Page et al., 2001) reveal a very similar trend to that of our BALQSO sample, with 4 significant detections above $`z=1.7`$ and only upper limits at lower redshifts, with similar detected levels of 850$`\mu \mathrm{m}`$ fluxes.

In Figure 1 bolometric luminosities (total far-infrared luminosity plus total UV/optical luminosity) are plotted against redshift for our BAL sample. Both the optical and submm luminosities have been K-corrected, assuming $`f_\nu \nu ^\alpha `$, with $`\alpha =0.5`$ for the optical correction, and greybody emission with $`\beta =1.5`$ and $`T=40K`$ for the submm correction. Infrared luminosity is integrated over the greybody from 1100$`\mu `$m to 60$`\mu `$m rest-frame, while optical luminosity is integrated over the powerlaw continuum from 2000Å to 5000Å. The plot also shows just the far-infrared luminosities overlaid (smaller symbols), showing that the bulk of the bolometric luminosity is generated through dust reprocessing, generating ULIRG-class luminosities for both submm detections. Note, however, that the deep SCUBA limits are insufficient to exclude the lower redshift BALQSOs from being ULIRG-class sources also. The mean submm signal from the five non-detections (calculated by weighting each measurement by its inverse variance) is 0.7$`\pm `$0.5 mJy, which translates into a far-infrared lumionsity of 8.6$`\times `$1011L. Note also that the negative K-correction in the observed 850$`\mu \mathrm{m}`$ flux results in a quasi-linear flux-luminosity relation, varying less than $`20\%`$ over the redshift range redshifts 1$``$3 (e.g. Blain, Barnard, & Chapman, 2003). Thus, while the negative K-correction eases the detection of distant submm sources, its small variation over the redshift range of interest here implies that the detections/upper limits would have been similar if the redshifts were interchanged.

Figure 2 presents a comparison of $`S_{850\mu \mathrm{m}}`$–to–R-band flux ratios for our BALQSO sample (filled circles) against those for a subsample of the Archibald, et al. (2001) radio galaxies for which R-band magnitudes were obtained (open squares). Both the optical and submm fluxes have been K-corrected, and hence represent rest-frame quantities. The errors reflect only uncertainties in the submm flux. The strong evolution of submm luminosities exhibited by these radio galaxies may be a key signature of the star formation history of massive elliptical galaxies which are believed to arise from gas-rich mergers. Thus, they may offer clues to test the evolution hypothesis for the BAL phenomenon. The most striking aspect of the submm–to–optical flux ratios in Figure 2 is that the BALQSOs are substantially less submm luminous than the radio galaxies (relative to their R-band flux). This is most likely due to increased obscuration of the radio galaxy optical nucleus, since the Archibald, et al. (2001) radio galaxies are steep spectrum sources and thus, according to AGN unification schemes, are probably being viewed at large inclinations with respect to their radio axes. Note also that that BALQSOs appear to fall into a similar distribution as the radio galaxies, with higher redshift objects possessing more firm submm detections than the lower redshift objects, for which only upper limits could be obtained. Indeed, this trend prevails in this entire sample of radio galaxies (Archibald, et al., 2001). One final and important point is that there does not appear to be a correlation between the submm and optical fluxes of our BALQSOs (or of the radio galaxies) and the high values of $`S_{850\mu \mathrm{m}}`$ relative to R-band flux is indicative of dust heating by starbursts rather than by the BALQSO nucleus. This is also found to be the case in several other large quasar samples (e.g. McMahon et al. 1999; Isaak et al. 2002 and see also discussion and references in Bertoldi & Cox 2002).

### 3.1. Selection effects

It is important to examine whether there are any selection biases in our BALQSO sample. The BALQSOs were drawn from quite inhomogeneous surveys; e.g. LBQS 0059-2735 was selected from objective prism plates (Hazard et al., 1987; Weymann, et al., 1991), while Hawaii 167 was found in a spectroscopic followup of K-band sources (Cowie, et al., 1994). Hence, there is no well defined selection function. One possible concern is that the highest redshift sources are intrinsically the most luminous, and thus a correlation between total luminosity and submm luminosity could be responsible for our bimodal detections between high and low redshifts. However, as Fig. 2 demonstrates, there is no evidence for a correlation between the submm and optical luminosities of our BALQSOs. Although dust obscuration may be affecting the submm–to–optical flux ratios, this lack of correlation is also readily seen in larger samples of optically bright quasars (e.g. Isaak et al. 2002 and see also references in Bertoldi & Cox 2002).

## 4. Discussion & Conclusions

This study suggests that the FIR luminosities, inferred dust masses and star formation rates of BALQSOs are comparable to those of ULIRG-class sources. Whilst this may at first seem consistent with evolutionary models for BALQSOs, two key aspects of these submm properties, namely the trend of submm flux with redshift and the lack of correlation between the submm and optical fluxes, are also seen in samples of non-BAL quasars and radio galaxies (McMahon et al., 1999; Archibald, et al., 2001; Isaak et al., 2002). If, as is proposed in evolutionary scenarios, BAL outflows are a rapid mass-loss phase triggered by a recent gas-rich merger or close interaction event involving vigourous star formation and associated enhanced dust emission, then it is difficult to understand why the submm properties of BALQSOs are similar to those of other AGN. While these dust properties are consistent with qualitative ideas about the role of star formation in the early universe and the evolution of the overall AGN population, they do not suggest that there is anything remarkable about the dust properties of BALQSOs that would be indicative of a connection between the BAL phenomenon and the presence of dust.

If the submm properties of our (albeit small) BALQSO sample are verified by larger samples, how do we then interpret the BAL phenomenon? It is now clear that standard orientation-based models are unable to account for all properties of BALQSOs; radio spectral index and radio axes measurements indicate that many BALQSOs simply cannot be viewed at large inclination angles (Becker et al., 2000). A possible model that can explain the radio properties of BALQSOs is one in which the BAL material resides in a poorly-collimated, weak radio jet (Kuncic, 1999). Such a model can also further explain weak, low-velocity absorption outflows in Seyferts, which are clearly low-luminosity counterparts to the BAL features in quasars. Again, such outflows are difficult to explain with evolution models in which the outflow is triggered by a merger; this is simply not the case for Seyferts with BAL-like outflows.

Owing to the small number of BALQSOs in our sample, our conclusions are tentative and clearly, larger samples are needed in order to verify our results. Ongoing intensive surveys, such as the Sloan Digital Sky Survey, are proving to be very successful in discovering quasars in a relatively unbiased fashion, with the 2dF Quasar Survey now cataloguing almost 17000 quasars (Boyle et al., 2000), providing an ideal basis for such studies.

## Acknowledgements

We thank the anonymous referee for comments that improved this paper. Furthermore, we thank the staff of the JCMT for their assistance with the SCUBA observations, and the weather for being so cooperative. GFL thanks the Australian Nuclear Science & Technology Organization (ANSTO) for financial support.