format: text

text:
# Statistical Properties of Collisionless Equal- and Unequal-Mass Merger Remnants of Disk Galaxies

## 1 Introduction

Detailed observations of individual giant elliptical galaxies have shown that they can be subdivided into two groups with respect to their structural properties (Bender, Döbereiner & Möllenhoff, 1988; Bender, 1988; Kormendy & Bender, 1996) \[BDM, and references therein\]. Faint giant ellipticals are isotropic rotators with small minor axis rotation and disky deviations of their isophotal contours from perfect ellipses. They might contain faint disks which contribute up to 30% to the total light in the galaxy. Therefore their disk-to-bulge ratios overlap with those of S0-galaxies (Rix & White, 1990; Scorza & Bender, 1995). Disky ellipticals have power-law inner density profiles (Lauer et al., 1995; Faber et al., 1997) and show little or no radio and X-ray emission (Bender et al., 1989). Boxy ellipticals, on the other hand, are in general more massive than disky ellipticals. They have box-shaped isophotes, and show flat cores. The kinematics of boxy ellipticals is generally more complex than that of disky ellipticals. They rotate slowly, are supposed to be supported by velocity anisotropy and have a large amount of minor axis rotation. A number of rotationally supported ellipticals also show boxy isophotes. These systems are either purely boxy (this is supposed to originate from tidal interactions with nearby massive companions) or show disky isophotes in the inner part and boxy isophpotes in the outer part (Nieto & Bender, 1989). Occasionally, boxy ellipticals have kinematically distinct cores (Franx & Illingworth, 1988; Jedrzejewski & Schechter, 1988; Bender, 1988). These cores inhabit flattened rapidly rotating disk- or torus-like components dominating the light in the central few hundred parsecs (Bender, 1990; Rix & White, 1992; Mehlert et al., 1998), but they contribute only a few percent to the total light of the galaxy. The fact that these cores are metal-enhanced shows that gas (at least in the inner regions) must have played an important role during their formation (Bender & Surma, 1992; Davies, Sadler & Peletier, 1993; Bender, 1996; Davies, 1996). Boxy ellipticals show stronger radio emission than the average and have high X-ray luminosities consistent with emission from hot gaseous halos (Beuing et al., 1999). The distinct physical properties of disky and boxy elliptical galaxies point to the fact that both types of ellipticals could have different formation histories. It has been argued by Kormendy & Bender (1996) and Faber et al. (1997) that the observed stellar disks and the high density power law centers in disky ellipticals are signatures of dissipation during their formation. In this sense they seem to continue the Hubble sequence from S0s to higher bulge-to-disk ratios (Kormendy & Bender, 1996). Boxy ellipticals with an even higher bulge-to-disk ratio show a stronger kinematical decoupling at their centers and no signature of a disk at all. Therefore the early type Hubble sequence from S0s to massive ellipticals might represent galaxies with a disk component embedded in more prominent spheroidal body.

The traditional view on the formation and evolution of giant elliptical galaxies is that they are very old stellar systems and formed very early at a redshift of more than two (Searle, Sargent & Bagnuolo, 1973). After an intensive initial star-formation phase they experienced very little mass evolution (Bruzual & Charlot, 1993). It has been argued by many authors that the stellar evolution of ellipticals is compatible with pure passive evolution models (Bower, Lucey & Ellis, 1992; Aragon-Salamanca et al., 1993; Bender, Ziegler & Bruzual, 1996; Ellis et al., 1997; Ziegler & Bender, 1997) or models with exponentially decaying star formation (Ziegler et al., 1999). Alternatively, hierarchical theories of galaxy formation predict that massive galaxies were assembled relatively late in many generations of mergers of disk-type galaxies or smaller subunits and mass accretion. It has been argued by Kauffmann (1996) and Kauffmann & Charlot (1998) that this merger scenario is consistent with observations of galaxies at different redshifts. The idea that elliptical galaxies can form from mergers of disk-galaxies has been originally proposed by Toomre & Toomre (1972). Thereafter the “merger hypothesis“ has been investigated in great details by many authors. For a recent review of current models of spheroid formation see Burkert & Naab (2003a, b). Negroponte & White (1983), Barnes (1988) and Hernquist (1992) performed the first fully self-consistent merger models of two equal-mass stellar disks embedded in dark matter halos. The remnants were slowly rotating, pressure supported and anisotropic systems and generally followed an $`r^{1/4}`$ surface density profile in the outer parts. However, due to phase space limitations (Carlberg, 1986) it was necessary to start with progenitors with massive central bulge components (Hernquist, 1993b) to fit the observed de Vaucouleurs profile also in the inner parts. These simulations showed that global properties of equal-mass merger remnants resemble those of ordinary slowly rotating massive elliptical galaxies. More detailed investigations of isophotal shapes of the merger remnants have shown that the same remnant can appear either disky or boxy when viewed from different directions (Hernquist, 1993b) with a trend pointing to boxy isophotes (Heyl Hernquist & Spergel, 1994; Steinmetz & Buchner, 1955). Barnes (1998) and Bendo & Barnes (2000) investigated a sample of disk-disk mergers with a mass ratio of 3:1 and found that the remnants are flattened and fast rotating in contrast to equal-mass mergers. Naab, Burkert & Hernquist (1999) investigated the photometric and kinematical properties of a prototypical 3:1 merger remnant in details and compared the results to observational data of disky elliptical galaxies. They found an excellent agreement and proposed that fast rotating disky elliptical galaxies can originate from pure collisionless 3:1 mergers while slowly rotating, pressure supported ellipticals form from equal-mass mergers of disk galaxies. In this paper we extend the analysis of Naab, Burkert & Hernquist (1999). A large number of 112 merger remnants resulting from a statistically unbiased sample of simulations of mergers between disk-galaxies with mass ratios of 1:1, 2:1, 3:1, and 4:1 is investigated. This large sample allows a much more thorough investigation of the statistical properties of merger remnants as all previous studies. After a short description of the simulation methods in section 2 we investigate the photometric and kinematic properties of the simulated remnants in section 3. We discuss the drawback of our results on the theory of the formation of elliptical galaxies in section 4.

## 2 The merger models

The spiral galaxies were constructed in dynamical equilibrium using the method described by Hernquist (1993). We used the following system of units: gravitational constant G=1, exponential scale length of the larger disk $`h=1`$ and mass of the larger disk $`M_d=1`$. Each galaxy consists of an exponential disk, a spherical, non-rotating bulge with mass $`M_b=1/3`$, a Hernquist density profile (Hernquist, 1990) with a scale length $`r_b=0.2h`$ and a spherical pseudo-isothermal halo with a mass $`M_d=5.8`$, cut-off radius $`r_c=10h`$ and core radius $`\gamma =1h`$.

We followed a sequence of mass ratios of the progenitor disks from $`\eta =1`$ to $`\eta =4`$ where $`\eta `$ is the mass of the more massive galaxy divided by the mass of the merger partner. The equal-mass mergers were calculated adopting in total 400000 particles with each galaxy consisting of 20000 bulge particles, 60000 disk particles, and 120000 halo particles. We decided to use twice as many halo particles than disk particles to reduce heating and instability effects in the disk components (Naab, Burkert & Hernquist, 1999). For the mergers with $`\eta =2,3,4`$ the parameters of the more massive galaxy were as described above. The low-mass companion contained a fraction of $`1/\eta `$ the mass and the number of particles in each component with a disk scale length of $`h=\sqrt{1/\eta }`$, as expected from the Tully-Fisher relation (Pierce & Tully, 1992).

The N-body simulations for the equal-mass mergers were performed by direct summation of the forces using the special purpose hardware GRAPE6 (Makino, Fukushige & Namura, 2003). With this highly efficient hardware one force calculation for 400000 particles takes approx. $`11`$ seconds. The mergers with mass ratios $`\eta =2,3,4`$ were followed using the newly developed treecode VINE (Wetzstein et al., 2003) in combination with the GRAPE5 (Kawai, Fukushige, Makino, & Taiji, 2000) hardware for which the code was optimized. VINE uses a binary tree in combination with the refined multipole acceptance criterion proposed by Warren & Salmon (1995). This criterion takes the mass distribution of every node into account. It enables the user to control the absolute force error which is introduced by the tree construction. We chose a value of 0.001 which guarantees that the error resulting from the tree is smaller than the intrinsic force error of the GRAPE5 hardware which is of the order of $`0.1\%`$. To compensate for the limited mass resolution of GRAPE5 we limit the maximum mass of a tree node to $`10^3`$ times the minimum particle mass. To further increase the speed of the calculation we only use GRAPE5 if forces for more than 50 particles are needed at the same time. Otherwise we compute the forces on the host computer as for small particle numbers the communication time between the board and the host computer exceeds the time for the force calculation on the host. With the parameters given above one force calculation with VINE and GRAPE5 for 400000 particles takes approx. $`12`$ seconds.

We used a gravitational Plummer-softening of $`ϵ=0.05`$ and a fixed leap-frog integration time step of $`\mathrm{\Delta }t=0.04`$. For the equal-mass mergers simulated with direct summation on GRAPE6 the total energy is conserved, VINE in combination with GRAPE5 conserves the total energy up to $`0.5\%`$.

The galaxies approached each other on nearly parabolic orbits with an initial separation of $`r_{sep}=30`$ length units and a pericenter distance of $`r_p=2`$ length units (same parameters as e.g. Hernquist (1992)). A study of orbits of merging dark matter halos in cosmological large scale simulations by Khochfar (2003) has shown that most of the merging halos are indeed on parabolic orbits. The inclinations of the two disks relative to the orbit plane were $`i_1`$ and $`i_2`$ with arguments of pericenter $`\omega _1`$ and $`\omega _2`$. In selecting unbiased initial parameters for the disk inclinations we followed the procedure described by Barnes (1998). For the spin vector of each disk we defined four different orientations pointing to every vertex of a regular tetrahedron. The initial orientations we used translate to the following set of angles: For the first galaxy $`i_1=(0,109,109,109),\omega _1=(0,60,180,60)`$. The second galaxy has $`i_2=(180,71,71,71),\omega _2=(0,30,30,90)`$. These parameters result in 16 initial configurations for equal-mass mergers and 16 more for every mass ratio $`\eta =2,3,4`$ where the initial orientations are interchanged. Following the simple hypothesis that the orientations of the merging disks are independent of each other and independent of their mutual orbital plane, every merger geometry has an equal probability to be realized (Barnes, 1998). The orbital parameters are listed in Table 1. In total we simulated 112 mergers. The total computing time for the entire sample was about 1600 hours wall clock time.

For all simulations the merger remnants were allowed to settle into equilibrium approximately 8 to 10 dynamical times after the merger was complete. Then their equilibrium state was analyzed.

## 3 Photometric and kinematical properties of the remnants

To compare our simulated merger remnants with observations we analyzed the remnants with respect to observed global photometric and kinematical properties of giant elliptical galaxies, e.g. surface density profiles, isophotal deviation from perfect ellipses, velocity dispersion, and major- and minor-axis rotation. Defining characteristic values for each projected remnant we followed as closely as possible the analysis described by BDM.

### 3.1 Isophotal shape

