text:
# The Effect of Neutron Star Rotation on the Properties of Thermonuclear X-ray Bursts

## 1 Introduction

Nearly-coherent millisecond brightness oscillations with frequencies between 270–620 Hz have been observed during thermonuclear X-ray bursts from 13 neutron star low-mass X-ray binaries (see Strohmayer & Bildsten 2003, for a review). Two of these systems also show persistent millisecond pulsations at the same frequencies in their non-burst emission (Chakrabarty et al., 2003; Strohmayer et al., 2003). It is therefore widely accepted that the oscillations result from brightness patterns that form on the surfaces of these rapidly-rotating neutron stars during X-ray bursts, thus probing two very different pieces of physics: the distribution of neutron star spin frequencies (e.g., White & Zhang, 1997; Bildsten, 1998; Chakrabarty et al., 2003), and how unstable nuclear burning proceeds on a neutron star’s surface (e.g., Nath, Strohmayer, & Swank 2002; Spitkovsky, Levin, & Ushormirsky 2002; Muno, Özel, & Chakrabarty 2003b).

In the sources that exhibit the millisecond burst oscillations, they are only detected from about half of all X-ray bursts. We previously showed that if a distinction is drawn between sources rotating at $`300`$ Hz and $`600`$ Hz based on the timing properties of their persistent emission (see van der Klis, 2000, for a review), then there is a physical difference in the properties of the bursts that exhibit oscillations (Muno et al., 2001). The fast, $``$600 Hz oscillations almost always occur during the strongest X-ray bursts, during which the photosphere of the neutron star is driven to a large radius by radiation pressure. In contrast, the slow, $`300`$ Hz oscillations occur with equal likelihood in bursts with and without photospheric radius expansion. Observational selection effects or differences in viewing angles cannot by themselves produce the correlations between burst properties and the presence of oscillations.

The proposed explanation for this difference has two parts (Franco, 2001; Muno et al., 2001): (1) that burst oscillations are almost exclusively detected when the accretion rates onto the neutron stars are relatively high ($`0.1\dot{M}_{\mathrm{Edd}}`$), even though the X-ray bursts themselves have also been observed at significantly lower accretion rates, and (2) that the properties of bursts at these high accretion rates are different in the fast and slow rotators. This explanation was motivated largely by observations of only two systems (4U 1728$``$34 and KS 1731$``$260), and still requires confirmation with a larger sample of sources. After seven years of operation, there is now sufficient data in the archive of observations taken with the Rossi X-ray Timing Explorer (RXTE) Proportional Counter Array (PCA; Jahoda et al., 1996) to examine this hypothesis for several sources. Therefore, in this paper we examine how the presence of millisecond oscillations and of photospheric radius expansion in thermonuclear X-ray bursts are related to the persistent accretion rates onto the neutron stars.

## 2 Observations and Results

We obtained public RXTE PCA observations of 12 out of 13 sources of millisecond burst oscillations.111We omitted MXB 1743$``$29 from our analysis because its persistent emission cannot be isolated from the several sources within 1 of the Galactic center. We searched for X-ray bursts using the using the algorithm

explained in Galloway et al. (in preparation). As of August 2003, we detected 394 bursts from these 12 sources.

### 2.1 Persistent Accretion Rate onto the Neutron Star

We first characterized the persistent accretion rates onto the neutron stars ($`\dot{M}`$) using X-ray color-color diagrams. Using data with 128 energy channels between 2–60 keV and 16 s time resolution, we defined soft and hard colors as the ratio of the background-subtracted detector counts in the (3.6–5.0)/(2.2–3.6) keV and the (8.6–18.0)/(5.0–8.6) keV energy bands, respectively. We used 64 s integrations to calculate the colors when the source intensity was above 100 counts s-1, and 256 s integrations otherwise. We corrected the raw colors to account for changes in the gain of the PCA using linear terms that were determined by assuming that the Crab nebula had a constant count rate in all four bands (Muno, Remillard, & Chakrabarty 2002a).

The resulting color-color diagrams are illustrated using the grey points in Figure 1. As $`\dot{M}`$ onto the neutron star increases, a source moves from the top-left to the bottom-right, roughly tracing a Z-shaped pattern (Hasinger & van der Klis, 1989; Muno et al., 2002a; Gierlinski & Done, 2002a).222Note that all of the burst oscillation sources are traditionally classified as so-called “atoll” sources. They are not to be confused with the so-called “Z” sources, which are more luminous and do not usually exhibit X-ray bursts (see van der Klis, 2000, for a review). The transient systems Aql X-1 and 4U 1608$``$522 were observed over the largest range of accretion rates, and therefore trace a full Z-shaped pattern in the color-color diagrams. The sources 4U 1702$``$429, 4U 1728$``$34, 4U 1636$``$536, and KS 1731$``$260 trace portions of the bottom and diagonal branches of the Z-track (see also Muno et al., 2000; Franco, 2001; van Straaten et al., 2001). For these six sources, we have parameterized the position on the color-color diagram using the arc-length ($`S_Z`$) along the solid curves in Figure 1 (see Dieters & van der Klis, 2000, for details of how this is defined). The remaining sources were observed on only a portion of the Z-track. The sources MXB 1659$``$298 (Wijnands et al., 2002), 4U 1916$``$053 (Boirin et al., 2000), SAX J1748.9$``$2021 (Kaaret et al., 2003), and SAX J1750.8-2900 (Kaaret et al., 2002) were observed almost extensively on the bottom portion of the Z-track. The persistent pulsars SAX J1808.4$``$3658 and XTE J1814-338 were only observed on the top (hard) branch of the Z-track, despite large variations in their luminosities. The color-color diagrams for these sources will be presented elsewhere (Galloway et al., in preparation).

In order to determine the accretion rate at the time the X-ray burst occurred, we calculated the hard and soft colors in the 256 s interval prior to the burst. These are indicated by the blue stars in Figure 1. Most of the conclusions that follow will be based on the six sources in which $`\dot{M}`$ varies over a wide range.

### 2.2 Presence of Millisecond Oscillations

We searched for millisecond oscillations during the X-ray bursts in data recorded with $`2^{13}`$ s (122 $`\mu `$s) time

resolution. We computed fast Fourier transforms of each 1 s interval of data for the first 15 s of the burst, and searched for signals within 5 Hz of the known spin frequencies. We considered a signal to be a detection if it had less than a 1% probability of occurring by chance given the 160 trial frequencies searched for each burst. A signal was considered significant if it passed any of three tests: (1) having a chance of $`<6\times 10^5`$ that it was produced by noise in a single trial, (2) persisting for two adjacent time and frequency bins with a chance probability of $`<(6\times 10^5)^{1/2}/6=1\times 10^3`$, or (3) occurring in the first second of a burst with a chance probability of $`<10^3`$. Millisecond oscillations were detected from 137 bursts in this manner; fewer than 5 ($`<3`$%, conservatively assuming the tests were independent trials) should be spurious detections of noise signals according to our selection criteria. Oscillations are detected in all of the bursts from the millisecond pulsars SAX J1808.4$``$3658 and XTE J1814$``$338 (Chakrabarty et al., 2003; Strohmayer et al., 2003), but are almost exclusively observed at high apparent accretion rates ($`S_Z2`$) in the remaining LMXBs. We have indicated the X-ray bursts with oscillations using filled blue stars in Figure 1.

We examined why the millisecond oscillations are not detected at low accretion rates by comparing the fractional rms amplitudes of the detected oscillations to the upper limits on the non-detections. The fractional rms amplitude $`A`$ is defined as

$$A=\left(\frac{P}{N}\right)^{1/2}\frac{N}{NB},$$
(1)
where $`P`$ is the power in the oscillations, $`N`$ is the total number of counts used to compute the power spectrum, and $`B`$ is the estimated background rate taken from the 16 s prior to the burst (Leahy et al., 1983). For the 6 sources in Figure 1, we plot in Figure 2 the amplitudes as a function of position on the Z-track. Each burst is recorded once in the figure. For detected oscillations, we plot the amplitude of the signal with the largest power. For upper limits, we plot the largest rms amplitude among the noise signals detected within 5 s after the start of the burst decay (most detected oscillations are observed at this time; see Muno et al. 2002b). Although the trend is not absolute, the upper limits on the non-detected signals at low $`\dot{M}`$ are typically a factor of two smaller than the amplitudes of the detected signals at higher $`\dot{M}`$. This indicates that the fractional amplitudes decline as the accretion rates onto the neutron stars decrease.

### 2.3 Photospheric Radius Expansion

We characterized the spectral properties of the X-ray bursts by extracting spectra of the emission at 0.25 s intervals. The data used provided at least 32 energy channels between 2–60 keV. We accounted for the background and non-burst emission using spectra from the 16 s prior to each burst (e.g. Kuulkers et al., 2002). The instrumental responses for the detector units that were active during each burst were generated using the FTOOL pcarsp. We modeled the spectra between 2–25 keV using a blackbody function absorbed by the interstellar medium. The interstellar hydrogen column density was held fixed at the mean value determined from a prior fit in which it was allowed to vary. The model provides an apparent temperature ($`T_{\mathrm{app}}`$) and a normalization equal to the square of the apparent radius ($`R_{\mathrm{app}}`$) of the burst emission surface.

As a proxy for the other properties of the bursts, we adopt the presence of photospheric radius expansion, using the criteria of Galloway et al. (2003). We considered a source to have exhibited radius expansion if (1) $`R_{\mathrm{app}}`$ reached a local maximum close to the time of peak flux, (2) $`R_{\mathrm{app}}`$ was observed to be lower shortly after the flux, with a significance of at least $`4\sigma `$, and (3) coincident with the maximum in $`R_{\mathrm{app}}`$, $`T_{\mathrm{app}}`$ also reached a local minimum. X-ray bursts with radius expansion that are identified with these criteria also tend to have the largest peak fluxes and fluences, and usually rise and decay most quickly (Strohmayer & Bildsten, 2003, Galloway et al., in preparation).

The X-ray bursts with radius expansion are indicated with red circles in Figure 1. Radius expansion bursts occur preferentially at low $`\dot{M}`$ ($`S_Z2`$) in the sources with $`300`$ Hz spin frequencies in Figure 1 (4U 1702$``$429 and 4U 1728$``$34), but not at high $`\dot{M}`$ ($`S_Z2`$). In contrast, radius expansion occurs predominantly at high $`\dot{M}`$ ($`S_Z2`$) in the sources with $`600`$ Hz spin frequencies in Figure 1 (4U 1636$``$536, 4U 1608$``$522, KS 1731$``$260, and Aql X-1). The corresponding trends for the persistent millisecond pulsars are unclear, because they only trace the top branch of the Z-track on the color-color diagram, and all the observed bursts exhibit oscillations. However, the bursts from the 401 Hz pulsar SAX J1808.4$``$3658 all exhibit radius expansion like the fast sources (Chakrabarty et al., 2003), while the bursts from from the 314 Hz pulsar XTE J1814–338 all lack radius expansion like the slow sources (Strohmayer et al., 2003). The trends are less clear for the remaining LMXBs because less data is available, but they are consistent with those in Figure 1 (Galloway et al. in preparation).

## 3 Discussion

