text:
# On the accuracy of the S/N estimates obtained with the exposure time calculator of the Wide Field Planetary Camera 2 on board the Hubble Space Telescope

## 1 Introduction

ETCs play an important role in modern instrument use as they allow observers to determine how to carry out specific investigations and, especially, to predict the amount of time these will require. Since the time needed for the various programmes is a very sensitive issue in the allocation process for most modern high visibility ground and space-based facilities, the accuracy of these simulators must be well understood both by the observers and the time allocation committees that must rely on their results for a fair and scientifically effective distribution of the available time. In this context, unfortunately, besides the documentation accompanying the software tools, there is practically no published information on the reliability of existing ETCs of imaging cameras.

The WFPC2 has been so far the principal instrument on board the HST and it is expected to be of extreme utility to image parallel fields even now that the Advanced Camera for Surveys (ACS) is installed on the HST. ETC software utilities are available on the internet site of the STScI which simulate analytically the photometry for a given target for each HST instrument. The accuracy of these programmes plays a fundamental role in the planning of observations, in particular when extremely deep imaging is required and whenever the performances of two different instruments have to be compared.

While performing simulations for an HST proposal for the WFPC2 and the ACS, in which high accuracy was needed in order to evaluate the limiting magnitudes for deep observations of a globular cluster, we found substantial differences between the WFPC2 ETC results and real photometry obtained on archival images. We found similar differences also in archival non crowded fields, so that we decided to analyse the problem by directly comparing the ETC predictions with our photometry in various circumstances and here we show the results and the way in which they depend on field crowding. We also compare our photometry with the result of the recently published “ETC++” software (Bernstein, 2001), whose calculations are based on statistical analysis tools and take into account the real effects of the pixel size and charge diffusion.

## 2 The WFPC2 ETC: comparison with real point–source photometry

The WFPC2 ETC computes the expected SNR of a point source from its input parameters, namely: the magnitude of the star in a given spectral band, the spectral type, the filter to use, the channel of the detector (PC1 or WF2, WF3, WF4), the analogue to digital gain, the position of the star on the pixel (centre or corner), the exposure time of the whole observation (i.e. the sum of all the exposure frames) and the sky coordinates of the target (Biretta, 1996). As of late, the option to manually select a specific value of the sky brightness has been added (Biretta, 2001).

First, the programme computes the source count rate, assuming a blackbody spectrum, if the user has not specified it, and multiplies it by the response curves of the detector and filter. Then, the programme takes into account the various noise sources, including photon noise, read noise, dark noise and sky noise. The latter depends on the target position on the sky, with the sky brightening by about one magnitude from the ecliptic pole to the ecliptic plane. The programme uses the values from Table 6.4 in the WFPC2 Instrument Handbook (Biretta & Heyer, 2001) to compute the sky count rate per pixel and hence its photon noise. The contribution of the total noise to the photometry of a star depends upon the number of pixels in the point spread function (PSF) and how these pixels are weighted during data reduction. The WFPC2 ETC assumes that the data reduction employs PSF fitting photometry, so that it weights the pixels in proportion to their intensity, which maximises the SNR. The multiple read errors for a “cosmic–ray splitted” (CR-split) image, i.e. an image composed by many shorter frames, is then computed for a set of default splitting values and the corresponding SNR is also given in the ETC result page.

In order to quantify the possible WFPC2 ETC deviations from real photometry, we performed accurate aperture photometry (using the DAOPhot package) on both crowded and non crowded archival fields. The average image used in our analysis was computed after aligning the individual frames in the dithering pattern and removing cosmic ray hits. A custom programme was used, which computes the offsets of the frames by measuring the mean displacement of the centroid of some reference stars. The task then registers all the images to the first one, creates a mask of the CR-contaminated pixels, by means of an iterative sigma clipping routine with respect to the median value of the corresponding pixels in all frames, and finally computes the mean image by averaging the corresponding pixel on all the images if not included in the CR-mask. The CR-corrupted pixel in the original un-shifted frames are also replaced by the value of the corresponding pixel in the mean image in order to allow us to perform photometry on both the combined and on the individual frames. Our instrumental magnitudes were then transformed to the Johnson/Cousin UBVRI system by following the prescription of (Holtzmann et al., 1995) — specifically, their equation 8, which also takes account of the colour correction by means of the coefficients in Table 7 therein — and by making an optimal choice of the aperture radius for each star, so as to minimise the associated photometric error.

For the crowded field, we have used the images of the Galactic globular cluster $`M\mathrm{\hspace{0.17em}4}`$ (HST proposal 5461), obtained in the F555W and F814W filters. These are deep images centered in a region at one core radius from the centre of a dense globular cluster and should be representative of the cases in which the field observed is filled with a multitude of very bright and saturated stars ($`V16`$), whose haloes overlap each other and cover a significant fraction of the frame (Figure 1).

Images of the field of ARP 2 taken from HST proposal 6701, also obtained through F555W and F814W filters, were used as representative of a sparsely filled region in which the field is populated with faint stars with no appreciable overlapping haloes (Figure 2).

### 2.1 How we measured the SNR

Aperture photometry was performed on both series of images (i.e. crowded and non crowded), with the following parameters. The flux of the object was sampled within an aperture of radius $`r_0`$, which is varied in steps of $`0.5`$ pixel. The background is sampled within an annulus drawn from an inner radius $`r_1=r_0+1`$ pixel to an outer radius $`r_2`$, with an annulus width which is varied from $`3`$ pixel up to $`20`$ pixel. As is discussed later in this section, an adjustable aperture radius and annulus size allow us to maximise the SNR, by limiting the noise generated by the contamination of the neighbouring objects. Moreover, the background is always estimated by taking the mode, rather than the mean or median, of the pixel distribution within the annulus. Appropriate aperture corrections were applied, which were directly measured from the most isolated non saturated stars in the field. A direct comparison with the encircled energy curve for the WFPC2 PSF (Biretta & Heyer, 2001) shows a perfect match, thus proving that the growth curves that we measured are reliable.

The DAOPhot task, used with the optimal aperture radius $`r_0`$ and the radii $`r_1`$ and $`r_2`$ for the sky annulus, gives the best estimate of both the magnitude and the associated error $`\sigma _m`$, from which we compute the SNR by using the equation:

$$\mathrm{SNR}^\mathrm{D}=\frac{1}{\mathrm{e}^{\sigma _m/1.08574}1}$$
(1)
that comes from inverting Pogson’s relation $`\mathrm{\Delta }m=2.5\mathrm{log}((F+\mathrm{\Delta }F)/F)`$, where the numerical constant $`1.08574`$ is equal to $`2.5\mathrm{ln}(10)`$. Hereafter, the acronym $`\mathrm{SNR}^\mathrm{D}`$ indicates the SNR estimated on the basis of the photometric error given by DAOPhot.

As an independent check, we have computed the SNR as indicated in equation 6.7 of the WFPC2 Handbook (Biretta & Heyer, 2001) which, in the practical case of observed quantities, becomes:

$$\mathrm{SNR}^\mathrm{H}=\frac{FG^{1/2}N^{1/2}}{\sqrt{F+(S+N_R^2/G)\pi r_0^2+N_S}}$$
(2)
where $`r_0`$ is the optimal aperture radius used by DAOPhot, $`N`$ is the number of frames combined together, $`N_R`$ the read-out noise (in units of electrons) of each specific CCD, $`S`$ the average background per pixel inside the annulus from $`r_1`$ to $`r_2`$, in units of DN, $`F`$ is the flux within the aperture of radius $`r_0`$ after subtraction of the background contribution $`S\pi r_0^2`$, in units of DN, and $`G`$ is the effective gain factor, i.e. the CCD gain times the number of frames averaged togheter. Finally, $`N_S`$ is a small (although non negligible) contribution to the error affecting the estimate of the background level which takes on the form:

$$N_S=\frac{(S+N_R^2/G)\pi r_0^4}{r_2^2r_1^2}$$
(3)
The computation of $`\mathrm{SNR}^\mathrm{H}`$ makes no use of the error estimate on the magnitude or flux provided by DAOPhot, so it is reassuring to find that $`\mathrm{SNR}^\mathrm{H}`$ is in excellent agreement with $`\mathrm{SNR}^\mathrm{D}`$. This, however, only happens if we use an adaptive choice for aperture radius and for the background annulus, as explained above. In fact, if we select a fixed radius and annulus size in a crowded environment, the contamination due to neighbouring stars alters the statistics of the sky within the annulus and we always find $`\mathrm{SNR}^\mathrm{H}>\mathrm{SNR}^\mathrm{D}`$. This is precisely the reason that made De Marchi et al. (1993) conclude that core aperture photometry, i.e. source and sky measurement conducted as close to the source as possible, as well as the use of the mode for the background are most advisable in crowded environments.

### 2.2 How the ETC expects the SNR to be measured

In light of the consistency between $`\mathrm{SNR}^\mathrm{D}`$ and $`\mathrm{SNR}^\mathrm{H}`$ and since the latter stems directly from equation 6.7 of the WFPC2 Handbook, on which the WFPC2 ETC is also based, we can now proceed and compare our measured $`\mathrm{SNR}^\mathrm{D}`$ with the ETC predictions. Before doing so, however, we must make sure that the way in which we measure the SNR (i.e. $`\mathrm{SNR}^\mathrm{D}`$) is consistent with the way in which the ETC software expects users to carry out the photometry. In fact, the latter assumes that the data reduction process employ PSF fitting photometry, i.e. that optimal weighting be assigned to each pixel in proportion to its intensity in the PSF. As discussed above, however, we have used aperture photometry to determine $`\mathrm{SNR}^\mathrm{D}`$. The WFPC2 ETC instructions would indeed offer a correction to apply to the ideal PSF fitting case $`\mathrm{SNR}^\mathrm{P}`$ (we call it “ETC optimal SNR”) in order to convert it to the equivalent SNR that would be obtained with canonical aperture photometry $`\mathrm{SNR}^\mathrm{A}`$ (“ETC aperture SNR”). Following the WFPC2 ETC instructions in (Biretta, 1996), we have:

$$\mathrm{SNR}^\mathrm{A}=\mathrm{SNR}^\mathrm{P}\frac{K}{r_0}$$
(4)
where $`K=0.11`$ for the PC camera and $`K=0.17`$ for the WF chips, particularly valid when the aperture radius is $`r_0>2.5`$ pixel for the PC and $`r_0>1.8`$ pixel in the case of the WF.

Since we determined $`\mathrm{SNR}^\mathrm{D}`$ by using aperture photometry, it would seem that we need to take into account the correction given by Equation 4. We show, however, that this correction is not necessary thanks to the adaptive method that we used for photometry. In Figure 3 we plot, for the PC chip, the measured $`\mathrm{SNR}^\mathrm{D}`$ against the prediction of the ETC for the aperture photometry case, i.e. $`\mathrm{SNR}^\mathrm{A}`$. We should like to clarify here how Figure 3 and others of the same type in the following were built. After having measured the calibrated magnitude of a star in the images, we folded the latter value through the WFPC2 ETC so as to calculate the estimated SNR for an object of that brightness and for the exposure time and CR-SPLIT pattern corresponding to those of the actual combined image. For this and all the other figures in this paper, unless otherwise specified, we used the “average sky” option for the sky brightness setting as allowed by the new WFPC2 ETC Version 3.0.