An artificial image of the remnant was created by binning the central 10 length units into $`128\times 128`$ pixels. This picture was smoothed with a Gaussian filter of standard deviation 1.5 pixels. The isophotes and their deviations from perfect ellipses were then determined using a data reduction package kindly provided by Ralf Bender.

To be confident that a once determined isophotal shape is characteristic for the remnant and does not change with time we investigated the time evolution of the ellipticity and the $`a_4`$ profile starting $`20`$ time units after the merger of the bulge components was complete and followed the evolution for the next $`50`$ time units. In intervals of four time units the luminous part of the remnant was transformed to the principal axes of its moment of inertia tensor. The tensor was evaluated using $`40\%`$ of the most tightly bound particles. The isophotal properties were then analyzed as seen along the minor axis

The characteristic ellipticity $`ϵ_{\mathrm{eff}}`$ for each projection was defined as the isophotal ellipticity at $`1.5r_e`$. Figure 1 shows the ellipticity profiles for the 3:1 merger remnant 10 in intervals of four time units starting at $`t=150`$. The ellipticity profile shows little evolution with time and so does $`ϵ_{\mathrm{eff}}`$ (small box). The behaviour of this simulation is characteristic for the whole set of simulations.

Following the definition of BDM for the global properties of observed giant elliptical galaxies, we determined for every projection the effective $`a_4`$-coefficient, $`a4_{\mathrm{eff}}`$, as the mean value of $`a_4`$ between $`0.25r_e`$ and $`1.0r_e`$, with $`r_e`$ being the projected spherical half-mass radius. In case of a strong peak in the $`a_4`$-distribution with an absolute value that is larger than the absolute mean value, we chose the peak value. Being characteristic for all merger remnants Figure 2 shows the $`a_4`$-profiles at different times for the 3:1 merger remnant with geometry 10. Figure 3 shows the corresponding $`a_4`$-profiles for the equal-mass remnant with geometry 2. Although the individual $`a_4`$-profiles differ in details at different time steps the global characteristic shape is conserved. This is reflected in the time evolution of $`a4_{\mathrm{eff}}`$ which was assigned to each remnant at the different time steps as shown in the small diagrams of Figure 2 and Figure 3. The change in $`a4_{\mathrm{eff}}`$ due to evolutionary effects (a real change in the isophotal shape and/or a change in three-dimensional shape leading to a different projection angle) is of the same order as the error due to the limited particle number. The bootstrap error bar is shown in the small diagrams of Figure 1 and Figure 2.

These measurements convinced us that the photometric properties of the inner regions (inside one effective radius) of the merger remnants evolve relatively fast into an equilibrium configuration. As the photometric properties do not significantly change with time they can be used for a further statistical analysis.

To investigate projection effects we determined for each simulation $`a4_{\mathrm{eff}}`$ and $`ϵ_{\mathrm{eff}}`$ for 500 random projections of every remnant. This resulted in 8000 projections for $`\eta =1`$ and 16000 projections for $`\eta =2,3,4`$, respectively.

The normalized histograms of projected ellipticities are shown in Figure 4. For equal-mass mergers the distribution rises from a small number of projections with zero ellipticity to a peak around $`ϵ_{\mathrm{eff}}=0.3`$ and then falls to zero at $`ϵ_{\mathrm{eff}}=0.7`$. 2:1 remnants show one peak around $`ϵ_{\mathrm{eff}}=0.25`$ and a strong peak at $`ϵ_{\mathrm{eff}}=0.45`$ falling to zero at $`ϵ_{\mathrm{eff}}=0.7`$. 3:1 and 4:1 remnants show a flat distribution for small ellipticities rising to a strong peak around $`ϵ_{\mathrm{eff}}=0.55`$ and $`ϵ_{\mathrm{eff}}=0.6`$, respectively. No projection leads to ellipticities larger than $`0.7`$.

The distribution function for the isophotal shapes of equal-mass remnants peaks at $`a4_{\mathrm{eff}}0.5`$ and declines rapidly for more negative values (Figure 5). In total 47% of the projected 1:1 remnants show disky isophotes distributed in a broader wing for positive $`a4_{\mathrm{eff}}`$. For higher mass ratios the distribution functions have a similar shape but are shifted to more positive values of $`a4_{\mathrm{eff}}`$. 2:1 remnants peak at around $`a4_{\mathrm{eff}}=0.5`$ and $`73\%`$ of the projected remnants show disky isophotes. The distributions for 3:1 and 4:1 remnants both peak around $`a4_{\mathrm{eff}}1`$. For these systems the fraction of disky projections is about $`89\%`$ and $`91\%`$, respectively.

In addition, the 500 projected values for $`a4_{\mathrm{eff}}`$ and $`ϵ_{\mathrm{eff}}`$ where used to calculate the two-dimensional probability density function for a given simulated remnant to be ”observed ” in the $`a4_{\mathrm{eff}}`$-$`ϵ_{\mathrm{eff}}`$ plane. We added up the probability densities assuming equal weights for every simulation geometry at a given mass ratio. Figure 6 shows the results for mergers with $`\eta =1,2,3`$, and $`4`$. The contours indicate the areas of 50% (dashed line), 70% (thin line) and 90% (thick line) probability to detect a merger remnant with the given properties. Observed data points from BDM are over-plotted. Filled boxes are boxy ellipticals with $`a4_{\mathrm{eff}}0`$ while open diamonds indicate disky ellipticals with $`a4_{\mathrm{eff}}>0`$. The errors were estimated applying the statistical bootstrapping method (Heyl Hernquist & Spergel, 1994). The area covered by 1:1 remnants with negative $`a4_{\mathrm{eff}}`$ is in good agreement with the observed data for boxy elliptical galaxies. In particular the observed trend for more boxy galaxies to have higher ellipticities is reproduced. The disky equal-mass remnants also follow the trend for observed disky ellipticals. Remnants with an $`a4_{\mathrm{eff}}`$ around zero can have slightly higher ellipticities than observed. Remnants with $`\eta =2,3`$, and $`4`$ predominantly populate the region of disky ellipticals and more disky ellipticals also tend to be more flattened. There is a trend for projections with $`a4_{\mathrm{eff}}1`$ to have slightly higher ellipticities than observed (or at a given ellipticity the remnants are not disky enough). Observed very disky ellipticals with relatively small ellipticities can not be reproduced.

In summary there is a clear trend for unequal-mass mergers to produce more disky remnants. Responsible for the disky appearance of the 3:1 and 4:1 remnants is the distribution of the particles of the massive disk. Figure 7 shows the different contributions from the small and the large progenitor galaxy and the resulting isophotal map of a characteristic 3:1 merger remnant. The particles originating from the small progenitor accumulate in a torus-like structure that has peanut-shaped or boxy isophotes. In contrast, the dominant luminous material from the larger progenitor still keeps its disk-like appearance. In combination, the contribution from the larger progenitor – since it is three times more massive – dominates the overall properties of the remnant. This result holds for all 3:1 and 4:1 merger remnants. The more massive disk component is not completely destroyed during the merger event and determines the overall structure of the remnant. For equal-mass mergers both disks are effected during the merger and they can loose the memory of their initial state.

The isophotal analysis also provides information about the radial change of the relative orientation of the major axes of the isophotes. In general, the amount of isophotal twist depends on the projection angle. This is demonstrated in Figure 8 for a more elongated and a nearly round projection of a 3:1 merger remnant. To get a quantitative measure for the isophotal twist we determined the relative position angle $`\mathrm{\Delta }\mathrm{\Phi }`$ between the isophote at 0.5 $`r_{\mathrm{eff}}`$ and 1.5 $`r_{\mathrm{eff}}`$ for every projection of the remnant. Figure 9 shows a comparison of isophotal twists for characteristic remnants with $`\eta =1,2,3`$, and $`4`$. The isophotal twist is in general larger for projections that appear nearly round. For ellipticities larger than $`ϵ_{\mathrm{eff}}0.4`$ the isophotal twist is $`\mathrm{\Delta }\mathrm{\Phi }20^{}`$. The distribution of $`\mathrm{\Delta }\mathrm{\Phi }`$ versus $`ϵ_{\mathrm{eff}}`$ for random projections of every remnant is consistent with observations of elliptical galaxies (Bender, Döbereiner & Möllenhoff, 1988). However, the simulated remnants can show larger isophotal twists at high ellipticities than observed galaxies. In addition there seems to exist no obvious correlation between the mass ratio of the galaxies and the amount of isophotal twist.

### 3.2 Kinematics

The central velocity dispersion $`\sigma _0`$ of every remnant was determined as the average projected velocity dispersion of the luminous particles inside a projected galactocentric distance of $`0.2r_e`$. We defined the characteristic rotational velocity along the major and the minor axis as the projected rotational velocity at $`1.5r_e`$ and $`0.5r_e`$, respectively. Figure 10 shows the time evolution for these velocity measurements for a 3:1 merger simulation which is characteristic for all simulations. The derived kinematical properties of the remnants stay nearly constant over a long time period and are therefore a good measurement of the intrinsic kinematics of the simulated remnants. Again, the bootstrap error is of the same order as the change of the measured values with time.

Like for the isophotal shape we computed the statistical kinematical properties of the simulated remnants and compared them with observational data of elliptical galaxies.

The normalized histograms for $`(v_{\mathrm{maj}}/\sigma _0)`$ are shown in Figure 11.

There is a clear trend for 1:1 mergers to produce slowly rotating ellipticals with $`(v_{\mathrm{maj}}/\sigma _0)<0.3`$. For 2:1 mergers the peak value is around $`(v_{\mathrm{maj}}/\sigma _0)=0.4`$. 3:1 and 4:1 merger remnants have peaks around $`(v_{\mathrm{maj}}/\sigma _0)=0.6`$ and $`0.8`$, respectively. There are indications for a weak secondary peak around $`(v_{\mathrm{maj}}/\sigma _0)=1`$ for 3:1 and 4:1 remnants.

Figure 12 shows the distribution function in the $`(v_{\mathrm{maj}}/\sigma _0)`$-$`ϵ_{\mathrm{eff}}`$ plane. The area of slowly rotating boxy ellipticals (filled boxes) is almost completely covered by remnants of 1:1 mergers while 2:1 mergers have rotational properties resembling faster rotating boxy and slowly rotating disky ellipticals. 3:1 and 4:1 remnants are clearly fast rotating, show high ellipticities and can be associated with fast rotating disky ellipticals. Mergers with an increasing mass ratio produce faster rotating ellipticals with higher ellipticities. All simulated data are in good agreement with observations.

However, there is a trend for the simulated remnants to have slightly larger maximum ellipticities than observed. Projections with high ellipticities ($`ϵ>0.5`$) show systematically smaller values for $`v_{\mathrm{maj}}/\sigma _0`$ than observed. The effect is strongest for equal-mass merger remnants and is discussed in detail below. This indicates that the simulated systems are more strongly supported by anisotropic velocity dispersions than the observed systems, despite significant rotation for 3:1 and 4:1 remnants. Observed ellipticals with $`(v_{\mathrm{maj}}/\sigma _0)>1.2`$ can not be reproduced as was already reported by Cretton, Naab, Rix, & Burkert (2001).

