text:
# New UBVRI colour distributions in E-type galaxies

## 1 Introduction

The ”classical” data on the large scale colour distributions of E-type galaxies relies on observations by Bender and Möllenhof (1987), Vigroux et al. (1988), Franx et al. (1989), Peletier et al. (1990), Goudfrooij et al. (1994), to quote only the papers discussing the 1-D profiles of colour against radius, as distinguished from studies of dust patterns. Most of these data were reconsidered by Michard (2000) (RM00), in an attempt to collect a significant sample of objects with a complete optical colour set, i.e. U$``$B, B$``$V, B$``$R and V$``$I in a coherent photometric system. This was adequate to confirm previous indications about the cause of colour gradients: these appear to be due essentially to population gradients within galaxies, with the dust playing no important role, except in galaxies with central intense dust patterns. Such objects are rather rare among the Es.

Similar to most spectral indices of stellar populations, the colours suffer from the well known age-metallicity degeneracy, and, except U$``$B or U$``$V, are not very sensitive to the two parameters. They are affected by dust, at least locally, or perhaps systematically in the central regions according to inferences based on a survey by Michard (1999) (RM99). On the other hand, they may be measured at lower surface brightnesses or larger radii than the line indices. They could therefore bring useful information to the study of fossil stellar populations, and further constraints upon models of the evolution of E galaxies. The present work aims to provide an enlarged sample of objects with complete colour data, extending farther in radius than in previous studies, and hopefully of improved accuracy.

In Paper I, we present the usual information about the observations and data reduction, and part of the results in tabular form. A larger set of results will be made available in electronic form. The frames, partly reduced, will be made available from the HYPERCAT database, Observatoire de Lyon.

In Paper II, under the assumption that the observed colour gradients reflect abundance variations along the radius, metallicity gradients will be computed from the present data, using new colour-metallicity calibrations derived from multi-population models for E-galaxies. These metallicity gradients allow an estimation of central and mean metallicities. Statistics of galaxies included in our sample indicate that mean metallicities are about solar, in agreement with the study by Trager et al. (2000) based on spectral indices.

Often used notations

* $`r`$ isophotal radius; $`r=(ab)^{1/2}`$ for an ellipse of semi-axis $`a`$ and $`b`$.
* $`\mathrm{\Delta }_{UB}`$ colour gradient in U$``$B; $`\mathrm{\Delta }_{UB}=\mathrm{d}(\mathrm{U}\mathrm{B})/\mathrm{d}(\mathrm{log}r)`$ and similar for other colours.
* diE, boE, unE: subclassification of E galaxies as disky, boxy and undetermined.

## 2 Observations

The observations were performed with the 120cm Newtonian telescope of the Observatoire de Haute-Provence, in three runs: April 1-11 2000, May 29-June 5 2000 and January 18-29 2001, noted below as run 1, 2 and 3. Tables 1 to 4 gives lists of the observed galaxies with some parameters relevant on the observing conditions. A CCD target Tek1024 is mounted in the camera, giving a field of view of 11.6x11.6 ′for a pixel size of 24 microns or 0.68 ″. The relatively large field is a favorable feature of this system for the observation of colour distributions in nearby galaxies. Less favourable are the rather poor seeing at the OHP, with the FWHM of star images usually in the 2-3″range, with values at 4 or more during periods of northern wind (mistral), and also the sky illumination by ever increasing urban lights.

The camera is unfortunately affected by the so-called ”red-halo” effect.

## 3 Data analysis

### 3.1 Outline of the operations

The analysis of the frames entails the following steps:

1. The usual corrections for offset, the ”flat-fielding”, and the interpolation of bad columns. As explained below we tried to improve the flat-fields by ad hoc ”superflats”.
2. The registration of the 5 frames in each passband to a common geometry, based on a set of measured coordinates for 6 to 12 stars. This is intended to simplify the derivation of colour maps if needed.
3. The preparation of each frame for measurement, involving a final attempt to measure and correct residual large scale background trends, corrections for parasitic objects, a treatment against cosmic rays peaks, and the calibration against the available results of aperture photometry (see Sect. 3.2.3).
4. Large errors in colour measurements may result from small differences between the widths of the PSFs of the two frames involved (see for instance Michard 1999, and the previous literature quoted therein). When the FWHMs of measured PSFs in the colour set for a given object differed by more than 10%, it was our practice to modify the frame PSFs and try to make them equal to that of the best frame of the group (see Sect. 3.2.4).
5. The isophotal analysis of the V frame was performed according to Carter (1978), as implemented in the Nice technique described in Michard and Marchal (1994) (MM94).
6. The correction procedures for the red halo effect in V$``$I, and eventually for the effects of different PSF far wings in other colours, were performed (see Sect. 3.2.5). The correction necessitates the crossed correlation of the V frame by the I PSF and conversely. This operation cancels out the errors in the colour distribution induced by the red halo, but degrades the resolution. A correction to the calibrations performed before the convolutions is needed.
7. Finally, colour measurements were performed along the previously found isophotal contours, and the average isophotal colours tabulated. Our routine at this stage involves corrections to the adopted values of the sky backgrounds, in order to eliminate the obvious effects of inaccuracies in these (see Sect. 3.2.6).
8. The V surface brightness and colours have been collected in ad hoc files, and corrected for galactic absorption (or reddening) and the K effect, according to the precepts and constants given in the Third Reference Catalogue of Bright Galaxies (RC3, de Vaucouleurs et al., 1991). The usual linear representation of colours against $`\mathrm{log}r`$ have been calculated in selected ranges, avoiding on the one hand the central regions affected by imperfect seeing corrections and/or by important dust patterns; and on the other the outer range presumably affected by poor background corrections and residual noise.

### 3.2 Detail of operations.

Some important details of the above summarized procedures will now be discussed.

#### 3.2.1 The background

The sky background of flat-fielded frames showed disappointingly large scale trends, specially in the U band. To improve upon this situation, it was tried to derive corrections by mapping the background of such frames, at least those which were not ”filled” by a large galaxy. Such maps were found to be correlated, although less so than expected, and their average used as a ”superflat”. The quality of the background was then generally improved: if not, the superflats were not used. Note that, in the observing runs of May 2000 and January 2001, a number of frames of ”empty” fields were obtained (sometimes in moonlight hours) to contribute to the derivation of superflats.

A final improvement was obtained by measuring, during the treatment of each frame, a number of background patches, and subtracting a linearly interpolated map of these, instead of a constant. The background residual large scale fluctuations were often measured at the three steps of the procedure, that is, after the application of the flats, of the superflats, and after the final treatment. Table 5 summarizes the results. It may be noted that the combination of flats and superflats left large errors in our first and second runs, specially in the U colour. The final background linear ”rectification” allowed quite significant improvements, as seen by comparing the upper and lower halves of the table.

If we consider an E galaxy observed under the typical conditions of the present series (see above for a tabulation of sky background values), the final residual fluctuations quoted here represent local errors of less than 0.1 mag near the isophote $`\mu _B=25`$. We will return later to the question of errors resulting from background uncertainties.

#### 3.2.2 Parasitic objects

In galaxy photometry it is necessary to remove parasitic objects, stars and galaxies, that overlap the measured object. In the present work we used concurrently the following techniques:

. replacement of pixels in a circle enclosing the ”parasite” by a circle symmetric about the center of the galaxy.

. replacement of the circle by another one chosen in a nearby area.

. marking of the pixels to be discarded in such a way that they are later left aside in the measuring programmes.

#### 3.2.3 Calibrations

The frames were calibrated by comparisons with the results of aperture photometry. Our first choice was to use the UBVRI data of Poulain (1988) (PP88), and Poulain and Nieto (1994) (PN94) that are available for 26 objects of our survey, and are in Cousins’s system, notably for R and I. In a few cases the data collected in the HYPERCAT catalogues were used. These contain both primary data (those with an independent photometric calibration), and secondary data (actually calibrated with part of the primary data). Only primary data were used, selected according to our previous experience or prejudices. A few completely missing calibrations in R and I were replaced by values calculated from the tight correlations of V$``$R and V$``$I with B$``$V, derived from Poulain’s data for E galaxies.

Although the I filter in the camera is of Gunn’s type, our photometry is transfered to Cousins’s system through the calibration. It is assumed that the difference of pass-bands has no significant effect on colour gradients.

#### 3.2.4 Equalization of PSF FWHM’s

Many studies of colour distribution in galaxies are affected by errors resulting from the difference in the PSFs of the two frames involved in a colour measurement. Franx et al. (1989), Peletier et al. (1990), Goudfrooij et al. (1994) calculated the radial range where the errors due to ”differential seeing” are larger than some accepted threshold, and discarded the corresponding colours. Our policy, for instance in RM99, was to correct for this effect by adjusting the two frames to have PSFs with a common FWHM.

In the present study we tried to equalize the PSFs of the 5 frames in a given colour set. This is feasible if the 5 frames are taken in rapid succession, so that the PSFs have similar widths. The frame with the narrowest PSF is selected, and we find by trial and error a narrow gaussian (or sum of two gaussians) which can restore another frame to the same quality, or rather the same FWHM, by deconvolution. The parameters at hand are the $`\sigma `$ of the gaussian and the number of iterations in the deconvolution. In Table 1 to 4 we give the original FWHM $`W_o`$ of each frame, and the improved $`W_i`$ after the procedure described here. Obviously this cannot lead to perfect results, and sometimes we find in our data the signature of ”differential seeing”, in the form of large colour variations, peaks or dips, within the seeing disk: these defects were edited out unless there was some good reason to suspect a genuine central colour anomaly, such as large dust patterns, or the jet of NGC4486.

The reader may notice that two observations of NGC4406 are listed in Table 2, one of May 30, the other of May 31 2000. The first one was taken through fog and with average seeing, while for the other the ”mistral” brought a clearer sky and very poor seeing. A special treatment was then applied: the central peak of the sharp images was ”grafted” on the corresponding regions of the unsharp but deeper images. This explains why the $`W_i`$ is so much narrower than the original $`W_o`$ for the frames of May 31.

#### 3.2.5  ”Red halo” and PSFs far wings.

Before the start of this survey, the CCD camera on the telescope used was known to be affected by the ”red halo”, an unfortunate property of thinned CCDs. The aureoles surrounding stellar images are obviously brighter and more extended in the I band than in B or V. Not only the red halo, but more generally the outermost wings of PSFs, were measured during our observing runs in 2000-1. The techniques and results are described in Michard (2001) (RM01). The choice of appropriate star fields allowed us to extend the measurements up to a radius of nearly 3 ′, and down to a level of about $`\mathrm{0.5.10}^6`$ of the central peak. Due to the red halo effect, PSF wings in I may be a factor of 3 brighter in an extended radius range than the V ones. Much smaller but still significant differences may also occur between the PSFs of various spectral bands, the V PSF wings always being fainter. The V PSF wings however, and all the others at the same time, were reinforced between our observing runs of spring 2000 and winter 2001, probably an effect of 10 months ageing of mirrors coatings. The final output of the measurements are average ”synthetic” PSFs in the format 512x512 pixels, or 5.8x5.8 ′, for each run and pass-band.

