text:
# Pseudorapidity Distribution of Charged Particles in d+Au Collisions at √𝑠__{𝑁𝑁}=200GeV

## Abstract

The measured pseudorapidity distribution of primary charged particles in minimum-bias $`\mathrm{d}+\mathrm{Au}`$ collisions at $`\sqrt{s_{_{NN}}}=200\mathrm{GeV}`$ is presented for the first time. This distribution falls off less rapidly in the gold direction as compared to the deuteron direction. The average value of the charged particle pseudorapidity density at midrapidity is $`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta _{\eta 0.6}=9.4\pm 0.7(\mathrm{syst})`$ and the integrated primary charged particle multiplicity in the measured region is 82 $`\pm `$ 6(syst). Estimates of the total charged particle production, based on extrapolations outside the measured pseudorapidity region, are also presented. The pseudorapidity distribution, normalized to the number of participants in $`\mathrm{d}+\mathrm{Au}`$ collisions, is compared to those of $`\mathrm{Au}+\mathrm{Au}`$ and $`\mathrm{p}+\overline{\mathrm{p}}`$ systems at the same energy. The $`\mathrm{d}+\mathrm{Au}`$ distribution is also compared to the predictions of the parton saturation model, as well as microscopic models.

The pseudorapidity distribution of charged particles in $`\mathrm{d}+\mathrm{Au}`$ collisions is important for understanding the evolution of the system created in more complicated $`\mathrm{Au}+\mathrm{Au}`$ collisions, and may provide a constraint on the initial state parton density. The energy and centrality dependence of pseudorapidity distributions in $`\mathrm{Au}+\mathrm{Au}`$ collisions measured at RHIC BacPRL2001 ; Bac2 ; Brah  are consistent with the approach based on the ideas of parton saturation KhPLB523 ; GriPR100  and semi-classical QCD MacPRD . Measurements of $`\mathrm{d}+\mathrm{Au}`$ collisions may be the ideal way to search for the onset of gluon saturation, since the system should be much simpler than the one studied in $`\mathrm{Au}+\mathrm{Au}`$ collisions. The results of $`\mathrm{d}+\mathrm{Au}`$ collisions are therefore crucial for testing the saturation approach Dima2003 .

In this Letter we present the first measurement of the minimum-bias pseudorapidity distribution of primary charged particles ($`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta `$) produced in collisions of deuterons with gold nuclei at a nucleon-nucleon center-of-mass energy, $`\sqrt{s_{_{NN}}}`$, of 200 GeV. The pseudorapidity, $`\eta `$, is defined as $`\eta =\mathrm{ln}[\mathrm{tan}(\theta /2)]`$, where $`\theta `$ is the emission angle relative to the direction of the deuteron beam. The data were obtained with the PHOBOS detector at the Relativistic Heavy Ion Collider at Brookhaven National Laboratory. The data were collected using the multiplicity array BacNim , covering $`\eta 5.4`$. The array consisted of a barrel of silicon detectors surrounding the beam pipe in the central rapidity region (“Octagon”), and six forward silicon counters, three on each side of the interaction point (“Rings”). The multiplicity array used in $`\mathrm{d}+\mathrm{Au}`$ collisions was the same as that for $`\mathrm{Au}+\mathrm{Au}`$ collisions at $`\sqrt{s_{_{NN}}}=200\mathrm{GeV}`$ BacPRL2003 . The detector setup also included two sets of 16 scintillator counters, covering $`3<\eta <4.5`$. These counters were used in the primary event trigger and in the offline event selection.

The results presented in this Letter were obtained by two analysis methods: a “hit-counting” method and an “analog” method. The details of the analysis procedure of the two methods leading to the measurement of $`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta `$ can be found in Ref. BacPRL2001 . The measured pseudorapidity distribution was corrected for particles which were absorbed or produced in the surrounding material and for feed-down products from weak decays of neutral strange particles. Uncertainties associated with these corrections, which we take as 20–50% of the corrections, dominate the systematic errors.

Due to the low multiplicity in $`\mathrm{d}+\mathrm{Au}`$, a new algorithm for collision vertex reconstruction was developed, using the hit position and energy deposited in the Octagon detector. Monte-Carlo studies using HIJING HIJING  and GEANT 3.21, as well as comparisons with the vertex information from reconstructed tracks, show a resolution of $`\sigma _{\mathrm{vtx}}(\mathrm{z})=1.4\mathrm{cm}`$ for the most peripheral events and $`\sigma _{\mathrm{vtx}}(\mathrm{z})=0.8\mathrm{cm}`$ for the most central events, $`z`$ being the distance along the beam axis.

In order to determine the minimum-bias pseudorapidity distribution for this analysis, it is necessary to correct for the trigger and vertex finding efficiencies. Since these efficiencies are functions of multiplicity (and therefore centrality), the most model-independent approach is to measure the pseudorapidity distribution in narrow bins of centrality and then integrate over centrality to produce a minimum-bias result. The alternative approach of multiplying by a global Monte-Carlo (MC) correction factor (reconstructed/true MC) would rely on our model to get the details of the shape and the centrality dependence right. To further reduce the model-dependence, we used five distinct measures for estimating the collision centrality based on the multiplicity (defined by the energy deposited in the silicon detectors) in different regions of pseudorapidity. The first of these centrality measures, $`\mathrm{E}_{\mathrm{Tot}}`$, used data from the full pseudorapidity coverage of PHOBOS, $`\eta 5.4`$. The second measure, $`\mathrm{E}_{\mathrm{Oct}}`$, used data from the more central region $`\eta 3`$. The next two measures, $`\mathrm{E}_{\mathrm{AuDir}}`$ and $`\mathrm{E}_{\mathrm{dDir}}`$, used data from the gold direction ($`5.4\eta 0.5`$) and deuteron direction ($`0.5\eta 5.4`$), respectively. The final centrality measure, $`\mathrm{E}_{\mathrm{Ring}}`$, used data far from the mid-rapidity region from the silicon ring counters at $`3.0\eta 5.4`$.

The centrality measures were calibrated using HIJING HIJING  and GEANT 3.21 simulations. In the case of $`\mathrm{E}_{\mathrm{Tot}}`$, $`\mathrm{E}_{\mathrm{Oct}}`$, $`\mathrm{E}_{\mathrm{AuDir}}`$ and $`\mathrm{E}_{\mathrm{dDir}}`$, the distributions in data and MC were very similar once the triggering and event selection were applied to the MC and an additional scaling factor near unity was included. For instance, for $`\mathrm{E}_{\mathrm{Tot}}`$ this scaling factor was 1.046. This allowed us to use the zero-bias “true” distribution in the MC to divide the data into equal 10% cross-section bins and to estimate the overall trigger and vertex-finding efficiency of $``$ 83%. Based on studies of reconstructed and true pseudorapidity distributions from different Monte-Carlo simulations, we include a global 5$`\%`$ systematic uncertainty on the pseudorapidity distribution in the reported errors. In the case of $`\mathrm{E}_{\mathrm{Ring}}`$, the MC and data distributions were sufficiently different that the procedure had to be modified. From the calculated efficiencies in each of the MC centrality bins, we could derive how many of the events in each “true” 10% bin would appear in the data. Cuts were then made on the $`\mathrm{E}_{\mathrm{Ring}}`$ distribution that gave this efficiency-corrected number of events.

Once the efficiency had been calculated and centrality cuts made, we extracted the average number of participating nucleons $`N_{part}`$ for each centrality bin using HIJING+GEANT. The measured minimum-bias pseudorapidity distribution and the estimated number of participants were obtained by averaging over all centralities for each centrality measure. This averaging was done as the minimum-bias distribution is less sensitive to centrality resolution and multiplicity bias effects than the individual centrality bins. The final distribution was obtained by averaging over the five minimum-bias distributions given by the different centrality measures. The total number of participant nucleons, $`N_{part}`$, was estimated to be $`8.1\pm 0.7(\mathrm{syst})`$. To cross check the final distribution obtained by using the silicon centrality measures, a further analysis was performed, only requiring at least one hit in one of the scintillator counter arrays. This analysis was consistent within the present systematic uncertainty.

Figure 1a shows the minimum-bias pseudorapidity distributions of primary charged particles measured for $`\mathrm{d}+\mathrm{Au}`$ collisions at $`\sqrt{s_{_{NN}}}=200\mathrm{GeV}`$ obtained with the five centrality measures described above. The distributions agree to within $`5\%`$ everywhere. Fig. 1b shows the final result obtained by averaging with equal weights over the five distributions in Fig. 1a. The systematic errors (90$`\%`$ C.L.) are shown as a gray band. The statistical errors are negligible. Fig. 1b clearly shows that the cross section for particle production is largest in the gold direction ($`\eta <0`$), and smallest in the deuteron direction ($`\eta >0`$). By fitting the gold-side using a Gaussian in the region $`3\eta 0.5`$, the peak centroid was found to be at $`\eta 1.9`$. This peak position is much closer to mid-rapidity than predicted by the parton saturation model Dima2003 , see Fig. 3a. It may be interesting to note that a double Gaussian decomposition of the Au+Au distribution at $`\sqrt{s_{_{NN}}}=200\mathrm{GeV}`$ BacPRL2003  for the $`06\%`$ most central collisions places the Gaussian centroids at $`\eta \pm 1.9`$, see Fig. 2. The significance of this observation is, however, not obvious, and may indeed be accidental.

For $`\mathrm{d}+\mathrm{Au}`$, the measured average value of the charged particle pseudorapidity density at midrapidity is $`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta _{\eta 0.6}=9.4\pm 0.7(\mathrm{syst})`$ and the integrated primary charged particle multiplicity in the measured region is $`\mathrm{N}_{\eta 5.4}^{\mathrm{ch}}=82\pm 6(\mathrm{syst})`$. Since $`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta `$ falls off to only about half of its maximum value at the most negative pseudorapidity of the detector acceptance, it is non-trivial to estimate the total charged particle multiplicity, $`\mathrm{N}_{\mathrm{ch}}`$, for this system. However, using a triple Gaussian fit to the data we obtain the solid curve in Fig. 1b, which corresponds to a total number of produced charged particles of $`\mathrm{N}_{\mathrm{tot}}^{\mathrm{ch}}=99`$ and the fit of the upper limit of the systematic errors (upper limit of the gray band) gives $`\mathrm{N}_{\mathrm{ch}}=110`$. We believe, however, that this represents an extreme upper limit for $`\mathrm{N}_{\mathrm{ch}}`$ since experimental data on p+A multiplicity distributions at lower energies wit1 ; wit2 ; wit3 ; wit4 ; wit5 , as well as various theoretical models discussed below, indicate that the fall-off at large negative pseudorapidities is more abrupt than indicated by this fit to the data. A more reasonable extrapolation, guided by the expectations based on lower energy p+A data and predictions of the AMPT model, yields an estimate of $`\mathrm{N}_{\mathrm{ch}}=5`$ in the extrapolated region. Therefore, we estimate that the total charged particle multiplicity in this reaction is $`\mathrm{N}_{\mathrm{ch}}=87_{07}^{+23}`$. It should be kept in mind, however, that even this estimate contains a significant contribution from the unmeasured region of pseudorapidity.

