text:
# On expansion parallax distances for planetary nebulae

## 1 Introduction

Distances to individual nebular objects, such as Planetary Nebulae (PNe) are important to quantitatively understand their structure and evolution, but notoriously difficult to determine. Typically uncertainties of a factor of two can be expected from so-called statistical methods, a review of which can be found in Terzian (1993).

Probably the best method for measuring the distance to an individual nebula is the ‘expansion parallax’ method. Here the nebular expansion in the sky as measured from images taken at different epochs is compared to the radial velocity as measured spectroscopically. This method became feasible with the advent of high resolution imaging, in the late 1980’s at radio wavelengths using interferometers (Masson 1986, 1989a, 1989b; Hajian et al. 1993, 1995; Kawamura & Masson 1996; Hajian & Terzian 1996; Christianto & Seaquist 1998) and nowadays also in the optical with the Hubble Space Telescope (HST, Reed et al. 1999; Palen et al. 2002; Li et al. 2002).

Clearly the simplest application of this method suffers from a number of drawbacks. One complication is the choice for the spectroscopic velocity, where different ions often give different values. The other is the assumption of spherical expansion. However, authors have been making corrections for the shape and aspherical expansion of the measured nebulae, and applying such sophisticated templates clearly has made the method more useful, in some cases reaching claimed errors as low as 10–20%.

Still, the application requires a reasonably close or particularly rapidly expanding nebula and to date has only been applied to a limited number of PNe. Table 1 list all cases with a well determined expansion parallax distance, leaving out the cases where only upper limits were found.

Among the corrections needed to accurately use the expansion parallax method there is one which to date has received little attention in the literature, even though it is quite basic. This is the fact that the expansion velocity as measured in the sky is a pattern velocity, whereas the spectroscopically measured velocity is a material velocity, and the two are generally not the same. This effect was touched upon by Marten et al. (1993) and Steffen et al. (1997), but has not been taken into account in any of the published expansion parallax distances.

The extreme case would be that of an R-type ionization front making its way through a stationary medium, and clearly here the expansion parallax method becomes useless, as has been mentioned by various authors. Luckily, this situation is thought to be rare in PNe. However, also in the more common cases of shock fronts or slower moving D-type ionization fronts, the two velocities will differ. Although typically not by much, the effect is systematic, not random, and it should therefore be taken into account when using the expansion parallax method. Especially as the measurements are becoming more accurate over time, this factor can no longer be neglected.

In this paper I calculate the magnitude of the discrepancy between the two velocities for shocks (Sect. 2), ionization front structures (Sect. 3) and so-called ‘shells’ (Sect. 4). In Sect. 5 I apply the derived corrections to the published results of the expansion parallax method and suggest strategies to optimize the expansion parallax method.

## 2 Shock waves

Since PNe are basically wind-driven bubbles, they contain shocks. Numerical hydrodynamic modelling has shown that in round or elliptical nebulae the inner bright rim is normally associated with the shock wave being driven into the surrounding material by the stellar wind. Key papers describing the radiation-hydrodynamic evolution of PNe are Marten & Schoenberner (1991); Mellema (1994, 1995), and a review is presented in Schönberner & Steffen (2003). The application to an individual PN is shown in Corradi et al. (2000). Whenever I mention numerical simulations in what follows, I refer to these papers.

The jumps of density and velocity across a shock are given by the Rankine-Hugoniot conditions. These are usually given for the reference frame of the shock, in which case the velocity jump is

$$\frac{v_0}{v_1}=\frac{(\gamma +1)^2}{(\gamma 1)^2+2},$$
(1)
which for infinite Mach number $``$ and an adiabatic index $`\gamma =5/3`$ gives the classical value of 4. The Mach number of the shock is given by

$$=\frac{|v_0|}{a_0},$$
(2)
where $`a_0`$ is the sound speed in the pre-shock gas, given by

$$a_0=\sqrt{\frac{\gamma kT_0}{\mu m_\mathrm{H}}},$$
(3)
with $`k`$ the Boltzmann constant, $`T_0`$ the gas temperature, $`\mu `$ the mean molecular weight, and $`m_\mathrm{H}`$ the mass of atomic hydrogen.

However, in the stellar frame, the shock has a velocity which we will call $`u_\mathrm{s}`$ and the pre- and post-shock velocities are given by $`u_{0,1}=v_{0,1}+u_\mathrm{s}`$, see Fig. 1.

The expansion parallax method measures $`u_s`$ from the angular expansion of the shock front, $`\dot{\theta }=u_s/D`$, and $`u_1`$ from the spectroscopy, and derives the distance $`D`$ from the ratio of the two, assuming $`u_s=u_1`$. From the expression above it is immediately clear that $`u_s`$ and $`u_1`$ are not equal. If the PN woud be expanding into complete vacuum the two velocities would be identical, but in the case where it is moving into a surrounding medium, it is more proper to view the expansion as a disturbance or wave travelling through the medium, raising its density and velocity. The speed at which the wave travels is not the same as the velocity to which the material is accelerated by the passage of the wave.

Let us consider the ratio $``$ of the two velocities, $`u_\mathrm{s}/u_1`$. The distance to the PN can then be expressed as $`u_1/\dot{\theta }`$. Using the relations between $`u_{0,1}`$ and $`v_{0,1}`$ together with Eqs. 1 and 2 (realizing that in our shock frame both $`v_0`$ and $`v_1`$ are negative) gives

$$=\frac{(\gamma +1)u_0+(\gamma +1)^2a_0}{(\gamma +1)u_0+2(^21)a_0}.$$
(4)
The limits of this ratio are $`(u_0+a_0)/u_0`$ for $`1`$ and $`(\gamma +1)/2`$ for $`\mathrm{}`$. This shows that only for isothermal ($`\gamma =1`$) hypersonic shocks the ratio tends to one for high Mach numbers. Choosing values for the pre-shock velocity $`u_0`$, pre-shock sound speed $`a_0`$, and shock Mach number $``$ gives a value for $``$. The material velocity $`u_1`$ can then be found from

$$u_1=u_0+2a_0\left(\frac{^21}{(\gamma +1)}\right).$$
(5)
Reversely, given values for $`u_0`$, $`u_1`$, and $`a_0`$, the shock’s Mach number can be found from

$$=\frac{(\gamma +1)(u_1u_0)+\sqrt{(\gamma +1)^2(u_0u_1)^2+16a_0^2}}{4a_0},$$
(6)
which can then be used to derive $``$ using Eq. 4.

In Fig. 2 I plot the ratio $``$ as function of $`u_1`$, for six pre-shock velocities $`u_0`$ in the range 1 to 25 km s-1, using a $`\gamma `$ value of 5/3 (monatomic ideal gas), and $`a_0=15`$ km s-1, the adiabatic sound speed for an ionized gas at $`10^4`$ K. Most PNe have high densities and relatively slow shocks, so that the shocks are expected to be isothermal. Using $`\gamma =1`$ is therefore appropriate, and in Fig. 3 I show $``$ as function of $`u_1`$ for this isothermal case, where I have used an isothermal sound speed of 11.7 km s-1, valid for an ionized gas of cosmic abundances at a temperature of $`10^4`$ K. As expected from the limits derived above, the ratio is largest for low velocities. I should note that the ratio does depend on the choice for the sound speed $`a_0`$. A value corresponding to an electron temperature of $`10^4`$ K is typical for PNe, but the temperature can be both higher and lower than this, ranging from 5000 to 15,000 K. Lower temperatures give lower values of $``$.

Looking at the data in Table 1, it is clear that for these PNe, we are never in the very high Mach number regime. Using Fig. 3 one can see that for the observed ranges of velocities and for $`\gamma =1`$, the typical ratio $``$ is between 1.3 and 1.5, although it can be as high as 1.8. Obviously this falls outside the formal errors of the method and becomes as important an effect as the geometric corrections applied for example by Li et al. (2002). Furthermore, the factor is always larger than 1, so it does not make sense to add a 20-30% extra uncertainty to the distances. Rather, the distances should be scaled up by 30% and an extra error of $`10`$% added to it. Figure 3 can be used to estimate the magnitude of this effect for individual PNe, which I will do in Sect. 5.

The conclusion then is that if the expansion parallaxes are measured from the shocked component of the PN, the distances should be multiplied by a factor $``$ given by Eq. 4, where uncertainties in the velocity of the material into which the shock is expanding adds a $``$10% error. I will come back to this point in Sect. 5.

## 3 Ionization fronts

Ionization fronts are either R- or D-type. The basic difference is the speed of the front which for R-type fronts is higher than twice the sound speed in the neutral gas, giving the gas no time to react to the presence of the ionization front, wheras for D-type it is lower than twice the sound speed in the neutral gas. The result is that D-type fronts actually consist of a combination of an ionization front and a preceding shock front. A good overview of basic ionization front theory can be found in Shu (1992).

The numerical simulations for the formation of PNe which studied the effects of ionization fronts in detail, showed that although R-type fronts occur, they never persist for long, and it would be highly unlikely to observe one in action. When they occur, they are of the so-called weak type, i.e. supersonic with respect to both the neutral and ionized flow regions.

Weak D-type fronts (moving subsonically with respect to both the ionized and neutral regions) are more commonly found in the models, and persist for a longer time, and could be present in a number of observed PNe, although it is generally hard to prove this. The simulations also show that the presence of a D-type front sets up a disturbance in the gas which persists long after the front itself has disappeared, creating the shells which I discuss in the next section.

Basic ionization front theory shows that the velocities $`v_1`$ and $`v_2`$ on either side of the front are related by

$$\frac{v_1}{v_2}=\frac{1}{2a_2^2}\left(a_1^2+v_1^2\pm \left[(a_1^2+v_1^2)^24a_2^2v_1^2\right]^{1/2}\right).$$
(7)
The indices 1 and 2 refer to the neutral and ionized sides, respectively. Just as in the previous section, $`v`$ refers to velocities in the frame of the discontinuity, and $`u`$ to velocities in the stellar frame. If the front moves with a velocity $`u_\mathrm{i}`$ in the stellar frame, the velocities on either side of it are $`u_{1,2}=v_{1,2}+u_\mathrm{i}`$. For the weak R-type fronts one has to use the $``$ sign in Eq. 7 and for the weak D-type fronts the $`+`$. The velocity $`v_1`$ has to be larger than $`v_\mathrm{R}=a_2+(a_2^2a_1^2)^{1/2}`$ for R-type fronts and smaller than $`v_\mathrm{D}=a_2(a_2^2a_1^2)^{1/2}`$ for D-type fronts.