We can see from Figure 3 that the prediction of the ETC for aperture photometry ($`\mathrm{SNR}^\mathrm{A}`$) are over-estimated for faint stars and under-estimated for bright objects with respect to the measured values for both the sparse and the crowded field. Figure 4 is the analogue of Figure 3 but here the reference is the ETC optimal SNR, $`\mathrm{SNR}^\mathrm{P}`$, i.e. without any correction for aperture photometry. As one can easily see, the ETC in this case always overestimates the value of the SNR with respect to the measured ones by up to $`100\%`$ for the fainter stars. As the right hand side axis shows, such a mismatch of the SNR corresponds to a time estimation error of the same amount (see Equation 7 ahead), i.e. the ETC appears to underestimate the exposure time actually needed to achieve a given SNR.

A closer look at Figure 4, however, reveals that the scatter of the representative points on the plot is smaller when our measurements are compared with SNRP than for SNRA and that the overall behaviour is closer to the ETC prediction at any magnitude. This is a consequence of our optimised aperture and annulus photometry closely approaching PSF fitting. In light of these results, in the following we ignore the correction for aperture photometry given by equation 4 and compare our measurements directly with the ETC optimal SNR, i.e. SNRP.

### 2.3 How the predicted SNR compares with the observed one

Figures 3 and 4 clearly witness the dependence of the actual SNR upon the level of field crowding and, at the same time, its independence of the filter used. In principle, one could question the validity of our latest assumption, i.e. that of ignoring the correction to be applied to the SNR measured with aperture photometry. In fact, in a crowded field, PSF fitting photometry is expected to give better results. We have, therefore, attempted a direct comparison between the predictions of the ETC and the results of PSF fitting photometry. Rather than carrying out the reduction ourselves, we have utilised one of the finest examples of photometric work carried out on these very M 4 data by Richer et al. (1997), who employed very accurate ALLFRAME photometry as described in detail in Ibata et al. (1999). In their paper, these authors measure the magnitude of each star from the individual frames in the dithering stack and compute the combined magnitudes as the weighted average of the corresponding fluxes, the error on them, $`\sigma _m^P`$, being related to the flux scatter amongst the frames.

In order to make a reliable comparison with our results, we have performed, in a similar way, optimised aperture photometry on the individual frames (i.e. the original, not yet aligned images, in which CR-hits had been removed as described above). The measured fluxes were averaged with a weight inversely proportional to the DAOPhot estimated uncertainty after rescaling for the flux ratio. Our final magnitude errors, $`\sigma _m^A`$, are thus derived from the standard deviation of the fluxes, divided by the square root of the number of images combined. Figure 5 displays the comparison between $`\sigma _m^P`$, $`\sigma _m^A`$ and the ETC prediction, showing that the two photometric uncertainties overlap each other, while the ETC largely overestimates the precision that can be attained with PSF fitting photometry, even by one of the most experienced teams.

Thus, in this crowded case, it is also apparent that the ETC deviations are independent of the photometric technique adopted. In sparse fields, where aperture photometry and PSF fitting are equally effective and reliable, Figures 3 and 4 already prove that the ETC predictions depart from the measured data, although by a smaller amount than that applicable in the crowded case. Finally, in Figure 6 we compare the predictions of the ETC with the actual measurements for both the PC and WF, to show that the behaviour of the ETC applies regardless of the channel.

The relevance of the above considerations becomes clear when one uses an ETC to simulate very deep observations, especially when a comparison between instruments, e.g. ACS/WFC and WFPC2, is required to compare the limiting magnitude in given exposure times. As experience shows, a star finding programme is able to detect a faint point source only when its brightest pixel is at least $`2`$ or $`3\sigma _{\mathrm{sky}}`$ above the sky background (where $`\sigma _{\mathrm{sky}}`$ is the standard deviation of the background), with a value of $`5`$ or more being the typical prerequisite in most faint photometry precision applications. If we plot the so called object detectability $`d`$, defined as:

$$d=\frac{\mathrm{Peak}\mathrm{Sky}}{\sigma _{\mathrm{sky}}}$$
(5)
as a function of the magnitude error $`\sigma _m`$, we obtain the graph in Figure 7. Here we notice that the detectability (which is practically independent of the filter and crowding) drops to the value of $`d=2.1`$ just when the magnitude error approaches $`0.5mag`$, which is usually considered the maximum allowed error in canonical photometric work. By relating the detectability $`d`$ with the ETC optimal SNR, $`\mathrm{SNR}^\mathrm{P}`$, as done in Figure 8, we see that $`d=2.1`$ corresponds to an ETC optimal SNR of $`3.0`$ for the non crowded case and to $`7.0`$ for the crowded case. This literally means that if we need to know the magnitude of the faintest detectable star in an observation of a stellar field with the WFPC2 we should query the ETC, setting “average sky”, for a SNR of $`7.0`$ and $`3.0`$, respectively in a crowded and in a sparse environment. It is normally assumed that a $`3\sigma `$ detection requires a SNR of 3, but in the case of the SNR provided by the WFPC2 ETC, this is only true for an isolated object.

## 3 Discussion and corrections

The direct consequence of what we have illustrated so far is that, if the ETC were used to plan observations of faint stars in a globular cluster like M 4 with the WFPC2, the predicted exposure time could be considerably underestimated. Conversely, the same predictions would be almost correct for a star of equal brightness in a sparse field. In the following we try and provide an empirical correction formula that can be applied to the SNR given by the WFPC2 ETC to compensate for the effects of crowding.

In order to understand the discrepancy between the expected and measured SNR and to clarify how to exactly account for the effects of crowding in the simulations, we artificially modified the background level and photon noise in the sparse field so as to reproduce the sky level and sky variance measured on the crowded field. In practice, we added to the sparse field a Gaussian noise with a mean equal to the difference in the sky level between the two fields and a variance equal to the quadratic difference of the sky variances between them. The SNR diagramme for the modified image (Figure 9) reveals that the locus of the modified sparse data points shifts towards and perfectly overlaps the crowded field locus. This tells us, as expected, that the increased background level resulting from crowding is responsible for the differences shown in Figures 3 and 4 between sparse and crowded fields.

It is, however, true that the ETC gives the SNR under the best possible sky conditions, which are rarely encountered, if ever, in real observations. Moreover, it is generally not expected of the ETC to take account of the position and brightness of all the stars in the field as would be necessary to simulate how crowding increases the background level. We have, therefore, manually set the ETC sky brightness to match the levels directly measured with the DAOPhot SKY task on the crowded image (i.e.. the mode of the levels distribution), hoping in this way to force the SNR simulated by the ETC to agree with our measurements. In fact, the results change only marginally, as shown in Figure 10, where $`\mathrm{SNR}^\mathrm{D}`$ and $`\mathrm{SNR}^\mathrm{P}`$ are plotted against the observed magnitude (Johnson $`V`$ in this case). The ETC simulation gets closer to the real data, but it does not still match them. Moreover, it seems as if a suitable value for the background cannot be found at all as shown in Figure 11, where one sees that the sky value that would force the ETC prediction to match $`\mathrm{SNR}^\mathrm{D}`$, changes significantly as a function of star brightness.

We must, thus, conclude that the treatment of the background is a major issue for the WFPC2 ETC, although that alone cannot explain the whole discrepancy. It goes without saying that we have verified and confirmed that the predictions of the ETC as concerns the count rates per pixel in the source and background are precise to within an accuracy of 10 %, as one would expect of a professional tool. We have also repeated all our tests on the individual frames, compared in turn with the predictions of the ETC for a case of CR-SPLIT=1. The result being the same, we can exclude an error in either the way in which we combined the data or in the way in which the ETC accounts for CR-SPLIT$`>1`$. The rest of the discrepancy, then, must be attributed to the way in which the noise is estimated, the signal being correct. A delicate issue could be, for instance, the value and operational definition of $`\sigma _{\mathrm{sky}}`$. We notice here that large variations in the value of $`\sigma _{\mathrm{sky}}`$ are possible, in the crowded environment, depending as to whether we measure it with the IRAF SKY task, which fits a Gaussian around the mode, or as the standard deviation that one obtains by manual analysis over the darkest regions of the background in the image. In fact the latter can be up to 3 times smaller than the former, and also 2 times smaller than the mean sigma as measured inside the photometric sky annulus around each star. Conversely, all these numbers turn out to be quite similar for the sparse field image.

To try and account for the possible sources of the residual error, we considered recent results published by Bernstein (2001), who uses Fourier analysis and Fisher information matrices to show to which extent the SNR of a point source depends on factors which normally are not considered in ETC programmes, such as pixel size, intra-pixel response function, extra-pixel charge diffusion and cosmic ray hits.

According to this work, a programme that does not take all these parameters into account may overestimate the SNR by up to a factor of 2. More precisely, whenever background limited point source photometry is involved, the key factor for the SNR calculation, namely the “effective area” $`A_{SN}`$ (see equation 12 in Bernstein 2001), strongly depends on the detector geometry, such as pixel size, under-sampling factor, intra-pixel response function and charge diffusion. The finite pixel size plays an important role, as even a Nyquist sampled pixel (i.e. one $`\lambda /2D`$ in size) causes a 13 % degradation in the SNR of a faint star and the same applies to extra-pixel charge diffusion.

In order to check whether these problems also affect the WFPC2, we configured Bernstein’s “ETC++” software to simulate WFPC2 point source photometry for the sparse field. The result is shown in Figure 12 where the measured SNR ($`\mathrm{SNR}^\mathrm{D}`$), the ETC optimal SNR ($`\mathrm{SNR}^\mathrm{P}`$) and the ETC++ SNR for aperture photometry are plotted against the stellar magnitude. The ETC++ gives a confidence level for its results as the value of the cumulative function of the stars distribution above the computed SNR. The ETC++ line in Figure 12 means that 50 % of the stars of any given magnitude should be above this line. The WFPC2 ETC does not give confidence levels, but we can assume that its SNR is computed as the mean of the SNR distribution at any given magnitude, i.e. at 50 % confidence level. If this is the case, Figure 12 indicates that the actual SNR is located in between the WFPC2 ETC and the ETC++ predictions, thus confirming the difficulty of any analytical ETC in reliably estimating the SNR.

Thus, a correction for the currently on-line WFPC2 ETC can only be empirical in nature. The following formula can be used to obtain a realistic estimate of the SNR:

$$\mathrm{SNR}^\mathrm{C}(60C+17)(e^{0.012\mathrm{SNR}^\mathrm{P}}1)+0.93\mathrm{SNR}^\mathrm{P}$$
(6)
where $`\mathrm{SNR}^\mathrm{P}`$ is the SNR estimated by the ETC without correction for aperture photometry and $`C`$ is a measure of the crowding, defined as the logarithm of the ratio between the total area of the chip and the number of pixel with value lower than the modal sky value plus one standard deviation. For example $`C`$ is equal to $`0.05`$ for our sparse field, whereas it grows to $`0.42`$ in the crowded case of M 4. For faint stars, e.g. for $`\mathrm{SNR}^\mathrm{P}20`$, this equation can be roughly approximated by the rule of thumb that the actual SNR is about $`1/2`$, or $`2/3`$, of the $`\mathrm{SNR}^\mathrm{P}`$, respectively for a crowded and non crowded environment. It should be noted that not even in an ideal case of zero crowding ($`C0`$) would the measured SNR match the prediction of the ETC, since there would still be a discrepancy of the same order of that found in the sparse case.

The advantage of this formula is that $`\mathrm{SNR}^\mathrm{C}=3`$ would now always imply a $`3\sigma `$ detection, regardless of the level of crowding in the image. The correction that we propose would allow an observer to accurately plan his observations and make the best use of the HST time. For the low SNR regime (e.g. $`\mathrm{SNR}^\mathrm{P}50`$), equation 6 can actually be rewritten to more explicitly show the effects of crowding on the exposure time:

$$\mathrm{t}^\mathrm{C}\mathrm{t}^\mathrm{P}\frac{\mathrm{SNR}^\mathrm{P}}{(60C+17)(e^{0.012\mathrm{SNR}^\mathrm{P}}1)+0.93\mathrm{SNR}^\mathrm{P}}$$
(7)
where $`\mathrm{t}^\mathrm{P}`$ is the exposure time predicted by the ETC to reach a certain SNR and $`\mathrm{t}^\mathrm{C}`$ is its actual value.

An example of how serious the underestimate of the exposure time can be when the ETC is not used with the above caveat in mind is given in Figures 13a and 13b for a crowded environment. There we show a simulation of the detectability of the white dwarf cooling sequence with the WFPC2 in NGC6397, the nearest globular cluster, through the filters F606W and F814W. We have adopted the theoretical WD cooling sequence of Prada Moroni et al. (2002) which provide a perfectly thin isochrone and have applied to it the colour and magnitude uncertainty that one obtains from the estimated SNR by inverting equation 1. Two cases are shown: one (a) as predicted by the WFPC2 ETC and one (b) for our corrected estimate of equation 6. The difference is outstanding, as the ETC predictions, taken at face value and ignoring the effects of crowding, would suggest that the sequence is not spread very much by photometric errors and its quasi-horizontal tail between $`m_{606}=29`$ and $`m_{606}=30`$ is clearly noticeable, whereas in our realistic simulation the sequence is widely spread and its lower part lies well below the detection limit.

The delicacy of the issue is immediately apparent when one considers that, based on the ETC estimates, one would deem that 120 orbits are sufficient to reliably secure the white dwarf cooling sequence in the colour–magnitude diagramme of NGC6397 down to $`m_{606}=30.5`$ and $`m_{814}=30`$, whereas, in fact, the correction shows that as many as 255 orbits would be needed to comfortably reach those limits with the WFPC2.

All of the above considerations are valid not only for the WFPC2, but also for any analytical ETC in general, especially when used to estimate the SNR of stars embedded in crowded environment or when the detector considerably under-samples the PSF, as suggested in Bernstein (2001). We should underline here, however, that this does not mean at all that the ETCs are unreliable nor that they are useless. One of the most important and practical reasons for having a standardised ETC is to allow the telescope time allocation committees to compare all the proposals on equal footing. In this respect, the ETC does not necessarily need to be accurate. Clearly, the better the detector’s cosmetics, intra-pixel response, charge diffusion and readout noise, the closer will the real photometry be to the ETC prediction. Thus, we expect, for example, a better behaviour of the ACS/WFC on-line simulator with respect to the WFPC2.

A non-analytical SNR calculator, which would simulate the whole observing session, including the dithering pattern, by numerically reproducing the real field (i.e. with the correct stellar positions and brightness, as imaged by a realistic model of the detector) and which uses the same photometric tools that will be adopted by the user (such as DAOPhot, ALLSTAR, and the like), would be, in our opinion, the best method to accurately predict the expected performances of any planned observing programme providing reliable results. Alternatively, at least for imaging ETCs which have very few configuration parameters and are fairly stable such as those in space telescopes, one should consider empirical modeling. One can take real results, such as we did in this paper, to calibrate an ETC which in turn interpolates between calibrations. In this way, the use of an empirical correction formula such as the one proposed here would guarantee a closer matching between simulations and real observations.

## 4 Summary and conclusions

The results of the WFPC2 exposure time calculator for point sources have been analysed by direct comparison with aperture and PSF photometry on real archival images. Significant deviations have been found between the ETC predictions and the actual photometry on the real data. Specifically, the analysis shows that the ETC deviations are i) independent of the filter, ii) independent of the choice of optimised aperture photometry or PFS fitting photometry, iii) independent of the PC or WF channel used, iv) strongly dependent upon the level of crowding in the field and that v) the ETC systematically overestimates the SNR, slightly for the bright sources and more seriously for faint sources close to the detection limit. Moreover, when data reduction follows the optimised aperture photometry method, the measured SNR will be as good as that obtained with PSF fitting and there is no need to apply the aperture photometry conversion suggested in the ETC documentation. An empirical correction formula is given to compute realistic SNR estimates, so as to assist observation planning when extremely faint sources have to be imaged, an example of which is presented. Manually increasing the value of the sky brightness in the simulator, so as to mimic the effects of crowding, shows that, although important, the background level is not the key parameter to explain the discrepancy, which is present even for data collected in rather sparse environments. Thus, it is not possible to correct the WFPC2 ETC predictions by just modifying the sky level. A comparison with a software tool developped by Bernstein (2001), whose predictions slightly underestimate the SNR at variance with the WFPC2 ETC, suggests that the effects of pixel size, charge diffusion and cosmic rays hits could be more important than previously thought.

It’s our pleasure to thank H. Ferguson, M. Stiavelli, S. Casertano, F. Massi, L. Pulone and R. Buonanno for helpful discussions. We are indebted to F. Valdes, the referee of this paper, for his useful comments and suggestions. G. Li Causi is particularly grateful to the ESO Director General’s Discretionary Fund for supporting his work. We also wish to thank Gary Bernstein for making his ETC++ software available to us.

FIGURE CAPTIONS:

Figure 1: Negative image of a crowded field in the globular cluster M4 obtained with the PC channel of the WFPC2 through the F555W filter (Richer et al., 1997).

Figure 2: Negative image of a sparse field obtained with the PC (Ibata et al., 1998).

Figure 3: The ratio between the SNR measured in crowded and sparse fields (SNRD in the text) and the WFPC2 ETC prediction for aperture photometry (SNRA in the text) is shown for the F555W and F814W filters.

Figure 4: The ratio between the SNR measured in crowded and sparse fields (SNRD in the text) and the WFPC2 ETC prediction for PSF-fitting photometry (SNRP in the text) is shown for the F555W and F814W filters. The right hand side axis applies to the low SNR regime ($`50`$) and indicates the amount of the time estimation error, i.e. the ratio between the actual exposure time (Equation 7 in the text) and that estimated by the ETC, for a given SNR in the abscissa.

Figure 5: Measured magnitude error from PSF-fitting photometry of Ibata et al. (1999) ($`\sigma _m^P`$) and from our optimized aperture photometry ($`\sigma _m^A`$), as a function of the magnitude, compared with the WFPC2 ETC prediction, for the crowded field and F555W filter.

Figure 6: The ratio between the measured SNR (SNRD) and the ETC optimal SNR (SNRP) is shown for the PC and the WF2 channels of the WFPC2, in the crowded field case.

Figure 7: Detectability $`d`$ versus measured magnitude error ($`\sigma _m`$). An uncertainty $`\sigma _m=0.5mag`$, usually the highest allowed in most photometric works, corresponds to a detectability $`d=2.1`$.

Figure 8: Detectability $`d`$ versus ETC optimal SNR (SNRP). A value of $`d=2.1`$ corresponds to a detection limit of $`SNR^P=3.0`$ or $`SNR^P=7.0`$ respectively for the sparse and the crowded case.

Figure 9: The ratio between the measured SNR (SNRD) and the ETC optimal SNR (SNRP) is shown for the crowded and sparse fields, before and after the artificial brightening of the sparse field background.

Figure 10: Comparison between the measured SNR (SNRD) and the ETC predictions for i) default low background, ii) default average background, iii) default high background and iv) actually measured background.

Figure 11: Values to enter in the User Specified Sky Background parameter of the WFPC2 ETC in order to force the ETC to match the measured SNR for crowded and non crowded fields.

Figure 12: Comparison between the prediction of the WFPC2 ETC v.3.0, the prediction of the ETC++ software and the measured SNR (SNRD), for the crowded and non crowded cases in the two filters (both ETCs were used here after setting the sky magnitude to the value measured in the real images).

Figure 13: Comparison between (a) the WFPC2 ETC predictions (SNRP) and (b) our correction of equation 6 (SNRC), in a simulation of a 120 HST orbits observation of the white dwarfs cooling sequence in NGC6397, in a colour–magnitude diagramme made through the filters F606W and F814W.


format: text

text:
# Comment: on possible time reversal symmetry breaking in normal state of ruthenate.

## Abstract

A spiral spin structure with nesting vector is expected to be the ground normal state of $`Sr_2RuO_4`$. It might be observed with use of the ARPES method with circularly polarised photons designed in .

PACS: 71.10.Fd, 71.27.+a, 71.10.Hf

The latest exciting results of A.Kaminski et al.  on the polarisational photoemission study of cuprates allows to propose another object with time reversal symmetry breaking for similar study. It seems to be interesting to apply the same ARPES method to study the polarisation phenomena in a single-layered ruthenate $`Sr_2RuO_4`$. It attracts an attention as a superconductor ($`T_c1.5K`$) with possible triplet type of pairing . The main argument in favour of such type of pairing is the Knight shift behaviour . Several models for a symmetry of the triplet superconducting (SC) order have been discussed . Experimental studies ( NQR relaxation rate , heat capacity, anisotropy of thermal conductivity and so on, see  and references therein) point to an existence of the node lines in the superconducting gap. It does not exclude the conventional d-wave pairing also.

Unconventional triplet SC order breaks a time reversal symmetry. But the symmetry breaking may occur in normal state also. Recently we suppose  that the normal state of $`Sr_2RuO_4`$ has a spiral spin structure. The mean field (MF) study shows that the energy of spiral state is lower than that of the para-, ferro- and antiferro-magnetic states if the spirality vector $`Q`$ coincides with a nesting vector $`Q=2\pi (\frac{1}{3},\frac{1}{3})`$ of $`\alpha ,\beta `$ bands of the ruthenate . This hypothesis is consistent with an incommensurate peak observed in inelastic neutron scattering at $`q(0.3,0.3,1)`$ in units $`(\frac{2\pi }{a},\frac{2\pi }{a},\frac{2\pi }{c})`$ . But a spiral structure means that the normal state is a state with a broken time-reversal symmetry. This is a state with the non-zero spin currents $`𝐉_{}=𝐉_{}𝐐`$ of opposite directions for two different spin polarisations perpendicular to the spin-rotation plane of the spiral state. It means that electrons with up (down) polarisation occupy preferentially the k-states with $`\mathrm{𝐤𝐐}<0`$ or $`\mathrm{𝐤𝐐}>0`$ correspondingly. Such asymmetry leads to the polarisation asymmetry of Fermi surfaces found in MF solutions in . It has been shown that in spiral state the contributions $`A_\sigma (k,\omega )`$ to the total one-particle spectral function $`A(k,\omega )`$ from different spin polarisations $`\sigma =`$ or $``$ display at $`\omega 0`$ the different sections of the Femri surfaces (see Fig.2 in ). So a strong dependence of the photoemission signal on a direction of a circularly polarised light should be observed. The 100% effect for the long ranged spiral order is expected to be greatly suppressed by any domain structure or any sort of a disorder. In particularly, it must be destroyed by the thermal cycling procedure, which removes a set of shadow Fermi surfaces from the ARPES map of $`Sr_2RuO_4`$ .

Work is supported by Russian Fund of Fundamental Research (Projects No. 00-03-32981 and No. 00-15-97334). Authors thanks V.Ya Krivnov for useful discussions and P.Fulde for possibility to work in Max Planck Institute for Physics of Complex Systems, Dresden.


text:
# A weak instability in an expanding universe?

## 1 Introduction

A general gravitational action can be arranged as a derivative expansion, where higher derivative terms contain more factors of the curvature tensor and its covariant derivatives. As such this could be viewed as the expansion of some underlying theory. But major problems arise if some truncation of the derivative expansion is treated as a fundamental theory. For example the inverse propagator of the high derivative theory is at best of the form