In a previous paper, we found that for neutron stars spinning at $`600`$ Hz, millisecond oscillations were observed preferentially in thermonuclear X-ray bursts with photospheric radius expansion, while for neutron stars spinning at $`300`$ Hz, oscillations were equally likely to be observed in bursts with and without radius expansion (Muno et al., 2001). The additional data that has entered the public RXTE archive in the past two years demonstrates that this initial distinction can be explained as a combination of two effects. First, in all of the sources oscillations are preferentially detected when the accretion rates onto the neutron stars are high. Second, the X-ray bursts that occur at high accretion rates in the fast rotators exhibit radius expansion, but do not in the slow rotators.

### 3.1 Millisecond Oscillations as a Function of $`\dot{M}`$

The detectability of millisecond burst oscillations is clearly not determined by the properties of the X-ray bursts that they appear in, because if one considers the entire sample of sources in Figure 1, bursts with and without photospheric radius expansion are equally likely to exhibit oscillations. The properties of the X-ray bursts are probably the best indicators of the conditions in the burning layer, while the low amplitudes of the oscillations (5–10% rms; Muno et al., 2002b) suggest that they are only a minor side-effect of the burning. Therefore, since the the properties of the bursts at low accretion rates are so dramatically different in the fast and slow rotators, it seems unlikely that there is some subtle change in the mechanism producing the oscillations that prevents them from being observed at low $`\dot{M}`$ from all of the sources.

Therefore, we suggest that the mechanism producing the millisecond oscillations always operates during X-ray bursts, but that their amplitudes are attenuated at low $`\dot{M}`$ by a mechanism external to the burning layer. The amount of attenuation required is quite small, since the fractional amplitudes of the detected signals are on average only a factor of 2 larger than the upper limits on the non-detections (Figure 2). The most likely source of attenuation is a corona of electrons that scatters photons from the surface of the neutron star (e.g. Brainerd & Lamb, 1987; Miller, 2000). A corona of optical depth $`\tau 3`$ would be sufficient to attenuate the burst oscillations by a factor of 2 (e.g., Miller, 2000). Such a corona of electrons is thought to produce the high-energy power-law tail that is often present in the X-ray spectra of LMXBs, by inverse-Compton scattering thermal photons form the neutron star (e.g., Barret et al., 2000; Gierlinski & Done, 2002b). This power-law tail is only observed at low $`\dot{M}`$, and is inferred to originate from a corona with temperature $`kT20`$ keV and optical depth $`\tau approx3`$ (Barret & Olive, 2002; Gierlinski & Done, 2002b; Maccarone & Coppi, 2003). Therefore, this corona could be the reason that millisecond burst oscillations are not observed at $`S_Z2`$ in Figure 1.333We note, however, that at high accretion rates the same authors infer the presence of a cooler ($`kT5`$ keV), more optically thick ($`\tau 5`$) Comptonizing corona. Our interpretation would require that at high $`\dot{M}`$ this corona is geometrically arranged so as not to intercept photons from the surface of the neutron star.

The two persistent millisecond pulsars provide the only exceptions to the above trend: all of the burst oscillations, and all of the X-ray bursts, are observed in the hard portion of the Z-track on the color-color diagram that corresponds to low $`\dot{M}`$. The main proposed difference between the millisecond pulsars and the other bursters is that the former have stronger magnetic fields (Chakrabarty et al., 2003). It is therefore plausible that oscillations are observed at low inferred $`\dot{M}`$ for the pulsars because magnetic effects either enhance the oscillation amplitudes, or supress a scattering corona.

### 3.2 X-ray Burst Properties as a Function of $`\dot{M}`$

The properties of X-ray bursts also are correlated with the accretion rate onto the neutron star, but in a manner that additionally depends on the neutron star rotation rate. We expect the burst properties to be determined by two factors (Fujimoto, Hanawa, & Miyaji 1981; Fushiki & Lamb, 1987; Bildsten, 2000). First, as $`\dot{M}`$ increases, the temperature at the burning layer also increases, and thus the column density of helium required to trigger a burst generally decreases. Therefore, if the accretion is spherically symmetric, X-ray bursts that occur at high $`\dot{M}`$ should be weaker.

Second, there is a competition between how quickly a sufficient column density is accumulated such that helium burning is unstable, and how quickly hydrogen in the accreted material can be stably fused into helium. The local accretion rate per unit area ($`\dot{m}`$) actually drives the competition, but the accretion is generally assumed to occur with spherical symmetry. At values of $`\dot{m}`$ thought to correspond to the lower end of those commonly observed from bursting LMXBs (equivalent to global accretion rates of $`0.01\dot{M}_{\mathrm{Edd}}<\dot{M}<0.05\dot{M}_{\mathrm{Edd}}`$), H is burned into He faster than it can be accreted, so the X-ray burst occurs from pure He fuel. As $`\dot{m}`$ increases ($`0.05\dot{M}_{\mathrm{Edd}}<\dot{M}<\dot{M}_{\mathrm{Edd}}`$), H is accreted faster than it can be burned into He, so the bursts occur from mixed H/He fuel.444We note that at $`\dot{M}<0.01\dot{M}_{\mathrm{Edd}}`$ hydrogen burning is unstable, which should also produced mixed H/He bursts. However, in this regime the expected recurrence times (30 h to 10 days; see Bildsten 2000; Narayan & Heyl 2003) are significantly longer than those observed (2–10 h; see Cornelisse et al. 2003). This regime is unlikely to apply to the current data. The X-ray burst properties change between these regimes because the relative amount of H and He in the fuel determines how rapidly the nuclear energy is released during the bursts. Helium burns via a strong triple-$`\alpha `$ process that releases energy quickly, so the low-$`\dot{m}`$ He bursts are more likely to exhibit radius expansion. In contrast, H serves to moderate the He burning at the start of the burst, and only burns through a slow $`rp`$-capture process onto the products of He burning at the end of a burst. Therefore, the high-$`\dot{m}`$ mixed H/He bursts should last longer, and be less likely to exhibit radius expansion.

As a result of these two effects, as the global accretion rate onto the neutron star increases, the X-ray bursts should become weaker and less likely to exhibit radius expansion. This is the case for the slow rotators, but in the fast rotators the bursts are more likely to exhibit radius expansion at high $`\dot{M}`$. One possible explanation for this is that the accretion is not spherically symmetric (Bildsten, 2000). In particular, if the local accretion rate ($`\dot{m}`$) decreases as the global rate ($`\dot{M}`$) increases in the fast rotators, then the high-$`\dot{m}`$ bursts with radius expansion could occur at low $`\dot{M}`$, and the low-$`\dot{m}`$ bursts without radius expansion could occur at high $`\dot{M}`$. In contrast, the change in burst properties in the slow rotators appears consistent with a local $`\dot{m}`$ that increases as the global $`\dot{M}`$ does. It is then possible that the rotation rate of the neutron star influences how accreted material spreads over its surface, either through a lower effective surface gravity or a stronger Coriolis force in the fast rotators.; (e.g. Bildsten, 2000; Spitkovsky et al., 2002).

Although we have determined observationally that the rotation rate of a neutron star influences how the properties of thermonuclear X-ray burst change with the accretion rate onto the neutron star, it is not clear what causes the observed correlations. Further progress should be made by studying this sample of sources to see how the burst time scales, peak fluxes, fluences, and recurrence times change with the accretion rates.

We thank R. Wijnands for providing the flux history of SAX J1808.4$``$3658, and J. Hartman for reducing the data on XTE 1814–338. We also thank the referee for helpful questions and suggestions. This work was supported by NASA under contract NAS 5-30612 and grant NAG 5-9184, and through a Hubble Fellowship to MPM from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555.


format: text

text:
# Magnetic excitations in La1.5Sr0.5CoO4

## Abstract

We report magnetization and neutron scattering measurements of the half-doped compound La1.5Sr0.5CoO4, which exhibits a checkerboard pattern of charge ordering below $``$ 800 K. In the antiferromagnetically-ordered phase below $``$ 40 K the spins are found to be canted in the ab plane. The spin excitation spectrum includes spin-wave excitations with a maximum energy of 16 meV, and diffuse magnetic modes at energies around 30 meV.

preprint: ATB-1

Layered transition metal oxides have been much studied in an attempt to understand the mechanism of high-Tc superconductivity in the cuprates. Recently there has been speculation that charge stripe correlations may play an important role, leading to investigations of compounds such as La2-xSrxNiO4 which exhibits well defined spin and charge stripe ordering tranquada-Nature-1995 .

Here we consider the structurally identical but less-studied family, La2-xSrxCoO4. Previous studies zaliznyak-2000 ; zaliznyak-2001  of the half-doped member La1.5Sr0.5CoO4 have concluded that below $`T_c`$ 800 K the holes introduced by Sr doping form a charge ordered phase, with a loosely correlated checkerboard arrangement of Co3+ and Co2+ ions. It has also been reported that very slightly incommensurate long range magnetic order occurs below $``$ 60 K moritomo-1997 ; zaliznyak-2000 , with a spin-freezing transition at $`T_s`$ 30 K zaliznyak-2000 .

We report here on polarized- and unpolarized-neutron scattering measurements and magnetometry studies of single crystals of La1.5Sr0.5CoO4. The crystals were grown by the floating-zone method in Oxford.

Magnetization data were collected using a superconducting quantum interference device (SQUID) magnetometer with the field parallel to the ab plane. DC field-cooled (FC) and zero-field-cooled (ZFC) data were recorded in a measuring field of 100 Oe.

Figure 1a shows the temperature variation of the in-plane FC and ZFC magnetizations, with the inset showing data up to 350 K. Our results agree with those of Moritomo et al. moritomo-1997 , with a broad maximum in the magnetization centred at about 60 K. However, we also see a definite splitting between the ZFC and FC data at low temperatures indicating a glassy ground-state, and an abrupt change in slope at 30 K. It is possible that this feature corresponds to a spin-reorientation transition as observed in La1.5Sr0.5NiO4 freeman .

To examine the magnetic order further we performed polarized-neutron diffraction on the triple-axis-spectrometer IN20 at the Institut-Laue-Langevin. We studied two rod-shaped crystals of La1.5Sr0.5CoO4 (taken from one original crystal) with a combined mass of 12 g.

Figure 1b shows the temperature variation of the amplitude of a magnetic Bragg peak. This data confirms that La1.5Sr0.5CoO4 is magnetically ordered below $``$ 40 K, but reveals a slow build-up of magnetic correlations with decreasing temperature below $``$ 60 K. It is likely that this gradual ordering explains the broad maximum in the magnetization data (fig. 1a).

Polarization analysis performed on the magnetic order peaks revealed that below 30 K the spins lie in the ab plane at an angle of $`12^{}`$ to the Co-O bonds. We also found tentative evidence that a small reorientation ($`5^{}`$) occurs above 30 K, which would correspond to the kink in fig. 1a. Such a spin reorientation should cause an anomaly in the temperature variation of the magnetic peak shown in fig. 1b but no anomaly is evident. This may be because of the small magnitude of the reorientation, or it may imply a different origin for the magnetization kink.