In Fig. 4 I show the ratio $``$, now defined as $`u_\mathrm{i}/u_2`$, against $`u_2`$ for weak R-fronts running into a medium with $`u_1=10`$ km s-1. For the sound speeds I chose $`a_1=1`$ km s-1 and $`a_2=10`$ km s-1, and tests show that the ratio $``$ is hardly sensitive to this choice. As can be expected this ratio is large for very fast moving fronts, and approaches 1 for the slower moving ones. Because of the large range of ratios, it would seem dangerous to apply any type of correction in case of R-type ionization fronts.

In Fig. 5 I show the same ratio $``$ against $`u_2`$ for weak D-fronts, a situation which is more likely to occur in real PNe. Normally the shock front which precedes the D-front will have accelerated the AGB wind to a higher velocity. The numerical simulations show this velocity to be typically around $`u_1=20`$ km s-1. I plot the ratio for three values of $`u_1`$. Interestingly, around the typically observed velocities, the correction is of the same order of magnitude as for the shock waves in Sect. 2.

The full solution for a D-type front would eliminate the choice for $`u_1`$. Shu (1992) showed how this can be done for a constant flux of ionizing photons impinging on a static and constant density environment. Attempts to do this for a stellar wind type environment have been only partly successful. Giuliani (1989) found a self-similar solution which requires the flux of ionizing photons to be time dependent and fall off as $`t^1`$, which does not apply to most PNe.

Masson (1986, 1989a, 1989b) applied a series of corrections to his determination of the distance to NGC 7027, BD+30 3639 and NGC 6591, which are partly related to the difference between the material velocity (which he referred to as the “bulk velocity”), and the (pattern) speed of the ionization front. The magnitudes of these corrections are actually similar to the ones found above, 1.1 to 1.2. Perhaps somewhat confusingly he applies these corrections to the measured angular expansion rate, so that the final figures he quotes for $`\dot{\theta }`$ are not actually the ones measured. He derives this correction from the fact that as the ionized shell expands, the density will go down, reducing the number of photons used up in recombinations, and allowing the ionization front to expand. This is a different approach to the one used above and does not use the jump conditions. However, the basic idea is still that the ionization front expands faster than the material flow of the gas, so the corrections are related. In not considering the jump conditions he implicitly assumes a weak R-type ionization front with a small density jump, and the correction corresponds therefore to the one from Fig. 4. The numerical simulations as well as the observed velocities show that this is not correct, although the correction factors come out at similar values. Note that Masson is the only author to actually apply any type of correction. In Sect. 5 I will comment some more on this.

The conclusion from this section is that it is more difficult to find a general solution for ionization fronts, which makes them less suitable than shock fronts for applying the expansion parallax method. However, the results in this section allow for corrections to be made, and below I will illustrate how these corrections could be applied.

## 4 Shells

Numerical simulations show that the density disturbance induced by the D-type ionization front will persist and can be identified with the observed attached shells around the bright core nebula, such as for example in NGC 3242. Corradi et al. (2003) compiled a list of PNe with haloes and shells, which can serve as a reference. We follow their nomenclature and refer to these structures as shells.

Analytical models for the expansion of such shells are available in Chevalier (1997) and Shu et al. (2002), even though the latter had a different application in mind. Using their results, it is possible to extract a similar ratio $``$ of pattern speed over material speed, where the pattern speed is now the movement of the edge of the shell. Since these models actually use the isothermal shock conditions, the results are basically identical to ones from Sect. 2. However, I list them separately, since shells are commonly observed and the Chevalier (1997) and Shu et al. (2002) results are the full solutions including the acceleration of the AGB wind after ionization.

For an AGB wind of constant mass loss and velocity, the results from Chevalier (1997) show that $``$ is between 1.10 and 1.34. I plot this ratio against $`u_1`$ in Fig. 6, together with the ratio from his solution including the effect of a stellar wind, which gives somewhat lower values. The ratio for an AGB wind density falling off steeper ($`\rho r^\alpha `$ with $`\alpha >2`$, Shu et al. 2002) are also shown.

The numerical simulations also show that when an attached shell forms, its velocity structure is that of a rarefaction wave, with a positive outward gradient. This means that the bright rim expanding into this will be moving into an area with a velocity lower than the original AGB wind velocity. The implication for the expansion parallax method is that when using the rims of PNe with attached shells a value of $`u_0`$ lower than 10 km s-1 is more appropriate when determining the correction factor from Sect. 2.

The conclusion from this section is that for attached shells a correction of around 20% is needed.

## 5 Implications for distances and the method

The results presented thus far show that for the published distances, the results should typically be scaled up with a factor $`1.3\pm 0.2`$, which interestingly enough eliminates the systematic discrepancy between the expansion parallax distances and statistical distances noted by Palen et al. (2002). In this section I will go through the list of PNe for which distances have been measured, consider corrections for each case, and give some suggestions on how to improve the usage of the expansion parallax method further. For calculating the corrections I use the values $`\gamma =1`$ and $`a_0=11.7`$ km s-1 for the cases of shocks (corresponding to Fig. 3), $`a_1=1`$ and $`a_2=10`$ km s-1 for ionization fronts (corresponding to Fig. 5). The first part of this section is intended not only to derive new distances, but also to illustrate how to derive and use the correction factors.

### 5.1 New distance estimates

Going through the list there are a number of PNe which contain bright rims and attached shells, and the rim can therefore be considered to be bounded by a shock. These are NGC 3242, NGC 6578, NGC 6884, NGC 7662, and IC 2448. Making conservative assumptions about the value of $`u_0`$ ($`13\pm 12`$ km s-1) and using the results from Sect. 2, the distances to these PNe should be increased by factors $`1.3\pm 0.2`$, $`1.45\pm 0.15`$, $`1.5\pm 0.15`$, $`1.5\pm 0.3`$, and $`1.5\pm 0.1`$ respectively. The uncertainties are due partly to the uncertainties in the reported spectroscopic velocities and partly to the uncertainty in the value of $`u_0`$. The latter could in principle be reduced if one assumes that the velocities in the shells follow the numerical models, in which case $`u_0`$ will be low ($`<10`$ km s-1). I list the new (conservative) values in Table 2.

NGC 6543 is a more complex nebula. The Chandra results show that it contains a wind-driven bubble, but the ionized material surrounding this bubble does not resemble a standard shell. This complicates the choice for $`u_0`$, but within reasonable limits the correction factor is $`1.55\pm 0.15`$.

Three PNe in Table 1 are candidates for the presence of an ionization front: BD+30 3639, NGC 7027, and Vy2-2.

Of these BD+30 3639 is a single shell, and relatively young PN with molecular emission around it, although the deeper exposures also seem to show an ionized component around the bright shell (Sahai & Trauger 1998). Kawamura & Masson (1996) used a detailed analysis of the expansion of this object in order to derive the distance, assuming we are observing an ionization front. Part of this analysis is a correction factor of 1.19 to account for the difference between the material and pattern velocities, as mentioned in Sect. 3. Using the results from Sect. 3 for a weak D-type front, I estimate a correction factor of $`1.2\pm 0.2`$ assuming a velocity of $`25\pm 5`$ km s-1 for the shocked neutral gas and a spectroscopic velocity of $`22\pm 4`$ km s-1 (as used by Kawamura & Masson 1996), confirming their result for the distance.

Li et al. used optical HST data, better spectroscopic data and a correction for the ellipticity of the nebula to arrive at a smaller distance. Using their preferred value for the spectroscopic velocity, 25.6 km s-1, I find a correction factor of $`1.1\pm 0.1`$ if we are seeing an ionization front running into a shocked component with a velocity of 25.6 – 30 km s-1, and a correction of $`1.3\pm 0.1`$ if we are seeing a shock. It therefore seems reasonable to increase the distance by 10% and to double the uncertainty on that number: $`1.3\pm 0.2`$ kpc, bringing it closer to the value from Kawamura & Masson (1996)

NGC 7027 is a more complex PN where an ionization front has to be present, given the large amounts of dust and molecules surrounding this object. The elliptical shell seen at radio wavelengths (and partly obscured at optical wavelengths) could be the ionization front and is analyzed as such by Masson (1986, 1989a). To correct for the discrepancy between spectroscopic and angular velocities, he uses a factor of 1.2 in the first and 1.15 in the second paper (note that the second paper quotes 1/1.15=0.83 as the correction factor). Further corrections for the decrease of the radio flux in the second paper actually largely cancel the effect of this factor.

Bains et al. (2003) analyzed high resolution optical long slit spectra of NGC 7027 and from the \[O III\] line derived an equatorial expansion velocity considerably lower than used by Masson (1989a), namely $`13\pm 1`$ km s-1. Cox et al. (2002) report K-band imaging and spectroscopy for NGC 7027. Their best fitting model has an equatorial expansion velocity in Br$`\gamma `$ of about 13 km s-1 (Cox & Huggins, private communication; this number is not given in the paper), consistent with the \[O III\] value. For the molecular emission (H2) the same authors find an equatorial velocity of $``$15 km s-1. Allowing all possible values for the velocity of the neutral material $`u_1`$ (13 – 23 km s-1), the correction factor would be $`1.4\pm 0.4`$. Assuming that the H2 expansion is indicative of the value of $`u_1`$ (13 – 17 km s-1), the ratio becomes $`1.2\pm 0.2`$. If the ionized shell is actually bounded by a shock, the results of Sect. 2 show that the correction factor would be $`1.75\pm 0.15`$ for all allowed values of the pre-shock velocity $`u_0`$.

Taking an angular expansion of $`4.84\pm 0.82`$ mas yr-1 (Masson 1989a, this includes the correction for flux variations $`0.47\pm 0.47`$), but not the correction factor of 0.83), combining this with an equatorial velocity of 13 km s-1, and applying a correction of $`1.2\pm 0.2`$, I arrive at a new distance of $`680\pm 170`$ pc. Note that this error is rather optimistic in view of the wide range of correction factors mentioned above.