Figure 2 shows the normalized minimum-bias $`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta `$ distribution for $`\mathrm{d}+\mathrm{Au}`$ in comparison to the most central (0–6%) $`\mathrm{Au}+\mathrm{Au}`$ collisions BacPRL2003  and non-single-diffractive $`\mathrm{p}+\overline{\mathrm{p}}`$ collisions from UA5 UA5  at the same energy, $`\sqrt{s_{_{NN}}}=200\mathrm{GeV}`$. To be consistent with our previous Au + Au pseudorapidity distribution publications, we have used a normalization factor of $`N_{part}/2`$ which is unity for $`\mathrm{p}+\overline{\mathrm{p}}`$. The distributions from different systems are normalized chiefly to enable a comparison of their shapes. The distribution seen in d+Au collisions at RHIC has features similar to those seen in p+A collisions at lower energies wit1 ; wit2 ; wit3 ; wit4 ; wit5 . Compared to $`\mathrm{p}+\overline{\mathrm{p}}`$ collisions, there is a significant increase in particle production in the gold fragmentation region and a reduction of particle production at positive pseudorapidity in the deuteron direction. The overall normalized production of particles is approximately the same as in proton-(anti-) proton collisions UA5 .

Figure 3a shows a comparison of the measured minimum-bias pseudorapidity distribution to the published predictions of the parton saturation model Dima2003  and the results from RQMD RQMD . It is evident that these models are inconsistent with the data. The parton saturation model overestimates the height of the gold-side peak, underestimates its width, and predicts the peak at $`\eta 3`$ rather than $`\eta =1.9`$ as in the data. On the deuteron side, the model underestimates the charged particle production. RQMD drastically underpredicts the particle production and has surprisingly little pseudorapidity asymmetry.

Figure 3b presents the comparison of the measured $`\mathrm{dN}_{\mathrm{ch}}/\mathrm{d}\eta `$ distribution to the predictions of HIJING HIJING  and the published AMPT calculation AMPT . The HIJING calculation (dash-dotted curve) reproduces the deuteron side and the peak of the gold-side, but fails to reproduce the tail in the gold direction ($`\eta <2.5`$). The AMPT model uses HIJING for the initial stage of the collision, adding a “popcorn mechanism” for baryon-antibaryon production AMPT ; AMPT1  and final-state interactions (FSI). The solid curve is the default, while the dashed curve excludes FSI. We see that both final-state interactions and the popcorn mechanism appear to broaden the gold-side peak, leading to a moderate increase of the particle multiplicity in the region $`\eta 3.5`$. The ratio of the model predictions to data is shown in Fig. 3c.

In summary, the pseudorapidity distribution of charged particles produced in $`\mathrm{d}+\mathrm{Au}`$ collisions at $`\sqrt{s_{_{NN}}}=200\mathrm{GeV}`$ has been measured. The distribution is rather broad and peaked in the gold direction similar to observations in p+A collisions at lower collision energies. The measured pseudorapidity distribution is compared with the predictions of the parton saturation model as well as microscopic models. While AMPT predictions fall close to the data, the saturation model overestimates the asymmetry between particle production in the gold and deuteron hemispheres.

We thank RHIC Operations for providing the variety of colliding systems. This work was partially supported by U.S. DOE grants DE-AC02-98CH10886, DE-FG02-93ER40802, DE-FC02-94ER40818, DE-FG02-94ER40865, DE-FG02-99ER41099, and W-31-109-ENG-38, US NSF grants 9603486, 9722606 and 0072204, Polish KBN grant 2-P03B-10323, and NSC of Taiwan contract NSC 89-2112-M-008-024.


format: text

text:
# UCLA/03/TEP/30 Detectability of the Sgr Dwarf Leading Tidal Stream with Auger, EUSO or OWL

( )

## Abstract

We point out that if heavy metastable particles composing the dark matter of our galaxy are responsible for the ultra-high energy cosmic rays (UHECR) then the leading tidal stream of the Sagittarius dwarf galaxy could be detected through UHECR. The signal would be an anisotropy in the UHECR flux smaller than the telltale anisotropy towards the galactic center that would first establish unstable dark matter as the origin of the UHECR.

Cosmic ray particles with energies above the Greisen-Zatsepin-Kuzmin cutoff  of about $`5\times 10^{19}`$ eV have been detected by a number of independent experiments over the last decade . The existence of these ultra-high energy cosmic rays (UHECR) presents us with a problem. Nucleons and photons with those energies have short attenuation lengths and could only come from distances of 100 Mpc or less, while plausible astrophysical sources for those energetic particles are much farther away. Recently, results from the HiRes Collaboration  brought the violation of the GZK cutoff into question and this issue will only be resolved conclusively by the Pierre Auger Observatory.

Among the solutions proposed for the origin of the UHECR is the decay of supermassive relic particles which may constitute all or a fraction of the dark matter . Such particles must be metastable, with lifetimes exceeding the age of the Universe. In this scenario the flux $`\varphi `$ of UHECR coming from a particular volume element in the Universe is proportional to the density of dark matter in this volume, and inversely proportional to the square of the distance $`r`$ between us and this particular volume element. Thus the flux is the integral of the column density over the solid angle

$$\varphi \frac{d^3r\rho (\stackrel{}{r})}{r^2}=𝑑\mathrm{\Omega }\rho (\stackrel{}{r})𝑑r.$$
(1)
Consequently, in decaying dark matter models the UHECR are produced mostly in our galaxy and the key test is the expected anisotropy in the arrival directions caused by the offset of the Sun from the center of the galaxy . In these models cosmic rays from dark matter dominate at energies above about 6$`\times 10^{19}`$ eV, and are mostly photons . Even for protons the effects of magnetic fields on UHECR trajectories can be neglected at these energies , thus the trajectories of these UHECR are straight lines and point towards the place of origin. Spherical halo models yield an anisotropy symmetric around the direction of the Galactic Center (GC) (while triaxial models, which will not be considered here for simplicity, give considerable angular deviations). The amplitude of this anisotropy is controlled by the halo core radius $`R_c`$, a parameter which appears in the halo density models, the simplest of which is the isothermal halo model

$$\rho _{\mathrm{halo}}=\rho _o\frac{R_c^2}{r^2+R_c^2}$$
(2)
where $`r^{}`$ is the radius from the GC.

The halo core radius $`R_c`$ is not well determined from astrophysical means and could in fact be measured through the anisotropy itself. In isothermal halo models, the GC asymmetry amplitude, $`A_{GC}(\theta )`$, defined as the ratio of the flux in the direction of the GC and the flux from the galactic anti-center (AC) within an integration cone of aperture $`\theta `$, depends strongly on $`R_c`$. For example in , for $`R_c=5`$kpc, the amplitudes are found to be (see Fig. 5 of ) $`A_{GC}(5^{})5`$ and $`A_{GC}(80^{})3`$, while for $`R_c=10`$ kpc, $`A_{GC}(5^{})3`$ and $`A_{GC}(80^{})2`$.

The observation of a GC-AC asymmetry would prove the decaying particle scenario for the origin of the UHECR. In this case, the leading tail of the Sagittarius Dwarf Galaxy would give rise to a different asymmetry in the UHECR fluxes. This is the new asymmetry we are considering here.

The Two Micron All Sky Survey (2MASS) and the Sloan Digital Sky Survey (SDSS)  have traced the tidal stream of the Sagittarius (Sgr) dwarf galaxy more than $`360^{}`$ around the sky. This galaxy is a dwarf spheroidal of roughly $`10^9M_{}`$. It is a satellite of the Milky Way Galaxy, located inside the Milky Way, $``$12 kpc behind the Galactic Center (GC) and $``$12 kpc below the Galactic Plane (GP) . There are two streams of matter that extend outwards from the main body of the Sgr galaxy and wrap around the GC. These streams, known as the leading and trailing tidal tails, are made of matter tidally pulled away from the Sgr galaxy. They are on the orbital plane of the Sgr galaxy, which crosses the Galactic plane near the position of the Solar system. The leading tidal tail results from stars (and dark matter particles \[see \]) that were originally between the Sgr dwarf center and the center of the Milky Way. These stars that are closer to the GC orbit faster than the Sgr dwarf center, they move at lower gravitational potential and have shorter orbital periods around the GC. Thus, these faster orbiting stars (and dark matter particles) stream out ahead of the main body of the Sgr dwarf and create the leading tail. The stars that give rise to the leading tail are stripped at roughly a distance of 10 kpc from the GC, on the opposite side of the GC from the Sun. This distance of 10 kpc is determined  by subtracting the current distance between the GC and the center of the Sgr main body, $`16`$ kpc, minus the Sgr tidal radius, $`6`$ kpc. This distance of 10 kpc determines the perigalacticon of the leading tidal tail, both where it was removed from the Sgr main body and again on this side (the solar side) of the GC. Similarly, the trailing tail is formed from stars that were on the opposite side of the dwarf when they were stripped, they are further out and thus move slower than the Sgr galaxy main body.

The trailing stream roughly extends in a doughnut shape in the Southern galactic hemisphere around the position of the Sun, at a distance of about 20 kpc from us. The leading tail goes around the GC in the Northern galactic hemisphere, turns around at a radius of about 35 kpc, then approaches the galactic plane from the Northern galactic hemisphere and crosses this plane passing through the Solar system (see Fig. 11 of ) in the general direction orthogonal to the Galactic plane. The streams carry dark matter as well as visible stars. Thus the leading stream brings additional dark matter into the vicinity of the Solar system. The effect of the leading Sgr stream for direct detection experiments, if the DM consists of weakly interacting massive particles, has recently been studied in  where the estimated range of DM density in the stream was found to be

