text:
# Recoiling DNA Molecule: Simulation & Experiment

## 1 Introduction

One of the most important characteristics of the DNA molecule is its high mechanical flexibility. Thanks to this, the entire genome of each living organism, which can reach several meters in length, fits inside each of its cells. The flexibility of the DNA molecule also plays a key role in all its cellular functions, such as replication, transcription and recombination. With advances in micromanipulation techniques during the last decade single molecule DNA experiments have become common, and are explored through many different approaches .

Once free in solution, single DNA molecules present entropic elastic behavior, well described, in both low and high stretching regimes, by the *worm-like chain* (WLC) model . The flexibility of the polymeric chain can be inferred through the ratio $`A/L`$, where $`A`$ is its *persistence length* and $`L`$ its *contour length*. Smaller $`A/L`$ ratios correspond to more flexible chains. For the DNA molecule, a complete chain can have contour length L varying from around $`2\mu m`$ to more than $`1m`$, while the persistence length $`A`$ is of the order of $`4050nm`$ or $`120150`$ base pairs.

The usual method to obtain $`A`$ employing an optical trap  is based on *force versus extension* measurements . One end of a DNA molecule is attached to a coverglass surface while the other end is tethered to a polystyrene microsphere *(Fig. 1)*. The optical trap can then be used to capture and move the microsphere, stretching the DNA molecule. Obtaining force versus extension data for a DNA molecule, however, requires a very precise characterization and calibration of the trap, allowing the determination of the trapping force as a function of the position of the trapped microsphere. Here we examine an alternate method to determine the persistence length $`A`$: once the DNA is stretched, one removes the trap and follows the motion of the microsphere in time, using videomicroscopy .

Using this method, Feingold  obtained a value of $`152nm`$ for the persistence length of DNA, approximately three times larger than the value obtained from force versus extension measurements. The proposed explanation for this discrepancy was that during the initial part of the recoil, when the DNA is more than $`80\%`$ stretched, the motion cannot be considered quasistatic, such that the force at each position is time dependent. Recently, Bohbot-Raviv *et al*  considered a complete nonequilibrium theory for the relaxation of highly stretched semiflexible polymers like DNA and found a better agreement with Feingold’s experiment.

We have performed recoil experiments in which the friction between the microsphere and the coverglass is carefully taken into account and found persistence lengths comparable to the values obtained from force versus extension measurements just using the standard quasistatic approximation, without the need of a nonequilibrium theory for the relaxation of the DNA. We also perform computer simulations of the recoil using the corresponding Langevin equation to assure the validity of the assumptions used in the data analysis.

The balance of this paper is organized as follows. In section 2 we describe the recoil method and the motion of the microsphere. Section 3 details the simulation of the recoil dynamics of the microsphere. The procedures adopted to obtain experimental recoil curves are described in section 4. In section 5 we describe the analysis of simulated and experimental recoils leading to our results. Conclusions and final remarks are presented in section 6.

## 2 Experimental Technique

In essence the recoil technique introduced by Shivashankar *et al*  is very simple. As in the force-based experiments, optical tweezers are used to find and trap a microsphere anchored to a coverglass surface through a single tethered DNA molecule. The trapped microsphere is moved, stretching the DNA molecule. The tweezers are then cut off suddenly, releasing the microsphere, which is dragged through the fluid by the recoiling DNA molecule. Following the recoil of the microsphere using videomicroscopy, we obtain a recoil curve $`R(t)`$, from which we can extract information about $`A`$ and $`L`$.

According to , the recoil movement can be considered *quasistatic* since the equilibration time for the DNA molecule is very short compared with the time scale for the motion of the microsphere. This means that, during relaxation, the DNA essentially passes through a sequence of equilibrium stretched states . Because the time resolution of our videomicroscopy experiments ($`(1/30)s`$) is large compared to the ratio $`m/\gamma `$ ($`10^7s`$) for the microsphere, inertial effects are negligible in the analysis of the recoils. As a result, the average motion of the microsphere can be described by

$$\gamma \frac{dR}{dt}=F_{dna}(r),$$
(1)
where $`R`$ is the average position of the microsphere, $`r`$ is the distance between the extremities of the DNA molecule, $`\gamma `$ is the Stokes friction coefficient and $`F_{dna}(r)`$ is the entropic elastic force for the DNA molecule given by the WLC model ,

$$F_{dna}(r)=\frac{k_bT}{A}\left[\frac{r}{L}+\frac{1}{4\left(1\frac{r}{L}\right)^2}\frac{1}{4}+\underset{i=2}{\overset{i7}{}}\alpha _i\left(\frac{r}{L}\right)^i\right],$$
(2)
where the final term represents the correction introduced by Bouchiat *et al* in  with $`\alpha _2=0.5164228`$, $`\alpha _3=2.737418`$, $`\alpha _4=16.07497`$, $`\alpha _5=38.87607`$, $`\alpha _6=39.49944`$ and $`\alpha _7=14.17718`$. The bouyancy force on the microsphere is $`5fN`$, negligible compared to $`F_{dna}`$, which is of the order of $`pN`$ during most of the recoil movement.

Far from any boundary, the Stokes friction coefficient between the microsphere and the surrounding medium is

$$\gamma =6\pi \eta a,$$
(3)
where $`\eta `$ is the viscosity of the medium and $`a`$ is the radius of the microsphere. There are, however, important considerations to be made regarding the use of Eq. (1). The wall-drag effect caused by the coverglass nearby changes the friction coefficient $`\gamma `$ of the microsphere . As a result, $`\gamma `$ can no longer be considered constant and isotropic, and must continuously be corrected as the microsphere moves. This correction depends on the height $`h`$ of the microsphere with respect to the coverglass below.

We can integrate Eq. (1) directly. However, we cannot obtain a simple explicit analytic expression for $`R(t)`$, forcing us to devise an alternative procedure to fit the experimental data. We chose to apply a discrete derivative to the experimental recoil curves, generating data that may be fit directly with Eq. (1). Since the microsphere undergoes Brownian motion while recoiling, the numerical derivative will, however, be noisy. We therefore smooth the recoil curve, minimizing the Brownian noise, prior to performing the numerical derivative. To test the validity and effects of the adopted procedure, we compare our result with the simulated the recoil dynamics of a tethered microsphere.

## 3 Simulation

To simulate the recoil dynamics of the microsphere we solve, numerically, the Langevin equation for the system, using $`F_{dna}`$ as an external force . In the non-inertial limit, we have

$$\frac{dR}{dt}=\frac{F_{dna}(r)}{\gamma }+\frac{(t)}{\gamma },$$
(4)
where $`(t)`$ is the Brownian force, which has the properties $`(t)=0`$ and $`(t)(t^{})=\mathrm{\Gamma }\delta (tt^{})`$, where $`\mathrm{\Gamma }=2\gamma k_BT`$.

The discrete form of Eq. (4) is

$$R(t+\mathrm{\Delta }t)=R(t)\frac{F_{dna}(r)\mathrm{\Delta }t}{\gamma }+\frac{\sqrt{\mathrm{\Gamma }\mathrm{\Delta }t}}{\gamma }\xi _n,$$
(5)
where $`\xi _n`$ is a Gaussian random variable with $`\xi _n=0`$ and $`\xi _n\xi _m=\delta _{nm}`$.

Until this point, we have treated the motion of the microsphere using a reference frame directed along the force axis of the DNA molecule. This one-dimensional reference frame is not appropriate to our purposes, since it does not contain the experimental reference frame which we want to simulate, in which the microsphere undergoes three-dimensional Brownian motion. Therefore, we extend Eq. (5) to three-dimensional space, obtaining

$$X_i(t+\mathrm{\Delta }t)=X_i(t)\frac{F_{dna}(r)\mathrm{\Delta }t}{\gamma }\frac{X_i(t)}{R}+\frac{\sqrt{\mathrm{\Gamma }\mathrm{\Delta }t}}{\gamma }\xi _{n,i},$$
(6)
with $`i=1,2,3`$ and

$$r=\sqrt{\underset{i=1}{\overset{i=3}{}}X_i^2}a=Ra,$$
(7)
where the set $`\{X_i\}`$ represents the coordinates of the center of mass of the microsphere while $`r`$ represents the relative distance between the two extremities of the DNA molecule. Since the DNA is attached to the surface of the microsphere, the two coordinate systems are connected through Eq. (7), where $`a`$ is the radius of the microsphere.

### 3.1 Correcting $`\gamma `$: The wall-drag effect

Use of a constant and isotropic Stokes friction coefficient $`\gamma `$ is justified if and only if the fluid velocity attains constant value $`v_0`$ far from the microsphere, in all directions. When an obstacle interferes with the fluid flow near the microsphere, such as the coverglass surface in the present case, this assumption is no longer valid and $`\gamma `$ depends on the direction of motion and on the distance from the microsphere to the obstacle . According to the analyses of Faxen , Brenner  and Goldman *et al*  in the early and mid 20th century, for a microsphere of radius $`a`$ whose center of mass is at a distance $`h`$ from a planar surface, the friction coefficient $`\gamma `$ is anisotropic, with

$$\gamma _{||}\frac{6\pi \eta a}{\left(1\frac{9}{16}\left(\frac{a}{h}\right)+\frac{1}{8}\left(\frac{a}{h}\right)^3\frac{45}{254}\left(\frac{a}{h}\right)^4\frac{1}{16}\left(\frac{a}{h}\right)^5\right)}$$
(8a)
for motion parallel to the planar surface *($`xy`$ plane)*, and
$`\gamma _{}=`$ $`\mathrm{\hspace{0.33em}6}\pi \eta a\times `$ $`[{\displaystyle \frac{4}{3}}\mathrm{sinh}\alpha {\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{n(n+1)}{(2n1)(2n+3)}}\times `$
$`[{\displaystyle \frac{2\mathrm{sinh}\left[(2n+1)\alpha \right]+(2n+1)\mathrm{sinh}\mathrm{\hspace{0.33em}2}\alpha }{4\mathrm{sinh}^2\left[(n+\frac{1}{2})\alpha \right](2n+1)^2\mathrm{sinh}^2\alpha }}1]],`$
for motion perpendicular to the planar surface *(z axis)*, where $`\alpha =\mathrm{cosh}^1(h/a)`$. The exact correction for $`\gamma _{}`$, shown in Eq. (LABEL:gammav-e), can be replaced by a far simpler approximation, given by
$$\gamma _{}6\pi \eta a\left(1+\frac{a}{ha}\right),$$
(8c)
with an average error $`5\%`$ for $`h/a>1`$. The loss of isotropy in $`\gamma `$ alters Eq. (6), which becomes

$$X_i(t+\mathrm{\Delta }t)=X_i(t)\frac{F_{dna}(r)\mathrm{\Delta }t}{\gamma _{||}}\frac{X_i(t)}{R}+\frac{\sqrt{\mathrm{\Gamma }_{||}\mathrm{\Delta }t}}{\gamma _{||}}\xi _{n,i}i=1,2,$$
(9a)
$$X_3(t+\mathrm{\Delta }t)=X_3(t)\frac{F_{dna}(r)\mathrm{\Delta }t}{\gamma _{}}\frac{X_3(t)}{R}+\frac{\sqrt{\mathrm{\Gamma }_{}\mathrm{\Delta }t}}{\gamma _{}}\xi _n,$$
(9b)
where Eq. (9a) generates the *x* and *y* components of the recoil curve and Eq. (9b) the *z* component. This equation, however, does not take into account the fact that the height $`h`$ of the center of the microsphere cannot be less than its radius $`a`$, in which case the microsphere is touching the coverglass. At this point, interactions between the microsphere and the coverglass can no longer be neglected and Eq. (8c) can no longer be safely used. We have nevertheless observed that, after the recoil, the microsphere exhibits Brownian motion in the plane of the coverglass for quite some time, the same behavior being observed for most of the free microspheres as well. This observation indicates that the microspheres remain hovering close to the coverglass for a period of time before sticking to it. To determine how close to the coverglass the microspheres were hovering, we prepared a sample containing only free microspheres under the same experimental conditions. After letting the sample rest on the microscope for more than 2 hours, we measured the microsphere diffusion coefficient for motion parallel to the coverglass (Fig. (2)). The measured value was $`(0.063\pm 0.001)\mu m^2/s`$, approximately $`2.75`$ times less than the theoretical free difusion coefficient, $`\frac{k_bT}{6\pi \eta a}=0.173\mu m^2/s`$, for $`a=1.42\mu m`$. Assuming this discrepancy is caused by the increased friction close to the coverglass and applying Eq. (8a) we obtained an average hovering distance of approximately $`1.47\mu m`$, indicating that the microspheres in fact were not touching the coverglass, their point of nearest approach remaining approximately $`0.05\mu m`$ above it.