Vy2-2 is a compact PN which is hard to categorize. Also here one may have to be aware of ionization fronts since this is a fairly young, low excitation PN. Assuming either an ionization front or an isothermal shock, the reported distance should be scaled up by $`1.3\pm 0.3`$

Two PNe in the list are hard to categorize in the standard wind-blown bubble plus photo-ionization scheme: NGC 6210 and NGC 6572. I therefore will not suggest any corrections for these, although shocks could very well be present. For NGC 6210 this would imply a 30% increase in the reported distance, but for NGC 6572 the correction could become very large (1.4–1.8) due to the low value of the reported spectroscopic velocity.

### 5.2 Improving the method

The correction factors derived in this paper depend on a number of parameters. As already indicated, due to the isothermal nature of the slow shocks in PNe, the preferred value for the adiabatic index $`\gamma `$ is 1. The sound speed in the ionized medium ($`a_0`$) I have taken to be 11.7 km s-1 in Fig. 3, which is the value for an electron temperature of $`10^4`$ K. If the PN is known to have a particularly high or low electron temperature the figure should be recalculated for the appropriate sound speed, since the values for $``$ unfortunately do depend on the choice of the sound speed. The appropriate sound speed to use here is the isothermal sound speed.

Figure 3 or its equivalent can then be used to find the correction factor for the range of values of $`u_0`$ and $`u_1`$ which seem reasonable. Observational data on $`u_0`$ is scarce, but if it cannot be further constrained, $`10\pm 10`$ km s-1 should cover most cases. For $`u_1`$ the question is which spectroscopic velocity to pick. Ultimately, for an individual PN this question can only be answered through detailed (photo-ionization or hydrodynamic) modelling of the PN, assuming it is actually possible to produce a unique model for it. The approach of Gesicki et al. (1998), who derive spatial velocity profiles using photo-ionization models, may be useful here. Generic hydrodynamic modelling may help establish what would be a good choice in general cases and I understand that a project to do this is under way (Schönberner, private communication). Until that time the choice for $`u_1`$ will introduce uncertainties in the method, as it has always done. It is important to realize that the jump conditions are valid just before and after the front, so the best choice would be a velocity as close as possible to the front.

Round or elliptical PNe with attached shells are particularly suitable for the expansion parallax method. Their inner bright rim is bounded by a shock, and the presence of a surrounding ionized shell shows that this is not an ionization front. For these almost round PNe, long slit spectroscopy can give a good indication of the spatial and velocity structure, since it will be quite similar in the directions perpendicular and parallel to the line of sight. This way a distinction can be made between the velocity of the rim and the shell (the latter can be both higher and lower than that of the rim). If possible it is best to choose species which lie closest to the outer edge of the rim so as to avoid any velocity variations within the rim. For these rims numerical models show that they are expanding into a region of low velocity, typically lower than that of the original AGB wind. Choosing a value for $`u_0`$ of around 5 km s-1 is a reasonable guess. In any case, the corrections as derived in Sect. 2 only weakly depend on this choice. Using the measured spectroscopic velocity, the ratio $``$ can be determined.

It would be interesting to derive an independent distance using the attached shells, if available. If the long slit spectroscopy provides a velocity for the shell, then the corrections from Sect. 4 can be used to find a distance. Models show that the velocities of the shells are not as strongly position dependent as those of the inner rims, something which is reflected in the only mild ellipticities found observationally. This may make the outer edges of shells actually better suited for distance determinations with the expansion parallax method.

It is best to avoid areas and features well away from where a shock is suspected. At least for a shock we can make a connection between the pattern and material velocities. Away from these discontinuities, the two will also differ, but there is no way to know how.

Bipolar PNe display a larger range of expansion velocities, which means that it is harder to correct for inclination effects. Palen et al. (2002) state that the method is unusable for extreme bipolars, which seems too pessimistic. However, a thorough understanding of the dynamics of the complete PN is essential to be able to apply the method to these PNe, and in most cases this information is not available. It would also be valuable to have more radiation-hydrodynamic modelling of such systems.

## 6 Conclusions

Pattern velocities and the material (or bulk) velocities in a gas are not necessarily the same. For discontinuities, such as shocks and ionization fronts, the relation between the two can easily be derived using the jump conditions across them. The pattern velocity is then found to be always higher than the material velocity. Since measuring the expansion of PNe in the sky is mostly done using sharp edges, which are associated with either shocks or ionization fronts, a correction should be applied before calculating the distance from the ratio of the two velocities. This correction is typically larger for velocities of the order one or two times the sound speed in the ionized material, which is actually what is measured in most PNe.

Not using this correction will systematically underestimate the distances to PNe. For the sample of PNe to which the expansion parallax method has been applied successfully, the correction factors are around 1.2 to 1.3 for both shocks and ionization fronts. Applying the corrections given in this paper should lead to improved distance determinations to PNe.

###### Acknowledgements.

I like to thank Pierre Cox and Patrick Huggins for providing me with a value for the equatorial expansion velocity of NGC 7027. My research has been made possible by a fellowship of the Royal Netherlands Academy of Arts and Sciences.


format: text

text:
# Complexation in Polyelectrolyte Solution with Divalent Surfactant

## 1 Introduction

Gene therapy represents a promising way for the treatment of both genetic and acquired diseases. The basic idea is to replace the sick gene with a healthy one or in some cases to add a new gene to get a resulting synthesis of a therapeutic protein. The process of getting the new gene to its target involves the crossing of several barriers, such as the cell membrane and the nuclear membrane. Since both the DNA and the cell membranes are negatively charged, the naked polynucleotides are electrostatically prevented from entering into the cell. Viral vectors such as retroviruses and adenoviruses are very efficient and able to target a wide range of cells . The gene transfection is accomplished by the use of a virus in which the native DNA has been replaced by the required DNA. Since the main purpose of the virus is to replicate itself, the new gene is successfully transfected to the cell nucleus by endocytosis or membrane fusion. This process suffers from the drawback that in cases where repeated treatment is needed, this might cause the immune system to react negatively due to the viral origin of the vector.

As a solution for this problem nonviral vectors have been developed. One of the approaches, pioneered by Felgner and Ringold  relies on the association between anionic nucleic acid and cationic lipid lipossomes. The process of association neutralizes the excess negative charge of the DNA and the DNA-lipossome penetrate into the cell by endocytosis. The efficiency of this method is low and they are toxic to the cell at the concentrated used.

Among the parameters needed for achieving efficient transfection is the requirement that the complex should have small size similar to a virus . This factor indeed limits the efficiency of the DNA-lipossomes complexes. In contrast, cationic surfactant have been shown to condense to DNA into discrete particles containing a single nucleic acid molecule  that can neutralize and revert the charge of the DNA - allowing the surfoplex to approach the membrane. Monovalent DNA-surfactant complexes exhibit a discrete first-order phase transition between elongated coil and collapsed globule  for a concentration of amphiphilic molecules well below the miscelar concentration. Despite this unique feature, monovalent detergent are poorly efficient in vitro in gene transfer . When the surfoplex approaches the membrane, the interaction between the surfactant molecules and the phospholipid bilayer overcomes the electrostatic attraction between the anionic nuclei acid and the cationic detergent. The surfactant incorporate into the phopholipid membrane resulting in the unfolding of the DNA that stays outside the cell .

Therefore, in order to have a stable surfoplex, one needs to enhance the electrostatic interaction. In this paper we present a model of DNA-amphiphilic solution where the surfactant is divalent. We find that in equilibrium, solution consists of complexes composed of DNA and associated counterions and amphiphiles. Even for an small amount of amphiphiles in the solution, the cooperative binding is found. Due to the high valence, the interaction between the nuclei acid and the detergent is expect to overcome the interaction between the surfactant molecules and the phospholipids in the membrane and stay associated to the DNA allowing it to transpass the membrane. The complex formed is more stable than the monovalent structure, forming a collapsed globule that we expect will transect with high efficiency .

## 2 The Model

Our model, illustrated in Fig. 1, consists of a solution of DNA strands of length $`L`$ and diameter $`a_p`$, divalent surfactant and monovalent salt. In aqueous solution, the polyions become ionized resulting in a negative charge $`Zq`$ distributed uniformly distributed with separation $`b=L/Z`$. The solvent, water, is modeled as a continuous medium of dielectric constant $`D`$. The ions of the salt are completely dissociated, forming an equal number of positive and negative ions. Similarly, the surfactants are assumed to be fully dissociated producing negative monovalent coions and polymeric chains with a divalent cationic head group. For simplicity, all the counterions and coions are treated as identical, independent of the molecules from which they were derived. The electrolytes are depicted as hard spheres of diameter $`a_c`$ and charge $`\pm q`$ and the surfactant is modeled as a polymer of $`s`$ monomers each one considered as a rigid sphere of diameter $`a_c`$ with the head monomer carrying a charge of $`+2q`$.

The interaction between the hydrophobic tails is short ranged and characterized by the parameter $`\chi `$. The density of DNA strands is $`\rho _p=N_p/V`$, the density of monovalent salt is $`\rho _s=N_s/V`$ and the density of divalent amphiphiles is $`\rho _a=N_a/V`$.

The strong electrostatic interactions between the polyions, the counterions, salt cations and surfactant heads leads to the formation of complexes, which in thermodynamic equilibrium will be made up of one polyion, $`n_c`$ monovalent counterions and $`n_a`$ divalent surfactant. We do not consider the effects of of polydispersity in the size of the complexes, since it does not affect the final result. Due to the association and to the charge conservation, there are only two free quantities and so,

$`\rho _c=(Zn_c)\rho _p+\rho _s\rho _{a+}=\rho _an_a\rho _p\rho _{}=\rho _s+2\rho _a`$ (1)
where $`\rho _c`$ is the density of free monovalent counterions, $`\rho _{a+}`$ is the density of free amphiphiles and $`\rho _{}`$ is the density of negative ion.

The objective of this theory is to determine the number of counterions $`n_c`$ and surfactants $`n_a`$ associated to each DNA strand. For this, we construct the Helmholtz free enery of the system and minimize it. The details about the model can be found elsewhere -. We give here the main steps. The relevant contributions for the Helmoltz free energy are two, the electrostatic and the entropic namely:

$$F=F_{el}+F_{ent}$$
(2)
In $`F_{el}`$ three types of interactions can be found: between the free ions and free surfactants ,$`F_{is}`$, between the complex, free ions and free surfactants ,$`F_{pis}`$, and between the complexes $`F_{pp}`$.