$$\rho _{\mathrm{stream}}=[2,200]\times 10^4M_{}/\mathrm{kpc}^3=[0.001,0.08]\mathrm{GeV}/\mathrm{cm}^3.$$
(3)
This density corresponds to a fraction $`\chi `$=(0.3-25)% of the local density of the isothermal Galactic halo for the usual estimate $`\rho _\mathrm{h}=0.3`$ GeV/cm3, and to a fraction $`\chi `$=(0.45-37.5)% of the local halo density for a lower estimate $`\rho _\mathrm{h}=0.2`$ GeV/cm3 (e.g. see the compilation of possible local halo density values for spherically symmetric halo models in ). This evaluation assumes a constant density along the stream .

The trailing stream may lead to some small increase of UHECR events coming from the Sgr dwarf plane, but this increase will be very small and would be more difficult to detect than the leading trail. Recall that the flux depends on the integral of the density over the line of sight (see Eq. (1)), and we can see along the axis of the leading trail but only transversally through the trailing tail.

We concentrate here on the leading tidal stream, which passes through the solar system. The width of the tidal stream in 2MASS M stars is estimated to be 4-8 kpc and the best-fit projection brings the center of the leading tail within 2 kpc of the Sun . We can see in Fig. 11 of  that the stream starts turning around at about (25-35) kpc from the Sun, but for smaller distances can be approximated by a cylinder. For simplicity here we model the leading stream as a cylinder of radius $`R`$=(2-4) kpc and height $`L`$=(25-35) kpc “above” the plane of the galaxy, i.e. on the side of the North Galactic Pole (NGP). We assume that the central axis of the cylinder passes through the position of the Sun. The NGP is located at $`\delta =+27^{}24^{}`$ of declination and $`\alpha =`$12h 49m of right ascension.

We then compare the expected UHECR flux $`\varphi _{}`$ in a solid angle of aperture $`\theta `$ in the direction of the NGP (i.e. along the axis of the cylinder), with the flux $`\varphi _{}`$ coming from an identical cone taken in a direction perpendicular to both the directions to the NGP and the GC (a cone with axis in the galactic plane). There are two directions with this property, one pointing to $`\delta =48^{}`$, $`\alpha =`$9h and the other pointing to $`\delta =+48^{}`$, $`\alpha =`$21h. We define the asymmetry amplitude

$$A_{\mathrm{stream}}(\theta )=\frac{\varphi _{}}{\varphi _{}}.$$
(4)
In a spherically symmetric halo without the Sgr stream there should be no difference in the fluxes from both directions, thus $`A_{\mathrm{stream}}=`$1. The presence of the stream leads to an asymmetry $`A_{\mathrm{stream}}>`$1.

Taking the evaluation of the density of the leading stream from  to be $`\rho _{\mathrm{stream}}=\chi \rho _\mathrm{h}`$, where $`\rho _\mathrm{h}`$ is the halo density at the location of the Sun, and assuming that the density is constant along the stream, we have computed $`A_{\mathrm{stream}}(\theta )`$, using the halo model in Eq. (2). Notice that $`A_{\mathrm{stream}}(\theta )`$ depends on the local overdensity $`\chi `$ due to the stream and not on the local halo density itself. In Fig. 1 we plot the column density as seen from the location of the Sun for the halo model in Eq. (2) with $`R_c`$= 5 kpc and for several of our leading stream models. The increase in column density in the directions close to the axis of the stream is clearly visible, except for the lowest values of $`\chi `$.

In Figs. 2 and 3 we plot the asymmetry $`A_{\mathrm{stream}}(\theta )`$ for different values of $`\chi `$ and two values of $`R_c=5`$ kpc and 10 kpc (Fig. 2 and 3, respectively). With $`R`$ and $`L`$ as the dimensions chosen for the cylinder representing the stream, $`A_{\mathrm{stream}}(\theta )`$ is constant for $`\theta <\mathrm{arctan}(R/L)`$, and decreases rapidly with increasing values of $`\theta `$, as can be seen in Figs. 2 and 3. The maximum asymmetry, $`A_{\mathrm{stream}}(\theta 6.5^{})=1.7`$, is obtained for the largest cylinder, lowest halo density (i.e. highest $`\chi `$) and small core radius, namely for $`R`$=4 kpc, $`L`$=35 kpc, $`\chi =`$ 0.375 and $`R_c`$=5 kpc. A larger core radius decreases the asymmetry to $`A_{\mathrm{stream}}(\theta 6.5^{})=1.57`$. Reducing $`L`$ to 25 kpc one obtains $`A_{\mathrm{stream}}(\theta 9^{})=1.5`$ for $`R_c=5`$ kpc (and 1.38 for $`R_c=10`$ kpc). Reducing the radius of the cylinder to $`R`$ = 3 kpc and reducing $`\chi `$ to $`\chi =0.25`$, with $`L`$=25 kpc the asymmetry amplitude becomes $`A_{\mathrm{stream}}(\theta 7^{})=1.36`$ (1.28) for $`R_c=`$ 5 (10) kpc. These are some of the largest expected asymmetries, while the asymmetries for the smallest values of $`\chi `$ are negligible. For example, $`A_{\mathrm{stream}}<1.005`$ for $`R`$=3 kpc, $`L`$=25 kpc and $`\chi =0.003`$.

In order to detect the asymmetry $`A_{\mathrm{stream}}(\theta )`$ above the three-sigma level we need a minimum number of events $`N`$ within the two cones of aperture $`\theta `$ that we are comparing,

$$\frac{N_{}N_{}}{\sqrt{N_{}+N_{}}}=\sqrt{N}\frac{(A_{\mathrm{stream}}1)}{\sqrt{A_{\mathrm{stream}}+1}}\mathrm{\hspace{0.17em}3},$$
(5)
which means

$$N\mathrm{\hspace{0.17em}9}\frac{(A_{\mathrm{stream}}+1)}{(A_{\mathrm{stream}}1)^2}.$$
(6)
For the maximum value of the asymmetry, $`A_{\mathrm{stream}}=1.7`$, we need at least $`N`$=50 events; for $`A_{\mathrm{stream}}=1.4`$, we need $`N`$= 135 events; for $`A_{\mathrm{stream}}=1.2`$, we need $`N`$= 500 events; for $`A_{\mathrm{stream}}=1.1`$, we need at least $`N`$= few 103 events; and for $`A_{\mathrm{stream}}<1.01`$, we would need $`N>`$ a few 105 events. These are the number of events within the cones considered. While the number of events increases as the solid angle of the cones increases, the stream asymmetry decreases for cones of aperture larger than $`\mathrm{arctan}(R/L)`$. Thus there is an optimum cone aperture which maximizes the reach of each observatory.

Let us see how large an overdensity $`\chi `$ the Pierre Auger Observatory, EUSO and OWL would be able to measure. The reach of each observatory is shown in Fig. 2 and Fig. 3 with thicker solid (red) lines. As can be seen in the figures, the integration cones for which the stream is best detectable through the asymmetry have an aperture $`\theta `$ of about $`15^{}`$ to $`20^{}`$. For larger cones, although the minimum asymmetries that can be detected decrease with increasing number of events $`N`$ within each cone, the predicted stream asymmetry decreases faster, so that the reach in $`\chi `$ does not improve.

To estimate the number of expected events, we need a value for the UHECR flux. At energies of $`6\times 10^{19}`$ eV the AGASA and HiRes cosmic ray fluxes differ by about a factor of 3 (see e.g. Fig. 1 in ). Since it is the AGASA data that do not show a GZK cutoff, we adopt the AGASA value for the flux, with an energy dependence of $`E^2`$, which is typical of decaying dark matter models . This gives us an isotropic UHECR flux integrated above $`6\times 10^{19}`$ eV of about 0.07 (km2 sr yr)-1. At higher energies the data points of AGASA, which are about one order of magnitude higher than those of HiRes, give an integrated isotropic flux of about 0.04 (km2 sr yr)-1 for $`E>1\times 10^{20}`$ eV.

The southern station of the Pierre Auger Observatory (see for example ) is already under construction since the year 2000 in the Province of Mendoza, Argentina, at a latitude of about $`35.2^{}`$ South. It will be completed by 2004. The northern station is proposed for Utah and would then enable continuous full-sky coverage. The southern array consists of 1600 particle detectors spread over 3000 km2 with fluorescence telescopes placed on the boundaries of the surface array. The expected angular resolution is less than 1 for all energies (0.5 at 1020eV). This angular resolution is very good to localize the flux from the stream. The limiting effective aperture for the full southern array is 7350 km2sr, including only zenith angles less than 60, which is the region of the sky where the reconstruction algorithm works best. Given that the latitude of the southern array is $`35.2^{}`$ South, this means that only objects at declination smaller than $`+24.8^{}`$ can be reconstructed well. However, the stream is centered at the NGP, which has a declination of $`+27.4^{}`$, and so it is outside the best detection region of Auger South.

Thus, the stream asymmetry could not be studied by the Southern Auger site unless events below a 60 zenith angle can be analyzed effectively. In this case the aperture would increase by about 50% , becoming 11250 km2sr, and there would be a total of 790 events per year with $`E6\times 10^{19}`$ eV. In the extreme case in which zenith angles as large as 85 could be reached, and we are not sure that this can be done, the number of events within the cones of solid angle $`\mathrm{\Omega }`$ we are studying would be a fraction $`f=\mathrm{\Omega }/(3.64\pi )=(1\mathrm{cos}\theta )/(1+\mathrm{cos}40^{})`$ of the total number of events observed, and only cones with aperture $`\theta <22.4^{}`$ would be completely included in the observable portion of the sky. In this case the minimum $`A_{\mathrm{stream}}(\theta )`$ observable by Auger South in 10 years as function of $`\theta `$ for $`\theta <22.4^{}`$ can be seen in Figs. 2 and 3 (thicker solid (red) line with the label Auger S). In particular, for $`\theta =20^{}`$, we have $`f=0.034`$, and after 10 years we would get $`N=268`$ events with $`E6\times 10^{19}`$ eV. This number of events would allow us to detect an asymmetry $`A_{\mathrm{stream}}(20^{})`$ 1.28. From Figs. 4 and 5, which plot the overdensity $`\chi `$ due to the stream as a function of the asymmetry within a cone of $`20^{}`$ aperture, we see that if Auger South could observe at zenith angles of up to 85, it could detect the stream if $`\chi `$ is at least $`\chi =0.3`$ (this value occurs for $`R`$ = 4 kpc, $`L`$ = 35 kpc and $`R_c=5`$ kpc).