The minor-axis kinematics of the simulated remnants was measured as the rotation velocity along the minor axis at $`0.5r_{\mathrm{eff}}`$. The amount of minor axis rotation was parametrized as $`(v_{min}/\sqrt{v_{\mathrm{maj}}^2+v_{min}^2})`$ (Binney, 1985). Minor axis rotation in elliptical galaxies, in addition to isophotal twist, has been suggested as a sign for a triaxial shape of the main body of elliptical galaxies (Wagner, Bender & Moellenhoff, 1988; Franx Illingworth & De Zeeuw, 1991). Indeed, 1:1 mergers show a significant amount of minor-axis rotation (Figure 13). This is consistent with their mostly triaxial shape and the measured isophotal twist. 3:1 and 4:1 remnants are more oblate and show weak minor axis rotation. However, the amount of isophotal twist is comparable to equal-mass merger remnants. A detailed analysis of the connection between intrinsic shape, isophotal shape, and internal kinematics is beyond the scope of this paper and needs further investigation.

The anisotropy parameter $`(v_{\mathrm{maj}}/\sigma _0)^{}`$ was defined as the ratio of the observed value of $`(v_{\mathrm{maj}}/\sigma _0)`$ and the theoretical value for an isotropic oblate rotator $`(v/\sigma )_{theo}=\sqrt{ϵ_{\mathrm{obs}}/(1ϵ_{\mathrm{obs}})}`$ with the observed ellipticity $`ϵ_{\mathrm{obs}}`$ (Binney, 1978). This parameter has been used by observers as a test whether a given galaxy with observed $`v_{\mathrm{maj}}`$, $`\sigma _0`$ and $`ϵ_{\mathrm{obs}}`$ is flattened by rotation ($`(v_{\mathrm{maj}}/\sigma _0)^{}0.7`$) or by velocity anisotropy ($`(v_{\mathrm{maj}}/\sigma _0)^{}<0.7`$) (Davies et al., 1983; Bender, 1988; Nieto, Capaccioli & Held, 1988; Scorza & Bender, 1995).

Figure 14 shows the normalized histograms for the $`(v_{\mathrm{maj}}/\sigma _0)^{}`$ values of the simulated remnants. The 1:1 remnants peak around $`(v_{\mathrm{maj}}/\sigma _0)^{}0.3`$ with a more prominent tail towards lower values. Only about 10% of the projections have $`(v_{\mathrm{maj}}/\sigma _0)^{}>0.6`$. In contrast 82% have $`(v_{\mathrm{maj}}/\sigma _0)^{}<0.4`$ and would clearly be interpreted as being supported by anisotropic velocity dispersion. 2:1 remnants peak at $`(v_{\mathrm{maj}}/\sigma _0)^{}0.5`$ with 41% showing $`(v_{\mathrm{maj}}/\sigma _0)^{}>0.6`$. 3:1 and 4:1 mergers with peaks at $`(v_{\mathrm{maj}}/\sigma _0)^{}0.6`$ and $`(v_{\mathrm{maj}}/\sigma _0)^{}0.7`$ are consistent with model predictions of oblate isotropic or slightly prolate rotators. The percentage of projections with $`(v_{\mathrm{maj}}/\sigma _0)^{}>0.6`$ is 64% and 83%, respectively. 3:1 and 4:1 remnants also have predominantly disky isophotes and cover the area populated by observed disky ellipticals in the $`log(v_{\mathrm{maj}}/\sigma _0)^{}`$ \- $`a4_{\mathrm{eff}}`$ diagram (Figure 15).

1:1 merger remnants show the most complex behaviour. A significant number of the projected remnants lie in a region of disky anisotropic systems where no elliptical galaxy is observed. Assuming that a projected remnant clearly fails to resemble an observed elliptical if it has $`a_4>0.3`$ and $`(v_{\mathrm{maj}}/\sigma _0)^{}<0.5`$, 28% of the 1:1 remnants would fall into the forbidden regime.

To understand the behaviour of equal-mass remnants in more detail we illustrate the influence of the initial orientation of the progenitor disks on the properties of the merger remnants. We defined mean values of the investigated quantities of all 500 projections for every merger geometry. Thereafter we divided the equal-mass remnants into four groups of geometries that produce remnants with almost distinct properties: slowly rotating boxy remnants with $`\overline{v_{\mathrm{maj}}/\sigma _0}<0.2`$ and $`\overline{a_4}<0.4`$ (group A; geometries 2,4,5, and 14), slowly rotating remnants with predominantly disky isophotes having $`\overline{v_{\mathrm{maj}}/\sigma _0}<0.2`$ and $`\overline{a_4}>0.25`$ (group B; geometries 1,8,9,11, and 13), modestly rotating remnants ($`\overline{v_{\mathrm{maj}}/\sigma _0}>0.2`$) that appear to be isotropic with $`\overline{(v_{\mathrm{maj}}/\sigma _0)^{}}>0.7`$ (group C; geometries 7, 12, and 15) and remnants with modest deviations in isophotal shape $`0.2<\overline{a_4}<0.2`$ (group D; geometries 3,6,12, and 16). The results are shown in Figure 16.

Remnants of group A are slowly rotating anisotropic systems with purely boxy isophotes. More than 90% of the projected remnants are consistent with observations of boxy anisotropic elliptical galaxies. This class of equal-mass merger remnants was sampled by Naab, Burkert & Hernquist (1999). Most projections of group B show disky isophotes and appear to be anisotropic. Around 60% of the projected remnants of this group fail to fall on the observed correlation. Furthermore the remnants of this group show too high ellipticities at low rotation velocities, a behaviour already described by Heyl Hernquist & Spergel (1996). The shaded regions in Figure 16 indicate the location of the highly elliptical projections ($`ϵ>0.4`$) . They dominate the disky anisotropic branch of ”not observed” merger remnants. It is possible that the very elongated projections of group B would be classified as S0 galaxies. Some of them, as NGC4550, show weak net rotation and are still very elongated. Geometry 1 produces a very axisymmetric remnant with two counterrotating populations of stars in as it is observed in NGC4550 (see e.g. Pfenniger (1999)). However, Rix, Franx, Fisher, & Illingworth (1992) argued that a merger origin for this particular galaxy is unlikely as it would have heated both disks too much.

Group C combines all rotating 1:1 remnants. They appear slightly disky or boxy depending on the orientation and are isotropic. It has to be noted that even the disky remnants of this group are consistent with observations and only 10% of the projections are failures. This group shows a clear connection to the initial conditions as the spins of the progenitor disks were almost aligned. Projections of group D can either have disky or boxy isophotes and are mostly anisotropic. 35% of the projections fall in the forbidden regime.

## 4 Conclusions & Discussion

We used a large set of N-body simulations of collisionless mergers of disk galaxies with mass ratios of $`\eta =1,2,3,`$ and $`4`$ to investigate for the first time statistical properties of disk-merger remnants. In contrast to previous studies where only projections along the principal axes were investigated we analyzed a large number of randomly projected merger remnants and compared the results in a consistent way to observations. Galaxies on the sky are also viewed along random lines of sight therefore this is the appropriate way of comparison.

We showed that the detailed isophotal and kinematical properties of the simulated merger remnants reach their equilibrium values relatively fast after the merger of the central parts of the galaxies is complete and thereafter do not significantly evolve with time. The high resolution of the simulations made it possible to keep the errors, especially for $`a_4`$, at a reasonably low value and enabled us to perform a statistical analysis.

The basic result is that the projected kinematical and photometric properties of remnants of major mergers of disk galaxies are in surprisingly good agreement with the observational data for elliptical galaxies. The mass ratio of the progenitor disks determines the global properties of the remnants. The influence of the initial orientation on the remnants is strongest for equal-mass mergers.

Purely boxy, anisotropic and slowly rotating remnants with a large amount of minor axis rotation can only be produced by equal-mass mergers with certain initial orientations. However, if the initial spins of the disks were almost aligned the remnants appear to be isotropic and disky or boxy, depending on the viewing angle. In total 28% of all projected 1:1 remnants show properties that are not observed at all. As some of these remnants have high ellipticities and very small $`(v_{\mathrm{maj}}/\sigma _0)^{}`$ they appear to be flattened by strongly anisotropic velocity dispersions. In addition they mostly show very disky isophotes. As we can not exclude that mergers with these geometries have occurred they should be observed, which is not the case. If this controversy is not solved in the future it might constitute a serious problem for the merger hypothesis. However, if observed, these elongated disky objects might have been classified as S0s and have therefore been excluded from the observed sample of bona fide ellipticals.

Based on the diverse properties of the 1:1 remnants they more likely resemble the class of intermediate mass giant ellipticals in the transition region from rotating disky to non-rotating boxy ellipticals. In this regime even boxy rotating galaxies have been observed (Nieto & Bender, 1989).

In contrast, 3:1 and 4:1 mergers form a more homogeneous group of remnants. They have preferentially disky isophotes, are fast rotating and show small minor axis rotation independent of the assumed projection. In general the properties are in very good agreement with observations of fast rotating disky ellipticals.

2:1 mergers have intermediate properties with boxy or disky isophotes depending on the projection and the orbital geometry of the merger. Globally, different projection angles and orbits do not change the fundamental properties of 1:1 mergers on the one side and 3:1 and 4:1 mergers on the other side.

In summary, many photometric and kinematical properties of low and intermediate mass giant elliptical galaxies can be understood as a sequence of major mergers of disk galaxies with varying mass ratios. The homogeneous group of massive boxy ellipticals populates the boxy anisotropic area of the $`(v_{\mathrm{maj}}/\sigma _0)^{}`$ \- $`a4_{\mathrm{eff}}`$ plane. Their distribution is inconsistent with the distribution of simulated 1:1 remnants. Therefore they are most likely not remnants of equal-mass mergers of disk galaxies. Boxy and mildly anisotropic remnants are only reproduced by a small subsample of initial conditions. There is no reason why other geometries should have been avoided. It is more likely that massive boxy ellipticals have formed by other processes like mergers of early type galaxies (Naab & Burkert, 2000) or multiple mergers in a group environment (Weil & Hernquist, 1996).

The simulations presented here were purely dissipationless, taking into account only the stellar and dark matter component of a galaxy. The importance of gas in galaxy-galaxy mergers and the detailed influence on the structure of elliptical galaxies is not fully understood up to now. Numerical simulations of galaxy mergers including gas by Barnes & Hernquist (1996) and Barnes (1998) have shown that the presence of gas can change the orbital structure and the shape of the merger remnants. As soon as the gas is driven to the center during the merger the mass concentration seems to be responsible for the destruction of stellar box orbits. This process makes it even more difficult to explain the formation of giant boxy ellipticals by binary disk mergers.