$$F_{el}=F_{is}+F_{pis}+F_{pp}$$
(3)
With the aid of the theory of Debye-Hückel-Bjerrum (DHBj) it is possible to find the electrostatic interaction between the complexes, ions and surfactant given by -:

$`\beta f^{pis}={\displaystyle \frac{\rho _pZ_c^2(a/L)}{T^{}(\kappa a)^2}}\left\{2\mathrm{ln}\left[\kappa aK_1(\kappa a)\right]I_0+{\displaystyle \frac{(\kappa a)^2}{2}}\right\}`$ (4)
with $`\beta =1/k_BT`$ and

$$I_0_0^{\kappa a}\frac{xK_0^2(x)}{K_1^2(x)}𝑑x$$
(5)
where $`\kappa `$ in $`(\kappa a)^2=4\pi \rho _1^{}/T^{}`$ is the inverse of the Debye screening length, $`\rho _1^{}=\rho _1a^3=\rho _ca^3+4\rho _{a+}a^3+\rho _{}a^3`$ is the reduced density and $`T^{}=Dk_BTa/q^2`$ is the reduced temperature. Furthermore, $`Z_c=Zn_c2n_a`$ is the valence of each complex and $`a=(a_c+a_p)/2`$ is the effective radius of the exclusion cylinder around each complex.

In the framework of the Debye-Hückel theory, the interaction between the free ions is given by:

$`\beta f^{is}={\displaystyle \frac{1}{4\pi a_c^3}}\left[\mathrm{ln}(1+\kappa a_c)\kappa a_c+{\displaystyle \frac{(\kappa a_c)^2}{2}}\right].`$ (6)
The electrostatic free energy interaction between two complexes for large separations is screened. The short-range of this interaction allows one to use a mean-field approximation resulting in:

$$\beta f^{pp}=\frac{2\pi a^3Z_c^2\rho _p^2\mathrm{exp}(2\kappa a)}{T^{}(\kappa a)^4K_1^2(\kappa a)}.$$
(7)
The calculation of the entropic contribution can be obtained with the aid of Flory  theory of mixing. The free energy is a sum of ideal free energies of various species, namely

$$\beta f^{ent}=\left[\rho _s\rho _s\mathrm{ln}\left(\frac{\varphi _s}{\zeta _s}\right)\right]$$
(8)
where $`s`$ represents the different species and $`\zeta _s`$ is the internal partition of the species $`s`$. In the case of the particles without structure, the internal partition function $`\zeta _{}=\zeta _c=\zeta _{a+}=1`$. The volume fraction $`\varphi _s`$ of the different species are:

$`\varphi _p`$ $`=`$ $`{\displaystyle \frac{\pi \rho _p^{}}{4(a/L)}}\left({\displaystyle \frac{a_p}{a}}\right)^2+{\displaystyle \frac{Z\pi \rho _p^{}}{6}}(s_am_a+m_c)\left({\displaystyle \frac{a_c}{a}}\right)^3`$
$`\varphi _c`$ $`=`$ $`\rho _c^{}{\displaystyle \frac{\pi }{6}}\left({\displaystyle \frac{a_c}{a}}\right)^3`$
$`\varphi _{a+}`$ $`=`$ $`{\displaystyle \frac{s_a\pi \rho _{a+}^{}}{6}}\left({\displaystyle \frac{a_c}{a}}\right)^3`$
$`\varphi _{}`$ $`=`$ $`{\displaystyle \frac{\pi \rho _{}^{}}{6}}\left({\displaystyle \frac{a_c}{a}}\right)^3.`$ (9)
Here we introduce the fractions counterions, $`m_c=n_a/Z`$, and surfactant, $`m_a=n_c/Z`$, associated to each DNA strand. Accounting for all these terms, the entropic contribution becomes

$`\beta f^{ent}`$ $`=`$ $`\rho _p\mathrm{ln}\left({\displaystyle \frac{\varphi _p(1+m_c+m_a)}{\zeta _{cl}(1+m_c+m_as_a)}}\right)\rho _p`$ (10)
$`+`$ $`\rho _c\mathrm{ln}\varphi _c\rho _c+\rho _{a_+}\mathrm{ln}{\displaystyle \frac{\varphi _{a_+}}{s_a}}\rho _{a_+}`$
$`+`$ $`\rho _{}\mathrm{ln}\varphi _{}\rho _{}`$
where the entropy for the free surfactants arises from the Flory theory for polymers  which is also the basis of the entropic contribution for the complex . The internal partition function of the complex, $`\zeta _{cl}`$, can be calculated by modeling the DNA by an one dimensional lattice with Z sites. If the number of associated ions to each site can be only zero or one, this problem becomes equivalent to finding the free energy of an one dimensional array with the three different states.

Due to the presence of two different valences, the exact solution of this model is not trivial, so we employ the mean-field Gibbs-Bogoliubov-Feymman inequality. The resulting partition function is given by

$`ln\zeta _c[m_c,m_a]`$ $`=`$ $`\xi K\left[{\displaystyle \frac{Z_c}{Z^2}}^21\right]+\beta \chi \left(Z1\right)m_a^2`$
$`+`$ $`Zm_c\mathrm{ln}m_c+Zm_a\mathrm{ln}m_a`$
$`+`$ $`Z\left(1m_cm_a\right)\mathrm{ln}\left(1m_cm_a\right)`$
where $`\xi \beta q^2/Da`$ is the Manning parameter, $`K=Z[\psi (Z)\psi (1)]Z+1`$, and $`\psi (n)`$ is the digamma function.

The equilibrium configuration of the system is found by the minimization of the Helmoltz free energy, leading to two equations, namely

$`{\displaystyle \frac{F}{m_c}}\delta m_c`$ $`=`$ $`0`$
$`{\displaystyle \frac{F}{m_a}}\delta m_a`$ $`=`$ $`0.`$ (12)
Solving this system of two equations, it is possible to obtain the values of $`m_c`$ and $`m_a`$.

## 3 Results and Conclusions

We define a ”surfoplex” to be a complex in which almost all of the DNA’s phosphate groups are neutralized by the associated surfactant molecules.As mentioned earlier, we are interested in the minimum amount of cationic surfactant needed to transform naked DNA into surfoplexes. To this effect, we study the dependence of the number of condensed surfactant molecules on the bulk concentration of surfactant $`\rho _a`$. In order to evaluate the relevance of hydrophobic interactions between the amphiphiles, the hydrophobic parameter was varied from 0 to $`\beta \chi =6`$. The effect of addition of high concentrations of salt to the system was analyzed by varying the amount of salt added to the system $`\rho _s`$.

Fig.$`2`$ and Fig.$`3`$ illustrate the surfactant and the counterion binding isotherm, $`m_a`$ and $`m_c`$ respectively, as a function of total amphiphilic concentration for two different hydrophobicity parameter $`\beta \chi =3.5,6`$ and salt concentrations $`\rho _s=10^3M,10^4M`$. At very low densities of divalent surfactant, the condensation is dominated by the monovalent ions, since the divalent amphiphilic molecules gain more entropy, and thus lower the total free energy, by staying free. However, as the concentration of divalent surfactant increases, the gain in electrostatic energy due to condensation and to the hydrophobicity of the amphiphilic molecules wins over the entropy.

We note that unlike the association with the ionic monovalent surfactants, which exhibits a large degree of cooperativity characterized by the sharp rise in the surfactant binding fraction, the replacement of the condensed monovalent counterions by the divalent surfactant proceeds more smoothly. This result could have been anticipated a priori. After the first amphiphile is associated, the condensation of additional molecules is energetically favored since the buildup of the hydrocarbon density in the vicinity of a polyion helps to exclude water and, thus, reduces the unfavorable hydrophobic energy of the alkyl tails. This effect competes with the electrostatic repulsion between the like-charged counterions that tends to inhibit it any further association. The repulsion is stronger for divalent ions than for monovalents and this explain why the cooperative effect is stronger in the later. In the case of the divalent amphiphilic molecules the surfoplex is formed for a density of surfactant around $`0.0005M`$ as illustrated in Fig. $`4`$ what is much lower than the one observed in the monovalent case -. For comparison with the case of monovalent surfactant illustrated in refs.  for producing Fig. $`4`$ we employ the hydrophobicity factor $`\beta \chi =4`$.

Fig. $`5`$ and Fig $`6`$ illustrates the effective charge of the complex. Due to the simplicity of our model that does not allow for the association of more than one surfactant molecule to each charged group along the DNA, there is no charge inversion. However, the decrease in the charge produced by the divalent surfactant is larger than the one observed when the amphiphilic molecules are monovalent for a similar model -.

In resume, we have presented a simple theory of DNA, for monovalent salt and divalent surfactant solutions. Our results should be of direct interest to researchers working on the design of improved gene delivery systems. In particular, we find that addition of cationic divalent surfactants leads to a cooperative binding. This binding happens far below the critical micell concentration and far below the concentration for this transition to happen in the presence of monovalent surfactant. Until now experimental attempts in employing surfactants as nonviral agent for transfection have been limited to the use of monovalent surfactants. The surfoplex produced with monovalent amphiphilic exhibit poor eficiency. Close to the membrane, most of the surfactants disassociate from the DNA and the transfection does not occurs . We propose the use of divalent cationic surfactants in the formation of the complex. Besides being more strongly connected to the charged groups of the DNA what will give more stability to the complex during the transfection, the amount of surfactant required is lower. Since the surfactant are toxic to the organism, this should reduce the risk of unnecessary medical complications.

Acknowledgments

This work was supported by the brazilian science agencies CNPq, FINEP, PRONEX and CTPETRO


text:
# More on “Atomic jumps in quasiperiodic Al72.6Ni10.5Co16.9 and related crystalline material”

## Abstract

We strongly disagree with a number of statements by Dolinšek et al. about phason dynamics in quasicrystals (QCs).

We are shocked by the methods used in the Reply of Dolinšek and Dubois  to a Comment  of ours to their paper, which unduly thwart our attempts of clarification.

(1) We would like to draw the attention of the reader to a whole series of repeated attempts to put two completely disparate phenomena under a common denominator (1a-c). They serve as a platform to qualify phason dynamics as not special.