The simulation routine was then altered to reproduce this effect, not letting $`h`$ drop below its minimum value of $`1.47\mu m`$ during the recoil, as observed experimentally.

The initial conditions of the microsphere for all simulated recoils are $`X_1(0)=\sqrt{\left(0.85L+a\right)^2h_0^2}`$, $`X_2(0)=0`$ and $`X_3(0)=h_0`$, where we have assumed that the DNA molecule is at least $`85\%`$ stretched, initially. The coordinate system is taken so as to place the direction of stretching in the *$`xz`$ plane*. The time increment, $`\mathrm{\Delta }t`$, was set to $`(1/30)s`$, matching videomicroscopy time resolution for better comparison between simulated and experimental curves (which also proved adequate for numerical convergence). A simulated recoil curve is shown in Fig. (3).

## 4 Experiment

### 4.1 Experimental Setup

The experiments are made on a Nikon *TE300 Eclipse* inverted optical microscope coupled to two CCD cameras and a piezo-driven stage (Fig. 4). We use a SDL 5422 near infrared laser *($`\lambda =832nm`$)* focused through an infinity corrected objective lens *(Nikon plan apo, DIC H, 100$`\times `$, 1.4 NA, oil immersion)* to create the optical tweezers. Imaging from CCD 1 is used to locate and move suitable microspheres into position while imaging from CCD 2, filtered from the near infrared light, is videorecorded for posterior analysis.

### 4.2 Sample Preparation

Samples are prepared using the method described in , with slight variations. We prepare a solution containing PBS $`(150mMNa^+)`$ pH 6 and microspheres *(Polybead Polystyrene $`2.8\mu m`$ microspheres, $`\rho =1.05g/ml,n=1.6`$)*. The density of microspheres in the solution is around $`450/\mu l`$. For every $`395\mu l`$ of solution, we add $`5\mu l`$ of DNA solution *(New England Biolabs $`\lambda `$ DNA $`553\mu g/ml`$)*, previously heated in a thermal bath at 62$`{}_{}{}^{o}C`$ for 5 minutes, mix gently and incubate together for 20 minutes. The solution is then pipeted into cells *($``$ $`4mm`$ radius, $`3mm`$ height)* on coverglasses and covered with cut pieces of coverglass to prevent excessive evaporation. At this pH condition, DNA molecules bind to the microspheres and to the coverglass, preferably by its extremities . After 24h of incubation at room temperature, the cells are opened and gently washed with PBS $`(150mMNa^+)`$ pH 7.4 to remove the excess of free beads and reset the pH for the DNA. The cells are covered again and are then ready for experiments.

### 4.3 DNA-Microsphere System Selection

With the samples placed on the microscope, the optical tweezers are used to locate a microsphere bound to a single DNA molecule. Although many microspheres are attached to DNA, some judgement is needed to choose a good DNA-microsphere system in order to avoid microspheres attached to more than one molecule, molecules bound to the coverglass at multiple locations and other possible undesirable configurations. Verifying if the selected microsphere presents radial symmetry with respect to stretching and recoiling around its apparent resting point is the minimal required test for a good candidate. Once a suitable DNA-microsphere pair is selected, the initial height $`h_0`$ for the recoil is set adjusting the focal distance of the tweezers. The microsphere is trapped and moved to stretch the DNA. The videorecording is started and the tweezers cut off, starting the recoil.

### 4.4 Data Extraction

After videorecording the recoil through CCD 2, the images are analyzed using the *Image 1.62*  software package. The *x* and *y* coordinates of the center of mass of the microsphere are extracted from each frame and combined with the time index of the recording to give $`X_1(t)`$ and $`X_2(t)`$.

### 4.5 Estimate for $`h(t)`$

As already mentioned, the appropriate correction to the friction coefficient for the microsphere depends on the distance $`h`$ of its center from the coverglass. We use the optical tweezers to set the initial height $`h_0`$ of the microsphere. During the recoil, the value of $`h`$ can be measured from the defocusing of the microsphere, a very time consuming procedure. We believe that, for the recoil experiments, sufficient accuracy can be achieved by using a reasonable estimate for $`h(t)`$ guided by computer simulations. To obtain $`h(t)`$ for the experimental recoils we used averaged simulated recoils to establish how $`h(t)`$, which corresponds to $`X_3(t)`$ in the simulated recoils, was related to the horizontal component of the recoil, $`X_1(t)`$. We found that, because of the anisotropy in $`\gamma `$, $`h(t)`$ plotted as a function of $`X_1(t)`$ does not correspond to a straight line trajectory pointing to the origin of the coordinate system, as would be expected if $`\gamma _{||}`$ and $`\gamma _{}`$ were identical. Nevertheless, we observed that the average $`h(t)`$ can be estimated, with good agreement with the simulations (Fig. (5)), by

$$h(t)=0.925+(h_00.925)\frac{X_1(t)}{X_1(t=0)},$$
(10)
until $`h(t)`$ reaches $`1.47\mu m`$, where the experimental evidence indicates that the microspheres interrupt, at least for a time longer than the duration of the recoil, their vertical motion. During the most relevant region of the recoil, the distance from the bottom of the microsphere to the coverglass is larger than $`0.5\mu m`$, a value much larger than any roughness of such a good optical quality coverglass used.

## 5 Data Analysis

Using the previously chosen coordinate system, in which the molecule is stretched along the *$`xz`$ plane*, the recoil is primarily along the *$`x`$ axis*. We may then use the $`x`$ component of Eq. (1), given by