A set of numerical experiments on model galaxies is also presented in RM01, to illustrate the effects of these far wings on the observed surface brightness and colour distributions. The most striking effect occurs for the gradient $`\mathrm{\Delta }_{VI}`$, which appears strongly positive, while it is negative according to the classical results of Bender and Möllenhof (1987), or Goudfrooij et al. (1994). More subtle effects are found for other colours, with relatively small but definite changes in gradients.

To correct for the consequences of the red halo, or other similar effects upon the colour distribution in the index C1$``$C2, frame C1 is convolved with the PSF of frame C2 and conversely. After this operation, the resulting images have been submitted to the same set of convolutions, one in the atmosphere plus instrument, the other in the computer: they lead therefore to correct colour distributions, but with a significant loss of resolution. As the convolutions attenuate the central regions of the galaxy, and much more so for the V frame convolved by the I PSF, the mean colours are biased: a correction to the calibrations performed before the convolutions is needed. This has been done by a comparison of simulated aperture photometry to the observed one, an operation also used to estimate the errors in calibration (see below).

Since the extended PSFs are found with limited accuracy, it is necessary to discuss the validity of the corresponding corrections obtained through crossed convolutions, the more so because of the obvious changes of the PSF far wings between run 1 and run 3. The mean values of the colour gradients for subsamples of E galaxies have been used for these checks, with the results of Table 6. For a subsample of 12 or more E galaxies the mean colour gradients and their dispersions cannot differ much, so that their values may be used as checks of the need for a correction and its eventual success. 111This assertion is questionable as pointed out by the referee. We therefore tested it by comparing the distributions of the gradients in subsamples of 15 to 25 objects, sorted out by NGC numbers, from the surveys of Peletier et al. 1990, Goudfrooij et al. 1994 and our discussion in RM00. This procures 8 tests, according to the number of colours in each of the sources. It turns out that the subsamples are statistically coincident in 4 cases. They differ by slightly more than the calculated errors in the others. The reference for these comparisons is the subsample in Michard (2000) (RM00), mostly a rediscussion of the ”classical” data by Peletier et al. (1990), Goudfrooij et al. (1994) and others.

Looking at the Table 6, it is clear that the red halo introduces enormous errors in the V$``$I gradients, but that the corrections are remarkably successful in restoring the agreement of the results with the accepted reference, both as regards the mean values and the dispersion. The same may be said about the V$``$R gradients. The wings of the V PSF were strongly reinforced between our run 3 and run 1 or 2, but much less so for the I and R PSF wings. As a result the red halo effect is less in V$``$I for the frames of run 3 and disappears in V$``$R.

The situation is less clear for the U$``$B distributions. Our mean uncorrected gradients are in good agreement with the ”classical” data, essentially from Peletier et al. (1990), as rediscussed in RM00. On the other hand, the U PSF wings are consistently above the V ones in all our runs, so that the true slopes of the U$``$B variations may be a bit smaller than the observed ones. This error might well be present in the classical observations. Incidentally, the data of Peletier et al. were obtained in U$``$R, and it is impossible to be sure that the far PSF wings of the used telescope were the same in both pass-bands!

Similar remarks might be made about our B$``$V data. Since the mean measured gradient is the same for our data of the year 2000 and the adopted reference (and as a good B PSF is not available) we take as correct this set of results. For our data of 2001, there is evidence that the wings of the B PSF were slightly above those of the V one, so that the B$``$V gradients might also be biased upwards.

It appears that a significant source of error in the measurement of the small colour gradients in E galaxies has hitherto been overlooked. It might be that small systematic errors, of the order of some 15-20%, are still present in the U$``$B or U$``$V gradients published here. Although such errors would not have significant astrophysical implications, control observations are planned.

#### 3.2.6 Measurements of isophotal colours

Our procedure uses Carter’s isophotal representation of the V frame. The successive contours at 0.1 mag intervals are fitted to each of the two frames to be compared, the mean surface brigtnesses calculated and the corresponding magnitudes and colours derived. A sliding mean smoothing is applied to the data for the outer contours and a graph of the colour against $`\mathrm{log}r`$ displayed. This might hopefully be linear, or nearly so, in the studied range. If it is not, it is our practice to introduce corrections to the provisional values of the sky background for one or both of the frames, in order to get rid of the ” breaks” in the colour-radius relation typically associated with a poor choice of the background constant. The reader is referred to the graphs published in Goudfrooij et al. (1994) for examples of such features. In Fig.1 we show a color profile with a rather important defect due to poor backgrounds, and its adopted correction.

The introduction of such ”aesthetic” corrections to the raw data might be criticized, since it assumes a regular behaviour of the colours at large $`\mathrm{log}r`$. This is however a reasonable hypothesis: the introduced corrections remain small, as shown by the statistics of Table 7. It should be noted that the mean sky background values derived for our large field frames are more precise than in previous works based on small field frames, where the sky was not reached at all. The problem lies in the presence of residual large-scale background fluctuations (see above): their effects are similar to those resulting from the poor evaluation of a constant background, and can be approximately corrected by the introduction of an ad hoc constant, or rather a set of constants, for the 5 frames in the colour set.

The linear fit was finally performed on a range of $`\mathrm{log}r`$ selected so as to avoid the central regions affected by known dust patterns or possible residual seeing-induced errors, and the outermost regions visibly affected by deviations from the expected straight line.

#### 3.2.7 Errors

The noise is not a significant source of error in this type of work, because averages can be performed upon thousands of pixels in the galaxy regions of low S/N ratio. Residual noise effects at large $`r`$ can be easily recognized in sample plots of the data (see Fig. 1, 2 , 3 , 4). The main sources of errors lie:

* in the calibration, giving a probable error $`\sigma _C`$, affecting equally all parts of a given object, but dependent upon the quality of the available aperture photometry in each colour.
* in the ”differential seeing”, i.e. the small differences in PSF between the frames in a colour set, also after the adjustments described above. This error $`\sigma _0`$ occurs only near the center of the objects, say in a diameter of twice the seeing FWHM (it would be a much larger range without the performed adjustments).There is no obvious reason for it to be colour dependent.
* in the sky background residual large-scale fluctuations, or equivalently, errors in the sky background estimates. This error $`\sigma _S`$ is strongly dependent upon colour and the studied radius, or rather the corresponding surface brightness.

These three components of the total error will now be considered, together with other relevant topics.

1. Errors in colours from calibration inaccuracies.

When a number of measurements are available for a given object, for instance 5 apertures in PP88 and PN94, a probable error of the resulting calibration is readily derived from the dispersion of these measured values about the results of simulated aperture photometry, calculated from our data, i.e. V magnitudes, isophotal parameters and colours. For the preferred calibrations with Poulain’s data, the computed error is often less than 0.01 in B$``$V or V$``$R but may rise to 0.02 in U$``$B or V$``$I. Still larger calibration errors, up to 0.04, have been estimated for objects with very scanty or uncertain aperture photometry. These probable errors apply to the zero point of the colour regression as shown in Table 9, and to all colour data from the same object.
2. Residual errors from PSF equalization.

The residual errors after this step in our data treatment can be objectively ascertained by comparing the ”central” B$``$R colours in our survey with the equivalent data from RM99, derived from high-resolution CFHT frames. For the present survey the ”central colours” are the integrated colours within a radius of 3 ″. From RM99, Table 6, we find the colours at the isophote of 1.5 ″, which are likely similar. The statistics of the difference B$``$R(new)–B$``$R(RM99) are for 31 objects in common: mean=0.01; $`\sigma =0.038`$. Assuming then that the errors are equal in the two surveys, the probable error associated with poor PSF adjustements is $`\sigma _0=0.027`$ in B$``$R. This source of error has no reason to vary significantly from one colour to another. It is independent of the error of calibration previously discussed.

The comparison between the two surveys is made possible because the same set of calibrations has been used, although the field of the CFHT frames was often not sufficient to use all calibrations apertures. A minor part of the above differences may come from this source. In RM99, the B$``$R of NGC2768 is quoted too red by 0.08 and that of NGC3610 too red by 0.10.
3. Errors from sky background inaccuracies

Given $`ϵ_A`$ and $`ϵ_C`$, the relative errors in the sky background evaluation for frames A and C, A and C beeing a pair among UBVRI, the magnitude error in the colour A$``$C may be expressed as $`\delta _{AC}=2.5\mathrm{log}(1+ϵ_A10^{0.4\mathrm{\Delta }\mu _A})+2.5\mathrm{log}(1+ϵ_C10^{0.4\mathrm{\Delta }\mu _C})`$. Here $`\mathrm{\Delta }\mu _A`$ is the magnitude contrast between the object and the sky in colour A. Using average values for the colours of E-galaxies and for the sky brightnesses, $`\mathrm{\Delta }\mu _A`$ for any colour may be expressed in terms of $`\mathrm{\Delta }\mu _V`$. Then, in the range of small $`\delta _{AC}`$, the expression reduces to $`\delta _{AC}=1.0857(K_Aϵ_AK_Cϵ_C)10^{0.4\mathrm{\Delta }\mu _V}`$ where $`K_A=10^{0.4\mathrm{\Delta }\mu _A}/10^{0.4\mathrm{\Delta }\mu _V}`$.

$`K_V=1`$ by definition; we find $`K_U=3.1`$ and $`K_I=1.55`$, while $`K_B`$ and $`K_R`$ are slightly below 1.

The $`ϵ_A`$,… are unknown, but it is feasible to get statistics of the linear combinations $`K_Aϵ_AK_Cϵ_C`$. Indeed, as explained above, we have adopted ad hoc corrections to provisional sky background values, in order to regularise the colour-$`\mathrm{log}r`$ relations. It is reasonable to assume that the errors left after these corrections are proportional to the adopted corrections themselves. We take $`\delta _{AC}=\eta K_{AC}=\eta (K_Aϵ_AK_Cϵ_C)`$, where $`\eta `$ is a small constant and the $`K_{AC}`$ may be derived from the statistics of the adopted corrections given in Table 7. In practice somewhat different statistics have been calculated, to take into account the fact that our corrections for two colours are not necessarily uncorrelated. The constant $`\eta `$, different for our observing runs of 2000 and 2001, was chosen so as to get a system of errors compatible with the appearance of the data and also with the errors found for the slopes of the colour-$`\mathrm{log}r`$ relation. The finally adopted errors $`\sigma _S`$ from sky background inaccuracies are given in Table 8. Note that this source of error is negligible for $`\mu _V=22`$ or smaller. Predicted errors are reduced for our run 3 as compared to the two others.
4. Total errors