(1a) Phason dynamics and regular atomic jumps in a B2-based phase would be the same physics since they are both leading to a broad (fast) quasielastic neutron scattering signal.

(1b) Phason dynamics and regular atomic jumps in a B2-based phase would be the same physics since phason sites in QCs and structural vacancies in B2-based phases would be both ”structural vacancies” according to an akward new definition, introduced ad hoc by the authors.

(1c) Phason dynamics and regular atomic jumps in a B2-based phase would be the same physics since both QCs and B2-based phases would be ”open structures”. In reality, there is no proof for the claim that a QC would be an ”open structure”. The authors “derive” it inductively from the wrong premises (1a)-(1b). But the fact that a number of sites are offered to an atom to jump to, does not forcedly imply that the structure is “open”, as the counterexample of interstitial sites clearly shows. (Furthermore, the fact that phason jumps in quasicrystals are assisted means that there is more to it than just the presence of a site available. Certain atoms might have to move temporarily out of the way to render the jump possible).

What characterizes if something is special, are the features that it has not in common with other things. Drawing a long list of indistinctive features shared with other things is not a valid procedure to prove that something is not special. Such an argument is based on an invalid inversion and/or negation of the logical implication.

The introduction of claim (1c) is a remarkable construction, that is based on the same leitmotiv of incorrect inversion and/or negation of the logical implication. It is easily recognized that claims (1a-b) are harsh and wrong. With claim (1c) the authors take us upwards an ad hoc deductive (implication) scheme, towards a less outspoken situation. It is much more difficult to discuss the claim (1c) and to prove conclusively that it is wrong on the ba-

sis of presently available data than claims (1a-b). The difficulty to come to terms with claim (1c) must permit doubts to linger on about the issues (1a)-(1b), by suggesting that they would follow rather logically from (1c), and that (1c) is plausible. In reality, as we mentioned, (1c) has been artificially postulated to yield (1a)-(1b): It is thus (1c) that has been “derived” from (1a)-(1b) rather than (1a)-(1b) from (1c). But we are not obliged to discuss the truth of postulate (1c), as it is totally obvious that (1a)-(1b) are wrong, independently of the issue if postulate (1c) is wrong or otherwise. It is also totally obvious that phason dynamics in a QC and regular atomic jumps in a B2-based phase are not the same physics, even if the signals they give rise to may look similar in certain selected aspects.

(2) In the repeated attempts to identify phason dynamics with conventional atomic hopping in a B2-based phase, each time contradicting evidence has been ignored.

(2a) Making use of SrTiO3 to depicture phason jumps in QC as nothing special is not appropriate because the corresponding high-temperature oxygen dynamics in SrTiO3 are not ultrafast as in QCs.

(2b) Similarly, in the B2-based phase, the jumps are indeed fast like in QCs, but here it is the geometry, temperature dependence, etc… that do not blend into the picture. The jump distances are not short as in QCs and instead the fast dynamics are a trivial consequence of a huge vacancy concentration, which has no counterpart in QCs.

In both cases it is tacitly denied that there is more to phason dynamics than just the one single aspect that is being picked. Phason dynamics are fast, there is a special geometry involved in the phason site that is an intrinsic part of the definition of the quasilattice itself, phasons have an unusual, special temperature dependence, there could exist correlated simultaneous jumps of several atoms, etc…

(3) The authors present the reader with a long, technical demonstration that the signal in AlCuNi is not magnetic. The motivation for this should be that I would have claimed that the NMR signal in AlNiCo is due to magnetism. The reader can verify in the text of my Comment that I never stated that the signal in AlNiCo is magnetic. I addressed AlPdMn which is well known to be magnetic for certain compositions, a fact the authors admit themselves. I clearly stated that I wanted to mention magnetism as a possibility in principle, in order to show that the interpretation of NMR data is not as straightforward as the authors try to make us believe (by not providing any discussion at all), and that the authors have already been misguided by their lack of rigor in the past, presenting a magnetic signal as evidence for phason dynamics.

(4) The purpose of my Comment was to indicate that the interpretation of the authors is wrong and that the authors have not bothered about eliminating alternative, more reasonable interpretations of their data. If my preference goes to an interpretation of their data in the form of some small non-phasonic atomic displacements, due to the presence of lattice distortions induced by structural defects, I pointed out that this does not mean that I claim any interpretation at all. The burden of the (perhaps presently inaccessible) correct interpretation and its proof remains entirely with the authors, and one cannot reverse this situation by picking one example in my objection that there are many other possible interpretations in principle, in order to mispresent my Comment as a wrong interpretation of the data. That the authors go into great efforts to explain that the signal in AlNiCo is not magnetic is perhaps fine, but not really needed. It distracts the attention of the reader from the essence of the discussion and it is already well established that AlNiCo is not magnetic.

(5) As the authors now admit, the attribution of the NMR signal in reference  was wrong. At the time I pointed out in a Comment  that the temperature behavior of their signal could not possibly agree with their interpretation in terms of phason dynamics. The Reply  thwarted the recognition of the fact that that conclusion was final. The authors postponed admitting that the signal was due to magnetism to another occasion. Nonetheless, in reference  the authors claimed again that the signal observed in reference  is a valid observation of phason dynamics and the Kalugin-Katz transition. And in the present reply they admit once again that this is wrong and that the AlPdMn data have been overinterpreted.

(6) We would like to draw the reader’s attention to the use the authors make of a paper by Yan et al.  who erroneously called phason sites vacancies. This was pointed out in a Comment  but clarification was thwarted by a statement in the Reply by Yan et al.  that they used a different definition of the compound word ”columnar vacancy”, in order to escape from the criticism that they introduced a different definition of the single word ”vacancy”. This non-conventional definition had not been given in the original paper and was just introduced a posteriori in the Reply. Such ad hoc methods are not viable and can only produce confusion. Dolinsek and Dubois built on that confusion twice: (a) They claimed that AlNiCo contains a lot of ”vacant sites” citing the paper by Yan et al., ignoring my Comment on it. (b) Now I have pointed out once again that phason sites are not vacancies, they take the same loophole of escape as Yan et al. and state that they use a different definition for the compound word ”structural vacancy”. Again, the use of a different definition is explicited a posteriori, and at variance with normal, legitimate expectations. The scientific convenience of such a definition has been discussed under (1).

(7) The authors say that they have been cautious by stating that the data are compatible with phason dynamics. This is perhaps not a direct claim, but nevertheless far too much of a claim. The caution does not serve to warn the reader but to thwart due criticism for trying to pass an unwarranted claim. Real caution would consist in stating that the data are also compatible with a plethora of (more reasonable) explanations. And there are enough other arguments in my Comment that show that the data do not correspond to phason dynamics.

(8) When I stated it is not proved that tunneling states are not phasons, I thought that perhaps, it might never be possible to prove that they are not phasons, due to the limitations of the techniques. This certainly does not imply that it would be reasonable to claim that tunneling states are phasons: There is more in being right than not being proved wrong. In fact, the authors do not mention that there is even less proof that tunneling states would be phasons. No argument is given to substantiate their claim, except the experimental difficulties one encounters to provide the absolute final proof of the contrary. The reader will certainly be surprized to learn that in reference  the authors state themselves that “tunneling are definitely not an intrinsic property of quasicrystals”, contradicting their own claims!

Let us add to that that the double wells in phason dynmaics are defined on Euclidean three-dimensional space, while those in tunneling states are defined on an abstract configuration space. The new interpretation proposed by Dolinsek and Dubois rests thus on a confusion between two very different realities covered by the same words ”double well”. We pointed out that phason dynamics in QCs explore a discrete set of jump vectors, in contrast to what is considered to happen in tunneling states. That was already an excellent indication that tunneling states are not phasons.

(9) That it is not possible to link the high-temperature neutron signals to their NMR data, cannot justify the high-profile claims of the authors. They try to save their interpretation of the NMR signals in an amended version by claiming that there could be a cross-over of dynamical regime from hopping to tunneling. The surprizing character of this claim should not disarm or confuse the reader: It cannot make away with our arguments based on the X-ray data. We pointed out that in AlNiCu there is clear evidence from X-ray diffraction for an order-disorder transition. Below the transition temperature, the jump dynamics are just frozen. X-ray diffraction would not give evidence for an order-disorder transition if the atomic motion just crossed over from over-the-barrier hopping to tunneling, because it would still detect the instantaneous positional disorder, created by the atomic motion.

Furthermore, it must be clear that the distances in AlNiCu are too large (and the transition temperature too high) to consider tunneling of the type proposed by the authors as a real possibility. A comparison with tunneling in SrTiO3 is thus totally inappropriate.

On the same footing, it is generally admitted that the freezing of phason hopping in QCs is responsible for the increase of diffuse scattering observed below a certain temperature (which is e.g. $``$ 600 C in AlPdMn). If there were a cross-over from phason hopping to phason tunneling, there would be no change in the X-ray data below the cross-over temperature.

(10) I already pointed out that in AlPdMn the assistance energies are not close to 0.6 eV. Until much more serious proof is given caution requires thus to consider the numerical similarities as coincidences.

(11) The authors deny me the use of the word ”phason”, although I come first in matters of having spelled out numerous times the limitations of this terminology. Moreover, the latter criticism could be extended to the whole literature on QCs - often on much more serious grounds  \- rather than just focusing it onto the subject of my work on phason dynamics, the significance of which extends far beyond a mere matter of terminology.

(12) The farfetched, invalid criticisms of this group on my work  make a long, damaging, disparate list. The motivations given for the present attacks contradict each other. The authors state to have selected on purpose a sample with an exceptional amount of structural vacancies to permit them to qualify phason dynamics as “not special”  or “not QC-specific”. But previously they have stated to be “surprized” as they “did not intend to question the excellent work on phason dynamics”. It is thus just that all this excellent work has been witlessly misaimed at forcing an open door. The tacit implications of this declaration of intention of the authors carry strong personal undertones. The same can be said about the gratuitous, subjective rating that phason dynamics are not special. But both issues can even not be touched upon without attracting further personal discredit.


text:
# SDSSJ115517.35+634622.0: A Newly Discovered Gravitationally Lensed Quasar

## 1 Introduction