$$\frac{dX_1}{dt}=\frac{F_{dna}(r)}{\gamma _{//}}\frac{X_1}{R},$$
(11)
to fit our data. The $`y`$ component of the recoil will fluctuate around 0, allowing us to write

$$R\sqrt{X_1^2+h^2}.$$
(12)
Analysis of the recoil curves (Fig. (6)) consists of three steps. First, the recoil curve is *smoothed* using the *Stineman* smoothing function from *KaleidaGraph*  software package, filtering out high-frequency Brownian noise. The smoothed recoil curve is then used to evaluate the finite difference $`\frac{\mathrm{\Delta }X_1}{\mathrm{\Delta }t}`$. The resulting velocity data are fit according to Eq. (11).

To check the validity of the results produced by this procedure, we compared the values of $`A`$ and $`L`$ obtained from the analysis of several simulated recoils with the known values of $`A`$ and $`L`$ used as parameters in each simulation. The process was repeated several times, with different values of $`A`$, $`L`$, $`h_0`$ and $`T`$ used in the simulated recoils. In most cases, there was good agreement between the values of $`A`$ and $`L`$ used in the simulation and the values obtained from the analysis of the recoils generated (Fig. (7)). We have, however, observed that in certain cases the analysis of a simulated recoil produced results quite different than those obtained from the analysis of other recoils generated with the same parameters. We believe that such occasional discrepancies occur when a slight deviation from the average recoil trajectory is enhanced when the recoil is smoothed and its derivative evaluated. To minimize this effect in the experimental analysis, we have always recorded, analysed and compared at least 4 recoils for each DNA-microsphere complex.

We analysed 21 experimental recoils from 5 different DNA-microsphere systems in 3 different samples. The average persistence length obtained from our data was $`A=(43\pm 5)nm`$, a value well within the range obtained from other techniques .

## 6 Conclusion

In recent years measurements of the persistence length of single DNA molecules have been made through several different methods, ranging from spectroscopic techniques to electron microscopy. Each of these approaches involves a particular set of assumptions, as pointed out by Wang *et al.* in , and is susceptible to different sources of error.

We have shown, based on simulated and experimental data, that the *recoil method*, in which we follow the recoil dynamics of a previously stretched DNA molecule through the motion of a tethered microsphere, can be used to determine the persistence length of single DNA molecules, provided that adequate correction of the friction for the microsphere is made. Although it is possible to perform the analysis of the recoil using only on the horizontal component of the movement, knowledge of the vertical dynamics of the microsphere is fundamental, since its vertical position is the critical parameter governing the friction correction. Both simulated and experimental data were treated using the quasistatic response of the DNA molecule, without considering the far from equilibrium dynamics proposed in . We observe a high degree of consistency between experiment and simulation.

We acknowledge helpful discussions with G. V. Shivashankar and M. Feingold. This work was supported by the Brazilian Agencies: Fundação de Amparo à Pesquisa do Estado de Minas Gerais (FAPEMIG), Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) and FINEP-PRONEX.


format: text

text:
# Calculation of Chiral Lagrangian Coefficients

We present a systematic way to combine the global color model and the instanton liquid model to calculate the chiral Lagrangian coefficients. Our numerical results are in agreement well with the experimental values.

PACS: 24.85.+p, 12.38.Lg, 11.15.Pg

Quantum chromodynamics (QCD) is successful to describe strong interaction in terms of basic force at high energy. However, the large coupling constant at low energy makes the perturbation theory break down. At low-energy scale, the freedoms of quarks and gluons realize themselves in hadronic level as mesons and baryons, which is characterized by chiral perturbation theory ($`\chi `$PT).$`^{\text{[1, 2]}}`$ The $`\chi `$PT can give a lot of satisfactory explanations for experimental data while the coefficients of the chiral expansion are input parameters and their values are determined from experimental observables. Since the present technique cannot derivate the chiral Lagrangian from QCD directly, a class of chiral invariant QCD based models are used in deriving the chiral Lagrangian which is characterized by the effective quark-quark interaction.

The global colour model (GCM), which has provided a lot of successful descriptions of the long distance properties of strong interaction and the QCD vacuum as well as hadronic phenomena at low energy, is based on a truncation of QCD Lagragian preserves global colour gauge transformation and permits a large $`N_c`$ expansion.$`^{\text{[3, 4, 5, 6]}}`$ The instanton liquid model (ILM) of the QCD vacuum is based on a semiclassical approximation, in which all gauge configurations are replaced by an ensemble of topological nontrivial fields, instantons and anti-instantons.$`^{\text{[7, 9]}}`$ Although it does not give rise to a long range confining force between quarks, the instanton vacuum has been shown to provide a good phenomenological description of many hadronic properties. Therefore it is interesting to calculate the chiral coefficients with the two coupling nonperturbative models. In fact, the method of the two coupling nonperturbative models to calculate physical quantities already exists in the literature and is successful in describing some hadronic phenomena.$`^{\text{[17]}}`$ In this Letter, we approximate the effective quark-quark interaction by a dilute instanton liquid ensemble and calculate the chiral coefficients from the GCM following Ref. , while the pioneer work for obtaining chiral hadronic phenomena through GCM was given by Roberts et al.$`^{\text{[10]}}`$

Here we describe briefly the main skeleton of the GCM. The global colour model, based upon an effective quark-quark interaction which is approximated by the dilute instanton liquid configurations, is defined by a truncation of QCD Lagrangian while maintaining all global symmetries.

The generation function for QCD in Euclidean space can be written as

$$Z[\overline{\eta },\eta ]=𝒟\overline{q}𝒟q𝒟A\mathrm{exp}\left\{S+d^4x(\overline{\eta }q+\overline{q}\eta )\right\},$$
(1)
where

$$S=d^4x\left\{\overline{q}\left[\gamma _\mu \left(_\mu ig\frac{\lambda ^a}{2}A_\mu ^a\right)+\psi \right]q+\frac{1}{4}G_{\mu \nu }^aG_{\mu \nu }^a\right\},$$

and $`G_{\mu \nu }^a=_\mu A_\nu ^a_\nu A_\mu ^a+gf^{abc}A_\mu ^bA_\nu ^c`$. We leave the gauge fixing term, the ghost field term and its integration measure to be understood. Here $`A_\mu ^a`$ are approximated by the fields of the dilute instanton liquid. In order to preserve the chiral invariance of the above action $`S`$ (not the full action), the quarks have to be coupled to a external field $`\psi (x)`$ in stead of current quark mass $`m`$, the field $`\psi (x)`$ transforms in a certain way under chiral transformations$`^{\text{[1, 2, 4]}}`$ (For detailed discussion one can see Ref.).

Here we introduce

$$e^{W[J]}=𝒟Ae^{{\scriptscriptstyle d^4x({\scriptscriptstyle \frac{1}{4}}G_{\mu \nu }^aG_{\mu \nu }^a+J_\mu ^aA_\mu ^a)}},$$
(2)
with $`J_\mu ^a(x)=\overline{q}(x)\gamma _\mu \frac{\lambda ^a}{2}q(x)`$ and rewrite the generation function as

$$Z[\psi ,\overline{\eta },\eta ]=𝒟q𝒟\overline{q}e^{{\scriptscriptstyle d^4x[\overline{q}(\gamma +\psi )q\overline{\eta }q\overline{q}\eta ]}}e^{W[J]},$$
(3)
then expand $`W[J]`$ in the current $`J_\mu ^a`$:

$$W[J]=\frac{1}{2}d^4xd^4yJ_\mu ^a(x)D_{\mu \nu }^{ab}(x,y)J_\nu ^b(y)+\frac{1}{3!}J_\mu ^aJ_\nu ^bJ_\rho ^cD_{\mu \nu \rho }^{abc}+\mathrm{}.$$
(4)
The GCM is defined by a truncation of the functional $`W[J]`$ in which the higher order $`n(3)`$-point functions are neglected, and only the gluon 2-point function $`D_{\mu \nu }^{ab}(x,y)`$ is retained. This model maintains all global symmetries of QCD and permits a $`1/N_c`$ expansion, however local SU(3) gauge invariance is lost by the truncation. For simplicity, we use a Feynman-like gauge $`D_{\mu \nu }^{ab}(xy)=\delta _{\mu \nu }\delta ^{ab}D(xy)`$. Using the functional integral approach, we introduce bilocal field integration variables $`^\theta (x,y)`$ and integrate out the quarks and gluons fields to change dynamical freedoms of fundamental QCD to hadronic variables. This method has an obvious advantage that the effective hadron field interactions retain knowledge of their subhadronic origin. The resulting function $`Z_{\mathrm{GCM}}`$ takes the following form:

$$Z_{\mathrm{GCM}}=𝒟e^{S[]},$$
(5)
where the action is given by

$$S[\overline{\eta },\eta ,]=\mathrm{TrLn}[G^1]+d^4xd^4y\frac{^\theta (x,y)^\theta (y,x)}{2g^2D(xy)}+d^4xd^4y\overline{\eta }(x)G(x,y;[^\theta ])\eta (y),$$
(6)
and $`G^1`$ is defined as

$$G^1(x,y)=(\gamma +\psi (x))\delta (xy)+\mathrm{\Lambda }^\theta ^\theta (x,y).$$
(7)
Here the quantity $`\mathrm{\Lambda }^\theta `$ arises from the Fierz reordering of the current-current interaction term in Eq.(4) and represents the direct product of Dirac, flavor SU(3) and color matrices.

The vacuum configurations are defined by minimizing the bilocal action: $`\frac{\delta S[]}{\delta }|_{\overline{\eta },\eta =0}=0,`$ following this mean field approximation, we obtain

$$_0^\theta (xy)=g^2D(xy)\mathrm{tr}[\mathrm{\Lambda }^\theta G_0(xy)].$$
(8)
These configurations provide self-energy dressing of the quarks through the definition $`\mathrm{\Sigma }(p)\mathrm{\Lambda }^\theta _0^\theta (p)=i\gamma p[A(p^2)1]+B(p^2)`$. In fact, Eq. (8) is the rain-bow Schwinger-Dyson equation. In this Letter, the Schwinger-Dyson functions $`A(p^2)`$ and $`B(p^2)`$ are approximated by the results of ILM through equating the two point Green functions of the corresponding quark in the two nonperturbative methods.$`^{\text{[17]}}`$

With the GCM generation function Eq.(5), it is straightforward to calculate the vacuum expectation value of any quark operator of the form $`O_n(\overline{q}_{i_1}\mathrm{\Lambda }_{i_1j_1}^{(1)}q_{j_1})(\overline{q}_{i_2}\mathrm{\Lambda }_{i_2j_2}^{(2)}q_{j_2})\mathrm{}(\overline{q}_{i_n}\mathrm{\Lambda }_{i_nj_n}^{(n)}q_{j_n})`$ in the mean field vacuum. Here $`\mathrm{\Lambda }^{(i)}`$ represents the operator in Dirac, flavor and color space.

Taking the appropriate number of derivatives with respect to the external sources $`\eta _i`$ and $`\overline{\eta _j}`$ of Eq.(5) and putting $`\eta _i=\overline{\eta }_j=0`$ , we obtain

$$0|:O_n:|0=()^n\underset{p}{}[()^p\mathrm{\Lambda }_{i_1j_1}^{(1)}\mathrm{}\mathrm{\Lambda }_{i_nj_n}^{(n)}(G_0)_{j_1i_{p(1)}}\mathrm{}(G_0)_{j_ni_{p(n)}}],$$
(9)
where $`p`$ represents the permutation of the $`n`$ indices. In fact, the above equation can only be correct with the assumption of vacuum saturate or vacuum dominance.$`^{\text{[22]}}`$ Furthermore, in calculation, we should subtract the perturbative part of the two point Green function for the quark,$`^{\text{[18]}}`$ which is here neglected due to the small current quark mass.

The quark condensate $`0|:\overline{q}(0)q(0):|0`$ can be inferred from Eq.(9),

$$0|:\overline{q}(0)q(0):|0=\frac{12}{16\pi ^2}_0^\mu 𝑑s\frac{B(s)}{sA^2(s)+B^2(s)}.$$
(10)
In the instanton liquid model, $`A(q)=1`$. Here $`\mu `$ is taken to be 1 $`GeV^2`$, the scale dependence of the quark condensate is thoroughly discussed in the previous paper.$`^{\text{[21]}}`$ The quark condensate is calculated for later use in the expressions for the Chiral coefficients.

We expand the bilocal fields about the saddle point,

$`^\theta (x,y)=_0^\theta (xy)+\widehat{}^\theta (x,y),`$ (11)
and then examine the effective interactions of the fluctuations, $`\widehat{}`$. The bilocal fields contain information about internal excitions of the $`q\overline{q}`$ in addition to center of mass motion which is associated with the usual local field variables. Thus, they can be projected onto a complete set of internal excitations $`\mathrm{\Gamma }_n^\theta `$ with the remaining center of mass degree freedom represented by the coefficients $`\varphi _n^\theta (P)`$. In fact, the coefficients $`\varphi _n^\theta (P)`$ can be identified with the meson fields. The bilocal fluctuations can be written as

$$\widehat{}^\theta (P,q)=\underset{n}{}\varphi _n^\theta (P)\mathrm{\Gamma }_n^\theta (P,q).$$
(12)
The functions $`\mathrm{\Gamma }_n^\theta `$ are in general eigenfunctions of the free kinetic operator of the bilocal fields. At the mass shell point, $`P^2=M_n^2`$, they satisfy the homogeneous Bethe-Salpeter equation in the ladder approximation for the given quantum number $`\theta `$ and modes n.$`^{\text{[6]}}`$

To obtain the chiral coefficients, we take a digression to discuss Goldstone theorem. Dynamical chiral symmetry breaking is associated with the occurrence of a massless Goldstone mode which can be related to the dynamical generated scalar self-energy of the quark propagator in the limit of vanishing current quark mass. By analysing the corresponding zero momentum behaviour of the axial-vector Ward identity in the chiral limit, we can obtain the Goldberger-Treiman relation for the quark-pseudoscalar vertex and then reduce the bound state function $`\mathrm{\Gamma }_n^\theta (P,q)`$ from ladder Bethe-Salpeter equation for the pseudo-scalar Goldstone mode with $`P^2=0`$ to the quark self-energy function.$`^{\text{[4]}}`$

Our main goal here is to obtain the effective action for the Goldstone mode, the higher mass fluctuations present in the bilocal fields are neglected. In fact, the higher mass fluctuations should have some effects if they have been integrated out. Based on the above discussion, we can write the full bilocal field of quark fields in the following form$`^{\text{[4]}}`$:

$$\mathrm{\Lambda }^\theta ^\theta =\mathrm{\Sigma }(xy)+B(xy)[U_5(\frac{x+y}{2})1],$$
(13)
where $`U_5(x)=P_RU(x)+P_LU^{}(x)`$ with $`P_{R,L}`$ the standard right-left projection operators. For the $`SU(3)`$ flavor case under consideration here the chiral field $`U`$ is defined as $`Ue^{i\lambda ^a\varphi ^a/f_\pi }`$ .

Next, we follow Ref.  to calculate the chiral coefficients and we consider only the real contribution to the effective action. The pioneer work for obtaining the chiral hadronic phenomenon through global color model was given by Roberts et al.$`^{\text{[10]}}`$ The restriction of the fluctuations to the Goldstone modes with $`UU^{}=1`$ entails that the second term of the action in Eq.(6) is independent of the filed $`U`$ and can be neglected. Thus, the real contribution to the action can be written as

$$𝒮Re[S]=\frac{1}{2}TrLn(G^1[G^1]^{}),$$
(14)
where

$$G^1=\gamma _xA(xy)+\psi (\frac{x+y}{2})\delta (xy)+B(xy)U_5(\frac{x+y}{2}).$$
(15)
Here $`A(xy)=\delta (xy)`$ in the ILM.

By expanding the logarithm and dropping the irrelevant constant terms, we obtain the chiral Lagrangian to the fourth order (in Euclidean space):

$`𝒮=`$ $`{\displaystyle }d^4x\{{\displaystyle \frac{f_\pi ^2}{4}}tr[_\mu U_\mu U^{}]{\displaystyle \frac{f_\pi ^2}{4}}[U\chi ^{}+\chi U^{}]`$ (16)
$`L_1(tr[_\mu U_\mu U^{}])^2L_2tr[_\mu U_\nu U^{}]tr[_\mu U_\nu U^{}]`$
$`L_3tr[_\mu U_\mu U^{}_\nu U_\nu U^{}]+L_5tr[_\mu U_\mu U^{}(U\chi ^{}+\chi U^{})]`$
$`L_8tr[\chi U^{}\chi U^{}+U\chi ^{}U\chi ^{}]\},`$
where $`\chi (x)=2\overline{q}q\psi (x)/f_\pi ^2`$ and the remaining trace is over flavor. The coefficients $`L_1`$, $`L_2`$ $`(L_2=2L_1)`$, $`L_3`$, $`L_5`$, $`L_8`$ and $`\pi `$ decay constant $`f_\pi `$ are functions of the scalar part of quark self energy and its derives, as the expressions are lengthy, one can see Ref. for details.

In performing the gradient expansion it is important to keep the $`x`$-dependence of the external field $`\psi (x)`$. After carrying out the gradient expansion to the fourth order, we employ the equation of motion which is obtained at the second order and depends on the external field $`\psi (x)`$. Finally we can identify $`\psi (x)`$ with the current quark mass matrix. Failure to keep the $`x`$-dependence of $`\psi `$ to the very end violates chiral invariance and will render some unphysical results for the low-energy coefficients.

From the above expressions, we can see that the coefficients $`L_4`$, $`L_6`$, $`L_7`$ vanish, which is consistent with the large $`N_c`$ QCD result. Some modifications may be obtained if we have integrated out higher mass states.

Instantons, as the solutions of the classical Yang-Mills equation of motion, give many successful descriptions of strong interaction phenomena, such as $`U(1)_A`$ problem, dynamical chiral symmetry breaking, tunneling the $`\theta `$ vacuum and so on. After the great achievement of ’t Hooft’s computing the functional determinant in the one-instanton background, a tempt was made to understand quantitatively the nonperturbative QCD vacuum configurations in terms of instanton ensembles.$`^{\text{[11, 12]}}`$ However, the notorious infrared problem due to large size instanton made many researchers no longer believe in the dilute instanton gas model. To circumvent the disaster, one may think the overcome that the vacuum gets more filled with instantons of increasing size for larger distances, at some scale there might be some repulsive interactions to stabilize the medium while the semiclassical treatment is still possible and the instantons are not much deformed through their interactions, thus form a dilute instanton liquid ensemble. In Ref., the authors used variational techniques and mean field approximation to deal with the statistical mechanics of the instanton liquid while stabilizing it by phenomenological hard core and repulsive interactions derived from ansatz of gauge field. Subsequently, the propagation of light quarks in the dilute instanton liquid ensemble was formulated quantitatively and refined due to the treatment of the nonzero mode propagation.$`^{\text{[15, 16]}}`$ Although it does not give rise to a long range confining force between quarks, the instanton vacuum has been shown to provide a good phenomenological description of many hadronic properties.$`^{\text{[7, 9]}}`$ Furthermore, lattice calculations support that the vacuum gluon configurations can be approximated by the dilute instanton liquid.$`^{\text{[20]}}`$ Its shortcoming will be that the vector self-energy function of quark propagator in the instanton ensembles $`A(q^2)=1`$ which is not consistent with the corresponding one in lattice QCD and other potential models.$`^{\text{[19, 8]}}`$ In addition, the $`q^21GeV^2`$ behaviour of the scalar self energy function $`B(q^2)`$ of quark propagator in the instanton medium falls as $`1/p^6`$ which disagrees with QCD. Thus, it is interesting to substitute the quark propagator which is dressed through propagating in the dilute instanton medium into the expressions for the $`\pi `$ decay constant and chiral coefficients. Our numerical results fit experimental values well. The quark propagator in the ILM can be written as

$$G(p)=[i\gamma p+M(p)]^1+O(N/VN_c),$$
(17)
where

$`M(p)=\sqrt{{\displaystyle \frac{N}{2VN_c}}}{\displaystyle \frac{p^2\phi ^2(p)}{\sqrt{\frac{d^4k}{(2\pi )^4}k^2\phi ^4(k)}}},`$
$`\phi (k)=\pi \rho ^2{\displaystyle \frac{}{z}}(I_0(z)K_0(z)I_1(z)K_1(z))|_{z=k\rho /2},`$ (18)
where $`I`$ and $`K`$ are the modified Bessel functions, $`N/V(200MeV)^4`$ and $`\rho (600MeV)^1`$.

In conlusion, we have provided a systematic way to combine the GCM with the ILM.$`^{\text{[17]}}`$ Then based on the nonperturbative quark propagator derived in the ILM, we calculate chiral coefficients which are expressed in terms of the scalar part of quark self energy and its derives.

Our results are compatible with the experimental data. Our results are $`f_\pi =86MeV`$, $`\overline{q}q=(217MeV)^3`$, $`L_1=0.8\times 10^3`$, $`L_2=1.6\times 10^3`$, $`L_3=3.9\times 10^3`$, $`L_5=2.0\times 10^3`$, $`L_8=1.0\times 10^3`$. The experimental values are $`L_1=0.9\pm 0.3\times 10^3`$, $`L_2=1.7\pm 0.7\times 10^3`$, $`L_3=4.4\pm 2.5\times 10^3`$, $`L_5=2.2\pm 0.5\times 10^3`$, $`L_8=1.1\pm 0.3\times 10^3`$

Following the GCM, we truncate the QCD Lagragian to include only two point gluon Green function which is approximated by the dilute instanton liquid ensemble and describe the mesons as fluctuations about the saddle point of the effective action. The chiral coefficients obtained in this way are consistent with the large $`N_c`$ QCD, and embody their underlying theory through the dependence on the quark self-energy and its derives. Some modification will be made if we have integrated out the high mass fluctuations instead of neglecting them.

To obtain the nonperturbative quark propagator, we have to resolve the rain-bow Schwinger-Dyson equation. However, the form of the gluon propagator $`g^2G(p)`$ in the infrared region is unknown, one often uses model forms as input in the previous studies of the rainbow Schwinger-Dyson equation.$`^{\text{[3, 6, 5]}}`$ Obviously, this dressing comprises the notation of constituent quarks by providing a mass $`m(p^2)=B(p^2)/A(p^2)`$ with $`(m(0)=418MeV)`$ in the ILM, which is corresponding to dynamical symmetry breaking. Here we use the dilute instanton liquid to approximate the QCD vacuum configurations, the quark propagator is dressed through propagating in the dilute instanton liquid ensemble and contains valuable nonperturbative information which manifest itself in the coefficients. Although ILM has shortcoming as mentioned above, it does give satisfactory chiral coefficients and can be extended to other systems. Furthermore, the lattice calculations support that the dilute instanton liquid ensemble is a good approximation of vacuum configurations.


text:
# Spatial distribution of galactic halos and their merger histories

## 1 Introduction

The advent of large wide-field redshift surveys of galaxies, such as the Two-Degree Field (2dF, Colless et al., 2001) and the Sloan Digital Sky Survey (SDSS, York et al., 2000), will allow detailed studies of clustering of galaxies as a function of their environment and internal properties. Indeed, hierarchical growth of structure via gravitational instability is thought to play a dominant role in shaping both the large-scale galaxy clustering and internal properties of galaxies such as luminosity, colors, and morphology. This close connection implies that studies of the spatial distribution of galaxies as a function of their internal properties and environment should provide us valuable insights into the process of galaxy formation. Previous observational studies and the first results from the 2dF and SDSS have shown that clustering strength depends on morphology (e.g., Hermit et al., 1996; Guzzo et al., 1997), luminosity (see, e.g., Hamilton 1988; Benoist et al. 1996; Norberg 2001), and colors (e.g., Zehavi, 2002) of galaxies. Greatly enhanced clustering of super-luminous IR-galaxies (Bouchet et al., 1993) and morphological (Dressler, 1980; Postman & Geller, 1984; Whitmore et al., 1993; Biviano et al., 2002) and color (e.g., Butcher & Oemler, 1978; Margoniner et al., 2001) segregation in clusters of galaxies indicate dependence of clustering on the merging history and large-scale environment.

In this paper we use a novel statistical tool, the mark correlation functions (hereafter MCFs), to study clustering of galactic halos as a function of their properties and environment in a high-resolution numerical simulation of the $`\mathrm{\Lambda }`$CDM cosmology. Mark correlation functions (Stoyan, 1984; Stoyan & Stoyan, 1994) have been introduced into astrophysics only recently (Beisbart & Kerscher, 2000), although some aspects of marked point processes were discussed by Peebles (1980). The mark statistics can be used to quantify the differences in the spatial distributions of various galaxy samples (similarly to the usual two-point correlation function) and, at the same time, to study the interplay between the spatial clustering and the distribution of galaxy properties (marks). In this respect, the MCFs are the natural extension of the spatial correlation functions for studies where it is advantageous to consider mark and spatial distributions simultaneously. Variants of the MCFs can be used to study both continuous (e.g., luminosity or angular momentum) and discrete distributions (e.g., color classes or morphological types) of galactic properties. This makes them valuable for quantitative studies of luminosity and morphological segregation of galaxies as well as dependence of spatial distribution of galaxies on various events in their evolutionary history (e.g., time since the last major merger), which can be used as marks. Indeed, the mark correlation statistics quickly proved to be a very useful tool for identification of physical processes that shape the observed spatial distribution of galaxies (Szapudi et al. 2000, see Beisbart et al. 2002 for a recent review).

The first step towards the use of clustering as a probe of processes shaping properties of galaxies is a good theoretical understanding of how these processes affect spatial distribution of galaxies. During the last several years, thanks to continuously improving spatial and mass resolutions of numerical simulations and development of semi-analytic models of galaxy formation, there was a significant progress in the theoretical understanding of galaxy clustering evolution and bias (e.g., Bagla, 1998; Jing, 1998; Kauffmann et al., 1999; Katz et al., 1999; Colín et al., 1999; Kravtsov & Klypin, 1999; Pearce et al., 1999; Schmalzing et al., 1999). The complete information about the internal properties and evolution of galactic halos, usually available in theoretical analysis, allows one to study in detail the interplay between different evolutionary processes and spatial distribution of objects. In the present paper we use mark correlation functions to study clustering of galaxy-size dark matter halos and its dependence on the halo properties (e.g., circular velocity) and evolution history in a high-resolution simulation of the currently favored flat Cold Dark Matter (CDM) model with cosmological constant (see § 2).

The paper is organized as follows. We discuss briefly the numerical simulation in Sect. 2. In Sect. 3 we introduce and explain the properties of mark correlation functions. In Sect. 4 we present analysis of the spatial distribution of dark matter halos using mark correlation functions. In Sect. 5 we discuss results and their implications and summarize our conclusions.

## 2 Numerical simulations

We used the Adaptive Refinement Tree (ART) code (Kravtsov et al., 1997) to simulate the evolution of collisionless DM in the currently favored $`\mathrm{\Lambda }`$CDM model ($`\mathrm{\Omega }_\mathrm{m}=1\mathrm{\Omega }_\mathrm{\Lambda }=0.3`$; $`H_0=h\times 100`$km s-1 Mpc$`{}_{}{}^{1}=70`$km s-1 Mpc-1; $`\sigma _8=1.0`$). The age of the Universe in this cosmology is 13.5 Gyrs and normalization is in accordance with the four year COBE DMR observations (Bunn & White, 1997) as well as the observed abundance of galaxy clusters (e.g., Pierpaoli et al., 2001; Ikebe et al., 2002).

The numerical simulation of the $`\mathrm{\Lambda }`$CDM model followed the evolution of $`256^31.67\times 10^7`$ particles in a periodic $`60h^1\mathrm{Mpc}`$ box. The particle mass is thus $`1.1\times 10^9h^1\mathrm{M}_{}`$. The ART code reaches high force resolution by refining all high-density regions with an automated refinement algorithm. The refinements are recursive: the refined regions can also be refined, each subsequent refinement having half of the previous level’s cell size. This creates an hierarchy of refinement meshes of different resolution covering regions of interest. The criterion for refinement is local overdensity of particles: in the simulation presented in this paper the code refined an individual cell only if the density of particles (smoothed with the cloud-in-cell scheme) was higher than $`n_{th}=5`$ particles. Therefore, all regions with overdensity higher than $`\delta =n_{th}2^{3L}/\overline{n}`$, where $`\overline{n}`$ is the average number density of particles in the cube, were refined to the refinement level $`L`$. For the simulation presented here, $`\overline{n}`$ is $`1/8`$. The peak formal dynamic range reached by the code on the highest refinement level is $`32,768`$, which corresponds to the smallest grid cell of $`1.83h^1\mathrm{kpc}`$; the actual force resolution is approximately a factor of two lower. The simulation that we analyze here has been used in Colín et al. (1999), and we refer the reader to this paper for further numerical details.

Identification of DM halos in the very high-density environments (e.g., inside groups and clusters) is a challenging problem. The goal of this study is to investigate spatial correlations of halo populations as closely related to the observed galaxy population as possible. This requires identification of both isolated halos and satellite halos orbiting within the virial radii of larger systems. The problems associated with halo identification within high-density regions are discussed in Klypin et al. (1999). In this study we use a halo finding algorithm called Bound Density Maxima (BDM). The main idea of the BDM algorithm is to find positions of local maxima in the density field smoothed at a certain scale and to apply physically motivated criteria to test whether the identified site corresponds to a gravitationally bound halo111The detailed description of the algorithm is given in (Klypin et al., 1999) and (Colín et al., 1999).. It is based on the ideas of the DENMAX halo finder (Bertschinger & Gelb, 1991), in the sense that the BDM makes sure that the density peaks are gravitationally bound and estimates parameters of the halos after removing unbound particles. The algorithm identifies both isolated halos and subhalos located in the virial regions of more massive halos. The distribution of halos identified in this way can be compared to the distribution of galaxies directly, because the halo and galaxy catalogs include both isolated systems and objects within clusters and groups.

Even with algorithms tailored for identification of sub-halos, additional problems exist. Interacting halos exchange and loose mass; the total mass of a halo depends on its radius, which is difficult to define in a dense environment within virialized regions. We alleviate the latter problem by using the maximum circular velocity instead of the mass. In practice, maximum circular velocity is a rather stable quantity which changes little even when halos looses most of its mass and can serve therefore as a useful mass-related “tag” of a halo. Numerically, the maximum “circular velocity” ($`\sqrt{GM/R}`$), $`v_{\mathrm{circ}}`$ can be measured more accurately then mass. In addition, the maximum circular velocity can be more readily compared to observations than, for example, virial mass or mass within the tidal radius of the halo.

## 3 Mark correlation functions

In studying galaxy clustering with the mark correlation functions, we view galaxies as discrete points in space with marks describing their intrinsic physical properties. Thus, we consider a point set $`\{𝐱_i\}_{i=1}^N`$ and attach a mark $`m_i`$ to each point $`𝐱_i^3`$ resulting in the marked point set $`\{(𝐱_i,m_i)\}_{i=1}^N`$ (Stoyan, 1984; Stoyan & Stoyan, 1994). The marks, in turn, can be either continuous or discrete. In the following, we use the circular velocity as a continuous mark and merging/stripping events of halos as discrete marks. In a subsequent paper, we will apply MCFs to study various other marks, such as the spin parameter (continuous scalar mark) and the angular momentum (vector mark).

Let $`\varrho `$ be the mean number density of the points in space and $`\varrho ^M(m)\mathrm{d}m`$ the probability that the value of a mark on a point lies within the interval $`[m,m+\mathrm{d}m]`$. The mean mark is then $`\overline{m}=dm\varrho ^M(m)m`$, the mark variance is $`\sigma _M^2=dm\varrho ^M(m)(m\overline{m})^2`$. We assume that the joint probability $`\varrho ^{SM}(𝐱,m)`$ of finding a point at position $`𝐱`$ with mark $`M`$, splits into a space–independent mark probability and the mean density: $`\varrho ^M(m)\times \varrho `$. The spatial–mark product–density

$$\varrho _2^{SM}((𝐱_1,m_1),(𝐱_2,m_2))\mathrm{d}V_1\mathrm{d}m_1\mathrm{d}V_2\mathrm{d}m_2,$$
(1)
is the joint probability of finding a point at $`𝐱_1`$ with the mark $`m_1`$ and another point at $`𝐱_2`$ with the mark $`m_2`$. We obtain the spatial product density $`\varrho _2(𝐱_1,𝐱_2)`$ and the two–point correlation function $`\xi (r)`$ by marginalizing over the marks:

$$\begin{array}{c}\varrho ^2(1+\xi (r))=\varrho _2(𝐱_1,𝐱_2)=\hfill \\ \hfill =dm_1dm_2\varrho _2^{SM}((𝐱_1,m_1),(𝐱_2,m_2)),\end{array}$$
(2)
where $`\xi (r)`$ is the spatial two-point correlation function which depends only on the separation $`r=|𝐱_1𝐱_2|`$ of the points for a homogeneous and isotropic point set.

We define the conditional mark density:

$$\begin{array}{c}_2(m_1,m_2|𝐱_1,𝐱_2)=\hfill \\ \hfill =\{\begin{array}{cc}\frac{\varrho _2^{SM}((𝐱_1,m_1),(𝐱_2,m_2))}{\varrho _2(𝐱_1,𝐱_2)}\hfill & \text{ for }\varrho _2(𝐱_1,𝐱_2)0,\hfill \\ 0\hfill & \text{ else }.\hfill \end{array}\end{array}$$
(3)
For a stationary and isotropic point distribution $`_2(m_1,m_2|r)\mathrm{d}m_1\mathrm{d}m_2`$ is the probability of finding the marks $`m_1`$ and $`m_2`$ at two points located at $`𝐱_1`$ and $`𝐱_2`$, under the condition that they are separated by $`r=|𝐱_1𝐱_2|`$. Now the full mark product–density can be written as

$$\varrho _2^{SM}((𝐱_1,m_1),(𝐱_2,m_2))=_2(m_1,m_2|𝐱_1,𝐱_2)\varrho _2(𝐱_1,𝐱_2).$$
(4)
If there is no mark segregation in space, $`_2(m_1,m_2|r)`$ is independent of $`r`$, and $`_2(m_1,m_2|r)=\varrho ^M(m_1)\varrho ^M(m_2)`$.

Starting from these definitions, especially using the conditional mark density $`_2(m_1,m_2|r)`$, one may construct several mark correlation functions sensitive to different aspects of mark segregation (Beisbart & Kerscher, 2000). The basic idea is to consider weighted conditional correlation functions describing the probability of finding two points at a separation $`r`$. For a positively defined weighting function $`f(m_1,m_2)0`$ we define the average over pairs with separation $`r`$:

$$f_\mathrm{P}(r)=dm_1dm_2f(m_1,m_2)_2(m_1,m_2|r).$$
(5)
$`f_\mathrm{P}(r)`$ is the expectation value of the weighting function $`f`$ (depending only on the marks), under the condition that we find a point pair with separation $`r`$. For a suitable defined integration measure Eq. (5) is also applicable to discrete marks. The definition (5) is very flexible, and allows us to investigate the correlations of both continuous and discrete marks.

In the following analysis, we calculated the mark correlation functions taking into account the periodicity of the simulation box. However, we obtain virtually identical results using the estimator without boundary corrections (see Appendix A of Beisbart & Kerscher (2000) for details).

### 3.1 Correlations of scalar marks

For scalar marks the following mark correlation functions have proven to be useful (Stoyan & Stoyan 1994, Beisbart & Kerscher 2000, Schlather 2001): The simplest weight to be used is the mean mark:

$$k_m(r)\frac{m_1+m_2_\mathrm{P}(r)}{2\overline{m}}.$$
(6)
It quantifies the deviation of the mean mark on pairs with separation $`r`$ from the overall mean mark $`\overline{m}`$. For example, $`k_m(r)>1`$ indicates mark segregation for point pairs with a separation $`r`$, specifically their mean mark is larger than the overall mark average.
Higher moments of marks like the mark fluctuations

$$\text{var}(r)\left(m_1m_1_\mathrm{P}(r)\right)^2_\mathrm{P}(r),$$
(7)
or the mark covariance (Cressie, 1991)

$`\text{cov}(r)`$ $`\left(m_1m_1_\mathrm{P}(r)\right)\left(m_2m_2_\mathrm{P}(r)\right)_\mathrm{P}(r)`$
$`=m_1m_2_\mathrm{P}(r)m_1_\mathrm{P}(r)m_2_\mathrm{P}(r),`$ (8)
may be used to quantify mark segregation. For example, a positive $`\text{cov}(r)`$ indicates that points with separation $`r`$ tend to have similar marks, whereas a negative $`\text{cov}(r)`$ indicates different marks.

### 3.2 Correlations of discrete marks

For discrete labels only combinations of indicator–functions are possible, and the integral degenerates into a sum over the labels. Supposing the marks of our objects belong to classes labeled with $`i,j,\mathrm{}`$, the conditional cross–correlation functions are given by

$$C_{ij}(r)\delta _{m_1i}\delta _{m_2j}+(1\delta _{ij})\delta _{m_2i}\delta _{m_1j}_\mathrm{P}(r),$$
(9)
with the Kronecker $`\delta _{m_1i}=1`$ for $`m_1=i`$ and zero otherwise (Stoyan & Stoyan, 1994). By construction

$$\underset{i}{}\underset{j}{}C_{i,j}(r)=1.$$
(10)
Mark segregation is indicated by $`C_{ij}2\varrho _i\varrho _j/\varrho ^2`$ for $`ij`$ and $`C_{ii}\varrho _i^2/\varrho ^2`$, where $`\varrho _i`$ denotes the number density of points with label $`i`$. The $`C_{ij}`$ are cross–correlation functions under the condition that two points are separated by a distance of $`r`$.

## 4 Results

For any study one needs to have a complete halo sample that is not affected by the numerical details of the halo finding procedure. We have tested the completeness of the halo samples using different parameters for the halo finder. For the given force and mass resolution the halo samples do not depend on the numerical parameters of the halo finder for halos with $`v_{\mathrm{circ}}\genfrac{}{}{0pt}{}{_>}{^{}}100`$ km/s (Gottlöber et al., 1999). In Fig. 1 we show the cumulative number of halos with a circular velocity larger than a value $`v_{\mathrm{circ}}`$, for redshifts $`z=0`$, $`z=1`$, and $`z=3`$, respectively. Assuming a minimum circular velocity of 100 km/s the samples are complete at $`z1`$ but we are missing a small fraction of halos with $`v_{\mathrm{circ}}<130`$ km/s at $`z=3`$.

### 4.1 Merging of halos

According to hierarchical structure formation halos formed early and grow during evolution due to accretion of matter and merging with other halos. In particular merger events are important because they are expected to lead to dramatic changes in the structure of dark matter halos and the galaxies they harbor. In–falling objects may damage or even destroy a stellar disk. The inflow of material may also serve as a source of fresh gas and therefore increase the star formation rate. At the same time, collisions between halos may result in shock heating of the gas, which would tend to delay or prevent star formation for some period of time.

Following Gottlöber et al. (2001) we identify major mergers as events when the mass of a halo grows by more than 25% during a time interval of about 0.5 Gyrs (approximate interval between simulation outputs). In the above paper we showed that for redshifts $`z<2`$ the merger rate can be fitted by a simple power law $`(1+z)^{3.0}`$. This merger rate evolution is in very good general agreement with observations (e.g., Le Fevre et al. (2000) measured a merger rate varying with redshift as $`(1+z)^{3.2\pm 0.6}`$). In addition, we found that evolution of the merger rate depends on the environment of the halo: halos that end up in clusters and groups by $`z=0`$ have a steeper evolution of merger rate and a higher rate of major mergers at early epochs compared to isolated “field” halos. This is because clusters and groups form in the regions that are overdense on large scales in which halos form and undergo the phase of active merging earlier than the overall field halo population.

For the $`z=0`$ sample of halos with $`v_{\mathrm{circ}}>100`$ km/s, about 32 % of halos had one major merger in the past and additional 19 % of halos had two or more major mergers. Now let us consider the distribution of epochs of the last major merger (relevant, for instance, for estimating a fraction of halos that could host old disks such as that of the Galaxy). We found that 55 % of the halos located in clusters at $`z=0`$ underwent a major merger after redshift $`z=4`$, but that corresponding fraction for the isolated halos is 43 %. In contrast, the fraction of isolated halos which underwent a major merger at a redshift $`z<1`$ is somewhat higher (19 %) than the corresponding fraction of halos in clusters (14 %); for $`z<0.4`$ (i.e., within the last $`5`$ Gyrs) these fractions are 8 % and 3.5 %, respectively. This reflects the fact that due to the high internal velocity dispersion of halos in clusters mergers are almost impossible. Since the merger rate of group halos is high compared to the overall merger rate of halos at all analyzed redshifts, the fraction of present-day group halos which merged after a fixed redshift is always higher than the fraction of isolated or cluster halos which merged after the same redshift.

Finally, let us consider one more effect. Due to the tidal interactions halos in dense regions (i.e., virialized regions of groups and clusters) tend to loose mass via tidal stripping. In order to take this effect into account, we follow the mass evolution of all halos and identify halos that lost more than 30% of their maximum (over their evolution) mass from $`z=1`$ to the present epoch. One would expect that galaxies hosted by such halos also lost their supply of fresh gas so that no star formation was possible in the recent past.

### 4.2 Spatial distribution and evolution history

In the preceding section we considered the overall fractions of halos in different environments and with certain merger history classes. This relatively straightforward analysis reveals existence of some environmental dependency of halo evolution histories. The goal of this section is to carry out a more quantitative analysis of how spatial clustering of different halo subsamples depends on evolution histories of their halos. As discussed in the previous subsection, major mergers (and tidal stripping) can be expected to result in dramatic changes in the properties of galaxies (i.e., morphology and color). One can expect, therefore, that the spatial distribution of halos that experienced a recent merger or stripping event is different from the distribution of the overall halo population. For example, Knebe & Müller (1999) found that massive halos undergoing mergers at present exhibit a much stronger bias with respect to the dark matter than relaxed halos do.

Figure 2 shows the two–point correlation function of all halos (a) with $`v_{\mathrm{circ}}>100`$km/s (solid line) compared with that of the subsamples of halos with different evolution histories. We divided the sample of all halos into four subsamples: halos which never (n) underwent a major merger in the past, halos which underwent a major merger before (e) and after (l) redshift $`z=1`$, and halos that lost mass since $`z=1`$ (s). Note, that the stripped halos constitute a separate sample, however every stripped halo belongs also per definition to one of the other subsample. In particular, a substantial part of “stripped” halos in clusters underwent a major merger before redshift $`z=1`$, i.e. they belong to the sample (e) of halos. They lost most of their mass later on due to interactions.

Figure 2 shows clearly that the subsample of halos which underwent the last major merger before redshift $`z=1`$ is more clustered than the sample of all halos. This is not surprising since most of such halos formed early in the regions of large-scale overdensity and ended up in groups and clusters by the present epoch. The stripped halos are even more biased with respect to the overall halo population. This is also due to the fact that halos that loose mass via stripping are located in the high-density regions where tidal forces are strong.

Let us now consider the spatial distribution of different halos using the mark correlation functions introduced in Sect. 3. To this end, we first split the total halo sample (i.e., sample of all halos with circular velocities $`v_{\mathrm{circ}}>100`$km/s) into two subsamples consisting of halos which experienced a major merger (sample m) and halos which never experienced a major merger (sample n), respectively. Figure 3 shows the conditional cross–correlation functions (eq. 9) of these samples at $`z=0`$. Positive $`C_{m,m}`$ on scales below 3$`h^1`$Mpc indicates that halos that experienced a major merger in their formation history are relatively overabundant in close pairs of halos, while halos without a major merger are underabundant. No significant cross–correlation between m and n exists. Qualitatively, this feature is independent from a lower cut in the circular velocity $`v_{\mathrm{circ}}`$, but the amplitude of the effect is reduced if we consider only the more massive halos with $`v_{\mathrm{circ}}>150`$km/s.

To investigate the significance of these deviations we performed a non–parametric Monte Carlo test (Besag & Diggle, 1977), similar to the one used in Kerscher et al. (2001) and Kerscher (1998). Our null hypothesis is “no mark correlation”. We simulate this null hypotheses by keeping the positions of the halos fixed and randomly re-assigning the marks (in this case the class assignments) to the halos. We repeat this to generate 99 realizations of this null hypothesis. The shaded areas in Fig. 3 are the one–$`\sigma `$ regions numerically determined from these samples. To quantify the significance we have to define a distance measure. Using $`M=4`$ equidistant radii $`r_i`$ in the range from 0.8$`h^1`$Mpc to 2.8$`h^1`$Mpc we define the “distance” of the $`k`$–th sample to the expected value for no mark correlation:

$$d^k=\frac{1}{M}\underset{i=1}{\overset{M}{}}\left(C_{m,m}(r_i)\frac{\varrho _m^2}{\varrho ^2}\right)^2,$$
(11)
with the overall number density $`\varrho `$ and the number density $`\varrho _m`$ of merged halos. Similarly we determine the distance $`d^{\mathrm{halo}}`$ of the original halo sample to the null hypothesis. Then we sort all the $`d`$’s and determine the position of $`d^{\mathrm{halo}}`$. In this case $`d^{\mathrm{halo}}`$ is the fifth largest distance, and we conclude that the original halo sample is incompatible with the null hypothesis ”no mark correlation” at a significance level of 95% (see the comments by Marriott 1978 concerning the significance level).

Figure 4 shows the conditional cross–correlations of the three subsamples which we considered above using the two–point correlation function (Fig. 2: the halos which never underwent a major merger in the past, sample n; early major merger at $`z>1`$, sample e; and late major merger at $`0z1`$, sample l).

As before, positive $`C_{\mathrm{e},\mathrm{e}}`$ at small separations indicates that for objects at distances less than about 2 $`h^1`$Mpc the halos with a major merger in their early formation history are relative over–abundant, at the expense of halos without a major merger, as deduced from the lowered $`C_{\mathrm{n},\mathrm{n}}`$. This signal is most prominent on scales below 3$`h^1`$Mpc but it extends out to 10$`h^1`$Mpc in agreement with the enhanced two-point correlation function of that subsample (Fig. 2). We interpret this as indication of a high number of early merged halos in clusters. The signal has high significance and it is not influenced by uncertainties in the normalization of the correlation function that may be caused by selection effects. Interestingly, the halos with a late major merger show no excess correlations but rather a lowered abundance on small scales, also manifested as the lower correlation function amplitude of that subsample.

The lower panel of Fig. 4 shows the cross-correlation of halos from different evolution classes. The over-abundance of pairs of never and early merged halos reflects the continuous accretion process onto high density regions. Infalling isolated halos from less dense regions accrete onto higher density regions with high velocity dispersions and, thus, low probability of merging. Therefore, type-n halos can survive in the high-density regions relatively long which explains the increasing of $`C_{\mathrm{n},\mathrm{e}}`$ towards small scales. The opposite is true for never $`n`$ and late $`l`$ merged halos. Type-l halos located predominantly in groups where mergers are more likely due to the lower velocity dispersions. The probability of accreting type-n halo (located close to an type-l halo) to experience a merger in high-density regions is therefore high. Many such halos will thus disappear (will become l-halos) resulting in suppression of $`C_{\mathrm{n},\mathrm{l}}`$ amplitude at small separations.

Note, that these features are qualitatively independent from a lower cut in the circular velocity $`v_{\mathrm{circ}}`$, but the amplitude and the spatial range is reduced if we consider more massive halos with $`v_{\mathrm{circ}}>120`$km/s or $`v_{\mathrm{circ}}>150`$km/s. This is due to the higher number of mergers within the low circular velocity halos.

Let us now consider halos which lost a substantial part of their mass due to tidal interactions. Figure 5 shows the conditional cross–correlation function of the halos using two classes: no stripping (sample: ns) and stripped (sample: s). Positive $`C_{\mathrm{s},\mathrm{s}}`$ at $`r<5h^1\mathrm{Mpc}`$ indicates that the number of stripped pairs is strongly enhanced at these separations, whereas the number of non–stripped pairs is reduced. This result is in accordance with the strongly enhanced correlation function shown in Fig. 2. In fact, we expect to find stripped halos only in the environment of clusters. The results for samples with a higher cut in the circular velocity $`v_{\mathrm{circ}}>150`$km/s are very similar.

### 4.3 Spatial distribution and circular velocity

As mentioned in Sec. 2, the mass and the maximum circular velocity of halos are tightly related. At the same time, the circular velocity can be determined more reliably in simulation as well as in observations, either through direct measurement using emission line width or rotation curve or via galaxy luminosity using the Tully-Fisher and the Faber-Jackson relations. Therefore, it is interesting to explore galaxy mark correlations with the maximum circular velocity as mark. This would mimic to some degree luminosity segregation effects in observed galaxy samples. Figure 6 shows the mark correlation functions of halos at $`z=0`$ with the circular velocity as mark. There is a strong signal in $`k_m(r)`$ at small separations but the signal is significant even out to 10$`h^1`$Mpc . This indicates that the mean circular velocity of pairs of halos with separations below $`10`$$`h^1`$Mpc is larger than the overall mean circular velocity $`\overline{v}_{\mathrm{circ}}`$ of the parent halo sample.

The negative signal of $`\text{cov}(r)`$ is confined to small separations ($`2`$$`h^1`$Mpc). This is because pairs at small separations are more frequently built from one halo with circular velocity larger than $`\overline{v}_{\mathrm{circ}}`$ and the other halo with circular velocity smaller than $`\overline{v}_{\mathrm{circ}}`$. Hence, this signal is dominated by pairs of small-mass subhalos and massive parent halos. The mark correlations results for both samples with $`v_{\mathrm{circ}}>100`$km/s and $`v_{\mathrm{circ}}>150`$km/s are shown. The mean mark correlation function, $`k_m(r)`$, exhibits a slightly stronger signal for the sample with a lower cut in circular velocity, but the signal $`\text{cov}(r)`$ is weaker for $`v_{\mathrm{circ}}>100`$km/s. The latter is due to the considerably larger number of isolated small-mass halos; i.e. in addition to pairs between large- and small-mass halos, for the $`v_{\mathrm{circ}}>100`$km/s sample there are many more small–small mass pairs.

Figure 7 shows the evolution of the conditional covariance $`\text{cov}(r)`$ of the circular velocity with redshift. The conditional covariance is negative at low redshifts out to scales of 2$`h^1`$Mpc as discussed above. At high redshifts ($`z=2`$ and $`z=3`$), significant positive amplitude of $`\text{cov}(r)`$ indicates that pairs with similar circular velocities are overabundant. At these redshifts the signal is significant out to the scale of $`5`$$`h^1`$Mpc due to the large number of smaller-mass progenitors of the present day halos.

## 5 Discussion and conclusions

In the previous section we used a novel statistical tool, the mark correlation functions, to study clustering of galactic halos as a function of their properties and environment in a high-resolution numerical simulation of the $`\mathrm{\Lambda }`$CDM cosmology. We applied MCFs using several types of continuous and discrete marks: maximum circular velocity, $`v_{\mathrm{circ}}`$, of halos (continuous), merger mark indicating whether halos experienced (m) or not (n) a major merger in their evolution history, a stripping mark (s) indicating whether the halo underwent a tidal stripping (i.e., mass loss) during its evolution (discrete marks). The halos which underwent major merger (m) are further classified by the epoch of the merger: late (l; $`z<1`$) and early (e; $`z>1`$) mergers. Our main results are as follows.

The two-point correlation amplitude is different for the halo subsamples with different marks. The halos that experienced an early major merger or mass loss (e and s) are clustered considerably more strongly than the overall halo population, while halos with late or no mergers have correlation function amplitude below that of the overall halo sample. This result indicates that halo clustering depends sensitively on the details of their evolution history. If existence of a major merger during halo evolution is related to the morphology of galaxies that halos host, the above result indicates that early type galaxies and galaxies in clusters and groups (hosted by halos that undergo tidal stripping) should be clustered more strongly than the late type galaxies and the overall galaxy population. Qualitatively, such trend exists in the observed galaxy samples (e.g., Hermit et al., 1996; Guzzo et al., 1997; Zehavi, 2002) implying that the morphology-dependent clustering may be largely due to the overall merger history of the galactic halos.

Using maximum circular velocity of halos as a continuous mark, we found that at $`z=0`$ the mean circular velocity of pairs of halos with separations $`10`$$`h^1`$Mpc is larger than the overall mean circular velocity $`\overline{v}_{\mathrm{circ}}`$ of the parent halo sample (manifested as significant enhancement of the mean mark at these separations; see Fig 6). Moreover, the negative mark covariance (Eq. 3.1) at small separations shows an enhanced abundance of pairs with halo circular velocities above and below the average circular velocity. This mean circular velocity enhancement increases steadily during the evolution of halos from $`z=3`$ to $`z=0`$. The mark covariance, $`cov(r)`$, has a more complicated behavior: it is negative at present, disappears at redshift $`z1`$ and becomes positive at higher redshifts due to the larger number of low circular velocity halos (the circular velocity function of halos steepens at high redshifts, see Fig. 1). Although the relation is not direct, the maximum circular velocity of halos should correlate well with the maximum circular velocity or velocity dispersion of the galaxies they host. The enhanced mean circular velocity in small-separation pairs should therefore correspond to the luminosity segregation or luminosity-dependent clustering in the observed galaxies. The luminosity dependence of galaxy clustering was recently convincingly detected in both 2dF (Norberg, 2001) and SDSS (Zehavi, 2002) galaxy surveys.

The mark correlation analysis indicate that galaxy-size halos ($`v_{\mathrm{circ}}>100`$km/s) which experienced a major merger in their evolution history are over-abundant in pairs with separations $`3`$ $`h^1`$Mpc with respect to the overall halo population, while halos which never experienced a major merger are under-abundant at these separations. We find no significant cross-correlation between these two halos classes. The overabundance of merger halos is due largely to the halos which experienced a major merger relatively early ($`z>1`$); halos with late ($`z<1`$) major merger are not over-abundant (this is also manifested in the low amplitude of their two-point correlation function relative to that of the overall halo population; see Figs. 2, 4). This result can be interpreted as correlation between the time since the last major merger and present-day environment of the halo (i.e., halos which underwent an early major merger tend to be located in clusters and groups). The significance of the results was estimated by a non–parametric Monte Carlo test which showed that the segregation and anti-segregation have significance of $`>95`$% in the distance range 0.8$`h^1`$Mpc to 2.8$`h^1`$Mpc . Similarly, the probability of finding stripped halos in pairs of separations $`5`$$`h^1`$Mpc is twice higher than the corresponding probability for the overall halo sample. Halos which experienced early major mergers and/or mass loss due to tidal stripping are likely to host early type galaxies. In this case, the above mark correlation results for DM halos indicate that morphological segregation of galaxies may be due to the specifics of the mass evolution histories and environment of their parent halos.

The analysis presented in this paper showed that MCFs provide powerful, yet algorithmically simple, quantitative measures of segregation in halo spatial distribution with respect to their properties (e.g., maximum circular velocity) and merger history (e.g., time since the last major merger). The mark correlation functions allow us to quantify the degree of segregation as a function of scale and can be used to quantify the differences in the spatial distributions of various galaxy samples (similarly to the usual two-point correlation function) and, at the same time, to study the interplay between the spatial clustering and the distribution of galaxy properties (marks). In this respect, the MCFs are a natural extension of the spatial correlation functions for studies where it is advantageous to consider mark and spatial distributions simultaneously. We believe that this will make the mark correlation functions very useful for analysis of spatial clustering and segregation as a function of various galaxy properties in current (SDSS and 2dF) and future (e.g., DEEP2) large redshift surveys.

###### Acknowledgements.

S.G. acknowledges support from Deutsche Akademie der Naturforscher Leopoldina with means of the Bundesministerium für Bildung und Forschung grant LPD 1996. M.K. would like to thank the people at the AIP for their hospitality on several occasions. M.K. was supported by the Sonderforschungsbereich 375-95 für Astro–Teilchenphysik der Deutschen Forschungsgemeinschaft. A.V.K. was partially supported by NASA through Hubble Fellowship grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555.


text:
# 1 Spacetime in Quantum Physics

## 1 Spacetime in Quantum Physics

The point of departure from classical to quantum physics is the algebra $`(T^{}Q)`$ of functions on the classical phase space $`T^{}Q`$. According to Dirac, quantisation can be achieved by replacing a function $`f`$ in this algebra by an operator $`\widehat{f}`$ and equating $`i\mathrm{}`$ times the Poisson bracket between functions to the commutator between the corresponding operators.

In classical physics, the functions $`f`$ commute, so $`(T^{}Q)`$ is a commutative algebra.

But the corresponding quantum algebra $`\widehat{}`$ is not commutative. Dynamics is on $`\widehat{}`$. So quantum physics is noncommutative dynamics.

A particular aspect of this dynamics is fuzzy phase space where we cannot localise points, and which has an attendent effective ultraviolet cutoff: The number of states in a phase space volume $`V`$ is infinite in classical physics and $`V/\mathrm{}^{2d}`$ in quantum physics when the phase space is of dimension $`2d`$. The emergence of this cutoff from quantisation is of particular importance for the program of fuzzy physics .

This brings us to the focus of our talk. In quantum physics, the commutative algebra of functions on phase space is deformed to a noncommutative algebra, leading to a “noncommutative phase space”. Such deformations, characteristic of quantum theory, are now appearing in different approaches to fundamental physics. The talk will focus on a few such selected approaches and their implications.

Before proceeding further, let us mention a few of these lines of thought leading to noncommutative geometry: 1) Noncommutative geometry has made its appearance as a method for regularising quantum field theory (qft) and in studies of deformation quantisation. This talk will more or less base itself on these aspects. 2) It has turned up in string physics as quantised D-branes. 3) Certain approaches to canonical gravity  have also used noncommutative geometry with great effectiveness.