$$k^2\underset{i}{}(k^2m_i^2)\mathrm{with}\mathrm{all}m_i^2>0.$$
(1)
Close to every other zero the $`k^2`$ will have the wrong sign. These ghost poles lead to loss of unitarity, the possibility of unstable solutions, etc. At a deeper level the $`m_i`$ probably have little to do with the spectrum of the underlying theory. Also a truncated derivative expansion will have solutions bearing no resemblance to the solutions of the true theory. Higher derivative gravity only makes sense as an effective field theory where the higher derivative terms are properly treated as perturbations . (For a contrary view see .)

Nevertheless we will show in this paper that higher derivative gravity in the case $`m_i^2>0`$ can serve as an instructive and suggestive framework for the investigation of the following question. Can there be a dynamical nonlinear coupling between physics at very large and very small time scales (UV-IR coupling) in a theory of gravity?

We will consider a class of solutions to higher derivative gravity theories of the FRW form, where the standard FRW time dependence for the scale factor can be recovered as a special case. This special case corresponds to particular initial conditions, and departures from this will yield oscillations with Planck scale frequencies in the scale factor of the FRW metric. It has been noted before  that rapid metric oscillations can emerge from higher derivative theories. We will find that the oscillations, even when averaged over, produce a departure from the FRW solution, where the fractional amount of this departure grows with time. We refer to this as a weak dynamical instability since the growth of these oscillations turns out to be logarithmic.

We claim that this is a generic property of higher derivative gravity, true at any order in the derivative expansion. We have reached this conclusion through a study of the field equations at high order in the derivative expansion. Since Planck frequencies are involved, derivatives of the oscillating component are not suppressed. The amplitude of the oscillation on the other hand is extremely small, at least at late times, if the impact of the oscillations on the equations is smaller than or of order the usual effects of cosmological expansion. We therefore consider a small amplitude expansion, where terms quadratic in the scale factor and its derivatives produce the leading nonlinear effects. It is these nonlinear effects that couple the oscillations and cosmological expansion together and give rise to the dynamical instability.

As we have said the new solutions can be arbitrarily close to the standard FRW solution, in which case there will be an arbitrarily long time where the effects of all higher derivative terms are much smaller than the Einstein curvature term. But eventually the effects of higher powers of curvature can become as important as the Einstein term, and the effect on the cosmological expansion can be significant.

One would naively think that the large and growing number of parameters in higher derivative actions would produce results that are very ambiguous. On the contrary we find that in the small amplitude expansion that the derivative expansion is extremely restrictive. Each order in the derivative expansion, at six derivatives and beyond, results in only two parameters, one of which multiplies both the leading and next-to-leading terms in the small amplitude expansion. And then quite surprisingly we find that the second parameter drops out of the quantities of interest. Thus a very particular pattern emerges in the derivative expansion that allows us to relate the amplitudes and energies of the oscillating modes to the matter sources.

We start with the gravitational action

$$S=M_{\mathrm{Pl}}^2d^4x\sqrt{g}(2\mathrm{\Lambda }+R+aR^2+bR_{\mu \nu }R^{\mu \nu }+cR_{\mu \nu \lambda \kappa }R^{\mu \nu \lambda \kappa }+\mathrm{})$$
(2)
where $`a`$, $`b`$, and $`c`$ are of order $`M_{\mathrm{Pl}}^2`$. In this work we will only consider a spatially flat metric,

$$ds^2=dt^2+a(t)^2\delta _{ij}dx^idx^j.$$
(3)
The field equations are

$$\widehat{G}_{\mu \nu }=T_{\mu \nu }+\mathrm{\Lambda }g_{\mu \nu }$$
(4)
where $`8\pi G`$ has been absorbed into $`T_{\mu \nu }`$ and $`\mathrm{\Lambda }`$. $`\widehat{G}_{\mu \nu }`$ is a power series in derivatives, where the two derivative terms make up the Einstein tensor. The matter described by $`T_{\mu \nu }`$ has the equation of state $`p(t)=w\rho (t)`$. We introduce a dimensionless Planckian time

$$\stackrel{~}{t}=\frac{t}{\sqrt{6a+2b+2c}}M_{\mathrm{Pl}}t,A(\stackrel{~}{t})2\mathrm{ln}(a(t)),$$
(5)
and in the following we drop the tilde. For the above metric, $`6a+2b+2c`$ is the only combination of the four-derivative parameters that appear in the field equations, and we must assume that it is positive. We absorb a further factor of $`6a+2b+2c`$ into the terms in (4), thus making $`T_{\mu \nu }`$ and $`\mathrm{\Lambda }`$ dimensionless.

In terms of this dimensionless time we will be considering oscillatory behavior for $`A(t)`$ with roughly unit frequencies and very small amplitude $`\epsilon `$. Thus $`^nA(t)/t^n\epsilon `$ and we consider a small amplitude expansion in powers of $`\epsilon `$. The linear terms are contained in $`\widehat{G}_{ii}`$.

$$\widehat{G}_{00}[A(t)]=𝒪(\epsilon ^2)$$
(6)
$$\widehat{G}_{ii}[A(t)]=\underset{n1}{}c_{2n}_t^{2n}A(t)+𝒪(\epsilon ^2)$$
(7)
The definitions above yield $`c_2=c_4=1`$. At 6 derivatives and above these linear terms originate from terms in the action containing two curvature tensors, or contractions thereof, and the appropriate number of covariant derivatives.

To describe our results at higher orders in $`\epsilon `$ we make the following definitions,

$$\widehat{G}_{00}[A(t)]\underset{n}{}c_{2n}𝒢_\rho ^{(2n)}$$
(8)
$$\widehat{G}_{ii}[A(t)]\underset{n}{}c_{2n}𝒢_p^{(2n)}$$
(9)
where $`𝒢_\rho ^{(2n)}`$ and $`𝒢_p^{(2n)}`$ contain $`2n`$ derivatives.

We display in the Appendix the results up to 12 derivatives, for $`𝒢_p^{(2n)}`$ up to $`𝒪(\epsilon ^2)`$ and $`𝒢_\rho ^{(2n)}`$ up to $`𝒪(\epsilon ^3)`$. As the number of derivatives grows the number of possible terms in the action that contribute grows very fast. Despite this we see that in addition to $`c_{2n}`$ only one new parameter $`\chi _{2n}`$ appears at each order in the total number of derivatives, to the order in $`\epsilon `$ displayed. (Actually at 2 and 4 derivatives there are no parameters given the definitions above and the terms at all orders in $`\epsilon `$ are displayed.)

These results were obtained without obtaining the general field equations for general metric. Since the $`t`$ dependence of the metric (3) appears as a general function of $`t`$ in $`g_{ii}`$, we can insert the metric (3) into the action and vary with respect to $`A(t)`$. This will give the $`ii`$ field equation, and thus the $`𝒢_p^{(2n)}`$, for this particular metric. The $`𝒢_\rho ^{(2n)}`$ can then be obtained through the use of the conservation equation (generalized Bianchi identity) $`D^\mu \widehat{G}_{\mu \nu }=0`$, which we satisfy to $`𝒪(\epsilon ^3)`$. In fact this completely determines $`𝒢_\rho ^{(2n)}`$ at order $`𝒪(\epsilon ^{m+1})`$ given $`𝒢_p^{(2n)}`$ at order $`𝒪(\epsilon ^m)`$. All the terms in the field equations, to any order in $`\epsilon `$ and for any given action, can be obtained in this way, but only the terms at the orders in $`\epsilon `$ relevant for our study are given in the Appendix.

## 2 Pure oscillations

We start with the oscillating vacuum solutions of the linearized equations, corresponding to the Lagrangian

$$=\frac{1}{2}A(t)(\underset{n1}{}c_{2n}_t^{2n})A(t).$$
(10)
For $`A(t)`$ of the form $`\mathrm{cos}(\omega t)`$ the allowed frequencies must satisfy

$$\underset{n1}{}(1)^nc_{2n}\omega ^{2n}=0.$$
(11)
When a negative or complex $`\omega ^2`$ solves (11) the theory is exponentially unstable, a case we will not consider further. For all roots, besides the vanishing one, to be real and positive it is necessary (but not sufficient) that all $`c_{2n}>0`$.

Given the set of allowed frequencies $`\omega _i`$ we can write the Lagrangian as

$$=\frac{1}{2}A(t)\left[\underset{i0}{}(_t^2+\omega _i^2)\right]A(t),$$
(12)
after absorbing a constant normalization into $`A(t)`$. We will assume that $`\omega _j>\omega _i`$ for $`j>i`$, and $`\omega _0=0`$. This can be rewritten as a sum of conventional free field terms by defining111 This generalizes the definitions in .

$$\varphi _i(t)=\frac{\left[\underset{ji}{}(_t^2+\omega _j^2)\right]A(t)}{\left|\underset{ji}{}(\omega _i^2\omega _j^2)\right|^{1/2}}.$$
(13)
Note that $`A(t)\mathrm{cos}(\omega _it)`$ corresponds to all $`\varphi _j(t)`$ vanishing except for $`\varphi _i(t)`$. The Lagrangian then becomes

$$=\frac{1}{2}\underset{i0}{}()^i\varphi _i(t)(_t^2+\omega _i^2)\varphi _i(t).$$
(14)
The $`()^i`$ factor means that every other term comes in with the wrong sign, and in particular this is the case for the mode with the lowest nonvanishing frequency $`\omega _1`$. The wrong sign in the Lagrangian means that the classical energy of these modes is negative. The energy density of the $`i`$th oscillating mode with amplitude $`\epsilon _i`$ is thus proportional to $`()^i\epsilon _i^2`$. We will refer to this as a gravitational energy. These minus signs appearing in the linearized theory cause problems in the attempt to quantize higher derivative gravity.

These vacuum oscillating solutions of the linearized theory are of course not solutions of the full nonlinear theory. In particular an oscillating $`A(t)`$ will generate non-oscillating terms in the field equations at $`𝒪(\epsilon ^2)`$, and a source $`T_{\mu \nu }`$ must contribute at this order for a solution to exist. We will include both a cosmological constant $`\mathrm{\Lambda }`$ and matter with equation of state $`p(t)=w\rho (t)`$. With $`A(t)`$ in the form

$$A(t)=\epsilon \mathrm{cos}(\omega t)+a_2\epsilon ^2\mathrm{cos}(2\omega t)+\mathrm{}$$
(15)
we expand the field equations to $`𝒪(\epsilon ^2)`$ using the results for $`𝒢_\rho ^{(2n)}`$ and $`𝒢_p^{(2n)}`$ in the Appendix. A solution would consist in relating $`\epsilon `$ and $`a_2`$ to the gravitational sources. Solutions at higher order in $`\epsilon `$ could be found by continuing the expansion of $`A(t)`$ in the obvious way.222Solutions at $`𝒪(\epsilon ^3)`$ were found in the $`2+4`$ derivative case in .

Despite the complexity of the expressions in the Appendix, we find very simple results for the oscillating $`A(t)`$, and a clear pattern emerges in the terms with 2 through to 12 derivatives. We assume that this pattern will continue at yet higher derivatives.

$$\epsilon ^2\frac{3}{8}(+𝒜\mathrm{cos}(2\omega t))=\rho (t)+\mathrm{\Lambda }$$
(16)
$$\epsilon 𝒜\mathrm{cos}(\omega t)+\epsilon ^2(\frac{3}{8}𝒜+𝒞\mathrm{cos}(2\omega t))=p(t)\mathrm{\Lambda }$$
(17)
$$𝒜=\underset{n1}{}()^nc_{2n}\omega ^{2n}$$
(18)
$$=\underset{n1}{}()^{n1}(2n1)c_{2n}\omega ^{2n}$$
(19)
$$𝒞=𝒞(a_2,c_{2n},\chi _{2n},\omega )$$
(20)
We see that $`\rho (t)`$ at $`𝒪(\epsilon ^2)`$ is completely determined by the same constants $`c_i`$ that appear at linear order; the additional $`\chi _{2n}`$ parameters have all dropped out of this equation. The terms linear in $`\epsilon `$ in the $`p(t)`$ equation must vanish separately, thus determining the allowed frequencies as before: $`\omega \{\omega _i\}`$ such that $`𝒜=0`$. Thus at $`𝒪(\epsilon ^2)`$ we find that $`\rho (t)`$ is a constant.