The present study will serve as the basis of a further detailed investigation on the influence of an additional dissipative component. Kormendy & Bender (1996) proposed a revised Hubble sequence with disky ellipticals representing the missing link between the Im-Spiral-S0 sequence and boxy ellipticals. They noted that gas infall into the equatorial plane with subsequent star formation could lead to a second disk-like subcomponent. Ellipticals with disks could appear disky when seen edge-on and boxy otherwise (Scorza & Bender, 1995). Our simulations indeed indicate that a disk-like substructure is responsible for producing disky isophotes (see Fig. 6). However, in the present case, the disk is the remnant of the more massive spiral which was not completely destroyed during the minor merger (see also Barnes (1998)). Naab & Burkert (2001) investigated line-of-sight velocity distributions of dissipationless merger remnants and found a velocity profile asymmetry that is opposite to the observed one. They concluded that this disagreement can be solved if ellipticals would contain a second disk-like substructure that most likely formed through gas accretion. The situation is however not completely clear, as another study by Bendo & Barnes (2000) found a good agreement of the observed asymmetries for some cases. Naab & Burkert (2001) have shown that extended gas disks can form as a result of a gas rich unequal-mass merger. Barnes (2002) recently presented a first detailed set of equal-and unequal-mass merger simulations of gas-rich galaxies resulting in the formation of extended gas disks. His simulations demonstrate a very complex dynamical evolution. It has to be investigated how the presence of gas changes the detailed properties of the remnants in detail. In addition, star formation during the merger will influence the stellar kinematics (see e.g. Springel (2000)) and the stellar populations by adding young and probably more metal rich stars. As metal enrichment provides a futher strong observational constraint on the formation history of early type galaxies (see e.g. Thomas, Maraston, & Bender (2002)) more simulations, including dissipation, star formation and chemical evolution will now be required in order to understand in detail the role of gas for the formation of elliptical galaxies by mergers.

We thank Ralf Bender and Hans-Walter Rix for helpful comments on the script and the referee, Joshua Barnes, for valuable suggestions that significantly improved the paper. All the simulations presented in this work have been performed on the GRAPE cluster at the Max-Planck-Insitute for Astronomy in Heidelberg (MPIA). Thorsten Naab is grateful to the MPIA for its support.


text:
# The role of the initial state in the collapse of a Bose condensed gas

## Abstract

The stability of a Bose condensed gas with a negative scattering length is studied with regard to the role played by the initial state. A trapped ideal gas ground state is shown to be unstable when the Gross Pitaevskii equation still predicts the existence of a stable state. A possible relation to a recent experimental study of the critical conditions for stability or collapse of a Bose Einstein condensate \[J. L. Roberts et al., Phys. Rev. Lett. 86, 4211 (2001)\] is discussed.

The possibility of manipulating interatomic forces using magnetic field tunable Feshbach resonances , has provided new opportunities for the study of collision phenomena in ultra cold gases. Recent experiments with 85Rb  have shown how one can tune the two body interaction, in a condensate from strongly repulsive to attractive over a wide range of different $`s`$ wave scattering lengths. This makes possible the study of non equilibrium phenomena under a variety of controlled initial conditions. One of the most striking experiments of this type is that in which a condensate, first formed with repulsive forces, is switched to the attractive case, and collapses. In an idealized sense collapse means an unrestricted contraction of a gas within a finite amount of time. As shown by Pitaevskii , the time dependent Gross Pitaevskii equation (GPE) predicts this remarkable phenomenon for a trapped Bose condensed gas with a negative scattering length in a harmonic trap as soon as certain conditions, on the interaction strength, and the initial condensate wave function, are fulfilled.

Very recently, the control of the interatomic force in 85Rb  has allowed an experimental determination of the critical conditions for the stability, with respect to collapse, of a condensate. In the course of these experiments, the interatomic force of the gas was tuned towards specific negative $`s`$ wave scattering lengths by means of magnetic fields, using a Feshbach resonance. Since the initial scattering length was chosen to be negligibly small, from the point of view of the GPE the prepared initial state can be interpreted as an ideal gas condensate wave function.

Since their first observation in experiments with 7Li , Bose condensed gases with a negative scattering length have been studied by numerous authors. Most of this work has used the GPE and focused on the stability conditions for a ground state solution (see, e.g.,  and references therein) or condensate formation . The experimental results in Ref.  indicate, however, that as soon as the condensate becomes unstable the GPE alone is not suitable for the description of the dynamics of the collapse. In particular the observed substantial loss of condensate atoms contradicts the number conservation inherent in the GPE. This failure can be expected because the collision term of the GPE does not account for collisional losses of condensate atoms due to non forward scattering , which should become relevant for the rapid dynamics occurring during the collapse.

More interestingly the gas was seen to be unstable even under conditions, where the GPE predicts the existence of a stable ground state. To clarify this aspect of the observations we have studied the stability of a Bose condensed gas with a negative scattering length using the GPE, taking into account the role of the initial state, prepared in the experiment in Ref. . We show that this initial ideal gas state leads to collapse condition that is different from that for the existence of a ground state of the condensate.

Throughout this article the $`s`$ wave scattering length of the interatomic force, $`a_0`$, is assumed to be negative. For a Bose condensed gas in a spherical trap, with the harmonic potential $`V_{\mathrm{trap}}=m\omega _{\mathrm{ho}}^2r^2/2`$, the natural units of time and length are the inverse trap frequency $`\omega _{\mathrm{ho}}^1`$ and the harmonic oscillator length $`a_{\mathrm{ho}}=\sqrt{\mathrm{}/m\omega _{\mathrm{ho}}}`$, respectively. Expressed through these units the GPE for an interacting Bose gas, consisting of $`N`$ atoms, assumes the well known form

$$i\frac{}{t}\psi =\frac{1}{2}\mathrm{\Delta }\psi +\frac{1}{2}r^2\psi 4\pi k|\psi |^2\psi .$$
(1)
Here the condensate wave function $`\psi `$ is normalized to unity, $`k`$ is given by $`N|a_0|/a_{\mathrm{ho}}`$ and the corresponding unit of the energy is then $`N\mathrm{}^2/ma_{\mathrm{ho}}^2`$. Equation (1) shows that $`k`$ is the only free adjustable parameter of the GPE for a spherical trap.

On first sight, the form of Eq. (1) might suggest that the critical conditions for collapse are governed solely by the parameter $`k`$. With respect to the existence of a stable ground state of a condensate this is true. In fact, the theoretical studies in Ref.  confirm the condition for stability

$$k<k_{\mathrm{cr}}=0.575.$$
(2)
When the initial condensate wave function deviates from the ground state, however, a collapse can occur even if the actual parameter $`k`$ is well below the above standard value.

As proposed by Pitaevskii  the conditions for stability or collapse can be visualized by means of two fundamental estimates provided by the GPE. The first one is the uncertainty relation between the mean square of the radius of a condensate, $`r^2=𝑑𝐫r^2|\psi |^2`$, and its kinetic energy, $`T=𝑑𝐫|\psi |^2/2`$, which reads

$$r^2T9/8,$$
(3)
where the choice of units corresponds to Eq. (1). It is worth noting that once $`r^2`$ vanishes the condensate wave function has collapsed. The kinetic energy then becomes infinite. The second estimate consists in an inequality between the mean field energy, $`E_{\mathrm{mf}}=4\pi k𝑑𝐫|\psi |^4/2`$, and the kinetic energy,

$$𝑑𝐫|\psi |^4/2\beta T^{3/2},$$
(4)
where $`\beta =0.0575`$ is a universal constant. Equations (3) and (4) are preserved for all times by the GPE.

As an additional constant of motion, the total energy,

$$E=T+V+E_{\mathrm{mf}},$$
(5)
contains the kinetic and the mean field part as well as the trap energy $`V=𝑑𝐫V_{\mathrm{trap}}|\psi |^2`$. In the spherical case the mean square of the radius of a condensate and the trap energy are related through

$$V=r^2/2.$$
(6)
The inequality given by Eq. (3) then provides an estimate for the trap energy Eq. (6) through the kinetic energy. In addition, according to Eq. (4), the mean field contribution to the total energy Eq. (5) can be estimated in terms of $`T`$. This in turn implies the general relation

$$ET+\frac{9}{16}T^14\pi \beta kT^{3/2}f(T,k),$$
(7)
which is valid for all times.

In Fig. 1 the lower bound Eq. (7) is illustrated for the parameters $`k=0.5`$ (Fig. 1a), $`k=0.53`$ (Fig. 1b) and $`k=0.575`$ (Fig. 1c) and an initial ideal gas ground state. For this Gaussian initial condensate wave function the total energy and the initial kinetic energy are given by $`E=E_{\mathrm{in}}=3/2k/\sqrt{2\pi }`$ and $`T_{\mathrm{in}}=3/4`$, respectively. In Fig. 1a the function $`f(T,k)`$ exhibits a local minimum which encloses the point $`(T_{\mathrm{in}},E)`$. As the total energy is a constant of motion and $`T(t)`$ varies continuously in time the kinetic energy is caught in the valley for all times. As pointed out by Pitaevskii , the inequality given by Eq. (3) then predicts the mean square of the radius of the condensate to be positive, i.e. $`r^29/8T>0`$, and a total collapse, i.e. $`r^2=0`$, can never occur. As can be seen in Fig. 1c, very close to the standard critical interaction strength Eq. (2) the function $`f(T,k)`$ becomes monotonous. Hence, a ground state of the GPE ceases to exist as soon as the local minimum of $`f(T,k)`$ vanishes.

In Fig. 1b, however, the total energy is above the local maximum of $`f(T,k)`$ and at least the lower bound in Eq. (7) allows the kinetic energy to access arbitrarily high values. In this case a total collapse of the condensate wave function cannot be excluded by the inequality Eq. (3).

To clarify this situation, Fig. 2 shows the results of numerical simulations of the GPE in Eq. (1) for $`k=0.52`$ (Fig. 2a) and $`k=0.53`$ (Fig. 2b), with an initial ideal gas ground state condensate wave function. In Fig. 2a, the simulation predicts the condensate to be stable with respect to collapse within the considered period of time. For the present parameter $`k=0.52`$ the total energy is already slightly above the local maximum of the function $`f(T,k)`$ and a rigorous assertion concerning stability or collapse cannot be obtained directly from Pitaevskii’s estimate.

In Fig. 2b, however, the condensate wave function exhibits a collapse although the actual parameter of the interatomic force of $`k=0.53`$ is still sufficiently weak for a stable ground state of the GPE to exist. To visualize the contraction of the condensate wave function, in Fig. 3 a typical density $`|\psi (r,t)|^2`$, close to the time of the collapse, is compared with the initial Gaussian density. The corresponding wave function is obtained from the simulation in Fig. 2b. As can be seen on the logarithmic scale, the predicted condensate density, close to the collapse, deviates from a Gaussian and exhibits a sharp peak at its center.

The GPE thus predicts the preparation of an initial ideal gas ground state condensate wave function to lead to a pronounced reduction of the critical parameter of the interatomic force to

$$k_{\mathrm{cr}}0.53$$
(8)
in comparison with the standard critical value for the existence of a stable ground state. This critical value for an initial ideal gas ground state has been reported also in Ref.  as an intermediate result in the numerical study of stable ground states of the GPE. A numerical simulation alone, however, cannot predict stability of a condensate on arbitrarily long time scales. The required lower estimate of $`k_{\mathrm{cr}}`$ is thus provided by the analysis illustrated in Fig. 1 and yields

$$k_{\mathrm{cr}}>0.5.$$
(9)
The sharply defined critical parameter of