Under certain conditions, the gravitational potential of a foreground galaxy can produce multiple images of a background quasar. These rare events are of considerable interest due their utility in a number of astrophysical investigations, including modeling the mass distribution of the lens galaxy (Kochanek, 1995), measuring the Hubble constant through time delay measurements (Refsdal, 1964; Kundic et al., 1997; Koopmans & Fassnacht, 1999), and constraining cosmological models through the use of lensing statistics (Turner, 1990; Fukugita et al., 1990; Chae et al., 2002). Many of these applications would benefit from a sample of gravitational lenses defined by a clearly understood selection procedure, so as to allow for the correction of selection effects.

The Sloan Digital Sky Survey (York et al., 2000, henceforth SDSS) dataset provides an excellent opportunity to discover gravitationally lensed quasars. The $`10^4`$ deg2 of five-band photometry and the spectroscopic sample of $`10^5`$ quasars should allow for the compilation of a lensed quasar sample which is both large and statistically well-defined in comparison to existing optical samples. For further discussion on this subject, the reader is referred to Pindor et al. (2003, henceforth P03). In this work, we report the discovery of a previously unknown gravitationally lensed quasar, found during our ongoing search for such objects in the SDSS dataset.

## 2 Observations

### 2.1 Selection as a Lens Candidate

Five-band photometry of SDSS J115517.35+634622.0 (henceforth SDSS J1155+6346) was obtained in the course of normal SDSS imaging111The SDSS photometric designation is (run/rerun/camCol/field/id) 2078/22/1/84/169 on 2001 January 26. The object was classified as a galaxy by the automated star-galaxy separator and is bright enough ($`r_{\mathrm{PETRO}}=16.8`$) that it was targeted for spectroscopy by SDSS galaxy target selection (Strauss et al., 2002). SDSS spectroscopy222The SDSS spectroscopic designation is (mjd/plate/fibre) 52316/598/87 on 2002 February 11 instead revealed the object to be a quasar at a redshift of $`z=2.89`$. Having been identified as a quasar, the object was initially selected as a lens candidate by a method similar to that used to select the first two SDSS-discovered gravitational lenses (Inada et al. 2003a, Inada et al. 2003b) as well as for having apparent Mg II absorption at a redshift of $`z=0.72`$. It finally received follow-up spectroscopy when observed as part of a sample which passed the selection criteria described in P03.

Briefly, both of the aforementioned selection algorithms work by searching among spectroscopically confirmed quasars for objects that are extended333We use the word “extended” to mean “not consistent with being a point source” and the word “resolved” to mean “displaying two distinct photometric peaks” and that can be decomposed into two components of similar colours. The algorithm of Inada et al. uses photometric model likelihoods measured by the SDSS photometric pipeline (Lupton et al., 2001, henceforth PHOTO) to identify extended objects and decomposes objects using flux moments. P03 uses photometric models consisting of one and two point spread functions (henceforth PSF) to identify extended objects and to decompose objects. Our two PSF photometric model of the SDSS atlas image predicted a component separation of $`1\stackrel{}{\mathrm{.}}87`$, a flux ratio of A/B = 5/4, and a position angle of $`87^{}`$.

The SDSS is a photometric and spectroscopic survey across 10,000 square degrees of the northern Galactic cap using the 2.5m SDSS telescope at Apache Point Observatory. SDSS imaging is carried out with a wide-field camera (Gunn et al., 1998) which makes nearly simultaneous observations of objects in five passbands: $`u`$ $`g`$ $`r`$ $`i`$ $`z`$. Together, the passbands cover the optical wavelengths from the atmospheric cut-off in the blue to the minimum detectable energy for the silicon CCDs in the red (Fukugita et al., 1996). Photometric calibration of the imaging survey is separately carried out by an automated 0.5m telescope which monitors a set of standard stars (Smith et al., 2002) while photometric data is being acquired(Hogg et al., 2001). The SDSS imaging camera also incorporates astrometric CCDs which provide astrometry of detected objects with an accuracy typically better than 0$`\stackrel{}{\mathrm{.}}1`$ (Pier et al., 2003). SDSS spectroscopy is carried out on the same telescope by two fiber-fed double spectrographs which produce spectra with a resolution ($`\lambda /\mathrm{\Delta }\lambda `$) of $``$ 2000 covering the wavelength range 3800–9200 Å. Together, the spectrographs have 640 fibers which are assigned based on previous SDSS imaging through an efficient tiling algorithm (Blanton et al., 2003). For more comprehensive documentation of the survey, readers are encouraged to consult Stoughton et al. (2002) and Abazajian et al. (2003).

### 2.2 Spectroscopic Observations

We observed the object on 2002 May 4, using the 6.5m Multiple Mirror Telescope (MMT). Spectra were obtained using the blue channel of the MMT spectrograph with the 300/mm grating and $`1^{\prime \prime }`$ slit. The slit was aligned along the separation axis of the pair, as predicted from SDSS imaging by a photometric model consisting of two point sources (see P03), so that the two components were observed simultaneously. The spectrum was obtained with a single 300s integration.

Although the seeing at the time of observation was slightly sub-arcsecond, the two objects were only marginally resolved in the focal plane due their small separation ($`1\stackrel{}{\mathrm{.}}8`$). Consequently, we implemented a deblending procedure to minimize cross-contamination between the extracted 1D spectra. We binned the 2D spectrum into 20 bins along the dispersion axis and integrated the flux in these bins to produce a series of spatial cross-sections. We then simultaneously fit a profile consisting of two Gaussians to each of these cross-sections, subject to the constraint that the separation of the two Gaussians be the same in each bin. Finally, we defined two apertures by assigning a fraction of the flux from each pixel to the first (second) aperture corresponding to the fractional contribution of the first (second) Gaussian to the model profile at that pixel. Hence, the model profiles can be thought of as relative weights which determine how much of the flux in a given pixel should be assigned to either of the two apertures. This deblending procedure is conceptually based upon the deblending algorithms used by PHOTO. Having been thus deblended, the spectra were then extracted using standard IRAF routines. Figure 1 shows the spectra of the two components. Note that the feature at $``$ 7600 Å is atmospheric. We cross-correlated the spectra and estimated the velocity difference between them to be $`100\pm 300`$ km s-1, consistent with zero.

We further examined the spectra for evidence of emission from a lens galaxy. The SDSS spectroscopic pipeline identifies in the SDSS spectrum a cross-correlation match to a series of galactic absorption features at a redshift of $`z=0.1756`$. In fact, this match is reported with a confidence level (99.73%) only marginally lower than the match to the quasar redshift (99.75%). This redshift identification is slightly suspect due to the fact that the purported calcium $`H`$ and $`K`$ lines appear in the Lyman-alpha forest. However, as shown in Figure 2, a number of the absorption features are clearly visible in both the SDSS and MMT spectra. Hence, we report a redshift measurement of $`z=0.1756`$ for the possible lens galaxy.

We also examined the absorption doublet seen at $``$ 4810 Å in both components. This doublet was identified in the SDSS spectrum as being Mg II in absorption at a redshift of $`z=0.72`$. Its location relative to the broad N V emission line might suggest that it is the N V doublet in absorption in the vicinity of the quasar. In Figure 3 we show the best-fit locations of the Mg II and N V doublet systems overlayed on the SDSS spectrum, which, although spatially unresolved, has better spectral resolution and signal-to-noise than our MMT spectrum. We conclude that we are in fact observing Mg II absorption at a much lower redshift. It is possible that this absorber is associated with a mass concentration which contributes to the lensing properties of the system, but it is impossible to test this hypothesis with the existing data.

We consider the spectra alone to be inconclusive with respect to whether this system was a lensed or binary quasar. The lensing hypothesis is supported by the fact that the spectra are similar in overall shape and have similarly peaked emission lines. Spectroscopic identification of a possible lens galaxy clearly also supports the lensing hypothesis. On the other hand, the equivalent widths of the emission lines in the fainter B component are about twice the widths in the brighter A component. This seems inconsistent with the spectra being two images of the same quasar. However, in the next section, we present further evidence that reconciles the observed spectral differences with the lensing hypothesis.

### 2.3 Near-IR Imaging

Having determined that SDSS J1155+6346 did indeed have two separate quasar images at the same redshift and with similar SED, we next obtained a near-IR image. On 2002 August 8, we imaged the object in the $`K^{}`$-band using the Near Infrared Camera (NIRC, Matthews & Soifer 1994 ) on the Keck I telescope. The observations consisted of a five-point dither pattern, integrating for 10s at each pointing. The data were flattened, sky subtracted, shifted, and stacked using the DIMSUM package in IRAF.

Figure 4 shows the $`K^{}`$-band of SDSS J1155+6346. It is evident from visually inspecting this image that component A is extended. It is our interpretation of this image that SDSS J1155+6346 is indeed a gravitational lens, and that the lensing galaxy is almost co-incident with the quasar image we have labelled A. To support this interpretation, we first constructed a model of the system consisting of only two point sources. We used an analytic PSF consisting of two Moffat functions, as recommended by Racine (1996). The half-width of our PSF model was determined by fitting to a nearby star, labelled P. Figure 5a shows the residuals of the best-fit two point source model. The peak residuals are 32% of the peak intensity. The highly negative (white) residuals near the center of component A are produced as the model attempts to account for the extended profile by over-estimating the central flux. Since the images are background-dominated, the measurement errors can be inferred from the observed sky variance, and our two point-source model has a reduced chi-square of 442. We subsequently constructed a model of the system as a combination of two point sources plus an extended component represented by a deVaucouleurs profile of the form:

$`I(x,y)=I_0\mathrm{exp}(7.67((x^2+(y/q)^2)^{1/2}/R_{\mathrm{eff}})^{1/4})`$ (1)
convolved with the PSF, where $`q`$ is the axis ratio ($`b/a`$). In general, our model should also allow the position angle of the galaxy on the sky to vary, but we found that the major and minor axes happen to be well-aligned with the pixel grid. Our best-fit model has a reduced chi-square of 6 (with 11014 degrees of freedom) and the peak residuals are 7% of the peak intensity. We interpret the obvious improvement observed for our two point sources plus deVaucouleurs model as definitively demonstrating the presence of a lens galaxy. Figure 4b shows the residuals left after subtracting our model from the data. The structure of the residuals near component B indicate that the inaccuracy of our PSF model probably dominates the uncertainties in the model. For this reason, estimation of the errors in the model parameters is not straightforward. We estimated the random errors through bootstrap re-sampling, but this does not adequately represent the systematic errors. Hence, we considered two further photometric models; one consisting of two analytic PSFs (as before) plus an exponential disc, and one consisting of two empirical PSFs, directly using the image of P as a PSF model, plus a deVaucouleurs profile (as before). The best-fit exponential disc model has a reduced chi-square of 21, and it differs from our initial model mainly in that it predicts roughly equal flux in the two quasar images and substantially less total flux in the lens galaxy. Figure 5c shows the residuals left after subtracting this model from the data. The substantially worse chi-square value indicates that the lens is an early-type galaxy, as is in fact indicated by visual inspection of the image. The best-fit empirical PSF model actually has a lower reduced chi-square value, 4.5, than our initial model. However, using the image of P directly requires by-hand masking of a fainter nearby source, as well as of the local gradient produced by the extended profile of the lens galaxy. For this reason, we retain our initial model as the favoured photometric model of the system. Figure 5d shows the residuals left after subtracting the empirical PSF model from the data. The main difference between the results of the empirical PSF model and our initial model is that the distance between the center of the lens galaxy and quasar image A \[henceforth $`\mathrm{\Delta }\theta _{AG}`$\] is about twice as great in the former as in the latter. Hence, the main systematic uncertainty appears to be in the relative position of the lens galaxy, and the extent of this uncertainty is indeed underestimated by our bootstrap error estimation. Our favoured model predicts $`\mathrm{\Delta }\theta _{AG}=0\stackrel{}{\mathrm{.}}15`$, but a more robust estimate can be obtained from average of $`\mathrm{\Delta }\theta _{AG}`$ as predicted by the analytic and empirical PSF models, whereby $`\mathrm{\Delta }\theta _{AG}=0\stackrel{}{\mathrm{.}}22\pm 0\stackrel{}{\mathrm{.}}07`$ (systematic). The fairly large reduced chi-square of our best-fit model indicates that further information (i.e. more reliable errors) could be obtained from the image if a more satisfactory PSF model could be found. Table 1 summarizes the best-fit parameters for all three models. The near coincidence of the lens galaxy and quasar A can lead to some confusion so, for clarity, let us explicitly state: when refering to, for instance in Figure 1, component A, we mean the unresolved combination of quasar A and the lens galaxy, which is brighter than the resolved image of quasar B. However, the results of our photometric modeling indicate that quasar image A itself is fainter than quasar image B.

In order to further verify the results of our photometric model, we constructed a difference spectrum by multiplying the spectrum of component A by the flux ratio predicted by our model and then subtracting the result from component B. Figure 6 shows this difference spectrum. The difference in C IV line appears to be consistent with zero, but significant residuals remain in the vicinity of Lyman-alpha. However, given that our sightline to quasar image A obviously passes quite close the center of the lens galaxy, it is not unreasonable to expect some difference in the quasar spectra due to absorption and/or micro-lensing. Further, if our identification of the lens galaxy redshift, $`z=0.1756`$, is correct, then the residuals are quite possibly due to the 4000 Å break which we would expect to observe in the vicinity of $``$ 4700 Å. We interpret the difference spectrum as broadly confirming the results of our photometric modelling, but as not having sufficient signal-to-noise to contribute to estimation of the lens model parameters.

We also matched SDSS J1155+6346 with the 2MASS Second Incremental Data Release Point Source Catalog and found it to have a $`K`$-band magnitude of 14.6$`\pm 0.1`$. We can use the 2MASS measurement to calibrate the total magnitude of the system, so that our photometric model of the NIRC image predicts $`K`$ magnitudes444We neglect the conversion between $`K`$ and $`K`$’ as it is a negligible correction relative to our photometric accuracy. of 14.9$`\pm 0.1`$, 16.6$`\pm 0.1`$ and 17.4$`\pm 0.1`$ for the galaxy, quasar B, and quasar A, respectively.

A measurement of the lens galaxy magnitude allows us to estimate its redshift, using the Faber-Jackson (1976) relation for lensing as presented by Rusin et al. (2003). We re-write Rusin et al.’s equation (A1) in the form

$`m_{obs}M_0=DM+2.5\gamma _{E+K}z_d1.25\gamma _{FJ}\mathrm{log}\mathrm{\Delta }\theta _{red}`$ (2)
where $`DM`$ is the distance modulus, $`\gamma _{E+K}`$ is a parameter subsuming evolution and spectral K-corrections, $`\gamma _{FJ}`$ is the slope of the Faber-Jackson relation in the chosen band, $`z_d`$ is the lens galaxy (deflector) redshift, and $`\mathrm{\Delta }\theta _{red}`$ is the reduced image separation, $`\mathrm{\Delta }\theta _{red}(\mathrm{\Delta }\theta /\mathrm{\Delta }\theta _{})(D_s/D_{ds})`$, where $`\mathrm{\Delta }\theta _{}`$ is the image separation produced by a singular isothermal sphere with a velocity dispersion of 225 km s-1, and $`D_s`$ and $`D_{ds}`$ are the angular diameter distances between the observer and source and deflector and source, respectively. The $`K`$-band values of $`M_0(=24.4+5\mathrm{l}\mathrm{o}\mathrm{g}h)`$, $`\gamma _{E+K}(=0.31)`$, and $`\gamma _{FJ}(=3.40)`$ are all given by Rusin et al., based upon HST photometry and spectrophotometric modelling of 29 early-type lens galaxies. We can then numerically solve Equation (2) to get $`z_{galaxy}=0.17\pm 0.03`$. Hence, our photometric redshift estimate is in excellent agreement with the tentative spectroscopic redshift presented in §2.2. If we accept the galaxy redshift to be $``$ 0.17, then an SIS model predicts a velocity dispersion of $``$ 190 km s-1 in order to reproduce the observed image separation. However, an SIS model also predicts a B/A flux ratio of $``$ 10:1 (for our favoured model, $``$ 5:1 for the galaxy position in the empirical PSF model), rather than the observed $``$ 2:1. Possible explanations for this discrepancy include some combination of a large internal or external shear, microlensing, or a significant softening of the inner slope of the mass profile. The discrepancy could also be explained if the distance between the galaxy and quasar image A is significantly greater than our models predict.

## 3 Conclusions

In the course of an ongoing search for strongly lensed quasars in the dataset produced by the SDSS, we have discovered an object which is almost certainly a previously unknown gravitationally lensed quasar. High angular resolution spectroscopy confirms two quasar images with an angular separation of $`1\stackrel{}{\mathrm{.}}8`$ having identical redshifts of $`z=2.89`$. We also identify the presence of absorption features corresponding to a lens galaxy at a redshift of $`z=0.1756`$. High resolution near-IR imaging clearly indicates the presence of a luminous, early-type galaxy between the quasar images. Our photometric model of the system implies an angular separation of $`1\stackrel{}{\mathrm{.}}832\pm 0.007`$ between the quasar images. We used 2MASS photometry to measure the lens galaxy magnitude as being $`K=14.9`$. This magnitude allows for a photometric estimate of the lens galaxy redshift as being $`z=0.17\pm 0.03`$, based on the Faber-Jackson relation for lensing. One of the quasar images is only $`0\stackrel{}{\mathrm{.}}22\pm 0\stackrel{}{\mathrm{.}}07`$ ($`0.1R_{\mathrm{eff}}`$) from the center of the lens galaxy. Hence, the optical depth to quasar micro-lensing (Wambsganss, 2001) for this image may be by quite high. Similarily, higher resolution multi-band imaging of the system could provide an interesting measurement of the differential reddening (Falco et al., 1999) produced by the central region of the prominent lens galaxy.

B.P. is supported by NASA grant NAG5-9274. D.J.E. is supported by National Science Foundation grant AST-0098577 and by an Alfred P. Sloan Research Fellowship. Part of the work reported here was done at the Institute of Geophysics and Planetary Physics, under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under contract No. W-7405-Eng-48. G.T.R. is supported by HST grant HST-GO-09472.01-A.

A portion of the results reported here made use of the Multiple Mirror Telescope Observatory, a facility operated jointly by the University of Arizona and the Smithsonian Institution.

Additional observations reported here were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California, and the National Aeronautics and Space Administration. The observatory was made possible by the generous financial support of the W. M. Keck Foundation.

This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation.

Funding for the creation and distribution of the SDSS Archive has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Aeronautics and Space Administration, the National Science Foundation, the U.S. Department of Energy, the Japanese Monbukagakusho, and the Max Planck Society. The SDSS Web site is http://www.sdss.org/.

The SDSS is managed by the Astrophysical Research Consortium (ARC) for the Participating Institutions. The Participating Institutions are The University of Chicago, Fermilab, the Institute for Advanced Study, the Japan Participation Group, The Johns Hopkins University, Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, University of Pittsburgh, Princeton University, the United States Naval Observatory, and the University of Washington.


text:
# ACS Observations of a Strongly Lensed Arc in a Field Elliptical

## 1 Introduction

The past decade has witnessed a dramatic growth in both the number of known gravitational lenses, and the astrophysical and cosmological applications for which they have been used. Strong lensing can provide valuable information about the detailed mass structure of intermediate-redshift galaxies (e.g., Rusin, Kochanek, & Keeton 2003). The frequency of lensing, as well as the time delays and image positions as a function of redshift, can be used to constrain the cosmological parameters (e.g., Im et al. 1997, Chae et al. 2002). Of course, lensing also provides the opportunity to study the properties of distant objects with the aid of strong magnification (e.g., Franx et al. 1997).

With its unsurpassed resolution at optical wavelengths, the Hubble Space Telescope (HST), has played a central role in the study and exploitation of many lensing systems. The installation of the Advanced Camera for Surveys (ACS) (Ford et al. 2002), with its higher resolution, greater sensitivity, and wider field, has made HST a yet more powerful tool for gravitational lensing studies. Here, we present our analysis of a strongly lensed arc serendipitously discovered around a field elliptical using ACS. We adopt a cosmology with $`(h,\mathrm{\Omega }_m,\mathrm{\Omega }_\mathrm{\Lambda })=(0.7,0.3,0.7)`$ throughout.

## 2 Observations and Data Reductions