In the $`p(t)`$ equation only the $`𝒞\mathrm{cos}(2\omega t)`$ term remains, which depends linearly on $`a_2`$ and the $`\chi _{2n}`$ parameters. The solution, to be compatible with $`p=w\rho `$, is thus completed by choosing $`a_2`$ so that $`𝒞=0`$. In this way the $`\chi _{2n}`$ parameters affect $`A(t)`$ only at order $`𝒪(\epsilon ^2)`$, through the value of $`a_2`$. The result is that the total pressure, $`p(t)\mathrm{\Lambda }`$, vanishes and thus

$$\rho (t)=\mathrm{\Lambda }/w.$$
(21)
Then to obtain a solution for the mode with frequency $`\omega _i`$, its amplitude $`\epsilon _i`$ is determined from (16) as

$$\epsilon _i^2\frac{3}{8}\underset{n1}{}()^{n1}(2n1)c_{2n}\omega _i^{2n}=\rho (1+w).$$
(22)
We need the value of $`=()^{n1}(2n1)c_{2n}\omega ^{2n}`$ for the various roots of $`𝒜=()^nc_{2n}\omega ^{2n}`$. It turns out that $``$ also has a set of positive roots that are interlaced with those of $`𝒜`$. The result is that the sign of $``$ alternates at each successive root of $`𝒜`$, with the implication (assuming $`w>1`$) that

$$\mathrm{sign}(\rho )=()^i$$
(23)
for the mode with frequency $`\omega _i`$. Thus the same modes that have negative gravitational energy also require a negative matter energy for their existence. In other words for these modes $`T_{00}+\mathrm{\Lambda }g_{00}<0`$ whereas the modes with positive gravitational energy have $`T_{00}+\mathrm{\Lambda }g_{00}>0`$. In both cases $`T_{ii}+\mathrm{\Lambda }g_{ii}=0`$. Thus in order not to violate positive energy conditions for matter, we may conclude that the set of oscillating solutions to the full theory does not include the modes with negative gravitational energy.

Note that the vanishing of $`p(t)\mathrm{\Lambda }`$ with $`\mathrm{\Lambda }0`$ requires $`w0`$. On the other hand if $`\mathrm{\Lambda }=0`$ then the same analysis would yield solutions for $`w=0`$. Thus pressureless matter can act as a source of the oscillations, with the amplitude of the oscillations determined by the energy density. In this case $`\rho `$ can be freely chosen, unlike the case with $`\mathrm{\Lambda }0`$.

Energy conservation for matter implies that $`\rho (t)`$ must oscillate in an oscillating background, since $`\rho (t)a(t)^{3(1+w)}`$ (where $`a(t)=e^{A(t)/2}`$ is the scale factor) must be a constant. Thus an oscillating component of $`\rho (t)`$ will occur at $`𝒪(\epsilon ^3)`$. In summary an oscillation amplitude $`\epsilon `$ of $`A(t)`$ corresponds to a source of $`𝒪(\epsilon ^2)`$ having an oscillating component of $`𝒪(\epsilon ^3)`$, a situation very different than the pure linearized theory. The implications of any of this for the development of a quantum theory will not be pursued here.

## 3 Oscillations and expansion

We turn now to consider the dynamics of these nonlinear oscillations in the context of the expanding universe. We will find that modes that were prohibited in the pure oscillating case become involved in a weak instability in the expanding case. This dynamics depends on the equation of state.

We will have to explore this numerically in theories truncated to various finite numbers of derivatives. We use the equation of state to eliminate $`\rho (t)`$ and $`p(t)`$ and obtain an equation for $`A(t)`$ alone. A single equation can be dealt with quite successfully numerically, even though it is nonlinear with high derivatives. Once $`A(t)`$ is obtained then the other equation determines $`\rho (t)`$. A good check on the accuracy of the numerical integration is provided by the quantity $`\rho (t)a(t)^{3(1+w)}`$, which should remain constant. We have confirmed this for all the cases studied.

We begin by truncating the equations at 4 derivatives and setting $`\mathrm{\Lambda }=0`$. We initially consider radiation domination $`w=1/3`$. Here it is easy enough to keep all the terms, to all orders in $`\epsilon `$. A range of initial conditions gives rise to the positive energy expanding solutions of interest to cosmology. In Fig. (1) we plot $`a(t)/a(t)_{\mathrm{FRW}}`$ and $`\rho (t)/\rho (t)_{\mathrm{FRW}}`$ for different ranges of $`t`$. For different initial conditions the size of the oscillating effect changes and the departure from the FRW solution can in particular be made much smaller, but the qualitative form of the plots remains the same. We will discuss the initial conditions further below; but for now we see that $`a(t)`$ increases faster, and $`\rho (t)`$ falls faster, than in the FRW case.

We can consider the departure from the FRW equation $`3H^2=\rho `$ as an indication of the importance of the higher derivative terms. This is seen in the plot of $`H^2(t)/H(t)_{\mathrm{FRW}}^2`$ and $`\rho (t)/\rho (t)_{\mathrm{FRW}}`$ in Fig. (2). We see the relative smoothness of $`\rho (t)`$ and the growing oscillations in $`H^2(t)/H(t)_{\mathrm{FRW}}^2`$, which implies the growing relative importance of the higher derivative terms. It is in this sense that the oscillations are growing and signaling a dynamical instability. We find that the dynamics is such that departure of $`\rho (t)/\rho (t)_{\mathrm{FRW}}`$ from unity grows like a power of $`\mathrm{ln}(t)`$, and thus the instability is weak and physically interesting.

Since the oscillations will not be directly observed it makes sense to consider quantities that are suitably time averaged, so as to smear out the oscillations. For example the apparent critical energy density of the universe will be given by $`\rho (t)_{\mathrm{crit}}3H(t)^2`$. It is then of interest to consider the quantity $`\rho (t)/\rho (t)_{\mathrm{crit}}`$, and we display this along with $`H(t)/H(t)_{\mathrm{FRW}}`$ in Fig. (3). We therefore effectively have $`\rho /\rho _{\mathrm{crit}}<1`$ in a flat universe, and thus the oscillations are playing the role of dark energy. But an accelerated expansion is not a result.

We now consider matter domination $`w=0`$. Here we find that the oscillations persist, but their relative importance does not grow. For example the size of the oscillations in the quantity $`H(t)^2/H(t)_{\mathrm{FRW}}^2`$ does not grow. This is reflected by the leveling off of $`\rho (t)/\rho (t)_{\mathrm{FRW}}`$ at some value determined by the size of the oscillations. Both the initial and final value of $`\rho (t)/\rho (t)_{\mathrm{FRW}}`$ is adjustable through the initial conditions. At the same time $`H(t)/H(t)_{\mathrm{FRW}}`$ tends to unity. The physically interesting quantity $`\rho (t)/\rho (t)_{\mathrm{crit}}`$ displays the least time dependence, and its essentially constant value, less than unity, is determined solely by the size of the oscillations. Thus oscillations that have grown during a radiation dominated era become “frozen in”, in relative importance, during the matter dominated era.

Pressureless matter $`w=0`$ turns out to be the dividing line between the growth of the oscillations, in relative terms, for $`w>0`$ and the fading away of the oscillations for $`w<0`$. In particular a positive cosmological constant eventually dominates the cosmological evolution and drives the effective $`w`$ negative. When this happens the oscillations fade away.

How natural it is to have the oscillations depends on how close to the Planck time the initial conditions are specified. The initial conditions may be characterized, for example, by how much the acceleration parameter deviates from the FRW value at that initial time. A given fractional deviation specified at one initial time produces much less oscillation than the same fractional deviation specified at an earlier initial time. Indeed at an initial time close to the Planck time $`t_{\mathrm{Pl}}`$, significant fine tuning is required to produce a solution close to the FRW one.

More appealing is to consider an initial time of order $`30t_{\mathrm{Pl}}`$. This was used to produce the $`w=1/3`$ results in Figs. (1,2,3), where a large fractional deviation, a factor of 10, in the acceleration parameter was used. For a small fractional deviation, say at the 10% level, the resulting oscillations are tiny, and they may not grow to be significant until cosmologically interesting times have elapsed. This is based on numerical analysis that indicates that the quantity $`\rho (t)/\rho (t)_{\mathrm{crit}}`$ in the radiation dominated era can be parameterized as $`1/(1+c\mathrm{ln}(t)^a)`$ for small $`c`$ depending on the initial conditions and $`a5`$. This is for times accessible to numerical analysis, whereas the era of radiation domination ends at $`t10^{28}t_{\mathrm{Pl}}`$, and so whether this approximation remains true throughout this era is uncertain. Clearly the main issue is to use the time of nucleosynthesis $`t10^{18}t_{\mathrm{Pl}}`$ to put a constraint on the deviations from standard FRW evolution.

We have studied the higher derivative equations in a similar way, using the results of the Appendix. This maintains energy conservation to $`𝒪(\epsilon ^3)`$, and we have verified the insensitivity to the $`\chi _{2n}`$ parameters. We will not display these numerical results, because they resemble the results from the four derivative equations. The new qualitative feature is the existence of more than one frequency mode of oscillation, where the relative size of the coefficients $`c_{2n}`$ determine how easy it is to excite the various modes simultaneously. The resulting oscillations can deviate quite dramatically from being pure sinusoidal or even periodic. But for the observable smeared out quantities, the results closely resemble those above.

We now look a little more closely at the origin of the result that $`\rho (t)<\rho (t)_{\mathrm{FRW}}`$. In the pure oscillation case the modes with negative gravitational energy required negative $`\rho (t)`$. Here it is clear that these modes are being excited and that they are responsible for driving $`\rho (t)<\rho (t)_{\mathrm{FRW}}`$. In fact in the 4 derivative case there is only one oscillating mode and it has negative gravitational energy. At higher derivatives (22) suggests that contributions to $`\rho (t)`$ alternate in sign depending on the order of derivatives. This is displayed in Fig. (4) in the case of the equations with up to 6 derivatives and $`c_6=1/5`$. The relative sizes of the contributions are consistent with (22) given that the lowest mode with $`\omega =\omega _1`$ dominates. The conclusion is that we are linking a dynamical instability in the expanding system to the growth of a negative component of the total gravitational energy.

As a side issue we mention that when the cosmological constant is negative there appear to be solutions with oscillations on both large and small scales , where the frequency of the large scale oscillation is proportional to $`\sqrt{\mathrm{\Lambda }}`$. Our results now allow us to determine the sign of $`\rho `$ needed to support these solutions at any order in the number of derivatives, and we find that it is typically negative. This is not surprising since the amplitude of the large time scale oscillations can be adjusted to zero with appropriate initial conditions; this returns us to the pure oscillation case we have treated above where we found negative $`\rho (t)`$ for some modes, in particular for the lowest frequency. This leaves open the possibility that with appropriate initial conditions it may be possible to excite the modes with positive energy in such a way as to have oscillations on large as well as small time scales. The practical interest in this is unclear.