$$k_{\mathrm{cr}}=0.461\pm 0.012\pm 0.054,$$
(10)
reported in Ref. , is indeed significantly smaller than the predicted standard value for the existence of a stable ground state in Eq. (2). This corroborates the predictions of the GPE in view of the destabilizing role of the particular initial state that was prepared in the experiment. Equations (8) and (9), however, are not strictly valid in the experimental case because the magnetic field was ramped linearly to its final value rather than abruptly. Interestingly, though, our result obtained from the GPE agrees with the observations within the experimental upper limit for $`k_{\mathrm{cr}}`$.

The analysis of stability or collapse of a dilute Bose condensed gas with a negative scattering length, performed in this work, can be applied to arbitrary initial condensate wave functions and trap geometries. The optimal estimate from Eq. (7) is obtained for a spherical trap. Related experiments with more general initial states should be feasible with the method described in Refs.  and could provide a further understanding of the range of validity of the GPE in non equilibrium situations.

###### Acknowledgements.

The author thanks Thomas Gasenzer, Samuel Morgan and Keith Burnett for stimulating discussions. This work was supported by the Alexander von Humboldt Foundation.


text:
# Results

## Results

Fig. 1 shows the model compounds used in the calculations. In a recent polymer study (Nakayama et al., 1999), experimentally observed differences in intensity for ester and carbonyl oxygen main lines could be explained directly from consideration of the main line intensities, where the main line intensity was estimated from the overlap obtained from the $`n=0`$ term in equation (1). This has the clear advantage for inherently large molecules, such as polymers or polypeptides, that main-line losses can be estimated without consideration of all the shake states to which intensity is distributed. For large molecules, the number of shake states goes beyond that which can readily be calculated at the configuration interaction level of theory. In a Gly-Gly -Gly tripeptide (Fig. 1a), the central glycine unit can be expected to display a shake spectrum which is similar to that of a residue within a longer polypeptide chain. To investigate effects arising from the truncations of the polypeptide chain, shake effects were calculated for all atoms in the model, including the atoms from the terminal carboxyl and amino functional groups.

There are six carbon atoms in the peptide model, showing very similar overlaps between the initial unionized state, and the core ionized ground state (cf. Table 1). This overlap is $`0.85`$ for all carbons, except for the carbon at the terminal carboxylic acid group, which has an overlap of $`0.87`$. The intensity loss from the main line is therefore close to ca. 30% for all carbon atoms in these models (cf. Table 1).

The oxygen shake contribution is generally the largest. In particular, the carbonyl oxygens display large shake effects, in agreement with the general behaviour for organic molecules. Here, strong shake is found for the two carbonyl oxygens in the peptide links, as well as for the double-bonded oxygen in the terminal carboxylate functional group. The overlap for these three oxygens is $`0.78`$ (cf. Table 1), corresponding to a total shake intensity loss of ca. 40% (cf. Table 1). The terminal hydroxyl oxygen has a slightly larger overlap of $`0.82`$, corresponding to the main line intensity loss of ca. 33%.

There are two nitrogens in the peptide links, and one in a terminal amino group. The shake effect for all these atoms is similar: the overlap is $`0.90`$ (Table 1) and the main line loss is ca. 20% (cf. Table 1). Nitrogen shows the smallest shake contribution, and the results agree well with data on a polyimide polymer (Nakayama et al., 1999), which also showed little nitrogen shake effect.

Sulphur shake effects have been calculated for three different common chemical environments (Fig. 1b-d and Table 2). The largest loss was found for a disulphide bridge (R-S-S-R, 25 %), and the smallest one for the thiol terminal group (R-SH, 21 %).

The initial models used for the sulphur calculations were smaller than the models used for the other atom types. This could, in principle, be a cause of concern, as these type of calculations are known to underestimate shake contributions in extended systems. In order to test the stability of the numerical results, additional calculations were performed on models with increasing chain lengths ( $`CH_3SH`$, $`CH_3CH_2SH`$, and $`CH_3CH_2CH_2SH`$ ). The results indicated negligible drift (not shown).

## Discussion

Radiation damage prevents the structural determination of single biomolecules and other non-repetitive structures (like cells) at high resolutions in classical electron or X-ray scattering experiments (Henderson, 1990; Henderson, 1995). Analysis of time-dependent components in damage formation suggests that the conventional damage barrier can be substantially extended at extreme dose rates and ultrashort exposure times (Neutze et al., 2000; Hajdu, 2000; Hajdu et al., 2000; Hajdu & Weckert, 2001; Ziaja et al., 2001). A quantitative description of damage formation and a detailed analysis of the ionisation dynamics of the sample are crucial for planning experiments at future X-ray lasers. Results described in this paper give the first assessments of shake contributions to the ionisation of carbon, nitrogen, oxygen and sulphur atoms within the most common chemical environments in a protein molecule. In previous applications, the calculated shakeup intensities were usually a few percentage units too high, so we judge the present results to provide upper bounds to the total shake intensities.

Data for noble gases show that shake-up and shake-off effects were of comparable magnitude (Svensson et al., 1988; Armen et al., 1985; Wark et al., 1991). Calculations on carbon, oxygen and nitrogen by Mohammedein et al., 1993 and El-shemi & Hassan, 1997 show that photoemission accompanied by double Auger decay has nearly zero probability for these light elements. Such triple ionisations can thus be neglected in the ionisation dynamics of C, O, N compounds. Sulphur, on the other hand, behaves differently, and the probability of multiple ionisation ($`n>2`$) after the ejection of a K-shell photoelectron is larger. The average charge left in sulphur after a single K-shell photoionisation is estimated to be about 4 (Mohammedein et al., 1993). For instance, for S: $`P(2)5\%`$, $`P(4)34\%`$, $`P(6)3\%`$ (Mohammedein et al., 1993; El-shemi & Hassan, 1997). This implies that for sulphur, multiple ionisation mechanisms should be taken into consideration in order to obtain more accurate predictions for the ionisation dynamics of proteins.

The energy of the photoelectrons, which emerge from shake-up processes is similar to the energy of photoelectrons released during outer-shell ionisation (Nakayama et al., 1999). In these calculations the deviation was within 10-40 eV. As a consequence, these electrons can leave a submicroscopic sample in the same manner as outer-shell photoelectrons. The inelastic mean free path of these electrons is of the order of a few hundred Ångstroms (Ashley, 1990 and references therein; Ziaja et al., 2001), and thus the damage by the departing electrons will only have to be considered in large samples where such electrons may become trapped and deposit energy.

Low energy shake-off electrons will behave similarly to Auger electrons, and may cause secondary ionisation in the sample. A detailed analysis of secondary electron cascades elicited by low energy electrons can be found in Ziaja et al., 2001.

## Materials and methods

The average total shake contributions were estimated from quantum chemical (QC) calculations, using the semi-empirical INDO/S-CI program ZINDO, developed by Zerner and co-workers (Ridley & Zerner, 1976; Bacon & Zerner, 1979; Zerner et al., 1980), with standard parameterisation. A single determinant description was used for the initial (unionized) state of the system, while the excited core-ionized states were obtained from Configuration-Interaction calculations with Single excitations (CIS). The equivalent core approximation was used to represent the core ionized species.

The calculations are based on the sudden approximation (Aoberg, 1967) valid for high energy photons, in which the intensity of a given shake line is proportional to the square of the overlap $`\mathrm{\Psi }_0|\mathrm{\Phi }_n`$ :

$$P_n=|\mathrm{\Psi }_0|\mathrm{\Phi }_n|^2,$$
(1)
where $`\mathrm{\Psi }_0`$ is the wave function of the N-1 remaining electrons in the neutral molecule, and $`\mathrm{\Phi }_n`$ is the (N-1)-electron wave function of the nth state of the ionized system. These intensities were calculated with the SHAKEINT (Lunell, 1987; Lunell & Keane, 1988) program package. This method has been successfully applied to a wide range of systems, including organic molecules (c.f. Lunell et al., 1978; Sjögren et al., 1992), C60 (Enkvist et al., 1995), polymers (PMDA-ODA polyimide) (Nakayama et al., 1999), and adsorbates on metals (CO adsorbed on Cu(100) surfaces) (Persson et al., 2000).

## Acknowledgements

We thank Svante Svensson for illuminating discussions and David van der Spoel for excellent comments. This work was sponsored by the Swedish Natural Science Research Council (NFR), the Swedish Research Council for Engineering Sciences (TFR), the EU-BIOTECH Programme and in part by the Polish Committee for Scientific Research with grants Nos. 2 P03B 04718, 2 P03B 05119 to B. Z. A. S. was supported by a STINT distinguished guest professorship. B. Z. is supported by the Wenner-Gren Foundations.

Figure 1. The model polypeptide for calculating the shake contributions for carbon, nitrogen, and oxygen atoms (Fig. 1a). Shake effects for sulphur were estimated in separate calculations for the three most common chemical environments of sulphur in proteins (Cys, Met, Cys-S-S-Cys; Figs. 1b-d).


text:
# The 𝜎-Meson Production in Excited Υ Decay Processes

## 1 Introduction

Whether the light $`\sigma `$-meson really exists or does not is an important problem in hadron physics. For many years, its existence had been neglected phenomenologically mainly due to the negative result of the analyses of $`I=0`$ $`S`$-wave $`\pi \pi `$ scattering phase shift.

In many $`\pi \pi `$-production experiments, a large event concentration or a bump structure in the spectra of $`\pi \pi `$ invariant mass $`m_{\pi \pi }`$ around 500 MeV had been observed, however, conventionally it was not regarded as $`\sigma `$-resonance, but as a mere $`\pi \pi `$-background, under influence of the so callled “universality argument.” In this argument, it is stated that because of the unitarity of $`S`$-matrix and of the analyticity of the amplitudes, the $`\pi \pi `$ production amplitude $``$ takes the form $`=\alpha (s)𝒯`$ ($`𝒯`$ being the $`\pi \pi `$ scattering amplitude), with a slowly varying real function $`\alpha (s)`$. The pole position of $`S`$-matrix is determined solely through the analysis of $`𝒯`$, which was believed to have no light $`\sigma `$-pole at that time.

Recently the data of $`\pi \pi `$-scattering phase shift had been reanalyzed by many groups including ours and the existence of light $`\sigma (450600)`$ was strongly suggested. The result of no $`\sigma `$-existence in the conventional analyses was pointed out to be due to the lack of consideration on the cancellation mechanism guaranteed by chiral symmetry, and shown to be not correct. Furthermore, we have pointed out that the “universality argument” should be revised, taking into account the quark physical picture of hadrons: The essential point is that the strong interaction is a residual interaction of QCD among all color-neutral bound states of quarks, anti-quarks and gluons, and accordingly the requirement of unitarity and analyticity must be made on the $`S`$-matrix elements with the right “quark physical bases,” including not only the stable $`\pi `$ meson but also the “stable” $`\sigma `$ meson (which is also considered to be a $`q\overline{q}`$ bound state). The universality argument (and the conventional application of final state interaction (FSI) theorem) with the $`S`$-matrix bases including only stable particles, is not correct: Since this $`\sigma `$-meson has, in principle, an independent production coupling generally with a strong phase from the $`2\pi `$ system.