We studied the spin excitation spectrum both with polarized- and unpolarized-neutron scattering. The latter was performed with the MAPS spectrometer at ISIS on a crystal of mass 35.5 g. The spin correlations were found to be highly two-dimensional (2D), showing no measureable inter-plane correlations for energies above 4 meV. For simplicity, therefore, we describe the scattering with respect to the 2D reciprocal lattice of the square CoO2 planes.

Figure 2 shows a cut through unpolarized neutron data in the ($`h,h`$) direction. The data were measured at $`T`$ = 9.5 K with an incident energy of 50 meV. The cut has been averaged over $`7E9`$ meV and $`0.215h0.295`$ in the $`(h,h)`$ direction. Lorentzian profiles were fitted to the two peaks, which correspond to spin-waves propagating parallel to the $`(h,h)`$ direction away from the $`(0.255,0.255)`$ magnetic zone centre.

By fitting (mainly) constant-energy cuts we determined the two-dimensional spin-wave dispersions in the ($`h,h`$) and ($`h,h`$) directions. These are shown in figure 3. The spin-wave velocity is found to be virtually isotropic in the ab plane. The spin-wave is seen to extend to $``$ 16 meV, in contrast to $``$ 50 meV in La1.5Sr0.5NiO4 freeman2 . This indicates a much weaker exchange interaction in the cobaltates than in the nickelates. A quantitative analysis of the spin-wave dispersion is in progress.

Another feature of the excitation spectrum is a pattern of diffuse scattering between $``$ 25 meV and 33 meV, centred close to positions such as $`(0,0)`$, $`(\pm 0.5,0)`$, $`(0,\pm 0.5)`$ and $`(\pm 0.5,\pm 0.5)`$. Figure 4 depicts a cut along $`(h,h)`$, with data averaged over $`25E29`$. Scattering is seen at 9.5 K but not at room temperature, which indicates that the diffuse features are not due to phonons. Polarised neutron measurements confirmed that the scattering is magnetic.

The origin of this diffuse magnetic inelastic scattering is not clear at present. One possibility is that the Co2+ spins have strong planar anisotropy, and that the diffuse feature corresponds to the excitation of the out-of-plane spin component. Another is that we are observing magnetic correlations among the Co3+ sites. These possibilities are currently being assessed by modelling of the spin excitations.

We would like to thank the Engineering and Physical Sciences Research Council of Great Britain for their help in funding this work.


text:
# Integrability of a disordered Heisenberg spin-1/2 chain

## I introduction

Random matrix theory has long been used to describe the spectra of complex systems, such as nuclei, molecules, and mesoscopic solids weide . Recently it has been used in the study of strongly correlated spin systems spin ; shepe98 . The statistical properties of the quantum energy spectrum are strongly influenced by the underlying classical dynamics. The level spacing distribution of a classical integrable system is Poissonian, $`P_P(s)=\mathrm{exp}(s)`$, while the level statistics of a chaotic system is given by the Wigner-Dyson distribution, $`P_{WD}(s)=(\pi s/2)\mathrm{exp}(\pi s^2/4)`$. The Wigner-Dyson distribution is obtained in random matrix theory and it reproduces the level repulsion of chaotic dynamics. These two distributions characterize, respectively, the localized and the metallic phase in the Anderson model of disordered systems. At the critical point between the two phases an intermediate level spacing statistics occurs inter .

The problem of localization and the statistical properties of the spectra for the case of just one particle has long been understood, but only recently has the problem of many-body systems been addressed spin ; shepe98 . When more than one particle is present in the system, the interaction between them has to be taken into account. The interplay between interaction and disorder is a challenging problem in today’s condensed matter physics and it can lead to new and unexpected effects us .

Here, we consider a one-dimensional Heisenberg spin-1/2 chain with defects and several excitations. A defect corresponds to the site where the energy splitting is different from all the others. It is obtained by applying a different magnetic field in the $`z`$ direction to the chosen site. In the absence of defects this homogeneous system is integrable and is solved with the Bethe ansatz bethe . Its level distribution is therefore Poissonian. As random on-site magnetic fields are turned on and their mean-square amplitude starts increasing, the system undergoes a transition and becomes chaotic. But by further increasing the mean-square amplitude, localization eventually takes place and the distribution becomes Poissonian again. We determine the crossover from integrability to quantum chaos in such disordered spin chain. In addition, we discuss that the cause for nonintegrability is not the randomness of the system, but the mere presence of defects. If only one defect is placed out of the edges of the chain and if the defect excess energy is of the order of the interaction strength, the system is also chaotic. The same sort of transition integrable-chaotic-integrable is obtained as the defect excess energy increases.

We consider only nearest neighbor interaction. The Hamiltonian describing the system is

$$H=\underset{n=1}{\overset{L}{}}(h_n\frac{J}{2}\delta _{n,1}\frac{J}{2}\delta _{n,L})\sigma _n^z+\underset{n}{\overset{L1}{}}\frac{J}{4}\stackrel{}{\sigma }_n.\stackrel{}{\sigma }_{n+1},$$
(1)
where $`\mathrm{}=1`$ and $`\stackrel{}{\sigma }`$ are Pauli matrices. There are $`L`$ sites. Each site $`n`$ is subjected to a magnetic field in the $`z`$ direction, giving the energy splitting $`h_n`$. The chain is ideal whenever all sites have the same energy splitting. A defect corresponds to the site whose energy splitting differs from the others.

For simplicity, we work with an isotropic chain, that is, the coupling constant $`J`$ for the diagonal Ising interaction $`\sigma _n^z\sigma _{n+1}^z`$ is equal to the coupling constant for the $`XY`$-type interaction $`\sigma _n^x\sigma _{n+1}^x+\sigma _n^y\sigma _{n+1}^y`$. This last term is responsible for propagating the excitation through the chain. A single-particle excitation corresponds to a spin pointing up.

In this model, the $`z`$ component of the total spin $`_{n=1}^LS_n^z`$ is conserved, so states with different number of excitations are not coupled. We therefore look at the level spacing distributions for sectors with the same number of excitations. Since we are interested in determining if the system is integrable or chaotic, we focus on the sector with the largest number of states, that is the sector with $`L/2`$ excitations. This is the region where chaos should set in first.

In a very large system the boundary conditions have no effects, but numerical calculations are limited to a finite number of sites. In a periodic (or closed) chain we found too many degenerate states, so we decided to work with a chain with free boundaries (or open chain). Both systems, closed or open, are known to be integrable in the absence of defects. They are solved with the Bethe ansatz method bethe . An open chain with defects only on the edges is also integrable edge . Here we choose an open chain with defects of values $`J/2`$ on the edges. Such values nullify the border effects.

We work with $`L=12`$ sites and 6 excitations, which gives us 924 states. A matrix of such size, 924x924, is sufficient to have good statistics, as illustrated in Fig. 1. In both plots we have the Poisson distribution (dot-dashed line) and the Wigner-Dyson distribution (long-dashed line). The spacings $`s`$ correspond to $`S/M`$, where $`M`$ is the mean level spacing and $`S`$ is the actual spacing. Both histograms are normalized to 1. The histogram at the top of the figure shows that the level spacings of the eigenvalues of a diagonal random matrix of such size is well described by the Poisson distribution. The one at the bottom of the figure shows that the level spacings of the eigenvalues of a random matrix of this size agrees very well with the Wigner-Dyson distribution. The random elements have a Gaussian distribution.

## II Randomly distributed defects

First we analyze the case of random magnetic fields along the $`z`$ direction. The energy splitting of each site is given by $`h_n=h+d_n`$, where $`d_n`$’s are uncorrelated random numbers with a Gaussian distribution: $`d_n=0`$ and $`d_nd_m=d^2\delta _{n,m}`$. When $`d=0`$ the system is integrable and a Poisson distribution is obtained, as the histogram at the top of Fig. 2 shows. As $`d`$ increases the system undergoes a transition and becomes chaotic, the Wigner-Dyson distribution is obtained, as can be seen from the histogram in the middle of Fig. 2. However, as we further increase $`d`$ and it becomes much larger than $`J`$, the system becomes localized. As expected, a Poisson distribution reappears (see the bottom of Fig. 2). Large $`d`$ corresponds to a random diagonal matrix with negligible off diagonal elements.

A more convenient way to analyze the evolution of the level spacing distributions with respect to the ratio $`d/J`$ is by using the parameter $`\eta =_0^{s_0}[P(s)P_{WD}(s)]𝑑s/_0^{s_0}[P_P(s)P_{WD}(s)]𝑑s`$, where $`s_0=0.4729\mathrm{}`$ is the intersection point of $`P_P(s)`$ and $`P_{WD}`$ shepe98 ; Shepelyansky . A regular system has $`\eta =1`$ and a chaotic system has $`\eta =0`$. The stars on Fig. 3 show the dependence of $`\eta `$ on $`d/J`$. There the transition integrable-chaotic-integrable is clear. The system is initially integrable, reaches its maximum chaoticity when $`d10^1J`$ and then becomes integrable again as $`d`$ becomes much larger than $`J`$. At this last step, the energy splitting of each site becomes very different from all the others and the excitations become localized.

Besides considering the site energies as random numbers, another way to introduce disorder is by taking the non-diagonal (hopping) matrix elements at random. We assume that the off diagonal elements also have a Gaussian distribution and the mean square is given by $`J_r`$. Contrary to the case discussed before, random coupling (circles on Fig. 3) leads to chaos even when the energy splittings of all sites are the same ($`d=0`$). Maximum chaoticity is again reached at $`d10^1J_r`$, but there is now just one transition, from chaos to integrability. Such transition takes a little longer to happen than in the previous case, though localization is also attained once the ratio $`d/J_r`$ becomes large.

## III One single defect

In the system treated here, chaos can be associated with randomness only when the coupling is random. In the case of a constant interaction strength, what is really responsible for the nonintegrability of the system is not the randomness of the energy splittings, but the mere existence of a defect. Let us consider again a constant coupling. The top of Fig. 4 shows that a Wigner-Dyson distribution can also be obtained when there is only one defect in the middle of the chain and the defect excess energy is of the order of the interaction strength. For this histogram all sites have the same energy splitting $`h_n=h`$, except site 6, which has $`h_6=h+J`$.

Here too we use the parameter $`\eta `$ to study the evolution of the level spacing distribution with respect to the ratio $`d_6/J`$. The bottom of Fig. 4 indicates a transitional behavior very similar to the case of constant coupling and random diagonal elements. The one-defect system is initially regular, but by increasing $`d_6`$ it becomes chaotic. The highest degree of chaoticity happens when $`d_6J`$. The level repulsion now settles in more slowly than in the previous situation of random magnetic fields. As $`d_6`$ is further increased, an excitation on site 6 will become site-localized. This means that when the ratio $`d_6/J`$ becomes very large there are two kinds of states in the system: states with a localized excitation on site 6 and states with no excitation on the defect. These two types of states are not coupled. This system becomes equivalent to two smaller and uncoupled ideal chains, and integrability is therefore recovered.

In the inset at the bottom of Fig. 4, we again consider one single defect, whose excess energy is equal to $`J`$. It shows how the parameter $`\eta `$ depends on the position of the defect and confirms that the Heisenberg spin-1/2 system is integrable when it has defects placed on its edges.

## IV Conclusion