In this paper we have been concerned with a dynamical instability of standard FRW expanding universes, where the matter energy density remains positive, but is reduced in association with the excitation of gravitation modes with negative energy. This is found to be a generic feature of higher derivative gravity with any number of derivatives. The phenomenological implications of this will have to be explored elsewhere, but the instability appears to be of physical interest due to its slow logarithmic growth. This instability is incompatible with a positive cosmological constant. This observation may be interest in a theory in which the effective 4 dimensional cosmological constant is a derived quantity, perhaps emerging from a higher dimensional theory . The instability may provide a dynamical reason for the 4D flat spacetime to be preferred. At the same time the oscillations and their gravitational coupling to matter may provide a mechanism by which positive vacuum energy is converted through particle production to matter energy.

## Appendix

We list here results for the quantities defined in (8-9), to the required orders in the small amplitude expansion, and up to twelve derivatives.

$$𝒢_{p}^{}{}_{}{}^{(2)}=(\frac{^2}{t^2}A(t))+\frac{3}{4}(\frac{}{t}A(t))^2$$

$$𝒢_{\rho }^{}{}_{}{}^{(2)}=\frac{3}{4}(\frac{}{t}A(t))^2$$

$$𝒢_{p}^{}{}_{}{}^{(4)}=(\frac{^4}{t^4}A(t))+\frac{9}{4}(\frac{^2}{t^2}A(t))^2+3(\frac{}{t}A(t))(\frac{^3}{t^3}A(t))+\frac{9}{4}(\frac{^2}{t^2}A(t))(\frac{}{t}A(t))^2$$

$$𝒢_{\rho }^{}{}_{}{}^{(4)}=\frac{3}{4}(\frac{^2}{t^2}A(t))^2+\frac{3}{2}(\frac{}{t}A(t))(\frac{^3}{t^3}A(t))+\frac{9}{4}(\frac{^2}{t^2}A(t))(\frac{}{t}A(t))^2$$

$`𝒢_{p}^{}{}_{}{}^{(6)}=(\frac{^6}{t^6}A(t))+(\frac{27}{4}+\chi _6)(\frac{^3}{t^3}A(t))^2+(\frac{21}{2}+\chi _6)(\frac{^2}{t^2}A(t))(\frac{^4}{t^4}A(t))`$
$`+\frac{9}{2}(\frac{}{t}A(t))(\frac{^5}{t^5}A(t))`$
$`𝒢_{\rho }^{}{}_{}{}^{(6)}=\frac{3}{4}(\frac{^3}{t^3}A(t))^2\frac{3}{2}(\frac{^4}{t^4}A(t))(\frac{^2}{t^2}A(t))+\frac{3}{2}(\frac{}{t}A(t))(\frac{^5}{t^5}A(t))`$
$`+\frac{9}{2}(\frac{^4}{t^4}A(t))(\frac{}{t}A(t))^2(9\frac{3}{2}\chi _6)(\frac{}{t}A(t))(\frac{^3}{t^3}A(t))(\frac{^2}{t^2}A(t))`$
$`(3+\frac{1}{2}\chi _6)(\frac{^2}{t^2}A(t))^3`$
$`𝒢_{p}^{}{}_{}{}^{(8)}=(\frac{^8}{t^8}A(t))+(\frac{89}{4}+6\chi _8)(\frac{^4}{t^4}A(t))^2+(36+9\chi _8)(\frac{^3}{t^3}A(t))(\frac{^5}{t^5}A(t))`$
$`+(19+3\chi _8)(\frac{^2}{t^2}A(t))(\frac{^6}{t^6}A(t))+6(\frac{^7}{t^7}A(t))(\frac{}{t}A(t))`$
$`𝒢_{\rho }^{}{}_{}{}^{(8)}=\frac{3}{4}(\frac{^4}{t^4}A(t))^2+\frac{3}{2}(\frac{^5}{t^5}A(t))(\frac{^3}{t^3}A(t))\frac{3}{2}(\frac{^6}{t^6}A(t))(\frac{^2}{t^2}A(t))`$
$`+\frac{3}{2}(\frac{^7}{t^7}A(t))(\frac{}{t}A(t))+\frac{27}{4}(\frac{}{t}A(t))^2(\frac{^6}{t^6}A(t))(\frac{69}{4}+\frac{9}{2}\chi _8)(\frac{^2}{t^2}A(t))^2(\frac{^4}{t^4}A(t))`$
$`(\frac{69}{4}\frac{9}{2}\chi _8)(\frac{}{t}A(t))(\frac{^2}{t^2}A(t))(\frac{^5}{t^5}A(t))`$
$`(\frac{69}{2}9\chi _8)(\frac{}{t}A(t))(\frac{^3}{t^3}A(t))(\frac{^4}{t^4}A(t))`$
$`𝒢_{p}^{}{}_{}{}^{(10)}=(\frac{^{10}}{t^{10}}A(t))+(\frac{223}{4}+\chi _{10})(\frac{^5}{t^5}A(t))^2+(\frac{191}{2}+\frac{8}{5}\chi _{10})(\frac{^4}{t^4}A(t))(\frac{^6}{t^6}A(t))`$
$`+(59+\frac{4}{5}\chi _{10})(\frac{^3}{t^3}A(t))(\frac{^7}{t^7}A(t))+(26+\frac{1}{5}\chi _{10})(\frac{^2}{t^2}A(t))(\frac{^8}{t^8}A(t))`$
$`+\frac{15}{2}(\frac{}{t}A(t))(\frac{^9}{t^9}A(t))`$
$`𝒢_{\rho }^{}{}_{}{}^{(10)}=\frac{3}{4}(\frac{^5}{t^5}A(t))^2\frac{3}{2}(\frac{^4}{t^4}A(t))(\frac{^6}{t^6}A(t))+\frac{3}{2}(\frac{^3}{t^3}A(t))(\frac{^7}{t^7}A(t))`$
$`\frac{3}{2}(\frac{^8}{t^8}A(t))(\frac{^2}{t^2}A(t))+\frac{3}{2}(\frac{}{t}A(t))(\frac{^9}{t^9}A(t))`$
$`(\frac{165}{2}\frac{3}{2}\chi _{10})(\frac{^4}{t^4}A(t))(\frac{}{t}A(t))(\frac{^5}{t^5}A(t))`$
$`(63\frac{9}{10}\chi _{10})(\frac{^3}{t^3}A(t))(\frac{^6}{t^6}A(t))(\frac{}{t}A(t))`$
$`(\frac{33}{2}\frac{3}{10}\chi _{10})(\frac{^4}{t^4}A(t))(\frac{^3}{t^3}A(t))^2`$
$`(\frac{93}{4}\frac{3}{10}\chi _{10})(\frac{}{t}A(t))(\frac{^7}{t^7}A(t))(\frac{^2}{t^2}A(t))+9(\frac{^8}{t^8}A(t))(\frac{}{t}A(t))^2`$
$`(33+\frac{3}{5}\chi _{10})(\frac{^4}{t^4}A(t))^2(\frac{^2}{t^2}A(t))(\frac{33}{2}+\frac{3}{10}\chi _{10})(\frac{^3}{t^3}A(t))(\frac{^5}{t^5}A(t))(\frac{^2}{t^2}A(t))`$
$`(\frac{93}{4}+\frac{3}{10}\chi _{10})(\frac{^2}{t^2}A(t))^2(\frac{^6}{t^6}A(t))`$
$`𝒢_{p}^{}{}_{}{}^{(12)}=(\frac{^{12}}{t^{12}}A(t))+(318+\chi _{12})(\frac{^6}{t^6}A(t))^2+(549+\frac{8}{5}\chi _{12})(\frac{^5}{t^5}A(t))(\frac{^7}{t^7}A(t))`$
$`+(\frac{1395}{4}+\frac{29}{35}\chi _{12})(\frac{^4}{t^4}A(t))(\frac{^8}{t^8}A(t))+(\frac{315}{2}+\frac{2}{7}\chi _{12})(\frac{^3}{t^3}A(t))(\frac{^9}{t^9}A(t))`$
$`+(48+\frac{2}{35}\chi _{12})(\frac{^2}{t^2}A(t))(\frac{^{10}}{t^{10}}A(t))+9(\frac{}{t}A(t))(\frac{^{11}}{t^{11}}A(t))`$
$`𝒢_{\rho }^{}{}_{}{}^{(12)}=\frac{3}{4}(\frac{^6}{t^6}A(t))^2+\frac{3}{2}(\frac{^5}{t^5}A(t))(\frac{^7}{t^7}A(t))\frac{3}{2}(\frac{^4}{t^4}A(t))(\frac{^8}{t^8}A(t))`$
$`+\frac{3}{2}(\frac{^3}{t^3}A(t))(\frac{^9}{t^9}A(t))\frac{3}{2}(\frac{^{10}}{t^{10}}A(t))(\frac{^2}{t^2}A(t))+\frac{3}{2}(\frac{}{t}A(t))(\frac{^{11}}{t^{11}}A(t))`$
$`(\frac{3825}{8}\frac{3}{2}\chi _{12})(\frac{}{t}A(t))(\frac{^5}{t^5}A(t))(\frac{^6}{t^6}A(t))`$
$`(\frac{207}{4}\frac{3}{35}\chi _{12})(\frac{}{t}A(t))(\frac{^9}{t^9}A(t))(\frac{^2}{t^2}A(t))`$
$`(\frac{315}{4}+\frac{6}{35}\chi _{12})(\frac{^3}{t^3}A(t))(\frac{^7}{t^7}A(t))(\frac{^2}{t^2}A(t))`$
$`(\frac{2115}{8}+\frac{51}{70}\chi _{12})(\frac{^4}{t^4}A(t))(\frac{^6}{t^6}A(t))(\frac{^2}{t^2}A(t))+\frac{45}{4}(\frac{^{10}}{t^{10}}A(t))(\frac{}{t}A(t))^2`$
$`(\frac{207}{4}+\frac{3}{35}\chi _{12})(\frac{^2}{t^2}A(t))^2(\frac{^8}{t^8}A(t))(\frac{855}{8}+\frac{27}{70}\chi _{12})(\frac{^5}{t^5}A(t))^2(\frac{^2}{t^2}A(t))`$
$`(\frac{855}{8}\frac{27}{70}\chi _{12})(\frac{^3}{t^3}A(t))(\frac{^5}{t^5}A(t))(\frac{^4}{t^4}A(t))`$
$`(\frac{315}{4}\frac{6}{35}\chi _{12})(\frac{^3}{t^3}A(t))^2(\frac{^6}{t^6}A(t))`$
$`(\frac{2745}{8}\frac{9}{10}\chi _{12})(\frac{}{t}A(t))(\frac{^7}{t^7}A(t))(\frac{^4}{t^4}A(t))`$
$`(\frac{729}{4}\frac{12}{35}\chi _{12})(\frac{}{t}A(t))(\frac{^8}{t^8}A(t))(\frac{^3}{t^3}A(t))(\frac{285}{8}+\frac{9}{70}\chi _{12})(\frac{^4}{t^4}A(t))^3`$
## Acknowledgments

This research was supported in part by the Natural Sciences and Engineering Research Council of Canada.


text:
# INTERACTION OF COUPLED HIGHER ORDER NONLINEAR SCHRÖDINGER EQUATION SOLITONS

## I Introduction

The subject of integrable model is very fascinating largely because of its innumerable symmetries and a special class of solutions known as soliton solutions. Only a few systems in the nature are known to be integrable. The coupled nonlinear Schrödinger equation (CNLS) and its higher order generalization namely coupled higher order nonlinear Schrödinger equation (CHNLS) are some of the examples of integrable equations which have direct relevance in the propagation of optical solitons in Kerr type nonlinear fibre. Coupled integrable systems have many important applications in photorefractive crystals 1  and also in all optical computations 2 . The cross phase modulation (CPM) phenomena in CNLS equation provides an interesting pulse shepherding effect to align the time arrival of the pulses 3 . The CPM phenomena along with group velocity dispersion (GVD) can also be utilized in compressing optical pulses at the one soliton level 4 . Recently the existence of multi-component solitons are experimentally established 5 .