The above estimated errors are independent and should be added quadratically. In the central region with $`r<6`$ ″, the total error will be $`\sigma _T=(\sigma _C^2+\sigma _0^2)^{1/2}`$. In the mean region with $`18<\mu _V<22`$ the total error equals the calibration error. Finally, in the outer range of $`\mu _V>23`$ one can use $`\sigma _T=(\sigma _C^2+\sigma _S^2)^{1/2}`$.
5. Errors from the corrections for red halo and other far wings effects

To our knowledge, the crossed convolutions used to correct for the red halo and similar effects do not give rise to random errors, but rather to systematic errors due to inaccuracies in the adopted PSFs. Such problems may be detected from the study of the distributions of the slopes of the colour-$`\mathrm{log}r`$ relations discussed below.
6. Errors in the slopes of the colour-$`\mathrm{log}r`$ relations

These have been calculated by two complementary methods. On the one hand, we can look for the correlations between the slopes derived here and those from the literature, notably the data collected and discussed in RM00. Assuming then that the errors are the same in both sources, we get an estimate of our slope errors, hopefully an upper limit. On the other hand we can consider the internal correlations between the slopes of the various colour-$`\mathrm{log}r`$ relations, specifically $`\mathrm{\Delta }_{UB}`$ and others with $`\mathrm{\Delta }_{BV}`$. Figure 5 and 6 shows the correlations of the U$``$B and V$``$I colour gradients with that in B$``$V. The coefficients of correlation are respectively 0.73 and 0.40. A weighted mean $`Gm4`$ of the slopes in the 4 colours has also been used as reference instead of $`\mathrm{\Delta }_{BV}`$ with analogous results. From the dispersions of such correlations the slope errors can be estimated, if the error for the reference $`\mathrm{d}BV/d\mathrm{log}r`$ or $`Gm4`$ is ”guessed”. The two techniques give results in very good agreement, the internal correlations indicating somewhat smaller errors.

The probable errors of the slope estimates are then 0.03 in U$``$B or U$``$V, 0.01 in B$``$V, 0.015 in V$``$R or B$``$R, 0.02 in V$``$I.

## 4 Observational results

The available data from the present study are:

1. Table 9, giving for each object the linear representation of the colour-$`\mathrm{log}r`$ relations, i.e. the selected inner and outer radii of the fit, the zero point colour with its probable error $`\sigma _C`$ and the slope. The probable errors for the slopes are given above. This table allows easy calculation of the colours at any radius, notably the effective radius $`r_e`$ and others as considered below.
2. In electronic form only, the tables giving for each galaxy, as a function of the radius $`r`$, the V magnitude $`\mu _V`$ and the colours U$``$B, B$``$V, V$``$R, V$``$I. These tables are presently available with the U$``$B indices as observed, and consistent with Table 9. They will be later made available after correction for suspected effects of PSF far wings (see discussion in 3.2.5).
3. Series of graphs from the above tables, showing the colours as a function of $`\mathrm{log}r`$ or of $`\mu _V`$. Examples of these graphs are shown here to illustrate a number of properties of the radius-colour relations.
4. Table 10 gives the ”central” colours according to several definitions: we have considered the colours integrated inside the area of radius $`r=3`$ ″, and the colours calculated for $`r=1.5`$ ″from the linear representations of Table 9. They should be nearly equal if these representations remain valid at small $`r`$, which is not always the case: see below for a description of typical deviations.
5. Table 11 collects colours measured at the outermost range of the available data expressed in V magnitude. This gross limit varies between $`\mu _V=23`$ (for NGC4472) and 24.5. It is controlled by the size of the object and the ”cleanliness” of the nearby field.

The SA0 galaxies NGC3115, 3607, 4550 and 5866 have been observed with the E-type sample. The corresponding results are given in the Tables, but they have been discarded from the discussion.

### 4.1 Statistical comparison with previous work

* Assuming that colours result only from population variations, ”perfect” correlations between zero point colours in different pass-bands are expected. We have compared colour-colour correlations for previous surveys and the present one. For B$``$I against B$``$V we find from Goudfrooij et al. (1994) a coefficient of correlation $`\rho =0.61`$ (41 objects) compared to $`\rho =0.94`$ from Table 9 (37 objects). For U$``$R against B$``$R the data of Peletier et al. 1990 (38 objects) give $`\rho =0.58`$, while the present data leads to $`\rho =0.94`$. The improvement is less if we compare our results with the discussion in RM00 where the calibration of the available photometry was reconsidered.
* The correlation between colour gradients should also be as good as allowed by errors of measurements. Again we compare the correlation coefficients between gradients from the literature and our results for the same colours. From Goudfrooij et al. (1994), considering the gradient $`\mathrm{\Delta }_{BI}`$ against $`\mathrm{\Delta }_{BV}`$, we find $`\rho =0.78`$ which is quite good. The result derived from Table 9 is still better, i.e. $`\rho =0.85`$. The improvement is more pronounced when the $`\mathrm{\Delta }_{UR}`$ against $`\mathrm{\Delta }_{BR}`$ correlation is calculated from Peletier et al. (1990): one obtains $`\rho =0.36`$ only, instead of $`\rho =0.88`$ with our data. 222The coefficients of correlation given here are of course different from the previously calculated ones in Sect. 3.2.7.6 , which refered to other colours with smaller gradients.
* Our colour measurements could be obtained at much lower surface brightness (or larger radii) than in previous work. The present data extend to $`23.2<\mu _V<25.2`$, with a median value near $`\mu _V=24.5`$ in all colours. The tables by Goudfrooij et al. (1994), available from CDS Strasbourg, are mostly limited to $`\mu _V=22.75`$ in B$``$V and $`\mu _V=22.25`$ in V$``$I (median values), due to the small fields of the frames, notably in the I band (and probably an ad hoc cut off). Their published graphs generally extend to the same radius in both colours. The tables from Peletier et al. (1990), again at the CDS, extend to $`\mu _V=23.2`$ in B$``$R and $`\mu _V=22.6`$ in U$``$R (mean values) due to a systenatic cut-off at 10% of the sky. The printed tables may be still more severely truncated.

In summary, our colour data generally extends 1.5 to 2 magnitudes deeper than in previous works, so that ”external colours” refering to the level $`\mu _V=24.5`$ whenever possible, are presented in Table 11 with realistic error estimates.

### 4.2 Description of isophotal colour profiles.

Most of the profiles relating colour to $`\mathrm{log}r`$, or equivalently to the surface brightness $`\mu _V`$, are regular, meaning that they deviate very little from a straight line in the range of abscissae relevant to the present data, roughly $`r>2`$ ″and $`\mu _V<24.5`$. In this case, the ”central colours” integrated within $`r<3`$ ″, and the colours calculated from the linear representation at the average radius of $`r=1.5`$, differ very little. Central colour according to these two definitions are given in Table 10: compare columns (5) and (6) for U$``$V, and (3) with (7) for V$``$R. Figure 2 gives an example of a set of regular colour profiles for NGC4473. Of course the linear colour–$`\mathrm{log}r`$ relation is only approximate and breaks down at small $`r`$ for all galaxies observed at high resolution. Carollo et al. (1997) obtained V$`i`$ maps of a number of E-galaxies from HST frames, disclosing minute colour structures in some cases. In RM99, the CFHT resolution proved sufficient for the detection and classification of B$``$R central ” red peaks”. These are smoothed out, at least partly, with the OHP seeing.

Non regular profiles have been observed in the following cases:

1. When an important dust pattern occurs near the center of an object, it produces a central red hump in the colour profile. This is the case for galaxies with the value of 3 for the dust pattern importance index (DPII) introduced in RM99, such as NGC2768, 4125 and 4374, but also for 5813 and 5831. The dust ring of NGC3607, type SA0, produces noticeable bumps in its colour profiles. Fig. 3 gives an example of the colour profiles of such a centrally dusty galaxy, i.e. NGC4125. In such cases, the integrated colours within $`r=3`$ ″are redder than the extrapolated colours at $`r=1.5`$: this ”extra reddening” is smaller in U$``$B than in other colours. The consideration of the central reddening in V$``$R, a colour nearly insensitive to age-metallicity variations but sensitive to dust, gives a possibility to correct other colours for dust effects. This has been done, in some applications, for galaxies of DPII 3.
2. A few objects show a central red hump in some of the colour profiles, specially U$``$B and also eventually B$``$V. This is the case of NGC3377, 3379, 3610(?) and 4494. Fig.4 illustrates the case of NGC3377. For such objects the ”extra reddening” defined above is near zero in V$``$R and V$``$I. It is therefore permitted to attribute it to a metallicity increase rather than to dust.
3. On the contrary, the central colours U$``$B, B$``$V may be bluer than the extrapolation of the linear portion of the profile. This occurs for objects with larger than average colour gradients. The best such case is NGC4636; others are 4283, 4478 and the SA0 4550. An ”extra blueing” in U$``$V is then apparent from the comparison of indices in columns (5)-(6) of Table 10. This is also marginally the case for such giant galaxies as NGC4406, 4472, 4649: according to RM99 the colour profile of such objects tend to flatten out near the center.
4. NGC4486 is remarkable in showing a central ”blue deep”, the U$``$B colour beeing 0.51 at $`r=0`$ as compared to 0.72 near $`r=6.5`$ ″. This central feature is probably somehow related to the famous non-thermal jet of this galaxy. The jet is of course very conspicuous in U$``$B , with a peak colour near -0.1. Needless to say, both the central ”blue deep” and the jet are affected by seeing (and our attempts to correct its effects).

### 4.3 Correlations of interest.

#### 4.3.1 Correlations between colour gradients.

The correlations between the various colour gradients have already been noted as useful tools to evaluate the probable errors in gradients. These correlations are displayed in Fig. 5 and 6.

The coefficient of correlation between $`\mathrm{\Delta }_{UB}`$ and $`\mathrm{\Delta }_{BV}`$ is 0.75; that between $`\mathrm{\Delta }_{VI}`$ and $`\mathrm{\Delta }_{BV}`$ is only 0.40. For $`\mathrm{\Delta }_{VR}`$ it falls to 0.20, because the errors are of the same order of magnitude as the V$``$R gradients. Imposing regression lines running through the origin, the relative slopes are $`<\mathrm{\Delta }_{UB}/\mathrm{\Delta }_{BV}=2.2>`$, $`<\mathrm{\Delta }_{VI}/\mathrm{\Delta }_{BV}=0.8>`$ and $`<\mathrm{\Delta }_{VR}/\mathrm{\Delta }_{BV}=0.25>`$. These relative slopes are in good agreement with the results in RM00, and its conclusion, i.e. the negligible influence of dust upon colour gradient, is confirmed.