## 2 Fuzzy physics and quantum field theory

In what follows, we will focus on fuzzy physics both as a means to regularise quantum field theories and in their relation to quantum spacetimes.

But as mentioned already, we will not talk about how they emerge from string physics.

## 3 Fuzzy spacetime as regulator

The original ideas for using fuzzy spaces to regulate qft’s are due to Madore . They concern quantising underlying spacetime itself, making it into a fuzzy spacetime.

We know since Planck and Bose that quantisation introduces a short distance cut-off, changing the number of states in a phase space volume $`V`$ from $`\mathrm{}`$ to $`V`$$`/`$$`\mathrm{}^{2d}`$.

Now qft’s on a manifold require regularisation. The usual nonperturbative regularisation involves lattice qft’s. The use of fuzzy spacetimes can be another. The latter has important advantages like maintaining symmetries and avoiding fermion doubling. The particular approach reported here involves many colleagues. Our representative papers are  and the work of Vaidya, Dolan et al. and Lopez et al. in . Related or overlapping work is due to  and .

We do Euclidean qft’s. Quantising $`S^4`$ would be of the greatest physical interest. But for now, we will focus on the simpler case of $`S^2`$.

## 4 Fuzzy $`^2`$, $`S^3`$ and $`S^2`$

### 4.1 Relation between $`^2`$, $`S^3`$ and $`S^2`$