In this paper we have considered some interesting features of CHNLS equation, which describes the propagation of the optical pulses of very short width, of the order of $`10^{15}`$sec. We have studied the novel inelastic collision of two solitons for an $`n`$-component CHNLS equation. The inelastic collisions among the three solitons are also studied. The three soliton interactions may be interpreted as a combination of two-soliton interactions occuring at three different points. This is in compliance with the three particle interactions in integrable systems. The energy exchange among the components of three solitons may occur at any of the three two-soliton interaction points. In this context, it is important to note that the phenomena of shape changing vis a vis inelastic collisions and their consequences have been observed for the CNLS equation by considering two soliton solutions of both two and three components and subsequently $`n`$ component generalization of two soliton solutions and n-complexes have been studied 2 . But the soliton solution being nonperturbative, the shape changing phenomena of CNLS equation does not ensure the same for CHNLS equation. Moreover, from integrability point of view study of the collisions of three solitons becomes indispensable and one needs to address this issue seperately.

The paper is organized in the following sequence. In section I the CHNLS equation is introduced and its $`N`$ soliton solution is given explicitly. In section III two-soliton interactions and its asymptotic analysis are obtained showing the novel inelastic collisions among the components of the solitons. Some interesting features of three soliton interactions, related to the integrability of the $`n`$-component CHNLS equation are addressed in section IV. Section V is the concluding one.

## II CHNLS equation and N-soliton solution

The $`n`$-component CHNLS equation incorporating the effects of Kerr type nonlinearity and stimulated Raman scattering may be written as

$$E_{iz}+iE_{i\tau \tau }+2i(\underset{j=1}{\overset{n}{}}E_j^{}E_j)E_i+\epsilon E_{i\tau \tau \tau }+3\epsilon (\underset{j=1}{\overset{n}{}}E_j^{}E_j)E_{i\tau }+3\epsilon (\underset{j=1}{\overset{n}{}}E_j^{}E_{j\tau })E_i=0$$
(1)
where $`E_i`$ is $`n`$-component electric field. $`z`$ and $`\tau `$ denote the direction of propagation and time variable respectively. The parameter $`\epsilon `$ is the ratio of the width of the spectra $`\mathrm{\Delta }\omega `$ to the carrier frequency $`\omega `$ such that $`\epsilon =\frac{\mathrm{\Delta }\omega }{\omega }<1`$.

A particular gauge equivalent form of (1) may be written as 6

$$q_{ix}+\epsilon q_{ittt}+3\epsilon \underset{j=1}{\overset{n}{}}(q_j^{}q_j)q_{it}+3\epsilon \underset{j=1}{\overset{n}{}}(q_j^{}q_{jt})q_i=0$$
(2)
where the transformation relations are given by

$`E_i(z,\tau )=q_i(x,t)e^{i(\frac{1}{27\epsilon ^2}x\frac{1}{3\epsilon }t)}`$
$`x=z`$
$`t=\tau {\displaystyle \frac{1}{3ϵ}}z`$
The transformed equation (2) is known as the coupled complex modified KdV equation (CCMKdV) whose $`N`$ soliton solutions for $`n`$-component field is known 6 . Notice that the CCMKdV equation in (2) is convenient particularly for the Lax representation of the system and consequently for applying the inverse scattering method to obtain $`N`$-soliton solution. However, as a result of the transformations (II), there is no change of the envelope function, but for a constant shift in velocity for all solitons. Consequently, the shape of the pulses for CCMKdV equation and CHNLS equation remains the same i.e., $`|E_i|=|q_i|`$ and the transformation will not affect the subsequent results of the two gauge equivalent systems.

The $`N`$\- soliton solution of the $`n`$ component CCMKdV equation may be written in the most compact form as 6 ,

$$q_i(x,t)=2\underset{j=1}{\overset{N}{}}(\text{B}\text{C}^1)_{ij}e^{i\lambda _j^{}t}$$
(3)
where, B and C are respectively $`n\times N`$ and $`N\times N`$ matrices whose explicit forms are given by

$$(\text{B})_{ij}=i𝒞_{n+1i}^{(j)}(0)e^{8i\epsilon \lambda _j^3xi\lambda _j^{}t}$$
(4)
and

$$(\text{C})_{ij}=\underset{p=1}{\overset{n}{}}\underset{k=1}{\overset{N}{}}\frac{𝒞_{n+1p}^{(j)}(0)𝒞_{n+1p}^{(k)}(0)}{\alpha _{n+1n+1}^{}(\lambda _j^{})\alpha _{n+1n+1}^{}(\lambda _k)}\frac{e^{i(\lambda _i^{}+\lambda _j^{}2\lambda _k)t+8i\epsilon (\lambda _k^3\lambda _j^3)x}}{(\lambda _k\lambda _j^{})(\lambda _k\lambda _i^{})}\delta _{ij}$$
(5)
In the above equation $`\alpha _{ij}`$ are the elements of the $`(n+1)\times (n+1)`$ scattering matrix and $`\lambda `$ is the spectral parameter. $`𝒞_{(n+1)p}^{(j)}`$ is related to the elements of the scattering matrix at the position of the simple poles, $`𝒞_{(n+1)p}^{(j)}=\alpha _{(n+1)p}(\lambda _j^{})|_{x=0}`$ and $``$ over $`\alpha `$ denotes the derivative with respect to $`\lambda `$.

The one soliton solution (1SS) for n-coupled system directly follows from (3,4,5) by considering $`N=1`$ and is given as

$$q_i(x,t)=\frac{𝒞_{n+1i}^{(1)}e^{\frac{R_{11}}{2}}e^{i\eta _{1I}}}{\mathrm{cosh}(\eta _{1R}+\frac{R_{11}}{2})}$$
(6)
where $`\eta _{1R}`$ and $`\eta _{1I}`$ respectively denote the real and imaginary parts of $`\eta _1=2i\lambda _1^{}t8i\epsilon \lambda _1^3x+i\frac{\pi }{2}`$ and $`e^{R_{11}}=\frac{\kappa _{11}}{\lambda _{11}^2}`$. The complex variables are defined as

$$\kappa _{11}=\underset{p=1}{\overset{n}{}}\left|\frac{𝒞_{(n+1)p}^{(1)}}{\alpha _{n+1n+1}^{}(\lambda _1)}\right|^2$$

and

$$\lambda _{11}=\lambda _1\lambda _1^{}$$

For convenience let us define $`\lambda _1=\frac{\lambda _{1R}+i\lambda _{1I}}{2}`$, where the subscripts $`R`$ and $`I`$ denote the real and imaginary parts and consequently $`\eta _{1R}`$ and $`\eta _{1I}`$ become

$$\eta _{1I}=\lambda _{1R}t+\epsilon \lambda _{1R}(3\lambda _{1I}^2\lambda _{1R}^2)x+\frac{\pi }{2}$$
(7)
$$\eta _{1R}=\lambda _{1I}[t\epsilon (\lambda _{1I}^23\lambda _{1R}^2)x]$$
(8)
From the argument of the $`\mathrm{cosh}`$ function of the 1SS (6) and (8), it is straight forward to identify the width of the soliton pulse as $`\mathrm{\Gamma }=|\lambda _{1I}|^1`$ and the soliton travels with a group velocity $`V_g=[ϵ(\lambda _{1I}^23\lambda _{1R}^2)]^1`$ in the positive $`x`$ direction when $`\lambda _{1I}^2>3\lambda _{1R}^2`$ and in the negative $`x`$ direction when $`\lambda _{1I}^2<3\lambda _{1R}^2`$.

## III  Two-Soliton interaction

The two soliton solution (2SS) can be obtained from (3,4,5) by putting $`N=2`$, the explicit form of the $`i^{ih}`$ component of 2SS being

$$q_i=\frac{2_{m=1}^2\left(𝒞_{n+1i}^{(m)}e^{\eta _m}+e^{\eta _1+\eta _2+\eta _m^{}+\delta _m^{(i)}}\right)}{Det[C]}$$
(9)
where,

$$e^{\delta _m^{(i)}}=\underset{r,s=1}{\overset{2}{}}𝒞_{n+1i}^{(r)}\frac{\kappa _{ms}}{\lambda _{ms}}\left(\frac{1}{\lambda _{mr}}\frac{1}{\lambda _{ms}}\right)$$
(10)
In (10) $`\lambda _{ij}=\lambda _i\lambda _j^{}`$ and $`\eta _m=2i\lambda _m^{}t8i\epsilon \lambda _m^3x+i\frac{\pi }{2}`$. $`Det[C]`$ is the determinant of $`2\times 2`$ dimensional matrix obtained from $`(\text{C})_{ij}`$ (5). The form of $`Det[C]`$ in terms of the spectral parameters emerges rather lengthy but it will be useful for the asymptotic analysis of the interacting solitons and will be seen later. In terms of the spectral parameters $`Det[C]`$ is given as

$$Det[C]=1+\underset{r,s=1}{\overset{2}{}}e^{\eta _r+\eta _s^{}+R_{rs}}+e^{2\eta _{1R}+2\eta _{2R}+𝒬}$$
(11)
with

$$e^{R_{rs}}=\frac{\kappa _{rs}}{\lambda _{rs}^2},$$
(12)
$$\kappa _{rs}=\underset{p=1}{\overset{n}{}}\frac{𝒞_{n+1p}^{(s)}𝒞_{n+1p}^{(r)}}{\alpha _{(n+1)(n+1)}^{}(\lambda _r^{})\alpha _{(n+1)(n+1)}^{}(\lambda _s)}$$
(13)
and $`e^𝒬=Det[\kappa ]Det[],[\kappa ]=(\frac{\kappa _{ij}}{\lambda _{ij}})`$ (no sum over $`i,j`$), $`[]=(\lambda _{ij})^1`$, $`[\kappa ]`$ and $`[]`$ being $`2\times 2`$ matrices. To understand the interaction in a more explicit manner we analyse the asymptotic limits of the 2SS (9). The asymptotic limit may be obtained by observing the 2SS when both the solitons are infinitely apart. This may be achieved by taking the limit $`\eta _{2R}\pm \mathrm{}`$. As a consequence $`i^{th}`$ component of the remaining soliton acquires the form

$$q_i^{(l\pm )}=\frac{\lambda _{lI}A_i^{(l\pm )}e^{i\eta _{lI}}}{\mathrm{cosh}(\eta _{lR}+\mathrm{\Phi }^{l\pm })}$$
(14)
where $`A_i^{(l\pm )}`$ and $`\mathrm{\Phi }^{(l\pm )}`$ may be interpreted as the amplitudes and phases of $`i^{th}`$ component of $`l^{th}`$ soliton respectively. Notice that in general $`A_i^{(l)}`$ are different from $`A_i^{(l+)}`$. This may occur due to exchange of energy between the solitons. The amplitude $`A_i^{(l)}`$ in terms of $`𝒞_{n+1i}^{(l)}`$ is given as

$$A_i^{(l)}=\frac{𝒞_{n+1i}^{(l)}}{\sqrt{\kappa _{ll}}}$$
(15)
The expression for $`A_i^{(l+)}`$ however, is more involved and may be written in terms of $`A_i^{(l)}`$ as

$$A_i^{(l+)}=A_i^{(l)}T_i^l$$
(16)
$`T_i^l`$ may be interpreted as the transition matrix and is defined as