The distribution of colour gradients for E-galaxies may be of interest. The following parameters are found: $`<\mathrm{\Delta }_{UB}>=0.152`$ with $`\sigma =0.044`$; $`<\mathrm{\Delta }_{BV}>=0.070`$ with $`\sigma =0.021`$; $`<\mathrm{\Delta }_{VR}>=0.018`$ with $`\sigma =0.012`$; $`<\mathrm{\Delta }_{VI}>=0.054`$ with $`\sigma =0.023`$. The dispersions are not much larger than the errors estimated above. The distributions are asymmetric: there are 4 objects with $`\mathrm{\Delta }_{BV}`$ larger than $`+1.8\sigma `$ above the mean, but none at less than the same deviation. These galaxies, with $`\mathrm{\Delta }_{BV}`$ clearly steeper than average, by about twice the estimated probable error of measurement, are NGC4283, 4478, 4564 and 4636, seemingly a random collection.

Remark: An attempt to sort the E-galaxies by flattening as measured in MM94, and to look for some relation to the gradients, lead to negative results. Similarly no significant difference was found between diE and other galaxies.

The SA0 NGC4550 has quite exceptional gradients in all colours, and an admixture of dust and relatively young stars could be invoked to explain its properties. This object also has very remarkable kinematics, as first described by Rubin et al. (1992); a model has been proposed by Rix et al. (1992). The few other S0s in the present sample are similar to Es with regard to their colour gradients.

#### 4.3.2 Colour-colour diagrams

Many colour-colour diagrams can be built from the present data. The indices used may be calculated from the linear representations in Table 9 at the effective radius $`r_e`$, at a near center position $`r_e/8`$, and at an intermediate position $`r_e/2`$. The system of $`r_e`$ here used is an average of estimates in MM94, Prugniel and Héraudeau, 1998, and the RC3. It is satisfactory that the various indices $`(\mathrm{B}\mathrm{V})_{r_e}`$, $`(\mathrm{B}\mathrm{V})_{r_e/2}`$, $`(\mathrm{B}\mathrm{V})_{r_e/8}`$, define a common diagram with the corresponding U$``$B. This may suggest that a common physical variable controls the variations inside an object, and the object to object changes, of the two colours. Such graphs readily show larger than average calibrations errors: for instance, the U$``$B colour of NGC3610 is clearly too red for its B$``$V.

We have also traced colour-colour diagrams for the ”central” colours, (Table 10), i.e. integrated in the radius $`r<3`$ ″. They are similar to those traced with the interpolated colours, but with larger dispersions: this is not surprising since the central colours suffer from larger errors (see 3.2.7).

Finally, one can trace colour-colour diagrams with the ”outermost colours” collected in Table 11. They are in fair agreement with the diagrams derived from Table 9, and extend these towards the blue. We show in Fig. 7 a composite colour-colour diagram U$``$B against B$``$V, using the colours ar $`r_e/8`$, $`r_e/2`$ and the outermost range. Similarly Fig. 8 displays the diagram of V$``$I against U$``$V.

#### 4.3.3 The U$``$V colour as a metallicity index for ellipticals ?

Burstein et al. (1988) showed a correlation between the central $`Mg_2`$ index and a global B$``$V colour measured in a large aperture (see also Bender et al. , 1993). This type of correlation is reconsidered here using the U$``$V colour which is much more sensitive to metallicity than B$``$V, and taking advantage of recent estimates of the $`Mg_2`$ index far from center.

Two correlations between the $`Mg_2`$ index and U$``$V are considered in Fig. 9 and 10. The first shows the relation between the two quantities near the galaxy center: the $`Mg_{2C}`$ index is taken from the tabulation by Faber et al. 1989. The $`(\mathrm{U}\mathrm{V})_C`$ index is the integrated colour in a circle of radius $`r=3`$ ″. The value for three galaxies with important central dust patterns have been corrected by reference to the central bump in (V$``$R), a colour sensitive to dust but less so to metallicity changes. The coefficient of correlation reaches 0.825. Taking $`(\mathrm{U}\mathrm{V})_C`$ as $`x`$ and $`Mg_{2C}`$ as $`y`$ we find the regression $`y=0.221\pm .027x0.061\pm .003`$.

The second, i.e. Fig. 10, displays the correlation between the U$``$V colour and the $`Mg_2`$ index at the effective radius $`r_e`$. $`Mg_{2re}`$ has been taken from Kobayashi and Arimoto (1999) (KA99). To increase the number of data points the $`Mg_b`$ gradients were introduced, using the linear relation $`Mg_b=15Mg_2`$ derived from the distribution of the values of both indices in the tables of KA99. The coefficient of correlation still reaches 0.72. Again, with the colour in $`x`$ and the $`Mg_{2re}`$ in $`y`$ we find $`y=0.283\pm .069x0.133\pm .007`$. The difference to the above correlation for central indices is barely significant.

The quality of these correlations proves that both indices are essentially controlled by the same physical variables, and leaves little room for the effects of diffuse dust upon the colours of E-galaxies.

## 5 Conclusions

New colour-radius relations have been derived for 36 E-type galaxies of the northern Local Supercluster, using UBVRI frames obtained with the 120cm telescope of the Observatoire de Haute-Provence. Four SA0, i.e. NGC3115, 3607, 4551 and 5866 were also observed.

We aimed to take advantage of the large field of the camera to observe the galaxies at larger radii than hitherto feasible, and thus improve the accuracy of colour gradients. The availability of the series of aperture photometry in PP88 and PN94 for most of the sample was also considered an asset towards a more coherent system of colours. It appears indeed that the colour calibrations are improved here compared to previous work, if this can be judged from the quality of correlations between zero point colours in various surveys (see Sect. 4.1)

Two steps in the reduction procedure were thought significant in improving the quality of colour profiles: the first was the adjustment of the FWHM of the PSFs in a given colour set of 5 frames to the best of the five. This allowed us to get significant colours much closer to the galaxy center than otherwise feasible. The second was a careful ”mapping” of the background of each frame, in order to lessen the background fluctuations remaining after the usual flat-fielding procedures. Both these precautions proved successful, and, as a result, the radial range of satisfactory colour measurements was greatly enlarged. Near the galaxy center, it proved feasible to obtain ”central colours”, i.e. colours integrated in the circle $`r=3\mathrm{}`$, in fair agreement with high resolution data (see Sect. 3.2.7.2 and Table 10).

On the other hand, colours could be obtained at much lower surface brightness (or larger radii) than in previous work. Our colour data extend to $`23.2<\mu _V<25.2`$, with a median value near $`\mu _V=24.5`$ in all colours. According to the comparisons in Sect. 4.1, this is 1.5 to 2 magnitudes deeper than in previous work. ”External colours”, refering to the level $`\mu _V=24.5`$ whenever possible, are published for the first time (see Table 11), and may be useful to give some indications about stellar populations at the outskirts of E-galaxies.

On the other hand, the ”red halo” effect of the camera was found to give enormous errors in V$``$I colours and gradients. These were corrected by a rigorous technique, and results in agreement with ”classical” data were obtained. Considering the V$``$I gradients, one is not happy however to introduce in their evaluation, corrections larger than the quantity to be measured! Besides this specific problem of the red halo of thinned CCD, the far wings of the PSFs have been proven in a recent paper (RM01) to have non negligible effects in the gradients of other colours, and also to vary with the age of mirror coatings. It is not impossible that the U$``$B or U$``$V gradients given here are overestimated by 15-20%, although they are in excellent statistical agreement with the well-known work of Peletier et al. (1990).

Various colour gradients against $`\mathrm{log}r`$ for a given object are well correlated, generally better than in previous work (see statistics in Sect. 4.1), which is interpreted as due to smaller measuring errors, notably in U$``$B. These improvements in accuracy did not bring out any obvious correlation between gradients and other galaxy properties. A few galaxies have exceptionally steep colour gradients (nearly at $`2\sigma `$) without sharing other properties.

Colour-colour relations can be built from the present data for several locations in galaxies, such as near center, various fractions of the effective radius $`r_e`$, or the ”outermost” measured range around $`\mu _V=24`$. All these diagrams overlap to form a single stripe with moderate scatter (except for one rather obvious calibration error?). These might prove useful to test theories of old stellar populations and of their host galaxies. Colour-colour diagrams based upon integrated colours have already been used for this purpose (Worthey, 1994).

The U$``$B or U$``$V colours correlate very well with the $`Mg_2`$ index, both near the galaxy center and at the effective radius $`r_e`$. This seems to rule out any large influence of diffuse dust in the colours and colour gradients in E-galaxies. This was considered likely by Witt et al. (1992) and discussed by Wise and Silva (1996) with inconclusive results. Previous arguments against such an influence were presented in RM00: they were based upon the relative average values of the gradients in various colours, and are reinforced in the present work, since the mean gradients are nearly unchanged, and their errors lessened.

It is well known that, for single-burst stellar populations, colours and line indices depend both on the metallicity and on the age of the system (Worthey 1994; Borges et al. 1995). However, E-galaxies are constituted by a population mix, having age and metallicity distributions which reflect their star formation histories. Therefore a colour-metallicity calibration requires the use of models able to provide those distributions and, consequently, the integrated colours along the galaxy lifetime. Such a calibration will be presented in Paper II.

###### Acknowledgements.

TPI acknowledges a Fapesp pos-doc fellowship n 97/13083-7.


format: text

text:
# Collapse arrest and soliton stabilization in nonlocal nonlinear media

## Abstract

We investigate the properties of localized waves in systems governed by nonlocal nonlinear Schrödinger type equations. We prove rigorously by bounding the Hamiltonian that nonlocality of the nonlinearity prevents collapse in, e.g., Bose-Einstein condensates and optical Kerr media in all physical dimensions. The nonlocal nonlinear response must be symmetric, but can be of completely arbitrary shape. We use variational techniques to find the soliton solutions and illustrate the stabilizing effect of nonlocality.

Collapse is a fundamental physical phenomenon well-known in the theory of waves in nonlinear media. It refers to the situation when strong contraction or self-focusing of a wave leads to a catastrophic increase or blow-up of its amplitude after a finite time or propagation distance (see  for reviews). Wave-collapse has been observed in plasma waves  (the famous Langmuir wave-collapse), electromagnetic waves or laser beams  (also called self-focusing), Bose-Einstein Condensates (BEC’s) or matter waves , and even capillary-gravity waves on deep water . The effect of collapse appears also in astrophysics, where the gravitational attraction plays the same role as the self-focusing of electromagnetic waves, tending to compress stars of sufficient mass, eventually leading to their collapse into a black hole .

Typically the contraction must be able to act freely in two or more dimensions to be strong enough to generate a collapse. Moreover, the so-called norm, which is the power for electromagnetic and plasma waves, the atom density for BEC’s, and the mass for stars, must be above a certain critical value for a collapse to occur. Most commonly the collapse has been discussed in the context of the nonlinear Schrödinger (NLS) equation, which is a universal model for dispersive (or diffractive) weakly nonlinear physical systems . The NLS equation models, e.g., all systems mentioned above, in which a wave-collapse has been predicted and verified experimentally.