Consider $`^2`$ with coordinates $`z=(z_1,z_2)`$. We have

$$S^3=<z:z{}_{}{}^{2}:=z_i^2=1>^2$$
(1)
and

$`S^2`$ $`=`$ $`<\stackrel{}{x}=z^{}\stackrel{}{\tau }z,zS^3>,`$ (2)
$`\stackrel{}{\tau }`$ $`=`$ $`Paulimatrices,`$ (3)
$`\stackrel{}{x}.\stackrel{}{x}`$ $`=`$ $`1.`$ (4)
$`S^2`$ is the quotient of $`S^3`$ by the $`U(1)`$ action $`zze^{i\theta }`$ (Hopf fibration).

The group $`SU(2)=\{g\}`$ acts on $`^2`$, $`S^3`$ and $`S^2`$:

$`zgz,x_iR_{ij}(g)x_j,R(g)=Rotationmatrixforg.`$ (5)
### 4.2 The fuzzy $`^2`$

Now we quantise $`^2`$ by the replacements

$`z_\alpha a_\alpha ,\overline{z_\alpha }a_\alpha ^{},`$ (6)
where $`a_\alpha `$ and $`a_\alpha ^{}`$ are harmonic oscillator annihilation and creation operators with the usual commutation relations. $`SU(2)`$ still acts with generators