$$T_i^1=\frac{(1\frac{𝒞_{n+1i}^{(2)}}{𝒞_{n+1i}^{(1)}}\mathrm{\Lambda }_2)}{\sqrt{1\mathrm{\Lambda }_1\mathrm{\Lambda }_2}}\sqrt{(\frac{\lambda _2^{}\lambda _1^{}}{\lambda _2\lambda _1})(\frac{\lambda _2^{}\lambda _1}{\lambda _2\lambda _1^{}})}$$
(17)
and

$$T_i^2=\frac{(1\frac{𝒞_{n+1i}^{(1)}}{𝒞_{n+1i}^{(2)}}\mathrm{\Lambda }_2)}{\sqrt{1\mathrm{\Lambda }_1\mathrm{\Lambda }_2}}\sqrt{(\frac{\lambda _2^{}\lambda _1^{}}{\lambda _2\lambda _1})(\frac{\lambda _2\lambda _1^{}}{\lambda _2^{}\lambda _1})}$$
(18)
where $`\mathrm{\Lambda }_1=\frac{\kappa _{12}\lambda _{11}}{\kappa _{11}\lambda _{12}}`$, $`\mathrm{\Lambda }_2=\frac{\kappa _{21}\lambda _{22}}{\kappa _{22}\lambda _{21}}`$. If $`|T_i^l|=1`$, the solitons pass through each other without being affected in their shapes and sizes. Otherwise, we will see that the solitons exchange energy at the time of interaction. The phase shift $`\mathrm{\Phi }^{(l)}`$ as a result of collision may be obtained from the following relation.

$$\mathrm{\Phi }^{(l)}=\mathrm{\Phi }^{(l+)}\mathrm{\Phi }^{(l)}$$
(19)
where $`\mathrm{\Phi }^{(l)}=\frac{1}{2}(R_{11})`$ and $`\mathrm{\Phi }^{(l+)}=\frac{1}{2}(𝒬R_{22})`$ and consequently the phase shift becomes $`\mathrm{\Phi }^{(l)}=\frac{1}{2}(𝒬R_{11}R_{22})`$. It is now interesting to note that for equal values of $`𝒞_{n+1i}^{(l)}`$ for each component of the solitons, $`|T_i^l|`$ becomes unity and consequently there is no exchange of energy among the components of each soliton. However, if the relative phases among $`𝒞_{n+1i}^{(l)}`$ are introduced the energy of each individual component no longer remains constant due to collisions and in that case $`|T_i^l|1`$. This analysis has been presented graphically by plotting two soliton interaction for a two component system with the $`\epsilon =0.1`$, $`\lambda _1=1.5i`$ and $`\lambda _1=.52i`$. One of the solitons is moving with group velocity $`1.74`$ and the other is moving with a group velocity $`3.07`$.

In figures(1) two soliton interaction is plotted with $`𝒞_{n+1i}^{(l)}=1`$ for each component, where the energy profile of each component of the solitons remains unchanged. In figures(2) two soliton interaction is plotted with nontrivial relative phases, $`𝒞_{n+\mathrm{1\; 1}}^{(1)}=𝒞_{n+\mathrm{1\; 2}}^{(1)}=𝒞_{n+\mathrm{1\; 2}}^{(2)}=1`$ and $`𝒞_{n+\mathrm{1\; 1}}^{(2)}=46(1i)`$ showing a considerable exchange of energy among the components. figure(3) shows that the overall intensity profile, represented by $`|q|`$ of the 2SS remains unchanged although there are appreciable energy exchange among the components of a soliton. It is an interesting result in the context of all optical computations 2  leading to the construction of the logical binary gates.

## IV Three-soliton interactions

The explicit form of the three soliton solution (3SS) can be obtained from (3,4,5) by putting $`N=3`$, the explicit form of the $`i^{ih}`$ component of 3SS being

$`q_p={\displaystyle \frac{2}{Det[C]}}\left[{\displaystyle \underset{i=1}{\overset{3}{}}}𝒞_{n+1p}^{(i)}e^{\eta _i}+{\displaystyle \underset{i,j,m=1_{ijm}}{\overset{3}{}}}(adj\stackrel{~}{\kappa }^j)_{mm}(adj\stackrel{~}{}^j)_{mm}e^{2\eta _{iR}+\eta _j}+{\displaystyle \underset{m=1}{\overset{3}{}}}Det[\stackrel{~}{\kappa }^m]Det[\stackrel{~}{}^m]e^{{\displaystyle \underset{{\scriptscriptstyle \genfrac{}{}{0pt}{}{i,m=1}{im}}}{\overset{3}{}}}2\eta _{iR}+\eta _m}\right]`$
where,

$`Det[C]=1{\displaystyle \underset{i,j=1}{\overset{3}{}}}e^{\eta _i+\eta _j^{}}e^{R_{ij}}{\displaystyle \underset{i=1}{\overset{3}{}}}{\displaystyle \underset{\genfrac{}{}{0pt}{}{k,j=1}{ijk}}{\overset{3}{}}}(adj\kappa )_{jk}(adj)_{jk}e^{2\eta _{iR}+\eta _k+\eta _j^{}}{\displaystyle \underset{\genfrac{}{}{0pt}{}{i,j,k=1}{ijk}}{\overset{3}{}}}(adj\kappa )_{kk}(adj)_{kk}e^{2\eta _{iR}+2\eta _{jR}}`$
$`+Det[\kappa ]Det[]e^{2\eta _{1_R}+2\eta _{2_R}+2\eta _{3_R}}`$
with $``$ is a $`3\times 3`$ dimensional matrix with elements $`(\lambda _{ij}^1)`$ and $`[\kappa ]`$ is also a $`3\times 3`$ dimensional matrix with elements $`\frac{\kappa _{ij}}{\lambda _{ij}}`$ (no sum over $`i,j`$). \[$`\stackrel{~}{\kappa }^j`$\] is the matrix \[$`\kappa `$\], where the $`j`$th row has been replaced by $`\genfrac{}{}{0pt}{}{\genfrac{}{}{0pt}{}{3}{}}{i=1}e_{ji}𝒞_{n+1p}^{(i)}`$ and $`\stackrel{~}{}^j`$ is the matrix where the $`j`$th row has been replaced by a unit row vector, $`\genfrac{}{}{0pt}{}{\genfrac{}{}{0pt}{}{3}{}}{i=1}e_{ji}`$. The phenomena of shape changing due to exchange of energy among the components of solitons becomes more interesting in the three-soliton interaction case. It is observed that the interactions take place at three different space time points with two solitons interact at each time. Thus the energy exchage among the components in three solitons interaction may occur at three different points. This is demonstrated in  figure(4) by the contour plot of the three solitons interaction with $`ϵ=0.1,\lambda _1=1.5i,\lambda _2=0.01+2i,\lambda _3=0.9+2i,𝒞_{n+\mathrm{1\; 1}}^{(1)}=𝒞_{n+\mathrm{1\; 1}}^{(2)}=𝒞_{n+\mathrm{1\; 1}}^{(3)}=10^3,𝒞_{n+\mathrm{1\; 2}}^{(2)}=𝒞_{n+\mathrm{1\; 2}}^{(3)}=1`$ and $`𝒞_{n+\mathrm{1\; 2}}^{(1)}=10^30.8i`$. Two of the solitons are moving in the positive $`x`$ direction with group velocities $`2.5`$ and $`6.36`$ and the third soliton is moving in the negative $`x`$ direction with group velocity $`1.74`$. The nature of three soliton interactions opens up another possibility of three solitons interaction where second and third solitons interact first instead of first and second solitons leaving the final configuration same in both the cases.

## V Conclusion

In conclusion, we have shown the novel shape changing phenomena asscoiated with two soliton interaction for $`n`$ component CHNLS equation. The three solitons interaction also exhibits energy exchange, but the exchange may occur at three points. This may lead to a more flexibility in constructing logical gates in all optical computing systems. The three solitons interaction also demonstartes an interesting feature conforming the exact integrability of the system  a la Zamolodchikov 7  and indicates the existence of Yang Baxter like relation, which will be publishedelsewhere.

Acknowledgments. S. G. and S. N. would like to thank CSIR Govt. of India for financial support under the project 03(0896)/99/EMR II. S. G. also gratefully acknowledge the hospitality of the Organizing Committee of the International Conference “Geometry, Integrability and Nonlinearity in Condensed Matter and Soft Condensed Matter Physics”, where the present results were first reported.


text:
# Quantum Computatinal Gates with Radiation Free Couplings

## Figure Captions

Fig.1: $`\mathrm{\Lambda }`$-shaped level configuration of the persistent-current normal-state Aharonov-Bohm qubit. 1,2,3 are the eigenstates of Hamiltonian (1) at $`\alpha =\pi /3`$ where 1 and 3 are the computational basis (qubit) registers $`|0`$, $`|1`$, whereas 2 is the qugate (control) register $`|c`$. Arrows show transitions between the degenerate states effected by the potential biases applied to the ring sites.

Fig.2: (a)A sketch of magnetically focused lines of magnetic field (arrows) of the superconducting fluxon with flux $`\mathrm{\Phi }_1=hc/2e`$ making one half of the normal-metal flux quantum, $`\mathrm{\Phi }_0=hc/e`$, and effecting the ring $`𝐑`$ with three normal islands into a $`\mathrm{\Lambda }`$-shaped configuration. The fluxon $`\mathrm{\Phi }_1`$ is trapped in the opening of superconducting foil ($`𝐒`$) and further compressed by a ferromagnetic crystal to fit into the interior of the ring. (b)Schematic of the multi-ring qubit arrangement with the sites (circles) on the surface of cylindrical wall $`𝐑`$ surrounding ferromagnetic cylinder ($`𝐅`$) which focuses lines of magnetic field in a cylindrical tube inside superconductor ($`𝐒`$).

Fig.3: A sketch of the electric field (shown by arrows) applied to the ring through the potential electrodes ($`𝐕`$) in direction perpendicular to the direction of magnetic field. $`𝐂`$ is a coupling loop providing the connection to the nondemolition-measuring setup of the persistent current.

Fig.4: Eigenenergies of the ring biased with a potential $`V_0`$. Eigenstates 1 and 3 are degenerate at $`V_0=0`$ where they form the qubit states. The commensurate situation, marked by arrows, appears at $`K=1`$ and at $`K=3`$ where it allows for the temporal (virtual) transition to a higher level 2 and back thus effecting the bit-flip transition (at $`K=1`$) and the Hadamard-like gate (at $`K=3`$).

Fig.5: Bit evolution at $`K=1`$. At point indicated by an arrow ($`t=t_1`$), the population of control register ($`|c`$) vanishes whereas the populations of the computational registers of qubit ($`|0`$) and ($`|1`$), interchange.

Fig.6: Bit evolution at $`K=3`$. At point indicated by an arrow ($`t=t_3`$), the population of control register ($`|c`$) vanishes whereas the computational basis of the qubit, originally in a state ($`|1`$, equally populates to states $`|0`$ and $`|1`$.

Fig.7: A sketch of the $`\mathrm{𝐂𝐍𝐎𝐓}`$ quantum gate. The ring No.1 which is pierced by a positive-$`\mathrm{\Phi }_0/2`$ flux (marked “+”), after the subtraction of the negative-$`\mathrm{\Phi }_0/2`$ flux (“-”), couples via the loop $`𝐂`$ to a quantum Hall bar (H) through which a fixed current (J) is fed. The bar generates a Hall voltage output effecting, through the potential electrodes $`𝐕`$, the flip transition in the ring No.2.

Fig.1

Fig.2

Fig.3

Fig.4

Fig.5

Fig.6

Fig.7