The collapse singularity is an artifact of the model and signals the limit of its validity. Close to the singularity, when the amplitude is extremely high and the temporal- and spatial scales are extremely short, new physical processes will come into play . A common effect is nonlinear dissipation, such as two-photon absorption of electromagnetic waves and inelastic two- and three-body recombination for matter waves, which efficiently absorbs the collapsing part of the wave. Thus collapse acts as an effective nonlinear loss mechanism, as is well-known in, e.g., Langmuir turbulence  and BEC’s . Effects, such as discreteness, saturation of the nonlinearity, and nonparaxiality, will also eliminate the possibility of a collapse (see ). In contrast, effects such as weak linear loss , temperature fluctuations , and spatial incoherence , cannot eliminate the collapse. In any case the collapse effect represents a strong mechanism for energy localization, which is important to study to understand the properties of a given physical system.

The inherent nonlocal character of the nonlinearity has attracted considerable interest as a means of eliminating collapse and stabilizing multidimensional solitary waves. Nonlocality appears naturally in optical systems with a thermal  or diffusive  type of nonlinearity. Nonlocality is also known to influence the propagation of electromagnetic waves in plasmas  and plays an important role in the theory of BEC’s, where it accounts for the finite-range many-body interaction .

In this work we consider NLS equations with a general nonlocal form of the nonlinearity. Turitsyn proved the absence of collapse for three particular shapes of the nonlocal nonlinear response . The analysis of the collapse conditions for general response functions is difficult and has been carried out only numerically . However, in many systems, such as BEC’s, one has no knowledge of the particular response function, and thus it is important to maintain its generality in the model. Here we prove rigorously that nonlocality eliminates collapse in all physical dimensions for arbitrary shapes of the nonlocal response, as long as it is symmetric.

We consider the evolution of a wave field $`u=u(\stackrel{}{r})=u(\stackrel{}{r},\tau )`$ described by the general nonlocal NLS equation

$$i\frac{u}{\tau }+^2uV(\stackrel{}{r})u+N(I)u=0,$$
(1)
where $`V`$=$`V(\stackrel{}{r})`$ is an external (linear) confining potential, $`I`$=$`I(\stackrel{}{r})`$=$`I(\stackrel{}{r},\tau )`$=$`|u|^2`$, $`\tau `$ is the evolution coordinate and $`\stackrel{}{r}`$=$`(r_1,r_2,r_3)`$ spans a D-dimensional ”transverse” coordinate space. The nonlinear term $`N`$=$`N(I)`$ is represented in the general nonlocal form

$$N(I)=R(\stackrel{}{r}^{}\stackrel{}{r})I(\stackrel{}{r}^{})𝑑\stackrel{}{r}^{},$$
(2)
where we assume the response function $`R(\stackrel{}{r})`$ to be real, localized, nonsingular, and symmetric. For localized waves Eq. (1) conserve the power (in optics) or number of atoms (for BEC) $`P`$ and the Hamiltonian $`H`$,

$$P=I𝑑\stackrel{}{r},H=[|u|^2+VI\frac{NI}{2}]𝑑\stackrel{}{r}.$$
(3)
In optics $`u`$ is the envelope of the electric field with intensity $`I`$ and $`V`$ represents a guiding structure (waveguide). Here Eq. (2) represents a general phenomenological model for self-focusing Kerr-like media, in which the change in the refractive index induced by an optical beam involves a transport process. This includes heat conduction in materials with a thermal nonlinearity  or diffusion of molecules or atoms in atomic vapours . A nonlocal response in the form (2) appears naturally due to many body interaction processes in the description of BEC’s , if the assumption of a zero-range interaction-potential is relaxed . For BEC’s with a negative scattering length Eq. (1) is the nonlocal Gross-Pitaevskii (G-P) equation for the collective wave function $`u`$ ($`\tau `$ is time), with $`I`$ representing the density of atoms and $`V`$ representing the magnetic trap.

In the limit when the response function is a delta-function, $`R(\stackrel{}{r})`$=$`\delta (|\stackrel{}{r}|)`$, the nonlinear response is local (see Fig. 1a) and simply given by

$$N(I)=I,$$
(4)
as in local optical Kerr media described by the standard NLS equation and in BEC’s described by the standard G-P equation. In this local limit multidimensional optical beams with a power higher than a certain critical value ($`P>11.69`$) will experience unbounded self-focusing and collapse after a finite propagation distance . It is also well-known that BEC’s would collapse when the total number of atoms is larger than a critical number .

With increasing width of the response function $`R(\stackrel{}{r})`$ the wave intensity in the vicinity of the point $`\stackrel{}{r}`$ also contributes to the nonlinear response at that point. In case of weak nonlocality, when $`R(\stackrel{}{r})`$ is much narrower than the width of the beam (see Fig. 1b), one can expand $`I(\stackrel{}{r}^{})`$ around $`\stackrel{}{r}^{}`$=$`\stackrel{}{r}`$ and obtain the simplified model

$$N(I)=I+\gamma ^2I,\gamma =\frac{1}{2}r^2R(r)𝑑\stackrel{}{r},$$
(5)
where the small positive definite parameter $`\gamma `$ measures the relative width of the nonlocal response. The diffusion type model (5) of the nonlocal nonlinearity is a model in its own right in plasma physics, where $`\gamma `$ can take any sign . It was also applied to BEC’s , nonlinear optics , and energy transfer in biomolecules . In weakly nonlocal media with $`N(I)=I+\gamma ^2I`$ it is straightforward to show that collapse cannot occur. This was first done for plasmas , and later for BEC’s .

In the limit of a strongly nonlocal response much broader than the characteristic width of the wave function (see Fig. 1d), one can instead expand the response function and obtain (to lowest order)

$$N(I)P(R_0+R_2r^2).$$
(6)
where $`R_0=R(\stackrel{}{0})`$ and $`R_2=\frac{1}{2}^2R(\stackrel{}{0})`$. The evolution of optical beams in such a strongly nonlocal medium was considered in . Since this relation is linear, the highly nonlinear effect of collapse cannot occur.

So in the two extreme limits of a weakly and highly nonlocal nonlinear response the collapse is prevented. For arbitrary degree of nonlocality it is difficult to prove anything rigorously. Just saying that the dynamics is described by either the weakly nonlocal model (5) or the linear oscillator model (6), which both have no collapse, is not enough. As it is well-known from studies of general NLS equations, the typical singularity is a so-called blow-up featuring the amplitude locally going to infinity on a broad background localized structure . Such a two-scale field distribution, which was also recently observed in BEC’s , is clearly described by neither of the two simple limiting systems.

The stabilizing effect of nonlocality of an arbitrary degree was proven by Turitsyn for three specific examples, including Coulomb interaction ($`R(\stackrel{}{r})`$=$`1/|\stackrel{}{r}|`$) . He bounded the Hamiltonian from below for fixed power, which proves that a collapse cannot occur and that the soliton solutions are stable in a Lyapunov sence. Here we consider the general case of arbitrarily shaped, nonsingular response functions and prove rigorously that the Hamiltonian is bounded from below in all dimensions.

Introducing the D-dimensional Fourier transform (denoted with a tilde) and its inverse

$$\stackrel{~}{I}(\stackrel{}{k})=I(\stackrel{}{r})\mathrm{e}^{i\stackrel{}{k}\stackrel{}{r}}𝑑\stackrel{}{r},I(\stackrel{}{r})=\frac{1}{(2\pi )^D}\stackrel{~}{I}(\stackrel{}{k})\mathrm{e}^{i\stackrel{}{k}\stackrel{}{r}}𝑑\stackrel{}{k},$$
(7)
it is straightforward to show that for $`N(I)`$ given by Eq. (2) the following relations hold

$`|I(\stackrel{}{k})|=\left|{\displaystyle I(\stackrel{}{r})\mathrm{e}^{i\stackrel{}{k}\stackrel{}{r}}𝑑\stackrel{}{r}}\right|{\displaystyle I𝑑\stackrel{}{r}}=P,`$ (8)
$`{\displaystyle NI𝑑\stackrel{}{r}}={\displaystyle \frac{1}{(2\pi )^D}}{\displaystyle \stackrel{~}{R}(\stackrel{}{k})|\stackrel{~}{I}(\stackrel{}{k})|^2𝑑\stackrel{}{k}}.`$ (9)
For response functions with a positive definite spectrum, $`\stackrel{~}{R}(\stackrel{}{k})0`$, the Hamiltonian is thus bounded as follows

$$Hu_2^2\frac{1}{2}R_0P^2,$$
(10)
where $`u_p^p|u|^p𝑑\stackrel{}{r}`$ and we have used the fact that the minimum of $`V(\stackrel{}{r})`$ is zero.

The inequality (10) is the main result of this letter. It shows that for all symmetric response functions with a positive definite Fourier spectrum and a finite value at the center, the Hamiltonian is bounded from below by the conserved quantity -$`\frac{1}{2}R_0P^2`$, or conversely, that the gradient norm $`u_2^2`$ is bounded from above by the conserved quantity $`H+\frac{1}{2}R_0P^2`$. This represents a rigorous proof that a collapse with the wave-amplitude locally going to infinity cannot occur in BEC’s, plasma, or optical Kerr media with a nonlocal nonlinear response, for any physically reasonable response function.

The stabilising effect of the nonlocality can be further illustrated by the properties of the stationary solutions of Eqs. (1). As a simple example we consider nonlocal optical bulk Kerr media with a Gaussian response

$$R(\stackrel{}{r}^{}\stackrel{}{r})=\left(\frac{1}{\pi \sigma ^2}\right)^{\frac{D}{2}}\mathrm{exp}\left(\frac{|\stackrel{}{r}^{}\stackrel{}{r}|^2}{\sigma ^2}\right).$$
(11)
The ground-state stationary solutions are then radially symmetric bell-shaped, nodeless solutions of the form $`u(\stackrel{}{r},z)=\varphi (r)\mathrm{exp}(i\lambda z)`$, where the profile $`\varphi (r)`$ is found from the Euler-Lagrange equations for the Lagrangian

$$L=\left[\lambda \varphi ^2+|\varphi |^2\frac{1}{2}N(\varphi ^2)\varphi ^2\right]𝑑\stackrel{}{r}.$$
(12)
To capture the main physical effects we use the approximate variational technique with a Gaussian trial profile $`\varphi (r)=\alpha \mathrm{exp}[(r/\beta )^2]`$, in view of the fact that the Gaussian profile is an exact solution in the strongly nonlocal limit with $`N(I)`$ given by the parabolic potential (6).

Inserting this ansatz into the Lagrangian (12), with $`N`$ given by the general expression (2), the Euler-Lagrange equations give the amplitude $`\alpha `$ and width $`\beta `$