$$L_\alpha =a^{}\sigma _\alpha a$$
(7)
which commute with the number operator:

$`N=a^{}a,[N,L_\alpha ]=0.`$ (8)
### 4.3 The fuzzy three-sphere

Consider

$`S_\alpha `$ $`=`$ $`a_\alpha {\displaystyle \frac{1}{\sqrt{N+1}}},`$
$`S{}_{\alpha }{}^{}^{}`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{N+1}}}a_\alpha {}_{}{}^{},`$
$`S{}_{\alpha }{}^{}{}_{}{}^{}S_\alpha `$ $`=`$ $`{\displaystyle \frac{N}{N+1}}1asN\mathrm{},`$
$`[S_\alpha ,S_\beta ^{}]`$ $``$ $`0asN\mathrm{}.`$ (9)
So $`S_\alpha `$ are normalised commuting vectors as $`N\mathrm{}`$ and the fuzzy three-sphere $`S{}_{F}{}^{}^3`$ is the algebra generated by $`S_\alpha `$ and $`S{}_{\beta }{}^{}^{}`$.

It is important to note that $`\frac{1}{(N+1)}`$ here plays the role of a quantised Planck’s constant.

This raises the following important questions:Is it possible that Planck’s “constant” is in reality an operator ? How can one experimentally test this possibilty?

The representation of the $`S{}_{F}{}^{}^3`$-algebra on the Fock space is irreducible.

### 4.4 The fuzzy two-sphere

Since $`[N,L_\alpha ]=0`$, we have that

$$[N,S{}_{}{}^{}\sigma _{\alpha }^{}S]=0.$$
(10)
So we can restrict the algebra of the fuzzy sphere $`S_{F}^{}{}_{}{}^{2}`$ generated by $`S{}_{}{}^{}\sigma _{\alpha }^{}S`$ to the finite- dimensional vector space spanned by the eigenvectors of $`N`$ with eigenvalue $`n`$. They are spanned by

$`|n_1,n_2>`$ $`:=`$ $`{\displaystyle \frac{1}{\sqrt{n_1!n_2!}}}(a{}_{1}{}^{}{}_{}{}^{})^{n_1}(a{}_{2}{}^{}{}_{}{}^{})^{n_2}|0>,`$
$`n_1+n_2`$ $`=`$ $`n,`$ (11)
and carry angular momentum $`J{}_{0}{}^{}=\frac{n}{2}`$. As this representation is irreducible, it follows that the fuzzy sphere algebra for angular momentum $`J_0`$ is the algebra of $`(2J_0+1)`$ x $`(2J_0+1)`$ matrices $`Mat(2J_0+1)`$. We denote its elements by $`M`$. It is just the vector space spanned by the tensor operators in the angular momentum $`n/2`$ representation.

Rotation acts on $`M`$ by $`MgMg^1`$. So the fuzzy sphere has angular momenta $`0,1,\mathrm{},2J_0`$ , cut off at $`2J_0`$. The algebra of functions on the two-sphere instead has all integral angular momenta upto $`\mathrm{}`$.

## 5 Scalar fields on $`S{}_{F}{}^{}^2`$

Scalar fields are power series in “coordinates” $`S{}_{}{}^{}\sigma _{\alpha }^{}S`$. So

$$ScalarField=(2J_0+1)dimensionalmatrix.$$
(12)
A scalar action, such as

$$𝒜=\frac{1}{2J_0+1}Tr([L_i,\varphi ]{}_{}{}^{}[L_i,\varphi ]+\varphi ^k)$$
(13)
can be quantised by functional integral methods. Renormalisation studies of such actions have been carried out in .

Gauge theories can be formulated on fuzzy spaces. For brevity, we will not enter into their discussion.

Summarising, we have the classical descent chain

$$S^3S^2.$$
(14)
It becomes after quantisation

$$_FS{}_{F}{}^{}{}_{}{}^{3}S{}_{F}{}^{}{}_{}{}^{2}.$$
(15)
The algebra dimension is $`\mathrm{}`$ for all except $`S{}_{F}{}^{}^2`$ for which it is $`(2J_0+1)^2`$.

## 6 On coherent states and star products

In quantum field theory, we calculate correlation functions like

$$<\varphi (\stackrel{}{n}{}_{1}{}^{})\varphi (\stackrel{}{n}{}_{2}{}^{})\mathrm{}\varphi (\stackrel{}{n}{}_{j}{}^{})>,\stackrel{}{n}{}_{j}{}^{}S^2.$$
(16)
To compare with such expressions, we have to know how to map our operators (finite-dimensional matrices) to functions. This is where coherent states and star products prove important. For us, they will be particularly important also when discussing issues of topology fluctuations. We now explain how star products can be defined on fuzzy spaces.

We start from the infinte dimensional Fock space associated with two oscillators and introduce the standard coherent states

$$z;\mathrm{}>=e^{z_\alpha a_\alpha {}_{}{}^{}\overline{z}_\alpha a_\alpha }|0>$$
(17)
where $`\mathrm{}`$ has been inserted in the state vector to indicate that it is associated with the Fock space. (It is omitted from the vacuum stae which will be commom to both the Fock space and its subspaces which will appear below.)

A theorem asserts that the diagonal matrix elements

$$<z;\mathrm{}\widehat{A}z;\mathrm{}>=A(z)$$
(18)
determines the operator $`\widehat{A}`$ (under suitable conditions on $`\widehat{A}`$): the map $`\widehat{A}A`$ of operators to functions is one-to-one. Both the Moyal and Perelomov star products follow from this result as we now indicate.

The star product $`AB`$ of two functions $`A`$ and $`B`$ is defined by

$$AB(z):=<z;\mathrm{}\widehat{A}\widehat{B}z;\mathrm{}>.$$
(19)
(It is not the Moyal star product, but equivalent to it.)

If $`\widehat{A}`$ is $`0`$ outside the subspace where $`N=n`$,then

$`\widehat{A}`$ $`=`$ $`P\widehat{A}P,`$
$`P`$ $`=`$ $`\text{Projector\hspace{0.17em} onto\hspace{0.17em} this\hspace{0.17em} subspace}.`$ (20)
The explicit expression for $`P`$ is

$$P=\underset{n_1+n_2=n}{}n_1,n_2><n_1,n_2$$
(21)
where we have used the definition (11).

The diagonal coherent state expectation value of $`P\widehat{A}P`$ is (upto a constant in the definition of the state)

$`A(z)`$ $`=`$ $`<z|\widehat{A}|z>,`$
$`|z>`$ $`:=`$ $`{\displaystyle \frac{1}{\sqrt{n!}}}{\displaystyle \underset{n_1+n_2=n}{}}(z_\alpha a_\alpha {}_{}{}^{})^n|0>.`$ (22)
The group $`SU(2)`$ with generators $`L_i`$ acts on these states according to

$`|z>`$ $``$ $`|gz>,`$
$`g`$ $``$ $`SU(2)`$ (23)
preserving $`|z_i|^2`$, so we can set it equal to $`1`$. This normalises these states.

The fuzzy sphere operators $`L_i`$ and $`\frac{1}{N+1}L_i`$ can be restricted to $`N=n`$. We can set them equal to $`0`$ on the subspace $`Nn`$, and treat them as the above $`\widehat{A}`$’s.

In this way, we have a star product on $`S{}_{F}{}^{}^2`$. There are explicit expressions for this star product.

## 7 Fluctuating topologies: A novel space-time uncertainty principle

We now pass to the 4-sphere $`S^4`$ and consider the fuzzy $`S^4`$, which is denoted by $`S{}_{F}{}^{}^4`$.( See in this connection.) We realise $`S^4`$ as a sphere in $`^5`$:

$$S^4=<x:=(x_1,x_2,..,x_5):|x_i|^2=1>.$$
(24)
For quantisation, it is necessary to introduce $`^4`$ with coordinates $`z=(z_1,z_2,..,z_4)`$. The unit sphere in $`^4`$ is $`S^7`$:

$$S^7=<z:|z_i|^2=1>.$$
(25)
We also introduce the five gamma matrices $`\gamma _\lambda `$:

$$\gamma _\lambda =standard\gamma matrices.$$
(26)
Then we can set

$$x_\lambda =z^{}\gamma _\lambda z,zS^7.$$

Now instead of $`2`$, we have $`4`$ sets of creation-aniihilation operators $`a_\alpha ,a_\alpha ^{}`$ ($`\alpha =1,2,3,4`$). $`S^4`$ will be thought of as a sphere in $`^5`$ and its quantised fuzzy version will hence be formulated.

Coherent states appropriate for $`S^4`$ are generalisations of those for $`S^2`$:

$$|z>:=\frac{1}{\sqrt{n!}}\underset{n_1+n_2+n_3+n_4=n}{}(z_\alpha a_\alpha {}_{}{}^{})^n|0>$$
(27)
where the sum on $`\alpha `$ is over $`4`$ values.

The quantised version of $`x_\lambda `$ is $`S{}_{}{}^{}\sigma _{\alpha }^{}S`$ with $`S`$’s being defined as before.

A simple calculation shows that

$$<z|S{}_{}{}^{}\sigma _{\alpha }^{}S|>=(\frac{n}{n+1})z^{}\gamma _\lambda z.$$
(28)
Or it gives $`S^4`$ as imbedded in $`^5`$.

But we now show that the fuzzy $`S{}_{F}{}^{}^4`$ emerges only approximately, having fluctuations of the order of $`\frac{1}{n}`$ into the six-dimensional fuzzy manifold $`P{}_{}{}^{3}_F`$, the fuzzy version of $`P^3`$.

To demonstrate this, consider the correlation function

$$<z|S{}_{}{}^{}\sigma _{\alpha }^{}SS{}_{}{}^{}\sigma _{\beta }^{}|z>.$$
(29)
A short calculation shows that it is

$$z^{}\gamma _\alpha zz^{}\gamma _\beta z+O(\frac{1}{n})z^{}\gamma _{\alpha \beta }z$$
(30)
where $`\gamma _{\alpha \beta }`$ is $`[\mathrm{\Gamma }_\alpha ,\mathrm{\Gamma }_\beta ]`$. Now $`z^{}\gamma _{\alpha \beta }z`$ is not a function on $`S^4`$. It is a function only on $`P^3`$. (It is invariant under the phase change of $`z:ze^{i\theta }z`$.)

By definition,

$$<z|S{}_{}{}^{}\sigma _{\alpha }^{}SS{}_{}{}^{}\sigma _{\beta }^{}|z>=(\frac{n}{n+1})z^{}\gamma _\alpha z(\frac{n}{n+1})z^{}\gamma _\beta z.$$
(31)
.

A measure of the lack of localisation of the star product on $`S^4`$ is possible to construct, but we will not discuss that here.

The behaviour of $`z^{}\gamma _\alpha z`$ under star product is generic for most manifolds, notable exceptions being $`P^N`$.

Thus we have an algebra of observables, such as that generated by $`z^{}\gamma _\alpha z`$ under \*, which is only approximately localised on a manifold $`M`$. It has fluctuations $`O(\frac{n}{n+1})`$ into another manifold $`M^{}M`$ which disappear in the classical limit $`\frac{1}{n}0`$.

## 8 More on quantum topologies

The fuzzy space $`S{}_{F}{}^{}^2`$ has representations in all dimensions.

In contrast the fuzzy space $`S{}_{F}{}^{}^4`$, or rather $`P{}_{F}{}^{}^3`$, has representations in symmetric products of its four-dimensional representation, namely in dimensions $`4,10,20,\mathrm{}`$.

So $`20\times 20`$ matrices can approximate either $`S^2`$ or $`P^3`$.

So what we see in these matrices depends on the operators we examine . If we work with $`S^2`$-coherent states and operators appropriate for them, then they will approximate $`S^2`$. But if instead we work with states and operators appropriate for $`P{}_{F}{}^{}^3`$, then these matrices will approximate that manifold.

In high dimensions, such approximate manifolds proliferate.

We can answer questions like: the probabilty of finding a $`P{}_{F}{}^{}^3`$\- localised state like

$`\rho `$ $`=`$ $`|z><z|,`$
$`z`$ $``$ $`P{}_{}{}^{3}.`$ (32)
($`|z>`$ being the coherent state for $`P{}_{F}{}^{}^3`$) in a $`P{}_{F}{}^{}^1`$-localised state

$`\stackrel{~}{\rho }`$ $`=`$ $`|\stackrel{~}{z}><\stackrel{~}{z}|,`$
$`\stackrel{~}{z}`$ $``$ $`P{}_{}{}^{1},thatisS^2,`$ (33)
($`|\stackrel{~}{z}>`$ being the coherent state for $`P{}_{F}{}^{}^1`$. The answer is

$$|Tr\rho ^{}\stackrel{~}{\rho }|^2.$$
(34)
We can see from this discussion that questions with fantasy about quantum spacetimes become accessible in fuzzy physics.

## 9 Causality and $`CPT`$ violation

Fuzzy models can be formulated for spacetimes with Minkowsky metric, not just for Euclidean spacetimes like $`S^2`$. They are natural models to quantise time, preserving many symmetries.

A popular model for quantising Minkowsky spacetime which has been carefully studied by Doplicher, Fredenhagen and Roberts  and others  is governed by the commutation relation

$`[x^\mu ,x^\nu ]`$ $`=`$ $`i\theta ^{\mu ,\nu },`$
$`\theta ^{\mu ,\nu }`$ $`=`$ $`\theta ^{\nu ,\mu }=arealconstant.`$ (35)
With spacetime fuzzy, the meaning of “spacelike separation” loses exact meaning and leads to causality and hence generically to $`CPT`$ violation.

A way to understand this is as follows. In $`1+1`$-dimensions, we can set

$$\theta ^{\mu ,\nu }=\theta ϵ^{\mu ,\nu }$$
(36)
where $`ϵ^{\mu ,\nu }`$ is the Levi-Civita symbol with $`ϵ^{0,1}=1`$. Then  the fuzzy version of a free field of mass $`M`$ has the representation

$$\mathrm{\Phi }=\psi (k)e^{ik_0x_0}e^{ik_1x_1}$$
(37)
where

$`k_0^2k_1^2`$ $`=`$ $`M^2,`$
$`\psi (k)`$ $`=`$ $`\delta (k^2M^2)a(k),k_0>0,`$ (38)
$`=`$ $`\delta (k^2M^2)a(k_0,k_1)^{},k_0<0,`$
$`a(k)`$ and $`a(k)^{}`$ being the usual annihilation and creation operators.

We can regard $`(x_0+ix_1)/\sqrt{2\theta }`$ and its adjoint as annihilation and creation operators. If $`|z>`$ denotes the corresponding coherent state, then the square of the fuzzy field localised approximately at $`(Rez`$,$`Imz)`$ is

$$<z|(\mathrm{\Phi })^2|z>:=(\overline{\mathrm{\Phi }})^2(z).$$
(39)
One checks that when $`z`$ and $`z^{}`$ correspond to spacelike points and $`\alpha `$ and $`\beta `$ are generic states,

$$<\alpha |[\overline{\mathrm{\Phi }})^2(z),(\overline{\mathrm{\Phi }})^2(z^{})]|\beta >e^{|zz^{}|^2}.$$
(40)
Such causality violations are generic in such models. A fundamental question, inspired by these models then is: How can we test them? . We can in principle do so using forward dispersion relations, but they are expected to be small and it is not clear if these dispersion relations can be evaluated with sufficient accuracy to establish significant limitations on possible causality violation.

## 10 Winding up

The following broad observations are suggested by the preceding discussions:

\- Models of spacetime based on noncommutative geometry are suggested by quantum physics itself, string physics and attempts to regularise quantum field theories.

\- They lead to strikingly novel spacetime models.

\- But lacking contact with experiments, these models for now remain metaphysical models, just as quantum gravity and string models.

## Acknowledgments

The new ideas reported here come from the joint efforts of a group of us, some of our representative papers being  and the work of Vaidya, Dolan et al. and Lopez et al. in . This work has been supported by DOE under contract number DE-FG02-85ER40231 and by the NSF-CONAcYT grant INT-9908763 which facilitated several mutual visits by the US and Mexican participants.