The object studied in this paper was found in a parallel field observed with the ACS Wide Field Camera (WFC) during NICMOS observations of the HST Ultra-Deep Field (GO program 9803, PI: Thompson) in 2003 September. It was observed for 9,9,18,27 orbits in the F435W,F606W, F775W,F850LP bandpasses, respectively, hereafter referred to as $`B_{435}`$, $`V_{606}`$, $`i_{775}`$, and $`z_{850}`$. The pointing pattern was a grid of nine positions, spaced by about one NICMOS field, or $``$30″. These data are public and were retrieved from the STScI archive where they received the standard CALACS on-the-fly reprocessing (OTFR) to the point of flatfielded images in units of electrons. Further processing, including cosmic ray rejection and final image combination, was performed with the “Apsis” data reduction software described in detail by Blakeslee et al. (2003a), with refinements as noted in Blakeslee et al. (2003b). Apsis now also performs automatic astrometric recalibration of the images; the resulting zero-point error is $``$0$`\stackrel{}{\mathrm{.}}`$1, dominated by systematic effects in the astrometric catalogues.

Two of the $`z_{850}`$ exposures taken during one of the orbits were adversely affected by a bright star that landed precisely at the edge of one of the CCD chips, resulting in charge bleeding into the overscan area used in the STScI CALACS processing. This led to a large erroneous slope in the fit to the overscan, making OTFR versions of these images unusable. We chose simply to omit these two images, which reduced the final $`z_{850}`$ signal-to-noise by less than 2%. The final exposure times in our summed images were then 20700 sec ($`B_{435}`$), 20679 sec ($`V_{606}`$), 41400 sec ($`i_{775}`$), and 59800 sec ($`z_{850}`$). We use photometric AB magnitude zero points in these four bands of 25.673, 26.486, 25.654, and 24.862 (Sirianni et al. 2003, in preparation) and correct for $`E(BV)`$ = 0.008 mag of extinction (Schlegel et al. 1998). Finally, we obtained follow-up spectroscopy in 2003 November with the Magellan Clay 6.5-m telescope using the LDSS-2 at a dispersion of 5.3 Å pix-1. The total exposure time was 7200 sec, and the data were reduced using standard techniques.

## 3 Analysis and Results

Figure 1 shows a composite color image constructed from a 0$`\stackrel{}{\mathrm{.}}`$8 portion of the $`B_{435}`$, $`V_{606}`$, $`z_{850}`$ images. A dramatic blue arc, subtending $``$120 is readily visible at a radius of 1$`\stackrel{}{\mathrm{.}}`$63 perpendicular to the major axis of a red early-type galaxy. An apparent counter-arc of the same color at a radius of 0$`\stackrel{}{\mathrm{.}}`$6 can be seen in the inset image. A number of other blue, tangentially oriented features (e.g., a faint arc to the south of the elliptical) may also be lensed. Several nearby galaxies have colors and/or photometric redshifts (see below) similar to the lensing galaxy and likely make up a group. Keeton et al. (2000) predict that 20–25% of galaxy lenses will be in groups massive enough to perturb the lensing potential. We note that this field has also been observed with Chandra (Rosati et al. 2002), and we find an upper limit to the X-ray flux of $`1\times \mathrm{\hspace{0.17em}10}^{15}`$ erg cm$`^2`$s-1 in the 0.5–2 keV band, indicating $`L_X<2\times \mathrm{\hspace{0.17em}10}^{42}`$ erg s-1 at the measured $`z=\mathrm{\hspace{0.17em}0.62}`$. The surrounding group must then be low-mass, $`2\times \mathrm{\hspace{0.17em}10}^{13}M_{}`$, assuming the $`L_XM`$ relation from Reiprich & Böhringer (2002), although it could still have a significant effect on the lensing.

For accurate photometric and structural measurements of the arc, it is necessary to first model and subtract the light from the lensing galaxy. We used custom 2-d isophote-fitting software, masking the arcs and other nearby sources before fitting. The fitted models resulted in small residuals and indicate an effective radius $`R_e1\stackrel{}{\mathrm{.}}6`$. Figure 2 illustrates the detailed knot structure that becomes apparent in the arcs following galaxy subtraction; the knot positions are given in the figure caption. A combination of the $`B_{435}`$ and $`V_{606}`$ images is used for this since the arcs are so blue relative to the galaxy. The main arc is brighter than the counter-arc by a factor of $`7.9\pm 0.5`$, depending on whether or not the concentration near the south end of the main arc is counted. We tabulate our key measurements on this lensing system in Table 1. The galaxy redshift was measured from the LDSS-2 spectra using template cross-correlation, and the arc photometric redshift was estimated using “BPZ” (Benítez 2000) with revised templates from Benítez et al. (2004).

If we estimate the Einstein radius $`\theta _E`$ of the system to be the projected distance of the primary arc and assume a singular isothermal sphere model (e.g., Schneider, et al. 1992) with the redshifts from Table 1, then we derive a lens velocity dispersion $`\sigma 305`$ km s-1. This changes by only $`\pm 3\%`$ if the source redshift is varied by $`0.4`$. Interestingly, we predict a very similar $`\sigma `$ from a fundamental plane (FP) analysis of the elliptical galaxy using the data in Table 1. We $`K`$-corrected our $`z_{850}`$ photometry to the SDSS $`r`$ band and derived a value for $`\sigma `$ using the coefficients and methods (including luminosity evolution of $``$0.85$`z`$ in the $`K`$-corrected $`r`$ band and conversion of the fitted $`R_e`$ to an effective circular radius) from the very large field elliptical study by Bernardi et al. (2003). The result is $`\sigma =295\pm 45`$ km s-1, where the error is based on the FP scatter.

We wish to investigate whether or not an elliptical galaxy of the observed luminosity can naturally produce the observed lensed features, including the counter-arc. We construct a model using a round source galaxy and the singular isothermal form of the elliptical effective lensing potential of Blandford & Kochanek (1987), parameterized by the velocity dispersion $`\sigma `$ and an ellipticity $`ϵ_p`$. We also include an external shear component, which may arise from a host group, the bright foreground galaxies seen in Figure 1, or misalignment between the luminous galaxy and dark matter halo (Keeton et al. 1997).

Because the system is under-constrained, we fix some model parameters to reasonable values. In particular, we set $`\sigma =295`$ km s-1 as expected from the FP and $`ϵ_p=0.12`$, or a third of the ellipticity of the light (Mellier et al. 1993). We then allow the position, orientation, and shear components to vary. We used an optimization procedure that minimizes the deviation of the delensed positions in the source plane while optimizing the gross coincidence of the model and data in the image plane. We did not include magnification weighting or detailed constraints from the internal structure of the arcs. Figure 3 shows that the model successfully reproduces the observed geometry of the bright arc and counter-arc using a single round source. The major axis of the mass ellipsoid is aligned to within 20° of the galaxy light, while the external shear axis is within 10°. We defer more detailed modeling of the lens mass distribution and internal structure of the arcs to a future work when we have better information on the source redshift and lens velocity dispersion.

## 4 Discussion and Summary

Early-type lensing galaxies constitute a small but unique mass-selected sample in generally low-density environments (Kochanek et al. 2000), lacking the biases inherent in luminosity selection. Kochanek et al. (2000) concluded that these galaxies have old stellar populations with formation redshifts $`z_f2`$ and lie on the same FP as cluster ellipticals of the same redshift. However, van de Ven et al. (2003) have analyzed the FP of 26 early-type lenses and find a range in colors and ages, with formation redshifts as low as $``$1 and slightly younger ages on average than cluster ellipticals. The new lens studied here adds to the mass-selected sample of field ellipticals, and the multi-band coverage allows for some constraint on age. The observed $`(V_{606}i_{775})`$ and $`(i_{775}z_{850})`$ colors transform fairly closely to rest-frame $`(UB)_z`$ and $`(BV)_z`$, respectively. We find $`(UB)_z=0.45`$, $`(BV)_z=0.85`$, with $`\pm `$0.03 mag uncertainty in the transformations. This implies an age $`\tau _L=4.4\pm 0.7`$ Gyr from the latest solar metallicity Bruzual & Charlot (2003) models, or formation at $`z_f2.0\pm 0.5`$. Thus, the population is clearly evolved, though marginally younger than derived for cluster ellipticals (e.g., van Dokkum & Franx 2001; Blakeslee et al. 2003b).

The excellent angular resolution and wide field of ACS on HST now makes it possible to conduct systematic optical searches for arcsecond-scale strongly lensed objects and obtain statistically meaningful results. The two Hubble Deep Fields observed with WFPC2 together produced only one probable lens candidate (Barkana et al. 1999), while the medium-deep survey (MDS) found 10 good candidates over about 130 WFPC2 fields (Ratnatunga, Griffiths & Ostrander 1999). However, the recent paper by Fassnacht et al. (2003) lists 6 probable lens candidates in ACS GOODS fields covering only about a quarter of the MDS area. The dramatic instance of lensing we report here lies in a field adjacent to the GOODS area and brings the number of likely galaxy lenses found with ACS to 7, or about 40% of the total found with HST in less than 15% of its lifetime. It seems likely that after a few more years, searches in ACS fields will provide the first robust optical determination of the rate of strong gravitational lensing by individual galaxies, and thus their mass function.

We have found and provided a first analysis of a clear case of strong gravitational lensing by a field, or small group, elliptical galaxy. We estimate $`z2.4`$ for the arc using our four-band ACS photometry and have spectroscopically measured $`z=0.62`$ for the lens galaxy. Several smaller galaxies with similar colors may form a group with the lens galaxy, although the Chandra X-ray flux limit indicates this must be fairly low-mass. We have successfully modeled the arc/counter-arc system as a single round source being lensed by an isothermal ellipsoid of velocity dispersion $`\sigma =295`$ km s-1, the value indicated by the fundamental plane, with a modest amount of external shear that could result from neighboring or foreground galaxies. A more thorough analysis including detailed constraints on the mass distribution and source morphology would benefit from measurements of the arc redshift and lens galaxy velocity dispersion. We expect these will be forthcoming, as will discoveries and analyses of many additional lensing systems with ACS.

ACS was developed under NASA contract NAS 5-32864, and this research has been supported by NASA grant NAG5-7697. L.I., V.M., and G.H. are supported by a Chilean Conicyt Fondap grant “Center for Astrophysics.” We thank our fellow ACS Team members and support staff, especially Dave Golimowski for retrieving the data and Jon McCann for writing the software used in making Figure 1.