We have shown that in a Heisenberg spin-1/2 chain, randomly distributed defects lead to the transition integrable-chaotic-integrable, according to the ratio $`d/J`$. When random off diagonal elements are also considered, there is only one transition from non-integrability to integrability. The transition integrable-chaotic-integrable is also observed when the coupling is constant and there is only one defect out of the edges of the chain. The level spacing distribution obtained in this case depends on how large the defect excess energy is in relation to the interaction strength.

Obtaining analytical solutions for disordered spin chains is not an easy task and in many cases is simply impossible. The algebraic version of the Bethe ansatz is often used to construct integrable Hamiltonians Zvy . The analysis of level spacing distributions should therefore be useful to identify which real or constructed Hamiltonians are indeed integrable.

Understanding under what conditions disordered spin systems become integrable is not just relevant for condensed matter physics, but also for quantum computation, since these systems are commonly used to model different proposals of quantum computers (see Ref. us  and references therein).

###### Acknowledgements.

We acknowledge support by the NSF through grant No. ITR-0085922 and would also like to thank C. O. Escobar for helpful discussions.


text:
# 1 Introduction

## 1 Introduction

Very recently the data from the Wilkinson Microwave Anisotropy Probe (WMAP) has allowed analyses of the cosmological parameters to a high degree of accuracy. These analyses also indicate unambiguously the existence of cold dark matter (CDM) and put sharp limits on it. At the same time over the past decade experiments for the direct detection of dark matter have made enormous progress with reliable limits emerging on the CDM component in direct laboratory experiments. Further, experiments are planned which in the future will be able to improve the sensitivities by several orders of magnitude. In this talk we will give a brief review of the recent developments in supersymmetric dark matter (For a sample of other recent reviews see Ref.). We will review the constraints on the analyses of dark matter from $`g_\mu 2`$ and from $`B_{s,d}^0\mu ^+\mu ^{}`$. We will also discuss the effects of nonuniversalities and the effects of the constraints of Yukawa coupling unification. One of the main focus of our analysis will be the study of dark matter on the hyperbolic branch (and focus point region which is a subpiece of the hyperbolic branch) and its implications for the discovery of supersymmetry. As is well known SUGRA models with R parity provide a candidate for supersymmetric dark matter. This is so because in SUGRA unified models one finds that over a large part of the parameter space the lightest supersymmetric particle (LSP) is the lightest neutralino which with R parity conservation becomes a candidate for cold dark matter (CDM). (An interesting alternate possibility discussed recently is that of axionic dark matter). In the simplest version of SUGRA models, mSUGRA, which is based on a flat Kähler potential the soft sector of the theory is parameterized by $`m_0,m_{\frac{1}{2}},A_0,\mathrm{tan}\beta `$, where $`m_0`$ is the universal scalar mass, $`m_{\frac{1}{2}}`$ is the universal gaugino mass, $`A_0`$ is universal trilinear coupling and $`\mathrm{tan}\beta =<H_2>/<H_1>`$ where $`H_2`$ gives mass to the up quark and $`H_1`$ gives mass to the down quark and the lepton. The minimal model can be extended by considering a curved Kähler manifold and also a curved gauge kinetic energy function. Specifically these allow one to include nonuniversalities in the Higgs sector, in the third generation sector and in the gaugino sector consistent with flavor changing neutral currents.

## 2 Constraints on dark matter analyses

There are a number of constraints that must be imposed in the theoretcal analyses of supersymmetric dark matter. These include the constraints from $`g_\mu 2`$, the flavor changing neutral current (FCNC) constraints, constraints from $`B_s^0\mu ^+\mu ^{}`$ limits, and constraints from the currents limits on the relic density for cold dark matter. Here we discuss some of these. We begin with a discussion of the $`g_\mu 2`$ constraint. The analysis of $`g_\mu 2`$ has been under scrutiny for a considerable period of time. In supersymmetry it is predicted that the effects of supersymmetric electroweak corrections to $`g_\mu 2`$ can be of the same size as the standard model electroweak effects. Furthermore, the sign of the supersymmetric electroweak correction is correlated with the sign of the Higgs mixing parameter $`\mu `$. It is also known that the effects of extra dimensions on $`g_\mu 2`$ are substantially smaller than the effects from supersymmetric loop corrections. Experimentally the situation is still somewhat unsettled due to the ambiguity in the errors of the hadronic corrections. A recent estimate of the difference between experiment and the standard model gives a $`3\sigma `$ effect. Further, the sign of the difference in these analyses is positive indicating a positive $`\mu `$. In Fig.1 we give a numerical analysis of neutralino-proton cross section $`\sigma _{\chi p}`$ which enters in the direct detection of dark matter. The results are taken from the analysis of Ref.. The analysis shows that the $`g_\mu 2`$ constraint is very strong and eliminates a large part of the parameters space consistent with the relic density constraints. The analysis also shows that a very substantial part of the parameter space consistent with $`g_\mu 2`$ and relic density constraints will be accessible to future dark matter detectors.

We turn now to the implications of a positive $`\mu `$ for Yukawa coupling unification and the implications of the Yukawa unification constraint on supersymmetric dark matter analyses. It is known that $`b\tau `$ unification prefers $`\mu <0`$. The above arises from the fact that a negative correction to the b-quark mass is desired for $`b\tau `$ unification and a negative correction is most readily manufactured if $`\mu `$ is negative. However, a closer scrutiny shows that the sign of the loop correction to the b quark mass is not rigidly tied to the sign of $`\mu `$. Thus the dominant correction to the b quark mass arises from the gluino exchange diagram and the sign of this correction is determined by the sign of $`\mu M_{\stackrel{~}{g}}`$ (see Ref. and the references quoted therein). At same time the supersymmetric correction to $`g_\mu 2`$ is governed mainly by the chargino exchange and the sign of that contribution is determined by the term $`\mu \stackrel{~}{m}_2`$. Thus it is possible to relax the rigid relationship between the $`\mu `$ sign and the sign of the loop correction to the b quark mass while maintaining the usual connection between the sign of $`\mu `$ and the sign of the supersymmetric loop correction to $`g_\mu 2`$. The solution to this relaxation is provided by nonuniversalities which allow one to switch the sign of $`M_{\stackrel{~}{g}}`$ relative to the sign of $`\stackrel{~}{m}_2`$. This switch in sign can be seen to arise group theoretically when the gaugino mass terms has nontrivial group transformations. Thus, for example, in $`SU(5)`$ the gaugino masses in general tranform like the $`(24\times 24)_s`$ representation of $`SU(5)`$. Now $`(24\times 24)_s=1+24+75+200`$ and for the $`24`$ plet on the right hand side one has that the $`SU(3)\times SU(2)\times U(1)`$ gaugino masses are in the ratio $`M_3:M_2:M_1=2:3:1`$ and there is a relative sign between the $`M_3`$ and the $`M_2`$. Similarly for $`SO(10)`$ the gaugino masses transform like the $`(45\times 45)_{sym}`$ representations of $`SO(10)`$ where $`(45\times 45)_{sym}=1+54+210+770`$. Here the $`54`$ plet representation on the right hand gives $`M_3:M_2:M_1=1:3/2:1`$ and one finds once again that $`M_3`$ and $`M_2`$ have opposite signs. The above possibilities allow for (see also) $`\mu >0,\mathrm{\Delta }a_\mu >0,\mathrm{\Delta }_b^{\stackrel{~}{g}}<0`$ allowing for $`b\tau `$ unification for a positive $`\mu `$. Detailed analyses of these can be found in where the implications of Yukawa unification on supersymmetric dark matter are also discussed.

Next we discuss the implications of the $`B_{s,d}^0l^+l^{}`$ constraint. In the Standard Model $`B(\overline{B}_s^0\mu ^+\mu ^{})=(3.1\pm 1.4)\times 10^9`$. However, this branching ratio may lie beyond the reach of RUNII which may reach a sensitivity of $`10^{(78)}`$. It turns out that in SUSY/SUGRA models the $`B(\overline{B}_s^0\mu ^+\mu ^{})\mathrm{tan}^6\beta `$ for large $`\mathrm{tan}\beta `$. This arises from the so called counterterm diagram. Because of the large $`\mathrm{tan}\beta `$ factor the branching ratio in SUSY models can get much larger than in the Standard Model. Specifically in mSUGRA the branching ratio can be as large as $`O(10^6)`$ and within reach of RUNII. Further, it was found that CP phases can provide an extra enhancement of $`10^2`$ in some cases (see Ibrahim and Nath in Ref.). Thus SUSY enhancement brightens the prospect for the observation of this decay at Fermilab. Further, the observation of $`B_{s,d}^0l^+l^{}`$ will be evidence for SUSY even before the sparticles are seen directly at colliders. In Fig. (2) the implications of the constraints of various limits on $`B_{s,d}^0l^+l^{}`$ branching ratio are given. Thus, for example, a branching ratio $`BR(B_{s,d}^0l^+l^{})=10^8`$ can probe the parameter space in $`m_{\frac{1}{2}}m_0`$ in the range up ot $`600`$ GeV-$`700`$ GeV for $`\mathrm{tan}\beta =50`$, $`A_0=`$ and $`\mu >0`$.

## 3 The Hyperbolic Branch, Supersymmetry and Dark Matter

It was shown quite sometime ago that the radiative breaking of supersymmetry has two branches: an ellipsoidal branch (EB) and a hyperbolic branch (HB). This can be exhibited easily by examining the relation for radiative breaking that determines $`\mu `$, i.e.,

$`C_1m_0^2+C_3m_{1/2}^2+C_2^{}A_0^2+\mathrm{\Delta }\mu _{loop}^2=\mu ^2+{\displaystyle \frac{1}{2}}M_Z^2`$ (1)
where $`\mathrm{\Delta }\mu _{loop}^2`$ arise from loop corrections to the effective potential. When $`\mathrm{tan}\beta `$ is small, the loop corrections are typically small, and furthermore $`C_1,C_2^{},C_3>0`$ and the variation of these parameters with the renormalization group scale are also small. In this case the radiative breaking of the electroweak symmetry is realized on the ellipsoidal branch and there is an upper limit to the soft parameters for a fixed value of $`\mu `$. However, for other regions of the parameter space specifically when $`\mathrm{tan}\beta `$ is large one finds that $`\mathrm{\Delta }\mu _{loop}^2`$ is large and furthermore the scale dependence of some of the co-efficients $`C_i`$ is rather large. Specifically in this case, if one chooses a scale $`Q_0`$ where the loop correction $`\mathrm{\Delta }\mu _{loop}^2`$ is minimized one finds that the co-efficient $`C_1`$ at $`Q_0`$ turns negative and the hyperbolic branch is realized. Thus on the hyperbolic branch $`m_0,m_{\frac{1}{2}}`$ can get very large for fixed $`\mu `$. A very interesting phenomenon of inversion occurs when $`M_i>>|\mu |`$. In this case an examination of the neutralino mass matrix and of the chargino mass matrix shows that $`\chi _1^0,\chi _2^0`$, $`\chi _1^\pm `$ are essentially degenerate with mass $`\mu `$. More specifically one finds that (for other analyses that explore the implications of WMAP data see )