For Figs. 4 and 5 we have chosen a cone aperture of $`20^{}`$ because this aperture is close to that one that maximizes the reach, although the latter depends somewhat on the observatory (see Figs. 2 and 3).

If the Northern Auger site would become a reality, the NGP and thus the stream would be in full view of this new observatory. The projected total effective area of 44000 km2sr  of the combined Southern and Northern Auger sites, with total 4$`\pi `$ sr sky coverage, would yield 3080 events with $`E6\times 10^{19}`$ eV. The number of events within one of the cones of solid angle $`\mathrm{\Omega }`$ we are studying is a fraction $`f=\mathrm{\Omega }/(4\pi )=(1\mathrm{cos}\theta )/2`$. For $`\theta =20^{}`$, this gives $`f=0.030`$. Thus in 10 years of combined Southern and Northern Auger data we would obtain $`N=924`$ events with $`E6\times 10^{19}`$ eV. This would allow the detection of asymmetries above the thicker solid (red) line labeled Auger S$`+`$N shown in Figs. 2 and 3. In particular we could detect $`A_{\mathrm{stream}}(20^{})`$ 1.14, which would translate into the detection of a stream with $`\chi 0.14`$ (see Figs. 4 and 5).

The Extreme Universe Space Observatory (EUSO) , a mission on-board the International Space Station, will image the air fluorescence of the extensive air showers produced by UHECR in the atmosphere. Its angular resolution is expected to be less than 2 and a conservative estimate of the duty cycle is 10%. It will operate for three years, with an effective trigger aperture which increases very fast with energy from $`1.0\times 10^4`$ km2 sr at $`E=4\times 10^{19}`$ eV to $`4.5\times 10^4`$ km2 sr at $`E1\times 10^{20}`$ eV (these effective apertures include the 10% duty cycle). The sky coverage is expected to be very close to complete, thanks to a suitable inclination of the EUSO detector axis with respect to the vertical and to the characteristics of the International Space Station orbit (the ISS orbit has a typical inclination of about $`52^{}`$ with respect to the ecliptic, a period of about a 16th of a day, and a precession of the line of nodes of about $`5^{}`$ per day ). Thus, within $`4\pi `$ sr of sky coverage, we expect about 1800 events per year with $`E1\times 10^{20}`$ eV. The number of events within the cones of solid angle $`\mathrm{\Omega }`$ we are studying is a fraction $`f=\mathrm{\Omega }/(4\pi )=(1\mathrm{cos}\theta )/2`$ ($`f=0.030`$ for $`\theta =20^{}`$). This means that in three years EUSO could detect asymmetries above the thicker solid (red) line labeled EUSO shown in Figs. 2 and 3. In particular, in three years of operation $`N=162`$ events would be collected with $`E1\times 10^{20}`$ eV for $`\theta =20^{}`$. Thus, EUSO could measure a minimum asymmetry $`A_{\mathrm{stream}}(20^{})=`$ 1.36, which is larger than the largest we obtained in our models at $`\theta =20^{}`$ (see Figs. 2 and 3).

The Orbiting Wide-angle Light-collectors (OWL)  mission would consist of a pair of satellites observing extended air showers in the atmosphere produced by UHECR with stereoscopic view. The angular resolution would be less than 1. It would have an instantaneous aperture rapidly increasing with energy, from $`4\times 10^4`$ km2 sr at $`E=3\times 10^{19}`$ eV to about $`2\times 10^6`$ km2 sr at $`E1\times 10^{20}`$ eV. Taking into account that OWL can only view the dark side of the Earth, and accounting for the effect of the Moon, clouds and man-made light, the duty cycle is about 10%. In fact the continuous effective aperture would be $`2.3\times 10^5`$ km2 sr for $`E1\times 10^{20}`$ eV  and about two orders of magnitude smaller at $`E=3\times 10^{19}`$ eV. Because of the steep rise of the aperture with energy, OWL would see more events at $`E1\times 10^{20}`$ eV than at $`E3\times 10^{19}`$ eV. With a flux of 0.04 (km2 sr s)-1, OWL would collect 9200 events per year. The number of events within the cones of solid angle $`\mathrm{\Omega }`$ we are studying is a fraction $`f=\mathrm{\Omega }/(4\pi )=(1\mathrm{cos}\theta )/2`$ of these 9200 events per year. In five years of operation the minimum asymmetry OWL could detect is shown as a function of $`\theta `$ in Figs. 2 and 3 with the thicker solid (red) line labeled OWL. In particular, OWL would detect 276 events per year in a cone of aperture $`\theta =20^{}`$. This corresponds to a detection of $`N=1380`$ events within the same cone in 5 years. Thus OWL could detect stream asymmetries above 1.12, i.e. could reach local overdensities due to the stream $`\chi 0.12`$. In terms of overdensity due to the stream, OWL’s five-year reach is better than the reach of the combined Auger South and North in 10 years.

In conclusion, EUSO would not be able to see the stream. The Southern Auger site alone cannot study the stream asymmetry unless events below a 60 zenith angle can be analyzed effectively. In the hypothetical case in which zenith angles as large as 85 could be reached, in 10 years Auger South alone could test local overdensities in the stream of $`\chi 0.3`$ (if $`R_c=5`$ kpc). This is possible only if $`\rho _\mathrm{h}0.3`$ GeV/cm3. For example, if $`\rho _\mathrm{h}=0.3`$ GeV/cm3, the minimum stream density testable by Auger South alone in 10 years would be $`\rho _{\mathrm{stream}}=`$ 0.09 GeV/cm3, while if $`\rho _\mathrm{h}=0.2`$ GeV/cm3, it would be $`\rho _{\mathrm{stream}}=`$ 0.06 GeV/cm3.

The best chance of detecting the Sgr leading tidal stream in UHECR would come either with 5 years of observation with OWL or 10 years of observation with the combined Northern and Southern Auger sites. In fact, the combined Northern and Southern Auger data, integrated for about 10 years, would allow us to detect local stream overdensities of $`\chi 0.14`$ (if $`R_c=5`$ kpc). This means $`\rho _{\mathrm{stream}}`$ 0.042 GeV/cm3 if $`\rho _\mathrm{h}=0.3`$ GeV/cm3, or $`\rho _{\mathrm{stream}}>`$ 0.028 GeV/cm3 if $`\rho _\mathrm{h}=0.2`$ GeV/cm3.

The OWL data, integrated for about 5 years, would allow us to detect local stream overdensitites of $`\chi 0.12`$ (if $`R_c=5`$ kpc). This means $`\rho _{\mathrm{stream}}`$ 0.036 GeV/cm3 if $`\rho _\mathrm{h}=0.3`$ GeV/cm3, or $`\rho _{\mathrm{stream}}>`$ 0.024 GeV/cm3 if $`\rho _\mathrm{h}=0.2`$ GeV/cm3.

The numbers given here are only indicative and we do not attempt to provide an estimate of their error. Our conclusions depend on the assumed UHECR fluxes which we took here to be at the level measured by AGASA. Clearly, the searches would become less sensitive for smaller fluxes. Our conclusions also depend on the halo model, in that the asymmetry due to the stream would decrease if the halo is flattened in the direction of the stream.

P. G. would like to thank Benjamin Stokes for helpful information on the sky coverage of UHECR experiments. We thank the referee for suggesting corrections and Dmitry Semikoz for usefull discussions. The work of G .G. was partially supported by DOE grant DE-FG03-91ER40662 and by NASA grant NAG5-13399.


text:
# Unintegrated parton distributions and particle production in hadronic collisionsPresented at the XXXIII International Symposium on Multiparticle Dynamics, 5-11 September 2003, Cracow, Poland

## 1 Introduction

The recent results from RHIC (see e.g. ) have attracted renewed interest in better understanding the dynamics of particle production, not only in nuclear collisions. Quite different approaches have been used to describe the particle spectra from the nuclear collisions . The model in Ref. with an educated guess for UGD describes surprisingly well the whole charged particle rapidity distribution by means of gluonic mechanisms only. Such a gluonic mechanism would lead to the identical production of positively and negatively charged hadrons. The recent results of the BRAHMS experiment  put into question the successful description of Ref.. In the light of this experiment, it becomes obvious that the large rapidity regions have more complicated flavour structure.

I discuss the relation between UGD’s in hadrons and the inclusive momentum distribution of particles produced in hadronic collisions. The results obtained with different UGD’s  are shown and compared. I present first results obtained with unintegrated parton distributions from the solution of the so-called CCFM equation.

## 2 Inclusive gluon production

At sufficiently high energy the cross section for inclusive gluon production in $`h_1+h_2g`$ can be written in terms of the UGD’s “in” both colliding hadrons

$$\frac{d\sigma }{dyd^2p_t}=\frac{16N_c}{N_c^21}\frac{1}{p_t^2}\alpha _s(\mathrm{\Omega }^2)_1(x_1,\kappa _1^2)_2(x_2,\kappa _2^2)\delta (\stackrel{}{\kappa }_1+\stackrel{}{\kappa }_2\stackrel{}{p}_t)d^2\kappa _1d^2\kappa _2.$$
(1)
Above $`_1`$ and $`_2`$ are UGD’s in hadron $`h_1`$ and $`h_2`$, respectively. The longitudinal momentum fractions are fixed by kinematics: $`x_{1/2}=\frac{p_t}{\sqrt{s}}\mathrm{exp}(\pm y)`$. The argument of the running coupling constant is taken as $`\mathrm{\Omega }^2=\mathrm{max}(\kappa _1^2,\kappa _2^2,p_t^2)`$.