$`\alpha ^2`$ $`=`$ $`(\lambda +D/\beta ^2)(2+2\sigma ^2/\beta ^2)^{D/2},`$ (13)
$`\beta ^2`$ $`=`$ $`[4D+\sqrt{(4D)^2+16\lambda \sigma ^2}]/(2\lambda ),`$ (14)
In Fig. 2 we show the power $`P_s`$=$`(\pi /2)^{D/2}\alpha ^2\beta ^D`$ and Hamiltonian of the stationary solutions in 2D. The dashed lines give the results of the weakly nonlocal approximation with $`N`$ given by Eq. (5), from which $`\alpha ^2=4\lambda `$ and $`\beta ^2=2/\lambda +2\sigma ^2`$ is found, resulting in the power

$$P_s=4\pi (1+\sigma ^2\lambda ),$$
(15)
where $`4\pi `$ is the ($`\lambda `$-independent) power of the Gaussian approximation to the soliton solution of the standard 2D NLS equation, recovered in the local limit $`\sigma `$=0.

In the 2D NLS equation the collapse is critical and the stationary solutions are ”marginally stable” with $`dP_s/d\lambda `$=0 . Typically perturbations act against the self-focusing, with several effects, such as non-paraxiality and saturability, completely eliminating collapse . This is also the case with nonlocality, as evidenced from Fig. 2 and the simplified expression (15), which shows that any finite width of the response function (non-zero value of $`\sigma `$) implies that $`dP_s/d\lambda `$ becomes positive definite. According to the (necessary) Vakhitov-Kolokolov (VK) criterion  the soliton solutions therefore (possibly) become linearly stable.

For small $`\lambda `$ the soliton width $`\beta `$ decreases as $`1/\lambda `$. Thus the assumption of weak nonlocality, i.e. that the soliton is much wider than the response function, applies only to sufficiently small values of $`\lambda `$ satisfying $`\lambda \sigma ^21`$, which is also clearly seen from Fig. 2. The accuracy of the assumption of weak nonlocality is further discussed in Ref.  in terms of modulational instability.

The 3D case shown in Fig. 3 is more interesting, because the nonlinearity is much stronger than in 2D. The collapse in the local 3D NLS equation is so-called supercritical . Again the soliton width $`\beta `$ decreases as $`1/\lambda `$, so a threshold width should exist, below which the nonlocality is not strong enough to stabilize the soliton. This is exactly what is observed in Fig. 3: For $`\lambda <\lambda ^{\mathrm{th}}`$ the solitons are still linearly unstable with $`dP_s/d\lambda <0`$, but above threshold the nonlocality is strong enough to bend the curve and make $`dP_s/d\lambda >0`$, i.e., the solitons become linearly stable according to the VK criterion. From the definition $`dP_s(\lambda ^{\mathrm{th}})/d\lambda `$=0 the variational results give $`\lambda ^{\mathrm{th}}=1/(2\sigma ^2)`$, corresponding to a threshold in the soliton power (dashed curve in Fig. 3) and width

$$P_s^{\mathrm{th}}=(5\pi )^{3/2}5\sigma /4,\beta ^{\mathrm{th}}=2\sigma ,$$
(16)
which are both proportional to the degree of nonlocality $`\sigma `$. Thus, sufficiently broad and high-power solitons are stable. In the Hamiltonian versus power diagram in Fig. 3 the lower (upper) branches correspond to stable (unstable) solutions while the threshold is represented by the cusp . This stable solution branch was recently found numerically in the context of BEC’s with a nonlocal negative scattering potential . It corresponds to high density, self bound states of the condensate.

In conclusion we studied the properties of localized wave packets in nonlocal NLS equations. We proved rigorously that nonlocality of an arbitrary shape eliminates collapse in all physical dimensions. We also demonstrated that multidimensional soliton solutions of the NLS equation may be stabilized by the nonlocality. This opens a new and interesting discussion as to what is actually observed in collapse experiments in nonlocal systems. It seems clear that it all comes down to oscillations between opposite extreme states and how strong and rapid they are. Such oscillations were recently found to occur in BEC’s through numerical and variational studies .

We acknowledge useful discussions with Yu.S. Kivshar, S.K. Turitsyn, and Yu.B. Gaididei. This work was supported by the Danish Technical Research Council (Talent Grant No. 26-00-0355) and the Danish Natural Sciences Foundation (Grant No. 9903273).


text:
# Colloid aggregation induced by oppositely charged polyions

## I Introduction

In metallurgy it has been known for several millenia that mechanical properties of many metals can be considerably improved by alloying two or several elements. Similar improvements of material properties of ”soft matter” can be achieved by mixing two or more macromolecular components in fluid solutions or dispersions. One of the oldest known examples is the stabilization of colloidal suspensions by the addition of adsorbing polymer, which was already used in ancient Egypt for the preparation of inks. On the contrary, flocculation of colloids can be driven by the addition of non-adsorbing polymer via the depletion mechanism . More recent is the realization of the key role of Coulombic interactions in determining and improving material or physiological properties of complex colloidal and biomolecular systems, and their exploitation in numerous technological and biomedical applications. Widely studied systems include suspensions of charged colloids and polyelectrolytes , polymer-clay nanocomposites  or DNA bundles and DNA-lipid complexes . For example, even small polyelectrolyte additives can have a substantial impact on the aggregation and kinetic stability of charged colloidal particles . On the other hand, just a small fraction of clay fillers doubles the tensile strength of a polymer matrix , while highly asymmetric mixtures of charged clay platelets are used to improve rheological properties of drilling fluids .

To design complex mixtures with optimal properties, it is essential to gain a quantitative understanding of correlations and effective interactions between various macromolecular species. This can be achieved using the well-established techniques of liquid-state theory . The PRISM (”polymer reference interaction site model”) theory, originally designed for the study of polymer solutions and melts , has been recently extended to investigate various charged and uncharged colloidal systems , as well as mixtures of neutral spherical colloids and polymers  and like-charged spherical colloids and polyelectrolytes .

This paper focuses on the structure and stability of dispersions of charged colloidal particles in the presence of oppositely charged polyions, which play the role of mesoscopic counterions. There is a considerable literature on the widely studied case of microscopic counterions and the resulting effective interactions between the highly charged colloidal particles . In particular, it is now well established that divalent counterions can lead to overscreening of the Coulombic interaction between spherical or rod-like colloidal particles, which results in a correlation-induced, short-range attraction between the latter . The case where the counterions are mesoscopic or macromolecular polyions has received much less attention so far, except for the simplest case of polyelectolyte chains confined between two infinite parallel plates carrying charges of opposite sign to the polyelectrolyte, where overcompensation of the charge on the platelets, leading to effective attraction between the latter, has been predicted .

Two systems involving mesoscopic counterions are examined in the present paper, namely a) mixtures of highly charged spherical colloids and oppositely charged polyelectrolyte chains and b) binary mixtures of oppositely charged platelets of widely different diameters. The key problem to be addressed is that of structural diagnostics of a possible aggregation of the larger particles induced by the smaller particles, over a range of physical conditions. This goal is achieved by determining pair correlation functions and partial structure factors within the framework of PRISM, which is adapted in this paper to the systems in hand.

## II Interaction site models

The two systems under investigation are aqueous dispersions or solutions, but in view of the mesoscopic scale of the particles, the solvent will be modelled as a structureless dielectric continuum providing a macroscopic permittivity $`\epsilon `$. Any microscopic counterions or ions from added electrolyte will be considered at the linear response (or Debye-Hückel) level, i.e., they will screen the electrostatic potential due to the interaction sites on the colloids or polyions on a scale given by the usual Debye screening length $`\lambda _D`$ . The underlying Born-Oppenheimer-like assumption entails that the charge distribution on the mesoscopic particles does not contribute to screening . Since sites with charge of opposite sign are involved, the site-site interaction must contain a short-range repulsion, which will be modelled by hard sphere (HS) cores of diameter $`d_i`$. The interaction potential between sites on particles of species $`i`$ and $`j`$, carrying the charges $`z_ie`$ and $`z_je`$, will hence be of the generic form:

$$u_{ij}(r)=\{\begin{array}{cc}\mathrm{},\hfill & r\frac{1}{2}(d_i+d_j),\hfill \\ \frac{z_iz_je^2}{\epsilon r}\mathrm{exp}(\kappa _Dr),\hfill & r>\frac{1}{2}(d_i+d_j),\hfill \end{array}$$
(1)
where $`\kappa _D=\lambda _D^1`$. On the polyelectrolyte, $`N_p`$ charged sites are distributed at regular intervals along the linear macromolecular backbone, each site being associated with a monomer or segment of diameter $`d_p`$. The clay platelets will be modelled as circular discs of radius $`R`$, carrying interaction sites uniformly distributed on the surface; the core diameter associated with these sites accounts for the finite thickness of the platelets. As regards the spherical colloids, two models may be considered. In model A, $`N_s`$ charged sites are distributed over the spherical surface of radius $`R`$. For a sufficiently large density $`N_s/(4\pi R^2)`$ of the surface sites, the hard cores of diameter $`d_s`$ associated with each site will ensure that the spheres are impenetrable. In model B, the impenetrability of the spherical colloids is ensured by an uncharged central ($`c`$) site, with a core radius identical to the particle radius, while the electrostatic interactions are carried by $`N_s`$ surface sites distributed on a spherical shell of radius $`Rd_s/2`$. Each surface site carries a charge equal to the total charge on the colloid divided by the number of sites. Invoking Gauss’ theorem, it might be argued that this distribution should be equivalent to a single site carrying the full colloidal charge at the centre of the sphere. This equivalence would hold within an exact theory, but since an approximate closure relation is to be used, as discussed in section III, more accurate results may be expected when interaction sites with lower charges are used, such that the Coulomb coupling is weaker.

Within the PRISM formalism to be discussed in the next section, all interaction sites of the same type (i.e., $`p`$, $`c`$, or $`s`$) are assumed to be equivalent, leading to a considerable simplification over the original RISM . This assumption is obviously exact for spherical particles with surface sites only (model A), but is an approximation for polyelectrolytes of finite length (since end effects are neglected), and even more so for platelets, where edge effects are not taken into consideration. However previous experience with neutral or charged rods  and platelets  shows that the basic PRISM assumption does lead to sensible results for the pair structure of these rigid particles.

## III Multicomponent PRISM

Consider a multicomponent system involving $`\nu `$ species of particles with number densities $`\rho _\alpha `$ ($`1\alpha \nu `$), each particle of a species involving $`n_\alpha `$ classes of exactly or approximately equivalent interaction sites. The total number of distinct classes of sites in the mixture is

$$\underset{\alpha =1}{\overset{\nu }{}}n_\alpha =N.$$
(2)
It proves convenient to assign an index $`i`$ ($`1iN`$) to each class of sites, and to order them such that classes $`1in_1`$ belong to particles of species $`1`$, classes $`n_1+1in_1+n_2`$ belong to species $`2`$, … , and classes $`Nn_\nu iN`$ belong to species $`\nu `$. Each class $`i`$ contains $`N_i`$ equivalent interaction sites, and if the positions of sites belonging to classes $`i`$ and $`j`$ are designated by $`𝐫_k𝐫_k^{(i)}`$, and $`𝐫_l𝐫_l^{(j)}`$ ($`1kN_i;\mathrm{\hspace{0.17em}\hspace{0.17em}1}lN_j`$), then one may define the usual intramolecular structure factors, or form factors

$$\omega _{ij}(𝐪)=\frac{1}{\sqrt{N_iN_j}}\underset{k=1}{\overset{N_i}{}}\underset{l=1}{\overset{N_j}{}}e^{i𝐪(𝐫_k𝐫_l)}$$
(3)
provided $`i`$ and $`j`$ sites belong to the same particle, while $`\omega _{ij}(𝐪)0`$ otherwise. The $`\omega _{ij}(𝐪)`$ are hence elements of a $`N`$ x $`N`$ box-diagonal matrix $`𝝎`$$`(𝐪)`$. In the present paper rigid platelets and spherical colloids are considered, so that

$$\omega _{ij}(q)=\delta _{ij}+\frac{1}{\sqrt{N_iN_j}}\underset{k=1}{\overset{N_i}{}}\underset{l=1}{\overset{N_j^{}}{}}\frac{\mathrm{sin}(qr_{kl})}{qr_{kl}},$$
(4)
where $`r_{kl}=|𝐫_k𝐫_l|`$, $`q=|𝐪|`$, and the prime signifies that the ”self” term $`k=l`$ is to be left out for the diagonal elements $`i=j`$. PRISM is based on the assumption that all direct correlation functions between sites on pairs of different particles are identical if these sites belong to the same class on each of the two particles. This leaves a total of $`N(N+1)/2`$ independent direct correlation functions $`c_{ij}(q)`$. Corresponding total correlation functions $`h_{ij}(q)`$ are defined by averaging over the $`N_iN_j`$ correlation functions between all pairs of equivalent sites belonging to classes $`i`$ and $`j`$ on two particles of the same or different species. The $`N`$ x $`N`$ matrices of $`c_{ij}(q)`$ or $`h_{ij}(q)`$ are symmetric, but it proves convenient to introduce auxiliary non-symmetric matrices $`𝐇(q)`$, $`𝐂(q)`$ and $`𝐖(q)`$, with elements

$`H_{ij}(q)`$ $`=`$ $`\rho _iN_iN_jh_{ij}(q),C_{ij}(q)=\rho _ic_{ij}(q),`$ (5)
$`W_{ij}(q)`$ $`=`$ $`\sqrt{N_iN_j}\omega _{ij}(q),`$ (6)
where $`\rho _i`$ is the number density of that species of particles to which sites of category $`i`$ belong. With these conventions, the set of $`N(N+1)/2`$ Ornstein-Zernike (OZ) relations can be cast in the compact form

$$𝐇(q)=𝐖(q)𝐂(q)[𝐖(q)+𝐇(q)]$$
(7)
which may be solved for $`𝐇(q)`$ according to:

$$𝐇(q)=\left([𝐈𝐖(q)𝐂(q)]^1𝐈\right)𝐖(q).$$
(8)
The set of $`N(N+1)/2`$ independent OZ relations must be supplemented by as many closure relations between each pair of total and direct correlation functions. For intersite distances $`r>d_{ij}=(d_i+d_j)/2`$, we adopt the partially linearized Laria-Wu-Chandler closure :

$`{\displaystyle \frac{h_{ij}(r)+1}{h_{ij}^0(r)+1}}\mathrm{exp}[h_{ij}^0(r)\omega _{ii}c_{ij}^0\omega _{jj}(r)]`$ (9)
$`=\{\begin{array}{cc}\mathrm{exp}[\chi _{ij}(r)],\hfill & \chi _{ij}(r)0\hfill \\ 1+\chi _{ij}(r),\hfill & \chi _{ij}(r)>0,\hfill \end{array}`$ (12)
with

$$\chi _{ij}(r)=h_{ij}(r)\omega _{ii}(\beta u_{ij}+c_{ij})\omega _{jj}(r),$$
(13)
where $`\beta =1/(k_BT)`$ is the inverse temperature. The superscript $`0`$ refers to the underlying reference model with hard core interactions only (i.e., without the screened Coulomb interaction in Eq. (1)), while the asterisk $``$ denotes a convolution product. The reference site-site correlation functions $`c_{ij}^0(r)`$ and $`h_{ij}^0(r)`$ are obtained by solving the set of OZ equations (7) with the atomic Percus-Yevick closure for hard spheres, i.e. $`h_{ij}^0(r)=1`$, $`rd_{ij}`$; $`c_{ij}^0(r)=0`$, $`r>d_{ij}`$ . The closure (9) takes the long-range asymptotics of the direct correlation functions $`c_{ij}(r)`$, in terms of the molecular hypernetted chain (HNC) approximation , into account, but prevents the exponential increase of the total correlation functions $`h_{ij}(r)`$ for distances corresponding to strong attraction between oppositely charged sites . By construction, the total correlation functions and their first derivatives are continuous when $`\chi _{ij}(r)=0`$.

The closed set of equations (8) - (13) are solved numerically by a standard iterative procedure to obtain pair distribution functions $`g_{ij}(r)=1+h_{ij}(r)`$, and partial structure factors:

$$S_{ij}(q)=\omega _{ij}(q)\delta _{ij}+\rho _ih_{ij}(q).$$
(14)
## IV Colloid-polyelectrolyte mixtures

In this binary mixture ($`\nu =2`$), species $`1`$, referred to with the index $`p`$, are semiflexible polyelectrolyte chains with $`N_p`$ equivalent segments carrying each a charge $`z_pe=Z_pe/N_p`$ (where $`Z_pe`$ is the total charge of each chain). The form factor of the semiflexible chain is :

$`\omega _p(q)`$ $`=`$ $`1+{\displaystyle \frac{2}{N_p}}{\displaystyle \underset{n=1}{\overset{N_p1}{}}}(N_pn)`$ (16)
$`\times \mathrm{exp}\left[{\displaystyle \frac{q^2}{3}}[nd_pl_pl_p^2(1e^{nd_p/l_p})]\right],`$
where $`d_p`$ is the segment length and $`l_p`$ is the persistence length. The spherical colloids of radius $`R_c`$ (species $`2`$), carry two categories of sites: $`N_s`$ rigidly arranged surface sites of diameter $`d_s`$ carrying an electric charge $`z_se=Z_ce/N_s`$ (where $`Z_ce`$ is the total charge of each colloidal particle), and one central site, labelled with the index $`c`$, which is non-interacting (model A), or the centre of a hard sphere repulsion of radius $`R_c`$ (model B). In this paper only calculations based on model A will be presented. The corresponding form factors are $`\omega _{cc}(q)=1`$,

$`\omega _{ss}(q)`$ $`=`$ $`1+N_s\left({\displaystyle \frac{\mathrm{sin}(qR_c)}{qR_c}}\right)^2,`$ (17)
$`\omega _{cs}(q)`$ $`=`$ $`\sqrt{N_s}{\displaystyle \frac{\mathrm{sin}(qR_c)}{qR_c}},`$ (18)
where the double sums in Eq. (4) have been replaced by double integrals corresponding to a continuous distribution of interaction sites on the surface of the sphere of radius $`R_c`$. Since this corresponds formally to the limit $`N_s\mathrm{}`$, the two form factors (17) and (18) are approximately related by $`\omega _{ss}(q)=\omega _{cs}^2(q)`$. In practical calculations $`N_s=60`$ was chosen. If no microscopic co- and counterions are present, charge neutrality relates the polyelectrolyte and colloid density via:

$$Z_c\rho _c+Z_p\rho _p=0.$$
(19)
In the absence of microscopic ions, $`\kappa _D=0`$, and the potentials $`u_{ij}(r)`$ in Eq. (1) reduce to bare Coulomb potentials. With three classes of sites ($`p`$, $`s`$ and $`c`$), the matrix $`𝐖(q)`$ is given by:

$$𝐖(q)=\left(\begin{array}{ccc}\omega _{pp}(q)\hfill & 0\hfill & 0\hfill \\ 0\hfill & \omega _{ss}(q)\hfill & \omega _{cs}(q)\hfill \\ 0\hfill & \omega _{cs}(q)\hfill & \omega _{cc}(q)\hfill \end{array}\right).$$
(20)
There are a priori 6 independent direct and total correlation functions $`c_{ij}(q)`$ and $`h_{ij}(q)`$ ($`i,j=p,s`$ or $`c`$). However, independently of any specific closure relation, the OZ equations lead directly to the following three relations between the total correlation functions:

$`h_{ss}(q)`$ $`=`$ $`\omega _{cs}^2(q)h_{cc}(q),`$ (21)
$`h_{cs}(q)`$ $`=`$ $`\omega _{cs}(q)h_{cc}(q),`$ (22)
$`h_{ps}(q)`$ $`=`$ $`\omega _{cs}(q)h_{pc}(q).`$ (23)
These relations reflect the fact that the ”free-rotation” approximation  for rigid molecules is in fact exact when the interaction sites are distributed on the surface of a sphere. Other decoupling schemes which have recently been put forward for polymeric systems  are also exact in that limit.

We have considered mixtures of colloids of radius $`R_c=150`$ nm, and of oppositely charged semiflexible polyelectrolyte chains of contour length $`L=600`$ nm and persistence length $`l_p=25`$ and $`50`$ nm in water at room temperature ($`e^2/(\epsilon k_BT)=0.714`$ nm).

The calculated pair distribution functions $`g_{pp}(r)`$ are shown in Fig. 1 for three values of $`Z_c`$. The interchain correlations are seen to increase dramatically with $`Z_c`$, and to be much more pronounced for the more flexible chains. Qualitatively similar results hold for a mixture involving smaller colloidal particles and shorter chains, as shown in Fig. 2. The reliability of PRISM is illustrated in the inset of that Figure, where the PRISM predictions for $`g_{pp}(r)`$ is compared to available Monte Carlo data  for a mixture of charged spheres and oppositely charged dimers. The unusually large amplitude of the first peak in $`g_{pp}(r)`$ is the signature of an enhanced polyelectrolyte density near the oppositely charged spherical colloids, reminiscent of the well-known ”counterion condensation” of small ions . This condensation may drive phase separation  which can be diagnosed by large fluctuations in the colloid density.

Examples of the colloid centre of mass structure factor $`S_{cc}(q)`$ are shown in Fig. 3. A central peak in the $`q0`$ limit is seen to build up as the colloid density $`\rho _c`$ increases. The enhancement of the density fluctuations, as characterized by $`S_{cc}(q=0)`$, is illustrated in Fig. 4, as a function of $`\rho _c`$ for several values of $`Z_c`$. Although no divergence occurs, the tendency towards segregation is clear. In a naive random phase approximation (RPA) picture, values of $`S_{cc}(0)>1`$ point to an effective long range attraction between the colloidal particles. The existence of an attraction between the colloidal particles, induced by the oppositely charged polyelectrolyte, can be more firmly established by a direct inversion of the pair distribution function $`g_{ss}(r)`$ of surface sites on the colloids.

In the limit $`\rho _c0`$, this effective potential between surface sites would coincide with the potential of mean force, $`v_{ss}(r)=k_BT\mathrm{log}g_{ss}(r)`$, but at finite concentration of colloids, an inversion based on an accurate closure relation must be used. In other words the inverse problem is solved, whereby the pair distribution function $`g_{ss}(r)`$, as calculated from the solution of the multicomponent PRISM equations, is known, and the effective potential $`v_{ss}(r)`$ is extracted from the one-component version of the theory, using the corresponding one-component OZ and closure relations. Note that the inversion does not require an iterative procedure. Similar inversions have recently been successfully used to determine effective interactions between the centres of mass of interacting polymer coils , and between the interaction sites of rodlike polyelectrolytes .

Examples of the effective potentials are shown Fig. 5. As the charge $`Z_ce`$ increases, an attractive well develops and deepens. Instability towards colloid clustering and genuine phase separation is expected to occur when the well-depth becomes of the order $`k_BT`$, so that yet stronger Coulomb interaction, i.e., more highly charged polyelectrolyte, would be needed to drive the phase separation. This is presently under investigation. A comparison with the bare Coulomb interaction (1) clearly demonstrates the enhanced screening effect of the polyelectrolyte as the charge increases. The same figure also shows the effective potential $`v_{cc}(r)`$ between the centres of the colloidal particles, obtained from an HNC inversion of the centre to centre pair distribution function $`g_{cc}(r)`$ for one value of $`Z_c`$.

## V Bidisperse clay suspensions

We now turn to binary mixtures of oppositely charged disc-shaped platelets of different radii, $`R_1`$ and $`R_2`$ suspended in water, in the presence of microscopic co- and counterions. $`N_1`$ and $`N_2`$ interaction sites are assumed to be uniformly distributed over the surfaces of the discs of species $`1`$ and $`2`$, and interact via the screened Coulomb potential (1); the hard core diameters $`d_1`$ and $`d_2`$ determine the effective thickness of the platelets. In the limit of a continuous distribution of sites, the form factors $`\omega _i(q)\omega _{ii}(q)`$ ($`1i2`$) are easily calculated to be :

$$\omega _i(q)=\frac{2N_i}{(qR_i)^2}\left(1\frac{J_1(2qR_i)}{qR_i}\right),$$
(24)
where $`J_1`$ denotes the cylindrical Bessel function of first order and all interaction sites are assumed to be equivalent, i.e., edge effects are neglected. This may be a more severe approximation for platelets than for spheres (where the equivalence is exact), or rigid or semiflexible polyelectrolyte chains, but our previous results for the one-component case show that it leads to good agreement with available simulation data . The correlation function matrices $`𝐇(q)`$ and $`𝐂(q)`$ are now $`2`$ x $`2`$ matrices, and W(q) is diagonal.

We have solved the PRISM equations for platelet parameters and physical conditions appropriate for a water base drilling mud, made up of large, negatively charged bentonite platelets, (index $`b`$), and much smaller, positively charged mixed metal hydroxide platelets (index $`m`$). The radii of the two species were taken to be $`R_1=R_b=10^3`$ nm and $`R_2=R_m=25`$ nm; the effective thicknesses were chosen to be $`d_b=d_m=10`$ nm, mimicking conditions where the platelets are not fully dispersed and form stacks. The strong screening in water base drilling muds is mainly due to the addition of NaOH (pH=12), and this is modelled by adopting an inverse screening length, $`\kappa _D=1`$ nm<sup>-1</sup>. A homogeneous interaction site density on both types of platelets is achieved by imposing $`N_b/N_m=R_b^2/R_m^2=1600`$; in practice $`N_b=30400`$ and $`N_m=19`$, while surface charge densities were fixed at $`\sigma _b=0.09`$ e nm<sup>-2</sup> and $`\sigma _m=0.9`$ e nm<sup>-2</sup>. Examples of bentonite pair distribution functions $`g_{bb}(r)`$ are shown in Fig. 6. At the two lower clay concentrations, the position of the first peak scales roughly like $`\rho _b^{1/3}`$, but at the highest density the peak shifts dramatically towards the origin, a behaviour indicative of a collapse of nearly parallel pairs of bentonite platelets, under the influence of a strong effective attraction, induced by the smaller platelets that appear to ”condense” on the surfaces of the bentonite discs. Note that all secondary minima and maxima have disappeared from $`g_{pp}(r)`$ at the highest density. The corresponding structure factors $`S_{bb}(q)`$ (minus the fixed intramolecular contribution) are pictured in Fig. 7, together with the results obtained from PRISM for the same bentonite densities, but ten times lower densities of the smaller platelets. The structure factors $`S_{bb}(q)`$ at the two lower densities show typical liquid-like structure, with a main peak at $`q=2\pi /D_b`$, where $`D_b\rho _b^{1/3}`$. Increasing $`\rho _b`$ or adding mixed metal hydroxide shifts the peak toward higher $`q`$-values, corresponding to a shorter distance between bentonite platelets. Increasing $`\rho _m`$ also leads to an enhancement of the small angle ($`q0`$) value of $`S_{bb}(q)`$; a clear central peak appears at intermediate densities, while at the highest densities the amplitude of the central peak appears to diverge, signalling strong concentration fluctuations indicative of platelet aggregation, or spinodal instability reminiscent of the behaviour observed in recent computer simulations of mixtures of oppositely charged spherical particles . A similar strong small $`q`$ upturn has been observed experimentally for monodisperse clay suspensions in the gel phase , polyelectrolyte gels  and mixtures , Poly(styrenesulfonate) ion exchange resins , and low ionic strength polyelectrolyte solutions , but the exact physical origin of these large-scale fluctuations for these systems is still under debate.

## VI Conclusion

We have illustrated the mechanism of mesoscopic particle aggregation induced by oppositely charged polyions on two examples of complex two-component dispersions. The partial ”condensation” of the small particles leads to a considerable reduction of the strong Coulomb repulsion between the larger colloidal particles, which may go as far as an effective short-range attraction between the latter. The attraction, and resulting aggregation give rise to a strong enhancement of the long wavelength concentration fluctuations of the large particles, which is well characterized by the small $`q`$ limit of the corresponding partial structure factor.

While effective attractions have been widely studied theoretically and by simulations in the case of microscopic counterions , PRISM has allowed us to extend such studies to cases where the ”counterions” are themselves complex objects (semiflexible polyelectrolyte chains or clay platelets). PRISM turns out to be an accurate and flexible tool to tackle highly asymmetric complex systems of charged particles of various shapes, which are beyond the search of present day simulation methodology. The PRISM integral equations provide detailed structural information, in the form of partial, site-site pair distribution functions and structure factors. Thermodynamic properties can then be determined via the energy route, or the compressibility route , but the structural information is not sufficient to give direct access to the equation of state via the viral route, and to the free energy necessary to determine phase diagrams, which requires thermodynamic integration. Work along these lines is in progress, as well as further calculations based on model B for the colloid-polyelectrolyte mixture.

###### Acknowledgements.

The authors are grateful to Edo Boek for helpful discussions.


text:
# On the nature of the Extremely Red Galaxies

## 1. Introduction

EROs proved to be a very difficult population to study. They are too faint both in the optical and in the near-IR to obtain large number of spectra with enough signal-to-noise ratio to allow for a secure classification. Even the most recent deep spectroscopic surveys (e.g., Cimatti et al, 2002) can only reach K$``$19.2 and even in this case the number of unidentified objects is similar to the number of objects of each class. We have therefore searched for a diagnostics that could be applied to large sample of EROs without large amounts of telescope time. We noticed (Pozzetti & Mannucci, 2000) that the optical and near-IR Spectral Energy Distribution (SED) of a dusty object is expected to be swallower than that of an old elliptical because in the latter case the red color is given by the sharp 4000Å break (see fig. 1). As a consequence, dusty starburst with z$`<`$2 are expected to have redder J-K colors than old ellipticals. The positions of the two populations in the J-K vs. I-K (or R-K) color diagram were estimated by a large number of synthetic galaxy models. The expected difference in the J-K color is about 0.3 mag, even if for some models the difference reduces to zero. The derived diagnostic correctly reproduces all the objects with known nature (see fig. 2). A good photometry is needed to obtain reliable results as no empty region is present between the two populations and the relative distance is not large.

We have applied this diagnostic to a complete sample of 57 EROs selected to have K$`<`$20.3 and R-K$`>`$5.3 by Thompson et al (2000) in the field of a radio galaxy (see Mannucci et al., 2001, for details). The color-color diagram is shown in the left panel of fig. 3, together with the separation line. According to the models and the observations of cold dwarfs, the objects with J-K$`<`$1.4 were classified as stars. Among the other objects, 21 (37% of the sample) can be classified as starburst, the same number as ellipticals, and 10 objects (17% of the sample) remain unclassified. Several objects fall near the separation line so that their classification can be affected by the uncertainties on the colours. The significance of the classification was estimated by the ratio of the distance of each point from the classification line divided by the radius of the error ellipse in the direction perpendicular to the line. These results are shown in the right panel of fig. 3. Eleven ellipticals and 2 starbursts have a distance from the line below 1 sigma, therefore their classification should be considered weak.

We have checked the effect of the uncertainties on the relative number of objects of the two populations. While the statistical errors are unlikely to affect the result very much, the systematic effect on the colors and on the position of the line can effect the ratio between ellipticals and dusty galaxies of about 35%.

We have studied the morphology of the objects to check the consistency of the classification. A correlation was found between the K-band morphology and the classification: as expected, the objects identified as ellipticals have regular morphology and small angular dimensions, while the starburst objects appear to be more extended and often elongated.

We have derived a crude estimate of the contribution of the detected starburst galaxies to the cosmic star formation density. The cosmic variance and uncertainties on the Star Formation Rate (SFR) of each galaxies and on the cosmic volume actually sampled make the resulting value vary uncertain. Nevertheless a value of 0.03 M/yr was derived for a $`\mathrm{\Lambda }`$ cosmology with H0=70, $`\mathrm{\Omega }_m=0.3`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$. This is of the same order of the value derived from the LBG at higher redshift. This is an indication that the ERO technique is not only a good method to detect high-redshift ellipticals, but also a promising way to select a population of starburst galaxies containing a significant fraction of star formation and having intermediate properties between the LBG and the SCUBA objects (see fig. 4). It is also interesting to note that the results on the relative number of objects and cosmic star formation history of this work are in a very good agreement with those of the following K20 spectroscopic survey (Cimatti et al., 2002) even if the techniques are radically different.

## References

Cimatti, A., et al., 2002, A&A, in press (astro-ph/0111527).

Mannucci, F., et al., 2001, MNRAS, in press (astro-ph/0110648)

Pozzetti, L., & Mannucci, F., 2000, MNRAS, 317, L17