$`M_{\chi _1^0}=\mu {\displaystyle \frac{M_Z^2}{2}}(1\mathrm{sin}2\beta )[{\displaystyle \frac{\mathrm{sin}^2\theta _W}{M_1\mu }}+{\displaystyle \frac{\mathrm{cos}^2\theta _W}{M_2\mu }}]`$
$`M_{\chi _2^0}=\mu +{\displaystyle \frac{M_Z^2}{2}}(1+\mathrm{sin}2\beta )[{\displaystyle \frac{\mathrm{sin}^2\theta _W}{M_1+\mu }}+{\displaystyle \frac{\mathrm{cos}^2\theta _W}{M_2+\mu }}]`$
$`M_{\chi _1^\pm }=\mu +{\displaystyle \frac{M_W^2\mathrm{cos}^2\beta }{\mu }}{\displaystyle \frac{M_W^2}{\mu }}{\displaystyle \frac{(M_2\mathrm{cos}\beta +\mu \mathrm{sin}\beta )^2}{(M_2^2\mu ^2)}}`$ (2)
In the inversion region the squarks and sleptons and the gluino may lie in the several TeV region and thus may not be easily observable at accelerators. Here the lightest particles are $`h^0,\chi _1^0,\chi _2^0,\chi _1^\pm `$ where $`m_{\chi _1^0}m_{\chi _1^\pm }m_{\chi _2^0}\mu `$. We note that the mass relations here are in gross violation of the scaling laws. Further, the quantities relevant for the observation of the lightest supersymmetric particles in this case are the mass differences $`\mathrm{\Delta }M^\pm =m_{\chi _1^\pm }m_{\chi _1^0}`$, and $`\mathrm{\Delta }M^0=m_{\chi _2^0}m_{\chi _1^0}`$. While $`m_{\chi _1^0}`$, $`m_{\chi _1^\pm }`$, $`m_{\chi _2^0}`$ may lie in the several hundred GeV region, $`\mathrm{\Delta }M^\pm `$ and $`\mathrm{\Delta }M^0`$ lies in the 1-10 GeV region. The smallness of the mass differences makes the observation of these particles rather difficult since the decay of the NLSP and also of the light chargino will result in rather soft particles which may not be observable at the LHC. Situations of this type have been discussed in other contexts in the literature. However, quite surprisingly the satisfaction of the relic density constraints can occur easily in the inversion region. This is so because as pointed out already that in the inversion region there is a near degeneracy of $`m_{\chi _1^0}`$, $`m_{\chi _1^\pm }`$, $`m_{\chi _2^0}`$ which implies a lot of coannihilation. Specifically one has coannihilations involving the processes

$`\chi _1^+\chi _1^{},\chi _1^0\chi _2^0u_i\overline{u}_i,d_i\overline{d}_i,W^+W^{}`$
$`\chi _1^0\chi _1^+,\chi _2^0\chi _1^+u_i\overline{d}_i,\overline{e}_i\nu _i,AW^+,ZW^+,W^+h`$ (3)
Because of coannihilation relic density constraints can be easily satisfied even though the squark and slepton masses may lie in the several TeV region. Detailed analyses show that in the HB region $`\sigma _{\chi _1^0p}`$ may still be accessible to dark matter experiment. Thus observation of dark matter may be the only means of observing SUSY effects if the inversion region of HB is realized in nature. We note in passing that the so called focus point region is a part of the hyperbolic branch and corresponds to low values of $`m_{\frac{1}{2}}`$. For further discussion see Ref..

## 4 Analysis of dark matter with WMAP constraints

In this section we discuss the implications of the WMAP constraints. We also discuss the implications for the detection of dark matter in dark matter detectors presently in operation as well as those that are planned for the future. We begin with the result of WMAP which gives for CDM

$$\mathrm{\Omega }_\chi h^2=0.1126_{0.018}^{+0.016}$$
(4)
In the theoretical analyses of the relic density coannihilation plays a central role. However, unlike the usual coannihilation phenomena where the particles that enter in the coannihilation process are the neutralino and the stau here the particles that coannihilate are $`\chi _1^0,\chi _2^0.\chi _1^\pm `$. Specifically the processes are listed in Eq.(3). In the inversion region the neutralino is essentially a Higgsino as opposed to being a Bino which is what happens over the most of the rest of the parameter space of SUGRA models. In Fig.(3) a numerical analysis of the region in the parameter space of $`m_0m_{\frac{1}{2}}`$ allowed under the WMAP relic density constraint with a $`2\sigma `$ error corridor is exhibited. One finds that the region consistent with the WMAP relic density constraint can indeed stretch to rather large values in ($`m_0`$, $`m_{\frac{1}{2}}`$) extending into several TeV in each direction. In Figs.(4) a further investigation of the parameter space consistent with the relic density constraints is carried out. Thus in Fig.(4a) an analysis of the parameter space in the $`m_0`$ and $`m_{\chi _1^0}`$ consistent with the relic density constraint is given. Here one finds that $`m_{\chi _1^0}`$ is sharply limited while $`m_0`$ gets large. The vertical patch at $`m_{\chi _1^0}`$ around a TeV is the inversion region where the neutralino is essentially a Higgsino. A similar phenomenon is visible in Fig.(4b). A plot of the neutralino-proton scalar (spin dependent) cross section $`\sigma _{\chi p}`$ as a function of $`m_{\chi _1^0}`$ is given in Fig.(4c)(Fig.(4d)). A comparison with Fig.(1) shows that a significant part of the inversion region will be accessible to future experiments on the direct detection of dark matter.

## 5 Conclusion

In this paper we have given a brief summary of some of the recent developments in supersymmetric dark matter. We have discussed the constraints of $`g_\mu 2`$ and of $`B_s^0\mu ^+\mu ^{}`$ on dark matter analyses. One of the most stringent constraints arises from the recent observation from WMAP which has measured the relic density for CDM to a high degree of accuracy. We discussed the allowed parameter space in mSUGRA satisfying the WMAP constraints. It was shown that quite surprisingly the allowed parameter space is quite large. Specifically one finds that a very significant region on the hyperbolic branch with ($`m_0`$, $`m_{\frac{1}{2}}`$) extending in several TeV still allows the satisfaction of the relic density constraints consistent with WMAP. The consistency with the WMAP data arises due to coannihilation of $`\chi _1^0,\chi _2^0,\chi _1^\pm `$. Further, one finds that the neutralino-proton cross section fall in range that may be accessible to dark matter detectors in the future. Thus if SUSY is realized deep on the hyperbolic branch, then direct observation of sparticles, aside from the light Higgs, may be difficult. However, degeneracy of $`\chi _1^0`$, $`\chi _2^0`$, $`\chi _1^\pm `$ would lead to significant coannihilation and satisfaction of relic density constraints and the direct detection of supersymmetric dark matter may still be possible. Finally, we note that in heterotic string models $`\mathrm{tan}\beta `$ is a determined quantity under the constraints of radiative breaking of the electroweak symmetry and thus dark matter analyses are more constrained in this framework. This constraint will be explored in further work. A similar situation may occur in models based on soft breaking in intersecting D branes.

Acknowledgments
This research was supported in part by NSF grant PHY-0139967


text:
# 1 Transverse momentum (left) and Feynman 𝑥 (right) dependencies of 𝑃_Λ.

SPIN EFFECTS AND ROLE OF CONSTITUENT QUARK SIZE

S.M. Troshin, N.E. Tyurin

Institute for High Energy Physics,
Protvino, Moscow Region, 142281, Russia

## Abstract

The nonperturbative mechanism for spin effects in inclusive production based on the constituent quark model is considered. The main role belongs to the orbital angular momentum and the polarization of the $`\overline{q}q`$–pairs in the internal structure of the constituent quarks.

It is well known fact that spin observables demonstrate rather complicated behaviour and it somehow is related to a scarcity of the experimental data in the field. However, it seems very important to find some general trends in the data even under such strained circumstances. One among the more or less stable patterns is the momentum transfer dependence of various spin observables in inclusive production. To be specific we will concentrate on a particular problem of $`\mathrm{\Lambda }`$–polarization. Experimentally, the situation is stable and clear. $`\mathrm{\Lambda }`$–polarization is negative and energy independent. It grows linearly with $`x_F`$ for $`p_{}>0.8`$ GeV/c, and for large values of the momentum transfer ($`0.8<p_{}<3.5`$ GeV/c) it is $`p_{}`$–independent . It is remarkable that both parameters $`A_N`$ and $`D_{NN}`$ show $`p_{}`$–dependence similar to polarization .

In perturbative QCD. a straightforward collinear factorization leads to very small values of $`P_\mathrm{\Lambda }`$ . pQCD modifications and in particular account for higher twists result in the dependence $`P_\mathrm{\Lambda }1/p_{}`$ . This behavior still does not correspond to the experimental trends. Account for $`k_{}`$–effects when the source of polarization is moved into the polarizing fragmentation functions also leads to falling $`P_\mathrm{\Lambda }k_{}/p_{}`$ at large $`p_{}`$ values . Potentially $`\mathrm{\Lambda }`$–polarization could become even a more serious problem that the nucleon spin problem. And in any case the both problems are interrelated.

The essential point here is that the vacuum at short distances is taken to be a perturbative one. However, polarization dynamics could have its roots hidden in the genuine nonperturbative sector of QCD. The models exploiting confinement and the chiral symmetry breaking have been proposed. Our model considerations  are based on the effective Lagrangian approach which in addition to the four–fermion interactions of the original NJL model includes the six–fermion $`U(1)_A`$–breaking term.

Chiral symmetry breaking generates quark masses:

$$m_U=m_u2g_40|\overline{u}u|02g_60|\overline{d}d|00|\overline{s}s|0.$$

In this approach massive quarks appear as quasiparticles, i.e. current quarks surrounded by a cloud of quark–antiquark pairs of different flavors. For example, for the $`U`$–quark the ratio

$$U|\overline{s}s|U/U|\overline{u}u+\overline{d}d+\overline{s}s|U$$

is estimated as $`0.10.5`$. The scale of spontaneous chiral symmetry breaking is

$$\mathrm{\Lambda }4\pi f_\pi 1\text{GeV/c}$$

and provides the momentum cutoff which determines a transition to the partonic picture. We consider nonperturbative hadron as consisting of constituent quarks located in the central part of the hadron and embedded into a quark condensate.

Respectively, spin of the constituent quark is given by the following ”spin balance equation”:

$$J_U=1/2=S_{u_v}+S_{\{\overline{q}q\}}+L_{\{\overline{q}q\}}=1/2+S_{\{\overline{q}q\}}+L_{\{\overline{q}q\}},$$
(1)
where $`L_{\{\overline{q}q\}}`$ is the orbital angular momentum of quark–antiquark pairs in the structure of the constituent quark. Its value can be estimated  with the use of the polarized DIS data:

$$L_{\{\overline{q}q\}}0.4$$
(2)
This means that the cloud quarks rotate coherently and significant part of the constituent quark spin is to be associated with the orbital angular momentum. In the model just this orbital motion of quark matter is the origin of asymmetries in inclusive processes. It is to be noted that the only effective degrees of freedom here are quasiparticles. The gluon degrees of freedom are overintegrated, and the six-fermion operator in the NJL Lagrangian simulates the effect of the gluon operator $`(\alpha _s/2\pi )G_{\mu \nu }^a\stackrel{~}{G}_a^{\mu \nu }`$ in QCD. It is also important to note the exact compensation between the spins of $`\overline{q}q`$-pairs and their orbital momenta:

$$L_{\{\overline{q}q\}}=S_{\{\overline{q}q\}},$$
(3)
which follows from Eq. (1).

Assumed picture of hadrons implies that overlapping and interaction of peripheral clouds and condensate excitation occur at the first stage of the collision. As a result massive virtual quarks appear in the overlapping region and some mean field is generated. Inclusive production of hyperon results from the two mechanisms: recombination of constituent quark with virtual massive strange quark (soft interactions) or from the constituent quark scattering in the mean field, its excitation and appearance of a strange quark as a result of decay of the parent constituent quark. The second mechanism is determined by interactions at the distances smaller than the constituent quark radius ($`r<R_Q1/\mathrm{\Lambda }_\chi `$) and is associated with hard interactions. Thus, we adopt a two-component picture of hadron (hyperon) production which incorporates interactions at long and short distances and it is the short distance dynamics which leads to production of polarized $`\mathrm{\Lambda }`$’s.

Polarization of a strange quark results from the multiple scattering of parent constituent quark $`Q`$ in the mean field where it gets polarized

$$𝒫_QI\frac{m_Qg^2}{\sqrt{s}}$$
(4)
and the polarization is nearly constant in the model since $`m_Qm_h/3`$ and $`I\sqrt{s}`$. The second crucial point is correlation between $`s`$–quark polarization and polarization of the parent quark $`Q`$. Indeed, the total orbital momentum of $`\overline{q}q`$–pairs in the constituent quark which has polarization $`𝒫_Q(x)`$ is

$$L_{\{\overline{q}q\}}^{𝒫_Q(x)}=𝒫_Q(x)L_{\{\overline{q}q\}},$$
(5)
where the value $`L_{\{\overline{q}q\}}`$ on the right hand side enters Eq. (1) written for the constituent quark with polarization $`+1`$. On the basis of Eq. (3) we suppose that there is a compensation between spin and orbital momentum of strange quarks inside the constituent quark

$$L_{s/Q}=J_{s/Q}=\alpha 𝒫_Q(x)L_{\{\overline{q}q\}},$$
(6)
where the parameter $`\alpha `$ determines the fraction of orbital momentum due to the strange quarks. Eq. (6) is quite similar to the conclusion made in the framework of the Lund model but has a different dynamical origin rooted in the mechanism of the spontaneous chiral symmetry breaking.

Final expression for the polarization is

$$P(s,x,p_{})=\mathrm{sin}[\alpha 𝒫_Q(x)L_{\{\overline{q}q\}}]\frac{R(s,x,p_{})}{1+R(s,x,p_{})}.$$
(7)
The function $`R`$ is the cross–section ratio of hard and soft processes. At $`p_{}>\mathrm{\Lambda }_\chi `$ the function $`R(s,x,p_{})1`$ and the polarization saturates

$$P(s,x,p_{})=\mathrm{sin}[\alpha 𝒫_Q(x)L_{\{\overline{q}q\}}].$$
(8)
Characteristic $`p_{}`$–dependence of $`\mathrm{\Lambda }`$–polarization follows from Eqs. (7) and (8): polarization is vanishing for $`p_{}<\mathrm{\Lambda }_\chi `$, it gets an increase in the region of $`p_{}\mathrm{\Lambda }_\chi `$ and polarization saturates and becomes $`p_{}`$–independent (flat) for $`p_{}>\mathrm{\Lambda }_\chi `$. The respective scale is the scale of the spontaneous chiral symmetry breaking $`\mathrm{\Lambda }_\chi 1`$ GeV/c. Such a behavior of polarization follows from the fact that constituent quarks themselves have slow (if at all) orbital motion and are in the $`S`$–state, but interactions with $`p_{}>\mathrm{\Lambda }_\chi `$ resolve the internal structure of constituent quark and feel the presence of internal orbital momenta inside this constituent quark.

To describe the data quantitatively we have to introduce some parameterization. The function

$$R(s,x,p_{})=C(x)\mathrm{exp}(p_{}/m)/(p_{}^2+\mathrm{\Lambda }_\chi ^2)^2$$
(9)
implies typical behavior of cross–sections for hard and soft processes. The form

$$𝒫_Q(x)=𝒫_Q^{max}x$$

is suggested by the ALEPH data . The value of Eq. (2) alongside with $`m=0.2`$ GeV and $`\alpha =0.8`$ provides a rather good fit to the experimental data (Fig. 1 and 2).

In the model the spin transfer parameter $`D_{NN}`$ is positive since $`P_\mathrm{\Lambda }`$ has the same sign as $`𝒫_Q`$. The model also predicts similarity of $`p_{}`$ dependencies for the different spin observables. The respective features were clearly seen in E-704 experiment .

Experimental prospects

In this approach the short distance interaction with $`p_{}\mathrm{\Lambda }_\chi `$ observes a coherent rotation of correlated $`\overline{q}q`$–pairs inside the constituent quark and not a gas of the free partons. The nonzero internal orbital momenta in the constituent quark means that there are significant multiparton correlations. The important point is what the origin of this orbital angular momentum is. The analogy with an anisotropic extension of the theory of superconductivity seems match well with the adopted picture for a constituent quark. An axis of anisotropy can be associated with the polarization vector of the valence quark located at the origin of the constituent quark.

It seems interesting to perform $`\mathrm{\Lambda }`$–polarization measurements at RHIC. When two polarized nucleons are available one could measure three–spin correlation parameters $`(n,n,n,0)`$ and $`(l,l,l,0)`$ in the processes

$$p_,+p_,=\mathrm{\Lambda }_,+X.$$
(10)
It would provide important data to study mechanisms of hyperon polarization.

Experimentally observed persistence and constancy of $`\mathrm{\Lambda }`$–hyperon polarization means that chiral symmetry is not restored in the region of energy and values of $`p_{}`$ where experimental measurements were performed. Otherwise we would not have any constituent quarks and should expect a vanishing polarization of $`\mathrm{\Lambda }`$. It is interesting to perform $`\mathrm{\Lambda }`$–polarization measurements at RHIC and the LHC. It would allow to make a direct check of perturbative QCD and allow to make a cross-check of the QCD background estimations based on perturbative calculations for the LHC. On the basis of the above model we expect significant $`P_\mathrm{\Lambda }`$ at RHIC energies. On the base of the model in one expects zero polarization in the region where QGP has formed, since chiral symmetry is restored and there is no room for quasiparticles such as constituent quarks. The absence or strong diminishing of transverse hyperon polarization can be used therefore as a signal of QGP formation in heavy-ion collisions. This prediction should also be valid for the models based on confinement, e.g. the Lund and Thomas precession model. In particular, the polarization of $`\mathrm{\Lambda }`$ in heavy–ion collisions in the model based on the Thomas precession was described in where nuclear effects were discussed as well. However, we do not expect a strong diminishing of the $`\mathrm{\Lambda }`$–polarization due to the nuclear effects: the available data show a weak $`A`$–dependence and are not sensitive to the type of the target. Thus, we could use a vanishing polarization of $`\mathrm{\Lambda }`$–hyperons in heavy ion collisions as a sole result of QGP formation provided the corresponding observable is non-zero in proton–proton collisions. The prediction based on this observation would be a decreasing behavior of polarization of $`\mathrm{\Lambda }`$ with the impact parameter in heavy-ion collisions in the region of energies and densities where QGP was produced: $`P_\mathrm{\Lambda }(b)0\text{at}b0`$, since the overlap is maximal at $`b=0`$. The value of the impact parameter can be controlled by the centrality in heavy–ion collisions. The experimental program should therefore include measurements of $`\mathrm{\Lambda }`$–polarization in $`pp`$–interactions first, and then if a significant polarization would be measured, the corresponding measurements could be a useful tool for the QGP detection. Such measurements seem to be experimentally feasible at RHIC and LHC provided it is supplemented with forward detectors.

One of the authors (S.T.) would like to thank Organizing Committee of the 10-th Workshop on High Energy Spin Physics (Nato Advanced Research Workshop) Spin-03, Dubna, September 16-20, 2003 for the valuable support and warm hospitality in Dubna.


text:
# Finite Precision Measurement Nullifies Euclid’s Postulates

(Department of Physics, Technion—Israel Institute of Technology, 32000 Haifa, Israel)

## Abstract

Following Meyer’s argument \[Phys. Rev. Lett. 83, 3751 (1999)\] the set of all directions in space is replaced by the dense subset of rational directions. The result conflicts with Euclidean geometry.

Meyer’s claim  that “finite precision measurement nullifies the Kochen-Specker theorem” (that is, makes it irrelevant to physics) and some of its generalizations  have caused considerable controversy that lasts until today . Meyer’s proposal was to replace the set of all directions in space by the dense subset of rational directions, arguing that a finite precision measurement cannot decide whether or not a number is rational.

Let us apply the same argument to ordinary geometry and consider only points with rational coordinates. Then the line $`x=y`$ and the unit circle $`x^2+y^2=1`$ are both dense but they do not intersect, in contradiction to Euclid’s postulates .

This work was supported by the Gerard Swope Fund.


text:
# The near-IR luminosity-metallicity relationship for dwarf irregular galaxies

## 1. The Project

Although early studies found a very well defined luminosity-metallicity (L-Z) relationship for dwarf irregular (dIrr) galaxies (Skillman, Kennicutt, & Hodge 1989; Richer & McCall 1995), more recent investigations seem to suggest a mild relationship with much scatter, or no relationship at all (Hidalgo-Gamez & Olofsson 1998; Hunter & Hoffman 1999; Skillman, Côté, & Miller 2003, hereafter SCM03). The main difficulties with the present status are that (a) all the studies use blue luminosities as tracers of the mass, (b) galaxies in different environments are used, (c) inhomogeneous data are used, and (d) the scatter in distance modulus is unknown. In order to overcome these limitations, we started a medium-term project whose aims are (i) collecting homogeneous samples of oxygen abundances, (ii) focusing on well-defined environments with reduced distance range, and (iii) collecting near-IR luminosities for all the objects. As a first step, we are collecting data for the three nearest groups of galaxies: the M 81, Sculptor, and Cen A groups at a distance of about $`3.5`$, $`2.5`$, and $`3.5`$ Mpc, respectively. The interaction rate in the M 81 and Cen A groups is higher than that in Sculptor, so our targets offer the chance to test the effect of the environment on the L-Z relationship.

## 2. First Results

Spectroscopy has been obtained at the Lick observatory with KAST@3.5m for the M 81 group, and at ESO with EFOSC2@3.6m for the Sculptor group. Near-IR imaging has been obtained at La Palma observatory with INGRID@WHT, and at ESO with SOFI@NTT for the M 81 and Sculptor groups, respectively. Observations of the Sculptor group are presented here. The reductions of the M 81 group are in progress, and observations of the Cen A group have been planned.