Here I shall not discuss the distributions of “produced” gluons, which can be found in . Instead I shall discuss what are typical values of $`x_1`$ and $`x_2`$ in the jet (particle) production. Average value $`<x_1>`$ and $`<x_2>`$, shown in Fig.1, only weakly depend on the model of UGD. For $`y`$ 0 at the RHIC energy W = 200 GeV one tests UGD’s at $`x_g`$ = 10-3 \- 10-2. When $`|y|`$ grows one tests more and more asymmetric (in $`x_1`$ and $`x_2`$) configurations. For large $`|y|`$ either $`x_1`$ is extremely small ($`x_1<`$ 10-4) and $`x_2`$ 1 or $`x_1`$ 1 and $`x_2`$ is extremely small ($`x_2<`$ 10-4). These are regions of gluon momentum fraction where the UGD’s is rather poorly known. The approximation used in obtaining UGD’s are valid certainly only for $`x<`$ 0.1. In order to extrapolate the gluon distribution to $`x_g`$ 1 I multiply the gluon distributions from the previous section by a factor $`(1x_g)^n`$, where n = 5-7.

## 3 From gluon to particle distributions

In Ref. it was assumed, based on the concept of local parton-hadron duality, that the rapidity distribution of particles is identical to the rapidity distribution of gluons. In the present approach I follow a different approach which makes use of phenomenological fragmentation functions (FF’s). For our present exploratory study it seems sufficient to assume $`\theta _h=\theta _g`$. This is equivalent to $`\eta _h=\eta _g=y_g`$, where $`\eta _h`$ and $`\eta _g`$ are hadron and gluon pseudorapitity, respectively. Then

$$y_g=\mathrm{arsinh}\left(\frac{m_{t,h}}{p_{t,h}}\mathrm{sinh}y_h\right),$$
(2)
where the transverse mass $`m_{t,h}=\sqrt{m_h^2+p_{t,h}^2}`$. In order to introduce phenomenological FF’s one has to define a new kinematical variable. In accord with $`e^+e^{}`$ and $`ep`$ collisions I define a quantity $`z`$ by the equation $`E_h=zE_g`$. This leads to the relation

$$p_{t,g}=\frac{p_{t,h}}{z}J(m_{t,h},y_h),$$
(3)
where $`J(m_{t,h},y_h)`$ is given in Ref.. Now we can write the single particle distribution in terms of the gluon distribution as follows

$`{\displaystyle \frac{d\sigma (\eta _h,p_{t,h})}{d\eta _hd^2p_{t,h}}}={\displaystyle 𝑑y_gd^2p_{t,g}𝑑zD_{gh}(z,\mu _D^2)}`$ (4)
$`\delta (y_g\eta _h)\delta ^2\left(\stackrel{}{p}_{t,h}{\displaystyle \frac{z\stackrel{}{p}_{t,g}}{J}}\right){\displaystyle \frac{d\sigma (y_g,p_{t,g})}{dy_gd^2p_{t,g}}}.`$
In the present calculation I shall use only LO FF’s from  or .

Let us analyze now how the results for pseudorapidity distributions depend on the choice of the UGD. In Fig.2 I compare pseudorapidity distribution of charged pions for different models of UGD’s. In this calculation FF from  has been used.

In contrast to Ref., where the whole pseudorapidity distribution, including fragmentation regions, has been well described in an approach similar to the one presented here, in the present approach pions produced from the fragmentation of gluons in the $`ggg`$ mechanism populate only midrapidity region, leaving room for other mechanisms in the fragmentation regions. These mechanisms involve quark/antiquark degrees of freedom or leading protons among others. This strongly suggests that the agreement of the result of the $`ggg`$ approach with the PHOBOS distributions  in Ref. in the true fragmentation region is rather due to approximations made in  than due to correctness of the reaction mechanism. In principle, this can be verified experimentally at RHIC by measuring the $`\pi ^+/\pi ^{}`$ ratio in proton-proton scattering as a function of (pseudo)rapidity in possibly broad range. The BRAHMS experiment can do it even with the existing apparatus.

In Fig.3 I compare the theoretical transverse momentum distributions of charged pions obtained with different gluon distributions with the UA1 collaboration data . The best agreement is obtained with the Karzeev-Levin gluon distribution. The distribution with the GBW model is much too steep in comparison to experimental data. This is probably due to neglecting QCD evolution in .

## 4 Unintegrated parton distributions from the solution of the CCFM equation

Many unintegrated gluon distributions in the literature are ad hoc parametrizations of different sets of experimental data rather than derived from QCD. Recently Jan Kwieciński and collaborators  have shown how to solve so-called CCFM equations by introducing unintegrated parton distributions in the space conjugated to the transverse momenta

$$f_q(x,\kappa ^2,\mu ^2)=\frac{1}{2\pi }\mathrm{exp}\left(i\stackrel{}{\kappa }\stackrel{}{b}\right)\stackrel{~}{f}_q(x,b,\mu ^2)d^2b.$$
(5)
The approach proposed is very convenient to introduce the nonperturbative effects like internal (nonperturbative) transverse momentum distributions of partons in nucleons. It seems reasonable to include the nonperturbative effects in the factorizable way

$$\stackrel{~}{f}_q(x,b,\mu ^2)=\stackrel{~}{f}_q^{CCFM}(x,b,\mu ^2)F_q^{np}(b).$$
(6)
In the following, for simplicity, I shall use a flavour and $`x`$-independent form factor

$$F_q^{np}(b)=F^{np}(b)=\mathrm{exp}\left(\frac{b^2}{4b_0^2}\right)$$
(7)
which describes the nonperturbative effects. In Fig.4 I show the distribution in pseudorapidity of charged pions calculated with the help of the CCFM parton distributions  and the Gaussian form factor (7) with $`b_0`$ = 0.5 GeV-1, adjusted to roughly describe the UA5 collaboration data. Now both gluon-gluon and (anti)quark-gluon and gluon-(anti)quark fussion processes can be included in one consistent framework.

As anticipated the missing up to now terms are more important in the fragmentation region, although its contribution in the central rapidity region is not negligible. More details concerning the calculation will be presented elsewhere .

## 5 Conclusions

I have calculated the inclusive distributions of gluons and associated charged $`\pi `$’s in the NN collisions through the $`ggg`$ mechanism in the $`k_t`$-factorization approach. The results for several UGD’s proposed recently have been compared. The results, especially $`p_{t,h}`$ distributions, obtained with different models of UGD’s differ considerably.

Contrary to a recent claim in Ref., I have found that the gluonic mechanism discussed does not describe the inclusive spectra of charged particles in the fragmentation region, i.e. in the region of large (pseudo)rapidities for any UGD from the literature. I think that at presently available energies the gluonic mechanism is not the only one. Some of the missing mechanisms have been estimated in the approach based on parton distributions originating from the solution of the CCFM equation. It was found that the dominant missing mechanisms are (anti)quark-gluon and gluon-(anti)quark fusion processes.

The existing UGD’s lead to the contributions which almost exhaust the strength at midrapidities and leave room for other mechanisms in the fragmentation regions. It seems that a measurement at RHIC of $`p_t`$ distributions of particles for different rapidities should be helpful to test unintegrated parton distributions.


text:
# Single Stranded Breaks Relax Intrinsic Curvature in DNA

## I Introduction

The intrinsic curvature in DNA was discovered about twenty years ago for regular repeats of $`\mathrm{A}_\mathrm{n}\mathrm{T}_\mathrm{m},\mathrm{with}\mathrm{n}+\mathrm{m}>3`$, called A-tracts Marini et al. (1982); Wu and Crothers (1984); Hagerman (1984). Several profound reviews exist of the large volume of related data accumulated during the last decades Diekmann (1987a); Hagerman (1990); Crothers et al. (1990); Crothers and Drak (1992); Olson and Zhurkin (1996); Crothers and Shakked (1999). Due to its outstanding role in genome functioning and a very controversial experimental behavior, the sequence dependent DNA bending continues to attract great attention. The hottest point of debates is the molecular mechanism of this phenomenon. Already in the seventies, it was realized that the DNA double helix is not a rigid cylinder, but is bendable depending upon direction and the base pair sequence Namoradze et al. (1977); Zhurkin et al. (1979). By postulating that stacking in ApA steps is intrinsically non-parallel the wedge model predicted that A-tracts repeated in phase with the helical screw should bend DNATrifonov and Sussman (1980). This proved true, but the model later appeared unsatisfactory. It turned out that X-ray DNA structures are often bent in random sequences whereas A-tracts are exceptionally straight Calladine et al. (1988); Maroun and Olson (1988); Dickerson et al. (1994). With the growing number of sequences checked, non-zero wedges had to be introduced for all base pair steps De Santis et al. (1990); Bolshoy et al. (1991); Liu and Beveridge (2001), nevertheless, experimental counter-examples were eventually found where bending could not result from accumulation of wedges Dlakić and Harrington (1996, 1998). Several theories tried to reveal other physical factors involved in DNA bending along with base stacking. The most popular junction model proposed by the Crother’s groupLevene and Crothers (1983); Wu and Crothers (1984) originated from an idea that a bend should occur when two different DNA forms are stacked Selsing et al. (1979). If poly(dA).poly(dT) double helix had a special B’ form as suggested by some data Alexeev et al. (1987) the helical axis should be kinked when an A-tract is interrupted by a random sequence. Drew and Travers Drew and Travers (1984, 1985) apparently were the first to notice that narrowing of both DNA grooves at the inner edge of a bend is a necessary and sufficient condition of bending. They, and later Burkhoff and Tullius Burkhoff and Tullius (1987), considered the preference of narrow and wide minor groove profiles by certain sequences as the possible original cause of this effect. Similar ideas were considered within the context of the junction modelSprous et al. (1999). In addition, the possible involvement of solvent counterions in DNA bending was discussed long agoLevene et al. (1986), and it has been extensively studied in the recent yearsYoung et al. (1997); Rouzina and Bloomfield (1998); Hud and Polak (2001). Detailed analyses of all these theories including comparison with experimental data can be found in the recent literature Olson and Zhurkin (1996); Crothers and Shakked (1999); Mazur and Kamashev (2002). Here we note only that, in spite of all efforts, the controversies remain and continue to accumulate Merling et al. (2003); Hud and Plavec (2003).