Thus, the production amplitude $``$ has “independent” properties from the scattering amplitude $`𝒯`$, except that the pole position of resonant particles should be universal through both amplitudes. Accordingly, the data of many $`\pi \pi `$-production experiments, which had been analyzed following the universality argument, must be reanalyzed independently from $`𝒯`$ by including the effect of $`\sigma `$-production. We parametrize $``$ as a sum of Breit-Wigner amplitudes for the relevant resonances (including $`\sigma `$-resonance) and of the direct $`2\pi `$ production amplitude, following the VMW method. In the VMW method the above mentioned universal pole singularity is explicitly taken into account, and this method is consistent with the unitarity of the total $`S`$-matrix (and also with the FSI theorem, applied correctly) with the above metioned right bases. The VMW method had already been applied to the processes of $`J/\psi \omega \pi \pi `$ decay, $`p\overline{p}3\pi ^0`$ annihilation by us, and essentially similar method to the process $`D^{}\pi ^{}\pi ^+\pi ^+`$ by E791. The large event concentration in low $`m_{\pi ^+\pi ^{}}`$ region in $`J/\psi \omega \pi \pi `$ and $`D^{}\pi ^{}\pi ^+\pi ^+`$ was satisfactorily explained by $`\sigma `$-production, while in $`p\overline{p}3\pi ^0`$ clear evidences for $`\sigma `$-meson prodution were shown. In this paper we apply this method to the analyses of the hadronic decays of excited $`\mathrm{{\rm Y}}`$, $`\mathrm{{\rm Y}}(2S,3S)`$, and $`\psi (1S,2S)`$.

## 2 Method of the analyses

We analyze the $`m_{\pi \pi }`$ spectra of the processes, $`\mathrm{{\rm Y}}(2S)\mathrm{{\rm Y}}(1S)\pi ^+\pi ^{}`$, $`\mathrm{{\rm Y}}(3S)\mathrm{{\rm Y}}(1S)\pi ^+\pi ^{}`$, $`\mathrm{{\rm Y}}(3S)\mathrm{{\rm Y}}(2S)\pi ^+\pi ^{}`$ and $`\psi (2S)J/\psi \pi ^+\pi ^{}`$ following the VMW method in one-channel form. The $``$ is given by a coherent sum of $`\sigma `$-Breit-Wigner amplitude and of direct $`2\pi `$-production amplitude as

$``$ $`=`$ $`{\displaystyle \frac{e^{i\theta _\sigma }r_\sigma }{m_\sigma ^2si\sqrt{s}\mathrm{\Gamma }_\sigma (s)}}+r_{2\pi }e^{i\theta _{2\pi }};`$ (2.1)
$`\mathrm{\Gamma }_\sigma (s)={\displaystyle \frac{g_\sigma ^2p_1(s)}{8\pi s}},p_1(s)=\sqrt{{\displaystyle \frac{s}{4}}m_\pi ^2},`$
where $`r_\sigma (r_{2\pi })`$ is the $`\sigma (2\pi )`$ production coupling and $`\theta _\sigma (\theta _{2\pi })`$ corresponds to the production phase. These parameters are process-dependent. The $`J/\psi \varphi \pi ^+\pi ^{},\varphi K^+K^{}`$ process is analyzed by VMW method in two-channel form, where $`f_0(980),f_0(1370),f_0(1500)`$, as well as $`\sigma `$, are included. All the relevant spectra are fitted by using common values of the parameters, the mass of $`m_\sigma `$ and the $`\sigma \pi \pi `$ coupling $`g_{\sigma \pi \pi }`$.

## 3 Results and Conclusion

The results of our analysis are shown in Fig. 1. The spectra of six different processes are reproduced well by the interference between the $`\sigma `$-Breit-Wigner amplitude with $`m_\sigma =526_{37}^{+48}\mathrm{MeV},\mathrm{\Gamma }_\sigma =301_{100}^{+145}\mathrm{MeV}`$, and the direct $`2\pi `$ production amplitude. The corresponding pole position is $`\sqrt{s_{\mathrm{pole}}}=(535_{36}^{+48})i(155_{53}^{+76})`$MeV. The total $`\chi ^2`$ is $`\chi ^2/(N_{\mathrm{data}}N_{\mathrm{param}})=86.5/(15037)=0.77`$. The contribution of $`\sigma `$ and of direct $`2\pi `$ production amplitudes to the spectra are given, respectively, by dot and dot-dashed lines in this figure. It is notable that the destructive interference between $`\sigma `$ amplitude and $`2\pi `$ amplitudes explains the suppression of the spectra in the threshold region of $`\mathrm{{\rm Y}}(2S1S)`$ and $`\psi (2S1S)`$ decays, while in $`\mathrm{{\rm Y}}(3S1S)`$ decay these two amplitudes interfere constructively, and the steep increase from the threshold is reproduced. These threshold behaviors of the production amplitudes are shown to be consistent with the restriction from chiral symmetry. It is especially interesting that the double peak structure, with the bottom around the $`\sigma `$-peak position, of the spectra in $`\mathrm{{\rm Y}}(3S1S)`$ decays is also reproduced well by the above interference between the direct $`2\pi `$ production amplitude with zero phase and the $`\sigma `$-production amplitude with a moving phase. That is, we are observing the very phase motion of the $`\sigma `$-Breit-Wigner formula through the variation of the spectra in this process.

The obtained values of masses and widths of $`\sigma `$-meson given above are almost consistent with the results obtained in our previous $`\pi \pi `$ phase shift analysis,
$`m_\sigma `$=535$``$675MeV and $`\mathrm{\Gamma }_\sigma `$=385$`\pm `$70MeV. These results give strong evidence for existence of the light $`\sigma `$(450–600).

The consistency of our result of analyses with the general constraints from chiral symmetry is discussed in the separate talk.

In this conference, there was raised a doubt on the $`\sigma `$-existence along the line of the universality argument;
i) on $`J/\psi \omega \pi \pi `$ process. The experimental cos$`\theta `$-distributions show the parabolic shape around cos$`\theta 0`$ in the range of $`m_{\pi \pi }`$, 300 through 700 MeV, coming from the interference between the $`S`$ and $`D`$ waves. However, the reverse of the direction of the parabolic shape, which is expected from the relative phase motion, to the background $`D`$-wave phase, of the $`\sigma `$-Breit-Wigner amplitude, was not observed. This criticism will become invalid, if we take into account the possible strong phase. For example, $`\theta _\sigma `$=$`90^{}`$.
ii) on $`pp`$-central collision, $`pppp\pi \pi `$. The result of the partial wave analysis by WA102 shows that the non-resonant $`\pi \pi `$-production from one pion exchange is enough to reproduce the $`S`$\- and $`D`$-wave components in low $`m_{\pi \pi }`$ region. We consider that it is necessary, in addition to this, to consider the effect of the $`\sigma `$-Breit-Wigner amplitude, since the relative phase between $`S`$\- and $`D`$-wave amplitudes shows a structure around $`m_{\pi \pi }500`$MeV region, implying the interference between the non-resonant $`r_{2\pi }e^{i\theta _{2\pi }}`$ and $`r_\sigma e^{i\theta _\sigma }`$ in our formula Eq. (2.1). By the way, the experimental relative phase is completely different from that in the $`\pi \pi `$-scattering. This fact clearly shows the universality-argument (where $`\alpha =/𝒯`$ is supposed to be real) is not correct.
iii) on $`\pi \pi `$-scattering. No rapid phase variation due to $`\sigma (600)`$ has been observed. We had mentioned in many occasions the reason why the $`\sigma `$-Breit-Wigner phase motion is not directly observed in $`\pi \pi `$-scattering: Because of the constraints from chiral symmetry the $`\sigma `$-amplitude must be strongly cancelled out by the non-resonant repulsive $`\pi \pi `$-amplitude in $`\pi \pi `$-scattering.


text:
# Synthesis, Crystal Structure and Magnetic Properties of the Linear-Chain Cobalt Oxide Sr5Pb3CoO12

## I Introduction

In order to reveal nature of correlated electrons in condensed matters and to construct advanced models for those, one dimensional (1D) electronic system as a basis of the models has been subjected for both experimental and theoretical investigations. In the experimental part, a variety of quantum phenomena has been indeed found for the quasi-1D compounds in past few decades; Tomonaga-Luttinger-type electric conductivity for BaVS3 PRB94MN , Peierls instability for CuGeO3 PRL93MH , charge- and spin-density waves for (TMTSF)2AsF6 and TTF-TCNQ PRB82KM ; JPSJ76SK , spinon-holon separation for SrCuO2 PRL96CK , Haldane gap for Ni(C2H8N2)2NO2(ClO4) PRL89AY , and superconductivity for (TMTSF)2PF6 JPL80DJ . These characteristic phenomena might reflect the nature of the correlated electrons, and intensive research on the materials played an important role to steadily advance understanding of the nature. We have recently been exploring novel quasi-1D compounds in order to find additional systems showing correlations among their magnetic, electronic transport properties, and crystal structure.

The linear-chain cobalt oxide Sr5Pb3CoO12 was recently discovered in a course of studies of quasi-1D magnetic materials and electrical-carrier-doped those PRB01KY ; SSC00KY ; JSSC99KY ; JSSC90JSK ; JSSC91TGNB . The polycrystalline sample of the compound was obtained by high-temperature solid-state reaction, and then subsequently investigated by x-ray diffraction, magnetic susceptibility, and thermogravimetric analysis (TGA) studies. Crystal structure was investigated in detail by powder-neutron-diffraction at 6 and 295 K. The compound was found isostructural to the 1D antiferromagnetic copper oxide Sr5Pb3CuO12 \[P6̄2m, $`a=`$ 10.1089(6) Å and $`c=`$ 3.5585(2) Å\], in which distorted CuO4 units are linearly connected by sharing those conner-oxygen PRB01KY ; JSSC90JSK ; JSSC91TGNB . Spin-singlet ground state, which potentially results from linear-alternation of quantum spins, was indeed suggested at the composition Sr5Pb2.6Bi0.4CuO12 PRB01KY . Although the title antiferromagnetic compound Sr5Pb3CoO12 was therefore expected to present one dimensionally anisotropic magnetism, however Curie-Weiss(CW)-type magnetism was eventually observed even at low temperature ($`>`$ 5 K). The rather isotropic magnetic property probably indicates a lack of 1D-magnetic uniformity, which may be affected by high degree of structural disorder in chains, as in the Bi-nondoped copper analogue PRB01KY . In this paper, the crystal structure and magnetic properties of the novel chain compound are reported and compared with those of the analogous copper oxides.

## II Experimental