Oxygen abundances are plotted vs. $`H`$ magnitudes in Fig. 1. Note that the observed Hii region of ESO 245-005 lies in the outskirts of the galaxy, so the abundance of the central star-forming regions might be closer to the trend defined by the other galaxies. The trend is rather well defined, although with some scatter. On the other hand, the existence of the relationship depends very much on the single point defined by AM 106-382, whose abundance has been found via the empirical method (Edmunds & Pagel 1984). More data are clearly needed, in particular for low-luminosity galaxies. The abundances of the objects from SCM03 agree with our determinations, reinforcing the correlation. The figure shows the advantage of working in the near-IR: the galaxies cover almost $`6`$ magnitudes in $`H`$, while the SCM03 sample spans only $`3`$ magnitudes in $`B`$. The impact of the unknown individual distance moduli is then reduced in the near-IR.

## References

Edmunds, M.G., & Pagel, B.E.J. 1984, MNRAS, 211, 507

Hidalgo-Gamez, A. M., & Olofsson, K. 1998, A&A, 334, 45

Hunter, D. A., & Hoffman, L. 1999, AJ, 117, 2789

Richer, M.G., McCall, M.L. 1995, ApJ, 445, 642

Skillman, E.D., Kennicutt, R.C., & Hodge, P.W. 1989, ApJ, 347, 875

Skillman, E.D., Côté, S., Miller, B.W. 2003, AJ, 125, 610


text:
# The buckling transition of 2D elastic honeycombs: Numerical simulation and Landau theory

## I Introduction

A two dimensional honeycomb structure formed by solid walls is the prototype of a cellular solid . These are materials widely used in applications due to their remarkable mechanical properties, for instance its capacity for energy absorption under impact, and its low weight. Energy absorption is related to plastic deformation under stress. But still the ideally elastic and perfectly uniform two dimensional honeycomb presents some not completely solved puzzles. Under compressive stress, it has a buckling transition in which some (or all) of their walls bend. This transition is reminiscent of the well known buckling transition of an elastic bar under compressive stress at its extremes. There has been some controversy on what the buckling mode of a regular honeycomb should be. On one hand, in their original work, Gibson and Ashby presented the results of an experiment using an elastomeric honeycomb, under what they called biaxial loading, in which they observed a non-trivial buckling pattern consistent with a symmetry breaking in which four original cells form the new repetitive motif of the material. In a posterior paper, Hutzler and Weaire performed numerical simulations and did not observe this pattern, but instead a buckling mode equivalent to that obtained under uniaxial stress. They argue that the pattern observed in Ref.  was a consequence of finite size effect, and the use of flat confining walls. Numerical results taking into account these effects did show the pattern observed by Gibson and Ashby. Very recently, Okomura, Ohno and Noguchi have studied the problem using a combination of a homogenization technique and finite elements numerical simulations. Their results do not agree with those of Hutzler and Weaire. Instead, they found buckling patterns that can be interpreted as a superposition of three individual buckling modes. They also found that whether one, two, or three of these modes are active depends on the degree of anisotropy of the externally applied strain.

In view of the aforementioned contradiction between  and , and considering that the techniques employed in both cases are quite different, an independent investigation to determine which of the two results is correct seems appropriate. In the first (numerical) part of this paper I will show that appropriately done numerical simulations using the technique used in  do not support the results claimed there, but instead those reported in . In the second (more theoretical) part I will show how the results obtained in the simulations are fully compatible with the predictions of a Laundau theory of second order phase transition applied to the buckling problem. This theory allows to obtain at once the buckled configuration of the system under a generic form of the macroscopically homogeneous applied deformation.

## II Numerical Simulation

I have simulated a two dimensional honeycomb through the technique used in Refs. , , namely by considering the honeycomb walls as one dimensional rods, and including stretching and bending energy as

$`E_{stretch}`$ $`=`$ $`{\displaystyle \frac{1}{2}}k_s{\displaystyle \left(\frac{dl}{dl_0}1\right)^2𝑑l_0}`$ (1)
$`E_{bend}`$ $`=`$ $`{\displaystyle \frac{1}{2}}k_b{\displaystyle c^2𝑑l_0}`$ (2)
where $`c`$ is the local curvature. To discretize these expressions I have used 7 intermediate points between any two neighbor vertices, but particular cases where checked using 18 points, to guarantee the absence of noticeable effects due to discretization. The only essential parameter of the model is $`k_b/L^2k_s`$, where $`L`$ is the length of the individual rods. This ratio is physically related to the fraction $`\mathrm{\Lambda }`$ of two dimensional space that is occupied by the rods , namely $`\mathrm{\Lambda }=4\sqrt{k_b/L^2k_s}`$. The simulations presented below were done at $`k_b/L^2k_s=4.5\times 10^4`$ (then $`\mathrm{\Lambda }=0.085`$), but additional checks indicate that the results are not qualitatively dependent on the precise value of the parameter. All previously obtained buckling patterns for perfect honeycombs can be accommodated within a $`2\times 4`$ unit cell. Then the elemental cell I simulated is precisely the $`2\times 4`$ cell shown in Fig. 1(a), with periodic boundary conditions. The simulation method consists in calculating the forces acting on all points of the discretized system, and updating their positions using a viscous dynamics. The control variable was the macroscopic strain, that can be changed varying the size and shape of the simulation box. Stresses in the system can be evaluated both by numerical differentiation of the total energy with respect to strain, and by direct summation in terms of the forces between particles. The equivalence of the two results allows to check for consistency and convergence of the simulation.

Before indicating the results obtained, it is clarifying to discuss qualitatively the behavior observed (see ). The buckling structures that appear are related to reaccommodation of the vertices of the honeycomb structure, in such a way that lines of vertices forming zig-zag chains shift relatively to neighbor chains, as qualitatively indicated in Fig. 1(b). There are three of these modes, that will be refereed to as the elementary modes of buckling. The patterns they generate will be called the uniaxial patterns. They are characterized by the unitary vectors shown in Fig. 1(e). Whether one, two, or the three elementary modes acquire non-zero amplitude at the buckling transition depends on the macroscopic strain applied. The qualitative pictures of Fig. 1(c) and (d) show the effect of combining two or three elementary modes. We see that the patterns of Fig. 1(b) and (c) are identical (except for the wall bending, that in this qualitative picture is not taken into account) to the patterns in Ref. . The pattern in Fig. 1(c) will be referred to as the Gibson-Ashby pattern. The pattern in Fig. 1(d), that combines the three elementary modes, preserves the hexagonal symmetry of the structure (this will be referred to as the symmetric pattern). A buckling mode of this symmetry has been observed and simulated in , but only in the case of plastic buckling. To my knowledge, the only prediction of this mode for a perfectly elastic honeycomb is contained in . It is remarkable that this pattern has lost the mirror symmetry plane of the original structure: it is a chiral pattern. The chirality can be defined as the sign of the product of the amplitudes of the three elementary modes from which the pattern is constructed.

In the numerical simulations, starting from the unstrained structure and isotropically compressing it, I have observed that the symmetric pattern is the stable configuration after buckling. However, it has to be mentioned that the transition from the unbuckled to the symmetrically buckled configuration, although continuous (then implying no surmounting of any energy barrier) may take a long computational time, and that other patterns (which correspond to saddles of the energy) can be transitory observed. For instance, I show in Fig. 2 the mechanical energy of the system as a function of the compressive strain $`s`$ in a single run increasing $`s`$ at a constant rate. At the buckling point (at $`s0.0037`$, note that I define $`s`$ as $`s=(s_x+s_y)/2`$), and in the particular run shown, the system seems to buckle to the Gibson-Ashby pattern shown in Fig 3(b). This is not a minimum but a saddle of the energy, however it lasts long enough, and the energy of this branch can be numerically followed (as seen in Fig. 2 with open triangles) up to some strain at which it is observed to transform to the real energy minimum which corresponds to the symmetric pattern. In fact, when preparing the system in the symmetric pattern for large strain, and upon reduction of the strain, we obtain the results indicated by full circles in Fig. 2, which correspond to the real energy minimum of the system. In addition, if the system is prepared in the uniaxial or in the Gibson-Ashby pattern at large strain, and strain is reduced (sufficiently rapidly for the configuration not to destabilize, and sufficiently slowly so as energy can be calculated accurately) we can follow the energy of these two configurations, as indicated by the stars, and the full triangles in Fig. 2.

The uniaxial and Gibson-Ashby patterns may become the stable buckling modes under appropriate non-isotropic loading. For instance, by compressing along the $`x`$ ($`y`$) direction, and keeping the perpendicular direction unstrained, I have observed the Gibson-Ashby (uniaxial) pattern to be stable after buckling. The stability of the Gibson Ashby pattern under particular uniaxial loading agrees again with the results in , but not with those claimed in . We will see now how all these results can be fully systematized within a theoretical framework.

## III Landau theory of the buckling transition

The results presented in the previous section are enough numerical input to construct the Landau theory of this peculiar second order transition. In fact, we see in Fig. 2 that the energies of all buckled configurations (whether the true minimum or the saddles) meet together in value and derivative at the buckling point, coinciding also at that point in value and derivative with the branch corresponding to the unbuckled system. This means that at the buckling point, in a generic parameter space, the state point of the system passes from a configuration with a single minimum (for the unbuckled state) to one with different minima and saddles in a continuous manner.

I will present a Landau description in which the order parameters are the three (small) amplitudes $`\varphi _1`$, $`\varphi _2`$, $`\varphi _3`$ of the elementary modes. For convenience, these three modes will be associated with the unitary vectors $`𝐯^1`$, $`𝐯^2`$, and $`𝐯^3`$ shown in Fig. 1(e). The free energy of the system (in the present case it actually corresponds simply to the elastic energy) must a scalar, and then it can only contain combinations of the amplitudes and the (eventually anisotropic) external strains that are invariant with respect to the symmetry operations of the lattice. Considering for the moment only the isotropically compressed case, we should look for invariant combinations of the amplitudes. Up to fourth order those available are:

$`\varphi _{1}^{}{}_{}{}^{2}+\varphi _{2}^{}{}_{}{}^{2}+\varphi _{3}^{}{}_{}{}^{2}`$ (3)
$`\left(\varphi _{1}^{}{}_{}{}^{2}+\varphi _{2}^{}{}_{}{}^{2}+\varphi _{3}^{}{}_{}{}^{2}\right)^2`$ (4)
$`\varphi _1^2\varphi _2^2+\varphi _2^2\varphi _3^2+\varphi _3^2\varphi _1^2`$ (5)
Then in this case the most general form of the free energy describing a second order transition is

$`F=\alpha (s_cs)\left(\varphi _1^2+\varphi _2^2+\varphi _3^2\right)+\beta \left(\varphi _1^2+\varphi _2^2+\varphi _3^2\right)^2+`$ (6)
$`\gamma \left(\varphi _1^2\varphi _2^2+\varphi _2^2\varphi _3^2+\varphi _3^2\varphi _1^2\right)`$ (7)
with numerical constants $`\alpha `$, $`\beta `$, and $`\gamma `$. This free energy has a single minimum at $`\varphi _p=0`$ ($`p=1,2,3`$) for $`s<s_c`$, representing the unbuckled state. For $`s>s_c`$, this expression has saddles (or minima) at the following values of the amplitudes:

$`\varphi _p^2`$ $`=`$ $`A,\varphi _q=0,\varphi _r=0`$ (8)
$`\varphi _p^2`$ $`=`$ $`B,\varphi _q=B,\varphi _r=0`$ (9)
$`\varphi _p^2`$ $`=`$ $`C,\varphi _q=C,\varphi _r=C`$ (10)
where $`p`$, $`q`$, $`r`$ are arbitrary permutations of 1, 2, 3, and $`A`$, $`B`$, $`C`$, are given by:

$`A`$ $`=`$ $`{\displaystyle \frac{(ss_c)\alpha }{2\beta }}`$ (11)
$`B`$ $`=`$ $`{\displaystyle \frac{(ss_c)\alpha }{4\beta \gamma }}`$ (12)
$`C`$ $`=`$ $`{\displaystyle \frac{(ss_c)\alpha }{6\beta 2\gamma }}`$ (13)
They correspond respectively to the uniaxial, Gibson-Ashby, and symmetric patterns. The corresponding values of the free energy are,

$`f_{uni}f_{unb}`$ $`=`$ $`{\displaystyle \frac{(s_cs)^2\alpha ^2}{4\beta }}`$ (14)
$`f_{GA}f_{unb}`$ $`=`$ $`{\displaystyle \frac{(s_cs)^2\alpha ^2}{4\beta \gamma }}`$ (15)
$`f_{sym}f_{unb}`$ $`=`$ $`{\displaystyle \frac{(s_cs)^2\alpha ^2}{4\beta 4\gamma /3}}`$ (16)
(I have subtracted the free energy of the unbuckled system $`f_{unb}`$, that is taken as zero in the Landau theory, but that should be included when comparing with the results of the numerical simulations). We see that the symmetric pattern is the minimum energy one for $`\gamma >0`$ (whereas the uniaxial pattern provides the absolute minimum if $`\gamma <0`$). Then, in order to agree with the simulation results, we will assume $`\gamma >0`$. The physical justification for the positive sign of $`\gamma `$ is not provided by the Landau theory, and should come from an explicit evaluation of the mechanical energy of the honeycomb. From the previous values of the free energy a parameter free relation can be obtained and compared with the numerical results. First note that he ratio $`\beta /\gamma `$ can be obtained for instance as

$$\frac{\gamma }{\beta }=3\left(\frac{f_{sym}f_{uni}}{f_{sym}f_{unb}}\right)$$
(17)
This in particular should be independent of the actual value of the compression $`s`$ (as long as $`s>s_c`$). From the numerical results it can be obtained that this ratio is in fact independent on $`s`$, and the numerical value is approximately $`\gamma /\beta 0.01`$ for the parameters of the simulation. Then note that to a very good approximation we have

$$\frac{f_{GA}f_{sym}}{f_{uni}f_{sym}}=\frac{1}{3}$$
(18)
This is a parameter free relation that has to be satisfied in our numerical simulations. From the data in Fig. 2, it can be in fact verified that this is very accurately satisfied. This is a strong evidence that the present Landau theory describes the physics of the present buckling transition.

In order to make the theory more complete, I want to consider now the possibility of non-isotropic external loading on the system. This means that instead of a single parameter $`s`$, we have now a generic (symmetric) strain tensor $`s_{ij}`$ ($`i,j=1,2`$) applied onto the system (the previously introduced isotropic compression $`s`$ is related to the trace of this tensor). This has to be introduced into the free energy in a symmetrically invariant form. To lowest order I will include it only in the second order term, which is the one that triggers the transition. Using the unitary vectors $`𝐯^1`$, $`𝐯^2`$, $`\mathrm{𝐯𝐮}^3`$ of the elementary modes, two different terms quadratic in the amplitudes can be written:

$`F^{(1)}`$ $``$ $`{\displaystyle \underset{i,j=1,2}{}}{\displaystyle \underset{p,q=1,2,3}{}}s_{ii}v_i^pv_i^qa_{pq}\varphi _p\varphi _q`$ (19)
$`F^{(2)}`$ $``$ $`{\displaystyle \underset{i,j=1,2}{}}{\displaystyle \underset{p,q=1,2,3}{}}s_{ij}v_i^pv_j^qb_{pq}\varphi _p\varphi _q`$ (20)
where $`a_{pq}`$ and $`b_{pq}`$ are arbitrary numeric matrices. However, these expressions have to be invariant under permutation of the elementary vectors (since this is a symmetry operation obtained by the mirror symmetry along the line containing the third vector) and sign change of any of the amplitudes (which is obtained by a particular spatial translation allowed by symmetry). Then it is obtained that both $`a_{pq}`$ and $`b_{pq}`$ matrices should be proportional to the identity, i.e,

$`F^{(1)}`$ $``$ $`{\displaystyle \underset{i,j=1,2}{}}{\displaystyle \underset{p=1,2,3}{}}s_{ii}v_i^pv_i^p\varphi _p^2`$ (21)
$`F^{(2)}`$ $``$ $`{\displaystyle \underset{i,j=1,2}{}}{\displaystyle \underset{p=1,2,3}{}}s_{ij}v_i^pv_j^p\varphi _p^2`$ (22)
Then $`F^{(1)}`$ becomes proportional to $`s_{ii}(\varphi _1^2+\varphi _2^2+\varphi _3^2)`$ and is the term considered in the isotropically compressed case.

To analyze the second contribution it can be more convenient to use the following definition of the three independent components of the strain tensor:

$`s`$ $`=`$ $`(s_{11}+s_{22})/2`$ (23)
$`s_2`$ $`=`$ $`(s_{11}s_{22})/2`$ (24)
$`s_3`$ $`=`$ $`s_{12}=s_{21}`$ (25)
which represent the applied deformation in a more physical way: $`s`$ represents an isotropic compression (we already used this), whereas $`s_2`$ and $`s_3`$ are the two independent shear modes, which are related by a $`\pi /4`$ rotation. In terms of these variables, and using explicitly the components of the unitary vectors we finally arrive to the following form of the free energy:

$`F`$ $`=`$ $`\alpha \left[(s_cs)\left(\varphi _1^2+\varphi _2^2+\varphi _3^2\right)\right]+`$ (26)
$`+`$ $`\delta \left[s_2\left(\varphi _1^2\varphi _2^2/2\varphi _3^2/2\right)+s_3{\displaystyle \frac{\sqrt{3}}{2}}\left(\varphi _3^2\varphi _2^2\right)\right]+`$ (27)
$`+`$ $`\beta \left(\varphi _1^2+\varphi _2^2+\varphi _3^2\right)^2\gamma (\varphi _1^2\varphi _2^2+\varphi _2^2\varphi _3^2+\varphi _3^2\varphi _1^2)`$ (28)
This is the final expression for the free energy close to the buckling transition. Minimizing it we can obtain the buckled state under any particular combination of the three independent strains $`s`$, $`s_2`$, and $`s_3`$. Note that as this expression has to describe also the unstrained system, and then $`\alpha s_c`$ and $`\delta `$ will be respectively proportional to bulk and shear modulus of the original structure, then typically $`\delta <<\alpha s_c`$.

I want to describe now the buckling mode map of the system, namely, what the amplitudes of the three elementary modes are for any choice of the strain tensor. First notice the following scaling of the free energy: If we consider the values of the three order parameters at the minimum of (28), namely $`\varphi _{min}^p`$, to be a function of $`s_cs`$, $`s_2`$, and $`s_3`$, then the following relation is satisfied:

$$\varphi _{min}^p(s_cs,s_2,s_3)=\lambda ^{1/2}\varphi _{min}^p(\lambda (s_cs),\lambda s_2,\lambda s_3)$$
(29)
This implies in particular that the borders between different regions in the parameters space $`s_cs`$, $`s_2`$, and $`s_3`$ are spanned by rays propagating from the origin.

I will analyze a couple of particular cases. First consider the case $`s_3=0`$, i.e., purely compressive strains along $`x`$ and $`y`$ (although non necessarily equal). I show in Fig. 4(a) the map of buckling modes in the $`s`$-$`s_2`$ plane for this case. The borders between different regions can be worked out analytically. All of them are straight lines emanating from the point $`s=s_c`$, $`s_2=0`$, as the previous argument indicates. The transition between unbuckled ($`\varphi _p0`$) and uniaxial pattern ($`\varphi _10`$) can be easily obtained setting $`\varphi _2=\varphi _3=0`$ in (28). The limit line is given by $`s_2=(ss_c)\alpha /\delta `$. The transition line between the unbuckled and the Gibson-Ashby pattern ($`\varphi _2=\varphi _30`$) is obtained along the same lines, as $`s_2=2(ss_c)\alpha /\delta `$. Increasing $`s`$ at $`s_2=0`$, the symmetric pattern appears at $`s=s_c`$, as we already know from the isotropically compressed case. The symmetric pattern looses its strict rotational symmetry for any $`s_20`$. However, it still has the three elementary model active as long as we are within the V-shaped region in Fig. 4(a). The limits of this region are given by

$$s_2^{(unisym)}=\frac{1}{3}\frac{\alpha \gamma }{\beta \delta }\frac{s_cs}{\left(1\frac{\gamma }{3\beta }\right)}$$
(30)
for the transition to the uniaxial pattern, and

$$s_2^{(GAsym)}=\frac{1}{6}\frac{\alpha \gamma }{\beta \delta }\frac{s_cs}{\left(1\frac{\gamma }{3\beta }\right)}$$
(31)
for the transition to the Gibson-Ashby pattern. Note the exact relation

$$s_2^{(GAsym)}/s_2^{(unisym)}=2$$
(32)
valid for any values of the parameters of the free energy.

The results of Fig. 4(a) are fully compatible with the results in  (in particular, relation (32) is very well satisfied). We note that for the present parameters the stability of the single elementary mode and Gibson-Ashby pattern that was obtain numerically under uniaxial compression is recovered.

As an additional example I show the map of buckling modes in the $`s_2`$-$`s_3`$ plane for some $`s>s_c`$ in Fig. 4(b). Note the nice symmetry of this pattern, which has one two or three elementary modes active depending on the particular choice of the applied strains $`s_2`$ and $`s_3`$ (remember that $`s_2`$ and $`s_3`$ are related by a rotation of $`\pi /4`$). Again, all borders between different sectors are straight lines. The analytical expression for the line separating sectors 1 and 1,2 is given by

$$s_3=\sqrt{3}s_2\left(1\frac{\gamma }{3\beta }\right)\frac{\alpha \gamma }{\sqrt{3}\beta \delta }(s_cs)$$
(33)
All other lines can be obtained then from symmetry.

To finish, we note that all transitions in the parameter space are continuous, namely, there are no jumps of the order parameters at any point, and there is no possibility of metastabilities either.

## IV Conclusions

The buckling mode of an elastic two dimensional honeycomb provides an example of non-trivial patterns with symmetry breaking appearing in a very simple mechanical system. Remarkably, for isotropic compression the symmetry breaking produces the appearance of a chiral ground state. This problem is also a realization of a second order transition that can be accurately modeled through a Landau theory constructed on the basis of the symmetry of the problem. The agreement between the Landau theory and the numerical simulation is seen to be very good.