Conformational calculations were always intensively used for testing various hypothesis concerning DNA bending Levene and Crothers (1983); Koo et al. (1990); Zhurkin et al. (1979); Ulyanov and Zhurkin (1984); Sanghani et al. (1996); von Kitzing and Diekmann (1987); Chuprina and Abagyan (1988); Zhurkin et al. (1991). The progress in methodology of molecular dynamics (MD) simulations of nucleic acids Cheatham and Kollman (2000) recently made possible realistic modeling of A-tract DNA fragments Young and Beveridge (1998); Sprous et al. (1999); Sherer et al. (1999); Strahs and Schlick (2000). In straightforward simulations of bending, the initial straight symmetrical DNA is allowed to bend spontaneously with no extra forces applied, that is due to generic atom-atom interactions Young and Beveridge (1998); Sprous et al. (1999); Mazur (2000, 2001a); Mazur and Kamashev (2002). Such computational experiments are very demanding because they involve relatively long DNA fragments and require multi-nanosecond (at least) trajectories to reveal reliable trends toward bent states. Nevertheless, for preselected A-tract sequences, it was possible to demonstrate spontaneous emergence of intrinsic curvature in the course of simulations as well as attracting character of bent states, which is necessary for static curvature. These simulations resulted in a new hypothesis of the physical origin of intrinsic curvature which we refer to as the compressed backbone theory (CBT)Mazur (2000, 2001a); Mazur and Kamashev (2002).

In a certain sense, the CBT continues the line of earlier views that attributed the origin of bending to the properties of the DNA backboneDrew and Travers (1984, 1985); Burkhoff and Tullius (1987); Sprous et al. (1999). It postulates that the intrinsic curvature in B-DNA results from a geometric mismatch between the length of the sugar-phosphate backbone and the base pair stacking. Starting from the first published drawing of the double helix Watson and Crick (1953), its two backbone strands are viewed as regular spiral traces that wind along the surface of the cylindrical B-DNA core formed by stacked base pairs. Because regular spirals represent the shortest lines that join two points on a cylindrical surface, and because the backbone strands are covalently linked to bases, an ideal B-DNA can exist only if the average backbone length either matches exactly or is shorter than that dictated by the base pair stacking Mazur (2000). The CBT assumes the opposite. It postulates that, under physiological conditions, the equilibrium specific length of the B-DNA backbone, considered as a restrained polymer attached to a cylindrical surface, is slightly longer that in the canonical B-form. In physics, similar mismatching relationship is often called geometric frustration Ediger (2000). The backbone tries to expand and "pushes" stacked bases while the stacking interactions oppose this. Eventually, the backbone finds a "compromise" by deviating from its regular spiral trace, which causes quasi-sinusoidal modulations of the DNA grooves. Concomitant base stacking perturbations cause local bends that accumulate to macroscopic static curvature when the period of these modulations corresponds to an integral number of helical turns. The CBT considers B-DNA double helices with any sequences and it may have a number of interesting consequences. It suggests, for example, that, even under normal temperature, the B-DNA double helix can occur in a glass-like state characterized by very long relaxation times Mazur and Kamashev (2002). The latter prediction agrees with some recent experimental data Song and Schurr (1990); Okonogi et al. (1999); Brauns et al. (2002). Under certain assumptions, one can even imagine the possibility of a glass-like aging of a single DNA molecule in vivo, with direct relation to the problems of biological aging.

One of the possible experimental tests of CBT is based upon the predicted effect of single-stranded breaks (nicks) upon the intrinsic curvature Mazur and Kamashev (2002). Reversible single chain breaks in DNA occur during diverse biochemical processes, including general and site-specific recombination, replication, and DNA repair. Therefore, structural perturbations due to nicks and single-stranded "gaps" have been earlier characterized in a number experimental studies Thomas (1956); Hays and Zimm (1970); Shore and Baldwin (1983); Hsieh and Griffith (1989); Pieters et al. (1989); Snowden-Ifft and Wemmer (1990); Aymami et al. (1990); Rentzeperis et al. (1993); Le Cam et al. (1994); Mills et al. (1994); Hagerman and Hagerman (1996); Furrer et al. (1997); Lane et al. (1997); Singh et al. (1997); Roll et al. (1998); Kozerski et al. (2001); Kuhn et al. (2002). According to all these data, nicks perturb the free DNA very little. The corresponding X-ray structuresAymami et al. (1990) are straight and similar to those of intact duplexes with analogous sequences. This agrees with early physicochemical tests Thomas (1956); Hays and Zimm (1970); Shore and Baldwin (1983) as well as NMR studies in solution Pieters et al. (1989); Snowden-Ifft and Wemmer (1990); Singh et al. (1997); Roll et al. (1998); Kozerski et al. (2001). The main detectable effect of a nick is local melting or fraying, especially at low ionic strength, mild denaturating conditions, or elevated temperatures. This increases the isotropic flexibility of DNA and reduces its gel mobility Le Cam et al. (1994); Mills et al. (1994); Furrer et al. (1997); Kuhn et al. (2002). In normal conditions the nick fraying is small and it almost undetectable below 10CMills et al. (1994); Kuhn et al. (2002).

Based upon these experimental data and the earlier views of the origin of intrinsic curvature one should expect that nicks in phased A-tract sequences can either increase bending or have zero effect. The increase may be expected if the curvature is caused by external forces, as in electrostatic counterion models, because the flexible hinge at a nick can allow these forces to increase the bend angle in their direction. A zero effect is expected for the wedge and junction models because nicks arguably do not affect the A-tract structures Aymami et al. (1990) and do not cause wedges detectable in gel migration assays Mills et al. (1994). In contrast, if the intrinsic curvature is really caused by the backbone compression, single-stranded breaks should relax it, which should increase the gel mobility of curved DNA fragments. This effect is opposite to that of the increased isotropic flexibility, moreover, the latter can be eliminated by reducing the temperature so that these two phenomena cannot be confused. In addition, because the degree of the compression should vary regularly between the inner and outer edges of the bent DNA, relaxation produced by a single nick is expected to change systematically according to its position with respect to the overall bend.

Here we present the results of experimental and simulation studies carried out according to the above plan. Two 35-mer sequences from our recent comparative investigationMazur and Kamashev (2002) were used for constructing two series of nicked DNA fragments. The first sequence contained an A-tract repeat motif that was preselected in MD simulations to reproducibly induce strong static curvature. The second sequence was “random”, but with the same base pair composition. Backbone breaks were introduced in one strand at different portions to span approximately one helical turn. The curvature was probed by PAGE as well as long-time MD simulations. In experiment, a noticeable increase in the gel mobility of the nicked A-tract repeat is observed depending upon the nick position. In contrast, no such effect is found for the random fragment. MD simulations fail to reproduce experimental effects qualitatively, but qualitatively they also confirm that single stranded breaks interfere with the A-tract induced curvature. The results qualitatively agree with CBT, but they cannot be accounted for by other models.

## II Materials and Methods

### II.1 Oligonucleotides and construction of 5’-labeled DNA probes

The double stranded DNA containing four A-tracts flanked poly(dA).poly(dT) termini were constructed by annealing the synthetic oligonucleotide B (the top line in Table 1) with its complementary (bottom) strand. Double stranded DNA fragments containing single stranded breaks (nicks) were constructed by annealing two shorter oligonucleotides (Table 1) with the same bottom strand. In the same way a series of nicked duplexes was constructed from the reference random sequence fragment S, with its base pair composition identical to that of B. The sequences of the oligonucleotides used for the construction of the DNA fragments are assembled in Table 1. The fragment codes in this table use the following mnemonics. The capital ’B’ and ’S’ stand for "bent" and "straight". The small letter prefix ’n’ stands for "nick". The numbers indicate the shift of the nick position with respect to the center of the sequence.

In all the constructs, one of the oligonucleotides annealed was labeled with T4 polynucleotide kinase and \[32-P\]-ATP. In most cases the label was attached to the 5’-end of the bottom (continuous) oligonucleotide whereas the two strands at the nick position were terminated with 3’OH and 5’OH groups. For some sequences, 5’-phosphorylated nicks were additionally constructed to check the effect of phosphate charges. In this case the continuous strand carried 5’-OH while \[32-P\] was introduced in the corresponding shorter partner. The annealing was carried out by incubating 300 nM of the unlabeled oligonucleotide(s) with 30 nM of end-labeled ones for 3 min at 80C in 20 mM Tris-HCl (pH 8.0), 150 mM NaCl, and then allowing them to cool slowly. To distinguish between the two nick series, in the text, we add small letter suffixes ’o’ and ’p’ to the codes given in Table 1. Thus the whole two nick series of the bent DNA are referred to as nBo and nBp, respectively, while the individual nicked fragments are denoted as nB-2o, nB-2p, etc.

### II.2 Gel mobility assays

The mobility of the DNA fragments was analyzed in 12% gels (acrylamide to bis-acrylamide, 29:1) buffered with 90 mM Tris-borate, 1 mM EDTA, pH 8.6. In order to check the effect of Mg2+ ions the buffer was supplemented with 10 mM MgCl2 in the absence of EDTA. Gels were pre-run under constant power until stabilization of the current. Labeled DNA in a buffer containing 20 mM Tris-HCl (pH 8.0), 50 mM NaCl, 4% Phicoll-400 and xylencianol was loaded onto the gel. The electrophoresis was performed under constant voltage and constant temperature of 4C. The dried gels were exposed to storage phosphor screens and visualized on a 400S PhosphorImager (Molecular Dynamics).

### II.3 Calculations

MD simulations were carried out for the central 35 base pairs of six nicked DNA fragments in Table 1 constructed from the A-tract repeat sequence (nB-4o, …, nB+6o). The corresponding calculations for the two mother fragments without nicks (B and S) were reported earlierMazur and Kamashev (2002). For consistency with the previous studiesMazur (2000, 2001a); Mazur and Kamashev (2002), we used the AMBER98 Cornell et al. (1995); Cheatham et al. (1999) force field and TIP3P water Jorgensen (1981). Trajectories generally started from the fiber canonical B-DNA conformationArnott and Hukins (1972) and were continued to 20 ns. In a few cases, dynamics was also run starting from the canonical A-DNA conformation. In total, these simulations sampled from more than 0.2 ms of all-atom dynamics for nicked and intact 35-mer DNA fragments with the same A-tract repeat sequence.