The samples of the cobalt oxide were prepared as follows. A mixture of pure ($``$ 99.9 %) and fine powders of SrCO3, PbO, and Co3O4 (Sr:Pb:Co = 5:3:1 molar ratio, $``$7 g) was placed in a dense-alumina crucible with a cap and heated in air at 750 C for 19 hours. The sample was quenched at room temperature and then ground. Thereafter, the powder was reheated at 850, 900, and 950 C in turn for 66 hours in total. The preheated powder was then molded into several pellets, followed by heating at 970 C for 73 hours in air. Subsequently, some of the sintered pellets were annealed at 500 C in mixed gas, 20 % oxygen in argon, at 100 MPa for 5 hours in a commercial apparatus (hot-isostatic-pressing system, developed by KOBE STEEL, LTD). To qualitatively analyze the polycrystalline samples, x-ray diffraction at room temperature with CuKα radiation was employed. The x-ray apparatus (RINT-2000 system, developed by RIGAKU, CO) was equipped with a graphite monochromator on the counter side and an auto-divergence-slit system. The oxygen content of the samples was measured in a commercial TGA apparatus (PYRIS 1, developed by PerkinElmer, Inc) by heating small amount of each sample (powder, $``$20 mg) in 3 % hydrogen in argon at a heating rate of 5 C per minute to 700 C and holding for 6 hours.

To investigate the crystal structure further, the sample annealed in the compressed gas was again studied by neutron diffraction. The neutron data at 6 and 295 K were obtained by the BT-1 high-resolution powder diffractometer at the NIST Center for Neutron Research. A Cu(311) monochromator was employed to produce a coherent neutron beam ($`\lambda =1.5401`$ Å) with 15, 20, and 7 collimators before and after the monochromator, and after the sample, respectively. The neutron diffraction patterns were measured between 8 and 160 degrees in diffraction angle at 0.05 degrees step. With the neutron profiles, crystal structure parameters were refined to a high degree of agreement by Rietveld calculations on the program GSAS LANLR90ACL . The total numbers of reflection and data point were 172 and 3039, respectively. Neutron scattering amplitude in the calculations were set 0.702, 0.940, 0.253, and 0.581 ($`\times `$10-12) cm for Sr, Pb, Co, and O, respectively LANLR90ACL .

The magnetic properties of the samples were studied by a commercial apparatus (MPMS system) between 5 and 390 K, developed by Quantum Design, Inc. The magnetic susceptibility data were collected at 50 kOe on cooling, and magnetization curves were recorded between -55 and 55 kOe after cooling the each sample at 5 and 150 K. All the pellets thus obtained were too electrically insulating to be subjected for further electronic transport measurements; beyond the 10 M$`\mathrm{\Omega }`$ limit of a conversional two-terminal tester at room temperature.

## III Results and Discussions

The powder x-ray patterns of the annealed sample and the as-made sample are presented in Figs.1a and 1b, respectively. Recalled the hexagonal unit cell of the structure of the analogous copper oxide Sr5Pb3CuO12 (P6̄2m, $`a=`$ 10.1089(6) Å and $`c=`$ 3.5585(2) Å) PRB01KY , a hexagonal unit cell was, at first, tested to qualitatively analyze the both patterns. Almost peaks were clearly found to be at the expected positions from the hexagonal symmetry and lattice parameters, as marked by each $`hkl`$ notation (Fig.1a), except several small peaks indicated by solid stars. The lattice parameters were refined to $`a=`$ 10.11(1) Å and $`c=`$ 3.567(1) Å by a least squares fitting for the annealed sample pattern (Fig.1a), and $`a=`$ 10.12(1) Å and $`c=`$ 3.558(1) Å for the other (Fig.1b). The both patterns indicate quality of the samples. Although the hexagonal model was eventually found reasonable to explain the almost x-ray peaks for the both cobalt oxides as well as the copper oxide, the star-marked minor peaks, however, remained to be uncharacterized, indicating possible presence of either small amount of impurities otherwise somewhat structural modulation. An attempt to prepare a much pure sample, which may not show the extra peaks, by means of optimizing the heating conditions and starting compositions has been made, however, that is unsuccessful thus far.

In order to measure the oxygen content of the cobalt oxide, the TGA study was made on the selected samples; the data are shown in Fig.2, where dotted and solid curves indicate weight change for the as-made and annealed samples, respectively. A clear weight loss was found above approximately 400 C in the both measurements. Based on the hypothesis that the observed weight loss is totally due to oxygen reduction, the oxygen quantity was then calculated to be 12.19 and 12.32 moles per the formula unit for the as-made and annealed samples, respectively. The high-oxygen pressure annealing does not appear to produce a significant increment of oxygen quantity in the as-made cobalt oxide. This fact is supported by the x-ray data; rather small changes (less than 0.25%) in lattice parameters was found after the annealing. As the samples continue to louse those weight after the oxygen reduction completed (null data secured the accuracy of the measurements), probably due to a volatility of lead, a small amount of extra loosing should be superimposed on the major steps, which may lead an over estimation somewhat into the oxygen-quantity calculations.

The annealed sample was subjected to neutron diffraction study at 6 and 295 K to obtain details of the local-chain structure and degree of oxygen nonstoichiometry. The analysis of the neutron data using the Rietveld technique fundamentally followed the way developed in the structural studies on the analogous copper oxide Sr5Pb3CuO12 PRB01KY ; JSSC90JSK ; JSSC91TGNB . As a result, a good achievement was successfully obtained by the Rietveld calculations on the neutron profiles (Fig.3), suggests the average crystal structural model (P6̄2m) is reliable with the cobalt oxide as well as the copper oxide PRB01KY ; JSSC90JSK ; JSSC91TGNB . In the main panel of Fig.3, the observed and calculated profiles are presented at the best quality (5–7 % levels in agreement factors). Shown below the profiles as a deference plot between those, quality of the refinement is indeed on a convincing level. The structure parameters at the best quality are listed in Table 1, and selected interatomic distances and bond angles were calculated from the parameters (Table 2). Lattice constants at 295 K of the hexagonal unit cell are $`a=`$ 10.1093(2) Å and $`c=`$ 3.562 51(9) Å, which essentially meet the x-ray parameters.

During the refinements, temporary fits using unfixed isotropic-atomic-displacement parameters and occupancy factors for O(4) and O(5) had lead only incredible results, probably due to too low levels of those occupancy factors less than 0.1. The situation was not improved even for the low temperature data. The isotropic atomic displacement parameters at the oxygen sites O(4) and O(5) were, therefore, fixed at a conceivable level, 2.0 Å2 (1.2 Å2 at 6 K); calculations thereafter appeared to improve the situation. Due to the same reasons, the occupancy factor of O(5) was constrained to be equal to that of O(4). The occupancy factors of O(2) and O(3), the normally occupied sites, were fixed to be fully occupied in the final calculation because oxygen vacancy was found within one standard deviation of 1.00. The possible mixing between Sr and Pb was tested in preliminary manner, however substantial degree of the mixing was not detected. The displacement parameter of Co is unusually large, probably reflecting the local structure disorder as in the copper oxide PRB01KY ; JSSC90JSK ; JSSC91TGNB . The occupancy factors in the partially occupied sites, O(1), O(4), and O(5), are slightly higher than the expected values 2/3, 1/12, and 1/12, respectively, from the oxygen-stoichiometric composition 12 moles per the formula unit. The oxygen quantity calculated from the present neutron parameters is 12.42 moles per the formula unit. Although the estimation from the neutron data matches that of the TGA data (12.32 moles) within 1 % level, we believe the oxygen quantity may be slightly overestimated because the occupancy factors of O(4) and O(5) are too low to be accurate. As the neutron diffraction study on the copper oxide did not help to figure out the probable local-structure disorders and oxygen nonstoichiometry PRB01KY ; JSSC90JSK ; JSSC91TGNB , further attentions on the present analysis, including calculations with anisotropically thermal parameters of metals, did not shed light on the problems, either. The oxygen composition of the cobalt oxide might be slightly superstoichiometric ($``$ 12.4), but not far from the stoichiometric Sr5Pb3CoO12. To determine the oxygen quantity conclusively, further investigations would be required after much high-pure sample becomes available.

Schematic structural views of Sr5Pb3CoO12 were drawn from the neutron data (295 K) in Figs.4a and 4b. It is clear that the cobalt compound has a linear-chain structure basis, in which cobalt–oxygen polyhedra form chains at inter-chain distance of 10.1 Å PRB01KY . A part of the chain along $`c`$axis is presented in Fig.4b. All of cobalt and oxygen sites are shown by the dotted circles smaller and larger, respectively, and a probable arrangement of those atoms is indicated by solid and filled circles. Intolerably short bond distances (Table 2) between cobalt and oxygen atoms were precluded from making the arrangement. The irregular arrangement of cobalt atoms, as shown in Fig.4b, and probable presence of local displacements of atoms may account for the observed unusually large atomic displacement parameters in the average-structure data such as 3.8 Å2 and 2.9 Å2 for Co and O(1), respectively. The $`c`$axis constant (3.56251(9) Å) reflects the average of Co to Co distances in the chains; almost are close to $``$3.56 Å, but few parts at $``$2.61Å may be somewhat involved. Since the oxygen-coordination environment of all cobalt atoms is not unique, a degree of magnetic uniformity of the present chain compound might be as low as that of the copper oxide Sr5Pb3CuO12, far rather than other well-studied spin-chain compounds PR64JCB ; PRB83JVCB ; PRB87DCJ ; PRB94TB .

The low-temperature structure was investigated at 6 K (inset of Fig.3) in the same manner, and the obtained structure parameters, calculated interatomic distances and bond angles are summarized in Tables 1 and 2 as well. The hexagonal unit cell slightly shrinks (0.21 % along $`a`$axis and 0.17 % along $`c`$axis) by the cooling, while the Co site splitting is much conspicuous than that at 295 K (+11.5 % in splitting distance) as shown in Figs.5a and 5b. Any traces of magnetic ordering were not detected in the diffraction profiles.

The magnetic-susceptibility data of the samples annealed (open circles) and as-made (closed circles) are shown in Fig.6; $`T`$ vs $`1/\chi `$ and $`T`$ vs $`\chi `$ are plotted in the main panel and the inset, respectively. The data of magnetic field dependence of the magnetization at 5 and 150 K are shown in Fig.7. Seen in the main panel of Fig.6, the inverse magnetic susceptibility is linearly dependent on temperature above approximately 130 K for both samples. The CW law was, then, applied to analyze the two sets of the high-temperature portion of the data. The formula to fit the data by a least-squares method, as indicated by the solid lines (main panel of Fig.6), was

$`\chi (T)={\displaystyle \frac{N\mu _{\mathrm{eff}}^2}{3k_\mathrm{B}(T\theta _\mathrm{W})}}+\chi _0,`$ (1)
where $`\theta _\mathrm{W}`$, $`\chi _0`$, $`\mu _{\mathrm{eff}}`$, $`k_\mathrm{B}`$, and $`N`$ were Weiss temperature, temperature-independent term, effective magnetic moment, Boltzmann constant, and Avogadro’s constant, respectively. In a preliminary fit, the temperature-independent term $`\chi _0`$ was estimated to be on the order of 10-6 emu/mol-Co, which was negligibly small, and then we decided not to employ the second term in Eq.1 hereafter. The subsequent fits yielded $`\mu _{\mathrm{eff}}`$ of 3.64 $`\mu _\mathrm{B}`$ per Co and $`\theta _\mathrm{W}`$ of -38.8 K for the annealed sample, and $`\mu _{\mathrm{eff}}`$ of 3.83 $`\mu _\mathrm{B}`$ per Co and $`\theta _\mathrm{W}`$ of -44.5 K for the as-made sample, those parameters were characteristic of antiferromagnetic interaction. As the oxygen stoichiometry is approximate at about 3% level (12 to 12.4 moles per the formula unit), the formal valence of oxygen also becomes approximate about +2 to +2.8, when the valence of lead is +4. Due to the rather large degree of uncertainty of the valence of cobalt at about 40% level, it is unable to determine the electronic configuration of cobalt solely from the magnetic susceptibility data. Further investigations, including spectroscopic studies, might help to make clear the issue.