Molecular dynamics simulations were carried out with the internal coordinate molecular dynamics (ICMD) method Mazur (1997, 2001b) adapted for DNA Mazur (1998a, 1999a) with the time step of 0.01 ps. In this approach, the DNA molecule has all bond length and almost all bond angles fixed at their standard values. The only variable bond angles are those centered at the sugar C1’,…,C4’, and O4’ atoms, which assures the flexibility of the furanose rings. In contrast, bases, thymine methyls, and phosphate groups move as articulated rigid bodies, with only rotations around single bonds allowed. The highest frequencies in thus obtained models are additionally balanced by increasing rotational inertia of the lightest rigid bodies as described earlierMazur (1998b, a). The possible physical effects of the above modifications have been discussed elsewhere Mazur (2001b); Mazur et al. (2001).

The so-called "minimal model" of B-DNA was used Mazur (1998a, 2001c). It includes only a partial hydration shell and treats counterion and long range solvation effects implicitly by reducing phosphate charges to -0.5 and applying linear scaling of Coulomb forces. This model produces B-DNA structures very close to experimental data and has no other bias towards bent or non-bent conformations except the base pair sequence. The interest to implicit solvation in DNA simulations is long standing and it continues in the literature Wang et al. (2002). Advantages as well as limitations of such approaches have been recently reviewed Cheatham and Kollman (2000).

The minimal model is not meant to be generally applicable and we have chosen it for the present studies because of the following considerations. We assume that effects it ignores, like sequence specific counterion binding, play some role, but are not critical for the A-tract induced curvature. We have shown earlier that, in these conditions, the A-tract repeat motif we are using has a bent state with very strong attracting properties that allow one to observe spontaneous transitions to stable bent conformations in the course of MD Mazur (2000); Mazur and Kamashev (2002). This property is rather exceptional and only due to it one could hope to obtain useful information from a comparison between nicked and intact DNA in MD. We also carried out several simulations for the mother 35-mer and 25-mer A-tract fragments by the Particle-Mesh Ewald methodEssmann et al. (1995) with full solvation and periodical boundaries. The corresponding trajectories were continued to 5-10 ns, but they showed much less demonstrative bending dynamics. The origin of this difficulty is currently unclear, but it agrees with recent reports by some other authors Cheatham and Young (2000-2001). Finally, the minimal model provides for substantial savings in computations. The saving factor is critical because DNA bending is detectable only in relatively long fragments and it is likely to involve very slow motions. In spite of a few reported analogous simulations using full-scale solvation Young and Beveridge (1998); Sprous et al. (1999) such calculations would be too costly for the volume of sampling reached here.

## III Results and Discussion

### III.1 Construction of DNA Fragments

The A-tract motif $`\mathrm{AAAATAG}`$ originally attracted our attention in MD simulations of the natural DNA fragment taken from the first curved DNA locus studied in vitroWu and Crothers (1984); Mazur (2001a). The central 35-mer A-tract repeat for the bent fragments in Table 1 was constructed by repeating this motif four times and it had to be inverted to make the two DNA termini symmetrical. Such inversion should not affect bending, Koo et al. (1986) but is essential for simulations because the 3’- and 5’-end A-tracts may represent qualitatively different boundaries. In repeated simulations with this and similar A-tract fragments, the static curvature emerged spontaneously and it became more evident as the chain length increased Mazur (2000); Mazur and Kamashev (2002). To obtain a reference non-A-tract DNA, we have re-shuffled manually base pairs of the A-tract repeat. We preferred this randomized sequence to commonly used GC-rich straight fragments in order to keep the base pair content identical and reduce the noise that could cause small variations in gel mobility. In order to amplify the PAGE resolution of similarly bent nicked DNA, the 35-mer fragments were extended to 71 bp by adding poly(dA).poly(dT) tails. The tails continue the mother 35-mer A-tract repeat sequence smoothly, which should reduce the possible perturbations that could affect the comparison of PAGE results with MD simulations that involved only the central 35-mer nicked fragments.

### III.2 Relative PAGE Mobilities

A representative PAGE plate of the random sequence series of non-phosphorylated nicks is shown in Fig. 1. All nicked DNA exhibit similar mobilities close to that of the reference straight fragment. A weak retardation effect distinguishable for some nick positions is at the limit of experimental accuracy. Qualitatively similar patterns, but with stronger retardation were earlier reported for nicked DNA fragments under elevated temperature and/or mild denaturating conditions Le Cam et al. (1994); Mills et al. (1994); Kuhn et al. (2002). We suppose, therefore, that some reduction of mobility in Fig. 1 results from the isotropic flexibility at nick positions that is strongly reduced at 4C. This small retardation should be taken into account in the interpretation of other tests below.

Fig. 2 shows the results obtained in the same conditions for the A-tract nicks. The pattern exhibited here is evidently more complex than that in Fig. 1. On the one hand, for the central nick positions, the mobilities are close to that of the mother fragment, which is qualitatively similar to Fig. 1. However, when the nick is moved farther from the center, the PAGE mobility grows and becomes higher than that of the mother curved fragment. It reaches maximum values for nB-4o and nB+6o, but for nB-6o it is reduced again. As a result, the mobilities of nB-6o and nB-4o appear close to those of nB+4o and nB+4o, respectively, and the overall pattern of nick bands in Fig. 2 appears skewed-sinusoidal, with the period corresponding to that of the double helix. The increase of PAGE mobility indicates reduced DNA curvature and this effect can well be due to the relaxation predicted by CBT. The amplitude of the reduction is very significant since it exceeds 30% of the curvature difference between the S and B fragments produced by four A-tracts. This mobility profile suggests that the backbone compression is maximal within A-tracts and minimal between them. Although this specific phasing is right opposite to our first guessMazur and Kamashev (2002) it is accounted for by CBT and not by other theories as we discuss further below. However, there are a few alternative interpretations that should be considered first.

The first evident factor that should be checked is the perturbed balance of phosphate charges. The two nick series displayed in Fig. 1 and Fig. 2 have one phosphate group less than the reference fragments B and S. Consequently, their mobilities may be uniformly reduced by around 0.7%, which would mean that the apparent reduction of curvature in the nicked A-tract fragments nB-4o and nB+6o is even stronger than it is seems from Fig. 2. This effect, however, is not distinguishable in Fig. 1 suggesting that it is too small and may be safely neglected. The influence of phosphate charges can be more subtle, however. The electrostatic models of DNA bending assume that it is the broken balance of phosphate repulsion at the opposite DNA sides that forces it to bend Hud and Polak (2001). If the A-tract curvature had an electrostatic origin, the phosphate "hole" at the nick could perturb it by adding a local bend in a direction that should rotate as the "hole" is moved along the DNA chain. Like that the integral curvature would be deviated and partially compensated in a way compatible with the results shown in Fig. 2.

To check this possibility, similar experiments were carried out for a series of analogous phosphorylated nicks. Their relative PAGE mobilities are compared with other results in Fig. 3. The mobility coefficients Qm used in this plate were computed as follows. The Qm value of the S fragment (see Table 1) is always the largest in the gel and it is arbitrarily assumed to equal 100. The Qm of fragment B is assumed to equal 0. For any given band its Qm is promotional to the distance from the fastest S band in the same gel and it is estimated relative to the distance between B and S. If a nicked fragment migrates slower than the B-band in the same gel its Qm is negative. In contrast, if it is faster, the corresponding Qm is positive below 100. Fig. 3 displays thus obtained Qm values for all DNA fragments used in the present study.

With free 5’-phosphates added to the nick sites, the phosphate "holes" are quenched or even inverted, which should drastically change the observed curvature modulations. The results in Fig. 3, however, show the opposite. The original skewed-sinusoidal dependence is persistent. At the same time, the phosphorylated nicks exhibit uniformly reduced mobilities, which agrees with earlier reportsLe Cam et al. (1994); Mills et al. (1994) and should be attributed to the increased strand fraying in nicks. The additional phosphate charge increases the inter-strand electrostatic repulsion at single stranded breaks, which should enhance the fraying effect as discussed in detail elsewhere Furrer et al. (1997). Whatever the reason of the uniformly reduced mobility, however, Fig. 3 clearly shows that the above hypothesis of electrostatic curvature compensation can be ruled out.

This figure also shows the profile of mobilities of the nBo series in the presence of Mg2+ ions. The Mg2+ ions are long known to affect the DNA curvature by increasing or reducing it depending upon the sequence Diekmann and Wang (1985); Diekmann (1987b). In our case the curvature was increased as judged from the absolute PAGE separation of fragments B and S. It is seen in Fig. 3 that the profile of the nick mobilities in the presence of Mg2+ ions remains qualitatively similar, but some nicks migrated slower than the mother B fragment. Interpretation of this effect depends upon the mechanism by which Mg2+ increase the A-tract curvature. If the additional bend is caused by the electrostatic attraction due to Mg2+ positioning between the phosphate groups in the narrowings of the major groove than the gain in the curvature may be larger when the DNA is nicked in the widening of the minor groove at the opposite side of the bend. This interpretation would agree with Fig. 3 because the negative Qm are observed for nick positions outside A-tracts in the probable widenings of the minor groove. It is also possible, however, that Mg2+ ions somehow increase the nick fraying especially when the nick position approaches the poly(dA).poly(dT) tails. In this case, the negative Qm values would have the same origin as for the phosphorylated nicks.

The results in Fig. 2 can also be interpreted with the wedge model of DNA bending by assuming that nicks have a special structure with a wedge in a fixed local direction. The overall curvature is reduced when this wedge direction is opposite to that of the initial bend and increased in the opposite orientation. Similarly to the "phosphate hole" model above, skewed-sinusoidal modulations of curvature should occur as the nick site is translated along DNA. Note that this wedge should have a very large angle since it compensates in nB-4o and nB+6o more than 30% of the curvature due to four phased A-tracts. A strong wedge like that should result in an increased curvature in nB+0o as well as in the nSo series in Fig. 1, which is not seen. This interpretation also would not agree with earlier studies of nick structures and gel mobilities Aymami et al. (1990); Mills et al. (1994); Furrer et al. (1997). However, earlier PAGE studies were carried out at higher temperatures where the isotropic flexibility of nicks could, in principle, dominate all other effects. Also, the non-linear dependence of the gel mobility upon the curvature may be such that the acceleration effect of the nick wedge in Fig. 2 is more pronounced than the expected retardation for nB+0o and the nSo series. Therefore, in order to check the above explanations, double nicked A-tract fragments were constructed, with the two sites separated by an odd number of helical half-turns so that the hypothetical wedges had opposite directions and should have compensated one another. The PAGE mobilities of double nicked fragments are compared in Fig. 4. These data evidently disagree with the foregoing wedge interpretation of mobility modulations. For both nB-4o and nB+6o, the addition of a second nick at 1.5 helical turn neither cancel nor even reduce the corresponding band shifts with respect to Fig. 2. As expected, analogous constructs for fragment S all show PAGE mobilities similar to that of the intact duplex.