Seen in Fig.6, a characteristic maximum expected from 1D antiferromagnetism in $`T`$ vs $`\chi `$ plot for either Bonner-Fisher or spin-gap models was not observed at all, although the chains were substantially antiferromagnetic PR64JCB ; PRB83JVCB ; PRB87DCJ ; PRB94TB . Below approximately 100 K, the inverse magnetic susceptibility gradually starts to apart from the CW line to go down on cooling, indicative of gradual growth of a ferromagnetic component in short range. To study the low-temperature portion further, the Brillouin function was employed to make comparison with the magnetization data at 5 K (Fig.7). The function to calculate the ideal $`M_S`$ vs $`H`$ curve at 5 K for free spins ($`S`$) was

$`M_S(\alpha )=NgS\mu _\mathrm{B}`$ (2)
$`\times \left({\displaystyle \frac{2S+1}{2S}}\mathrm{coth}{\displaystyle \frac{2S+1}{2S}}\alpha {\displaystyle \frac{1}{2S}}\mathrm{coth}{\displaystyle \frac{1}{2S}}\alpha \right),`$
where $`\alpha =gS\mu _\mathrm{B}H/k_\mathrm{B}T`$. The effective $`S`$ was estimated to be $``$1.39 for the annealed and $``$1.48 for the as-made samples from calculations using the formula $`\mu _{\mathrm{eff}}=2\sqrt{S(S+1)}`$ and the $`\mu _{\mathrm{eff}}`$ parameters, obtained by the CW fits. The computed curves for free spins at the both $`S`$ vales are shown in Fig.7. Due to the substantial antiferromagnetic interactions, the magnetization at 50 kOe for the cobalt oxides is reduced to approximately one third of that for the ideal-free-spin system. Although the absolute Weiss temperatures are much higher than the studied temperature 5 K, the magnetization depends on applied field with rather positive curvature than linear relationship, which may indicate somewhat an influence of the ferromagnetic spin coupling (Fig.7). These observations, including the deviation from the CW line on cooling, suggest that minor ferromagnetic interactions coexist with the antiferromagnetism, and might help to destroy the potential 1D-antiferromagnetic characters. Since the chain involves local-structure disorder and potentially be modulated, characters of all probable magnetic bonds were unable to be uniquely determined in this way. The observed $`\theta _\mathrm{W}`$, which is apparently negative, is probably due to a balance between the major antiferromagnetic and the possible minor ferromagnetic interactions. The observed ferromagnetic component, otherwise, results from minor magnetic impurities.

In summary, we have reported the average structure and the magnetic properties for the quasi-1D novel cobalt oxide. Due to the probable presence of local-structure disorders such as the local displacement of cobalt and ligand oxygens, and the irregular arrangement of cobalt-oxygen polyhedra in short range, degree of magnetic coherence should not be as high as expected to the low dimensional magnetism. The local structure is too complicated to be clearly observed by the neutron diffraction which probes a positional average of the structure. Further studies with local structural and magnetic probes, more sensitive to the microscopic magnetic environment, could clarify the probable coupling between the local-structure disorders and the rather unusual magnetism.

###### Acknowledgements.

This research was supported in part by the Multi-Core Project administrated by Ministry of Education, Culture, Sports, Science and Technology of Japan.


text:
# Thermodynamics with 2+1 and 3 Flavors of Improved Staggered Quarks Preprint FSU-CSIT-01-50. Talk given by U.M. Heller at “Statistical QCD”, August 26–30, 2001, Bielefeld, Germany.

## 1 INTRODUCTION

With the Relativistic Heavy Ion Collider (RHIC) now producing data, it has become even more important to understand the phase diagram of QCD at finite temperature, and to determine properties of the high temperature quark-gluon-plasma phase with confidence, i.e. with controlled lattice spacing errors.

It is fairly well established that QCD with two flavors of massless quarks has a second order finite temperature, chiral symmetry restoring phase transition. This transition is washed out as soon as the quarks become massive. QCD with three flavors of massless quarks has a first order finite temperature, chiral symmetry restoring phase transition, which is stable for small quark masses. Not well known is how large the quark masses can be before the phase transition turns second order and then into a crossover, both for degenerate quarks and especially for the physically relevant case of two light and one heavier strange quark.

In previous studies, the second question is particularly badly answered due to the flavor symmetry breaking in Kogut-Susskind quarks, usually used for this purpose: how can one study the influence of the strange quark when most of the (non-Goldstone) pions are heavier than the (Goldstone) kaon?

Adding a few terms to the conventional Kogut-Susskind action, namely three-link, five-link and seven-link staples and a third-neighbor coupling, removes all tree-level $`𝒪(a^2)`$ errors . This “Asqtad” action shows improved flavor and rotational symmetry , and, at least in quenched QCD, good scaling properties . This is illustrated in Fig. 2. There, and in other figures below, we plot results in units of $`r_1`$, a scale defined in terms of the static $`\mathrm{Q}\overline{\mathrm{Q}}`$ potential by $`r_1^2F_{\mathrm{Q}\overline{\mathrm{Q}}\mathrm{static}}(r_1)=1`$, which leads to $`r_10.35`$ fm .

The three-link “Naik” term insures a good dispersion relation and thereby helps decrease the lattice artifacts in energy density and pressure at small temporal lattice size $`N_t`$, as can be seen in Fig. 2.

Based on the zero temperature simulations of Ref.  and the estimate of $`T_c150170`$ MeV  we deduce that for $`N_t=810`$ the kaon will be heavier than the heaviest non-Goldstone pion at the finite temperature transition.

We have zero temperature results, in particular, the value of the (bare) strange quark mass, at fixed lattice spacing $`a0.13`$ fm and $`a0.2`$ fm. Since we want to keep the physical quark masses approximately constant, when we vary the temperature which, at fixed $`N_t=1/(aT)`$, means varying the gauge coupling $`\beta `$, we interpolated (extrapolated) between the values at the two lattice spacings.

Given a gauge coupling $`\beta `$ and quark mass $`am_q`$, for the three-flavor simulations, or quark masses $`am_s`$ and $`am_{u,d}`$, for $`2+1`$ flavors, we determined the tadpole factor $`u_0=(\text{Tr}U_p/3)^{1/4}`$ self-consistently in short zero temperature simulations on $`L^4`$ lattices.

To convert lattice scales to physical scales we interpolated $`a/r_1`$ with the form advocated by Allton

$`a/r_1=c_0f(g_0^2)\left[1+c_2g_0^2f^2(g_0^2)\right],f(g_0^2)=(b_0g_0^2)^{b_1/(2b_0^2)}\mathrm{exp}\left({\displaystyle \frac{1}{2b_0g_0^2}}\right).`$ (1)
For the (one-loop, tadpole) Symanzik improved gauge action the bare coupling $`g_0^2`$ is related to $`\beta `$ by $`g_0^2=10/\beta `$. $`b_{0,1}`$ are the universal first two coefficients of the beta-function for three massless flavors, and $`c_{0,2}`$ are determined from the measured values of $`a/r_1`$ at $`a0.13`$ fm and $`a0.2`$ fm.

The Hybrid Molecular Dynamics R-algorithm was used for all simulations, with two noise vectors for the $`2+1`$ flavor simulations.

## 2 RESULTS FOR THREE DEGENERATE FLAVORS

For our first simulation with $`m_qm_s`$ we used a linear interpolation of $`am_q(\beta )`$ in $`\beta `$. Since the range of interpolation/extrapolation is quite large, we subsequently used a linear interpolation of $`\mathrm{log}(am_q(\beta ))`$, compatible with the leading behavior of $`a`$ (with $`m_q`$ fixed) as function of $`\beta `$. In the region of the finite temperature transition or crossover for lattices with temporal extent $`N_t=6`$ and 8, the difference between the two interpolation schemes is small enough that it will not affect our results significantly.

We show in Fig. 4 the real part of the Polyakov line and in Fig. 4 the condensate as function of the temperature, with the physical scale obtained via $`r_1`$, for our simulations with three degenerate quarks of mass $`m_qm_s`$ and $`m_q0.6m_s`$. The Polyakov line shows a crossover from confined behavior at low temperature to a deconfined behavior at high temperature. The condensate decreases as the temperature is increased, but shows no sign of a transition or sharp crossover.

In Fig. 5 we show the chiral susceptibility, $`\chi _{\overline{\psi }\psi }=\overline{\psi }\psi /m_q`$ from the simulations with $`m_qm_s`$ for lattices with $`N_t=6`$ and 8. A peak is seen in both cases. But the peak height does not increase when we change the spatial lattice size from $`N_s=2N_t`$ to $`3N_t`$, indicating that we are observing a smooth crossover from a chirally symmetric phase at high temperature (large $`\beta `$) to a phase of broken chiral symmetry at low temperature. The location of the peak coincides approximately with the steep increase of the Polyakov line in Fig. 4. The chiral susceptibility for $`m_q0.6m_s`$ looks similar to Fig. 5, but in that case we have, so far, results only for one spatial size, $`N_s=2N_t`$.

## 3 RESULTS FOR 2+1 FLAVORS

We have also preliminary results from simulations with 2+1 dynamical flavors, with the heavier quark kept at the strange quark mass, and the lighter masses at $`m_{u,d}=0.6m_s`$ and $`0.4m_s`$. The real part of the Polyakov line and the light quark condensate from these simulations are shown in Fig. 6. Again we see a crossover behavior, particularly for the Polyakov line. For these simulations we do not have the data yet to compute the chiral susceptibility.

## 4 CONCLUSIONS

We have made first simulations to explore the finite temperature phase diagram with an improved staggered fermion action, the “Asqtad” action, which reduces flavor symmetry breaking so that all pions are lighter than the kaon already at larger lattice spacing, and which improves rotational symmetry and the dispersion relation, leading to diminished lattice artifacts in energy density and pressure.

For three flavor simulations down to quark mass $`m_q0.6m_s`$, and for 2+1 flavor simulations with the light quark masses down to $`m_{u,d}=0.4m_s`$ while the heavier quark is kept at the strange quark mass, we observed finite temperature crossover behavior, but so far no sign of a genuine phase transition. This result is compatible with other recent simulations  which found phase transitions only at quark masses lighter than those studied by us so far. The temperature at the crossover, with the scale set by $`r_1`$, $`T_c190`$ – 200 MeV, is a little higher than expected from previous determinations . This is presumably also due to our quark masses still being rather high.

This work is supported by the US National Science Foundation and Department of Energy and used computer resources at Florida State University (SP), NERSC, NPACI, FNAL, and the University of Utah (CHPC).