The experimental results presented above evidence that single stranded breaks really reduce the overall curvature of the A-tract repeat in remarkable qualitative agreement with CBT. Other theories of DNA bending cannot account for these data unless strong ad hoc assumptions are introduced like sequence-and-position dependent nick wedges, for example.

### III.3 DNA Curvature in Simulations

As regards the overall stability of dynamics and its closeness to the B form, all MD trajectories computed in the course of this studies were similar to our earlier simulations under similar conditions Mazur (1998a, 2001c, 2000, 2001a); Mazur and Kamashev (2002). Table 2 shows parameters of the final 1ns-average conformations. They all have remarkably similar average helicoidals corresponding to a typical B-DNA. For example, the average helical twist estimated from the best-fit B-DNA experimental values Kabsch et al. (1982) gives $`34.0\pm 0.2`$ and $`33.8\pm 0.2`$ for the A-tract fragment and the randomized sequence, respectively. The nicks did not cause readily visible perturbations and it was commonly difficult to distinguish them in snapshots. Statistical analysis of the local inter-base parameters at nick positions is summarized in Table 3. As expected, the broken DNA strand exhibits somewhat increased flexibility at the nick site. The most noticeable are much stronger deviations of Twist form its canonical values. A similar feature is smaller, but distinguishable for Tilt and Roll. On the other hand, the nanosecond time scale fluctuations characterized by the standard deviations are not very different with and without the nick. The last observation suggests that, as long as the bases at the nicked step remain in the stack, their rapid motions are restrained to the same degree as in the intact structure, and that the above strong deviations of some averages are due to slower global motions that occur in much longer time scales. Apart from these strong fluctuations, however, we did not observe any repetitive static perturbations that might be attributed to strand breaking. Moreover, strongly non-canonical helical parameters are sometimes encountered in intact steps as well. As it was for dynamics of the mother A-tract DNA fragmentMazur and Kamashev (2002), the rms deviations from the canonical structures seen in Table 2 are largely due to variation of the overall curvature.

Fig. 5 compares the overall DNA curvature in the last five 1ns averaged conformations from all trajectories mentioned in Table 2. It is seen that the mother A-tract fragment (B) reached the largest bending amplitudes non accessible to other fragments. The most stable and significant curvature among the nicks is obtained for nB-4o and nB+0o. Detailed in Fig. 6 is the dynamics of the bend direction. Its short time scale fluctuations are smaller when bending is strong, therefore, Fig. 6 also gives an estimate of the bending amplitude. Again we note that the mother A-tract fragment (B) is distinguished by the rapid establishment of curvature. However, convergence to a statically bent state is clear for nB-4o and nB+0o as well. Following to our previous approachMazur (2000, 2001a); Mazur and Kamashev (2002) an attempt has been made to confirm the attracting property of the bent state for nB+0o by running dynamics starting from the canonical A-DNA conformation. The corresponding trajectory was continued to 17 ns and it showed convergence of the same bending direction as in Fig. 6, but with a somewhat smaller amplitude.

Fig. 5 and Fig. 6 show that our MD simulations failed to reproduce a regular variation of the bending magnitude with the nick position as seen in the PAGE experiments above. Fragments nB-4o and nB+0o exhibit noticeable curvature comparable, but smaller than that of the intact duplex. According to experiment, the curvature in nB+0o is similar to that in the mother fragment, but for nB-4o it should be definitely smaller. Fragment nB+6o exhibited less stable bending in a similar direction whereas for nB-2o, nB+2o, and nB+4o the 20 ns duration of trajectories apparently was not sufficient for convergence. At the same time, our calculations confirm that, in qualitative agreement with experiment, single strand breaking interferes with bending and generally reduces it. It is remarkable in Fig. 5 that strong bends beyond 40 were observed only in directions roughly corresponding to that in the mother fragment, which confirms the strong attracting property of this bent state. Taking into account the unclear physical nature of the phenomenon being studied and the complexity of the model system we consider the partial qualitative agreement achieved here as very satisfactory.

In spite of this only limited agreement with experiment, we tried to use the ensemble of computed curved DNA conformations to clarify the origin of the effect produced by nicks upon the A-tract curvature. To this end, for each of the fragments B, nB-4o, nB+0o, and nB+6o a nanosecond of dynamics have been selected for a more detailed analysis. The corresponding 1ns average structures are marked by arrows in Fig. 5. They all correspond to the largest degree of curvature reached during these trajectories. For B and nB+6o these are the last nanosecond structures. In contrast, for nB-4o and nB+0o the 17th and the 19th nanoseconds were used, respectively. According to Fig. 5, these bends are more or less convergent suggesting that the corresponding dynamics should be similar and subtle differences caused by nicks might be detectable. As seen in Fig. 7, for fragments B, nB-4o, and nB+6o the averaged conformations are rather close to one another. The curved DNA axis is almost planar with two or three zones of increased bending found between the A-tracts. In contrast, for nB+0o the curved axis is strongly non-planar due to significant deviation of bend direction in its central zone from the other two. As a result, for B, nB-4o, and nB+6o the overall bend is close to the sum of local bends, but for nB+Oo it is much smaller. The minor groove profiles averaged over the corresponding 1 ns periods are compared in Fig. 8. Again one can note a similarity between B, nB-4o, and nB+6o versus a more complex profile for nB+Oo. The widenings of the minor groove are found between A-tracts and they roughly correspond to the zones of increased bending in Fig. 7. According to Fig. 7 and Fig. 8, a reduced overall bending amplitude in nB-4o and nB+6o with respect to the mother A-tract fragment is attributable to the central zone where the widening is also less significant. In contrast, the smaller curvature in nB+0o is mainly due to the lost phasing between the local bends.

Relaxation of A-tract curvature by single stranded breaks has been originally proposed as a conceptual test of CBTMazur and Kamashev (2002). The results reported here qualitatively agree with its predictions and strongly support this model versus its alternatives. At the same time, the nick position dependence as revealed by our PAGE tests is inverted with respect to to the original guessMazur and Kamashev (2002). We thought that the backbone compression is maximal in the widenings of the minor groove, therefore, nicks between A-tracts should relax bending stronger than nicks inside them. In contrast, the experiments reported above suggest that the backbone compression should increase inside A-tracts and reduce outside them. Because A-tracts are found at the internal edge of the curved double helix the last suggestion agrees with the simple physical intuition that says that the surface of a curved cylinder is compressed at its inner edge and stretched at the opposite side. In order to understand which backbone component can be responsible for the compression/frustration relationship postulated by CBT we tried to measure the specific backbone length by using pairs of similar atoms in neighboring residues.

For all atom pairs the resulting profiles of the backbone length appear very noisy and regular modulations can be seen only after averaging with a sliding window. Inspection of all atom pairs allows one to divide them in two groups. For atoms of phosphate groups and their neighbors, the measured backbone length usually is larger at the inner side of the bend and smaller at the opposite side, which may be viewed as a trend towards the A-from in the widenings of the major groove. However, these modulations are very weak, in agreement with the conventional view of the B-DNA backbone as being non-compressible Packer and Hunter (1998); Calladine and Drew (1992). In contrast, for a group of backbone atoms including O4’ and its close neighbors, strong regular oscillations of the backbone length could be revealed, with their phases being opposite to those in the first group. The average backbone profiles produced by inter-O4’ distances are shown in Fig. 9. They were computed for the 1 ns intervals of trajectories marked in Fig. 5. It is seen that the backbone length correlates with the minor groove profiles in Fig. 8, but does not repeat them. Thus measured specific backbone length reaches its local minima at the inner edge of the bent DNA cylinder and corresponds exactly to the nick positions in nB-4o and nB+6o. Similar oscillations are observed in X-ray structures of bent DNA available in NDBBerman et al. (1992) (AKM, unpublished). It may seem surprising that profiles in Fig. 9 exhibit no perturbations at strand breaks, but this results from smoothing. In the non-smoothed plots the single stranded breaks in nB+6o and nB+0o were distinguished by very small and very large distances, respectively, but they have been compensated by strong opposite deviations in the neighboring steps.

The foregoing analysis demonstrates that the DNA backbone behaves as a complex mechanical system rather than a chain of strings or rigid rods, and that different its components compress and stretch in different zones of the curved DNA. Nevertheless, the correlations revealed in Fig. 9 suggest that the backbone compression, if it really exists, should mainly affect sugar-sugar interactions and that they can be at the origin of the intrinsic backbone frustration. Even though in B-DNA consecutive O4’ atoms do not interact, the neighboring sugar rings make a number of direct contacts non-mediated by the solvent.

## IV Concluding Discussion

The results presented here provide new insights in the putative physical mechanism of intrinsic DNA curvature. This intriguing phenomenon seems to be exhaustively studied, but it still attracts great interest. The A-tract repeats and single stranded breaks both have been earlier recognized as elements specifically involved in DNA bending. To our best knowledge, however, nobody checked what could be their cumulative effect. Here it is found that nicks relax the intrinsic curvature induced by A-tract repeats in conditions where they do not affect the random sequence DNA. The relaxation effect depends upon the nick position with respect to the bend. It is very significant when the nick is at the inner edge of the bent DNA inside A-tracts, and is gradually reduced essentially to zero when the nick is moved to the outer surface of the bend. This behavior is in qualitative agreement with the CBT and it suggests that the backbone compression in the curved DNA is larger in the minor groove narrowings at the inner surface.

The MD simulations confirm that single stranded breaks interfere with, and tend to relax the A-tract induced curvature. For about a half of the DNA fragments studied, the duration of trajectories provided reproducible convergence to statically bent states in good agreement with experiment as well as earlier simulations. Inspection of the computed conformations reveals strong regular modulations of the local backbone length as measured by distances between certain sugar atoms, with compression reaching its maximum at the inner edge of the bend where nicks produce the strongest relaxation in experiments. These results provide additional support to the CBT and suggest that the frustration in the B-DNA backbone may result from interactions between consecutive sugar rings.